Tuesday, May 29, 2012

Google Earth and ocean depth contours

Been playing around in R for a while and then a bit in the GIS environment - within R of course. Do not know very much about GIS, but know what I want: Looking at GIS fisheries data across various scales (macro to global). Now we have myriads of websites that provide a view, on a 'data' on a scale. But herein lies the limitation, it is piecemeal information. Now what GIS software medium actually operates on all scales? Who else to call than Google Earth?
Here is an snapshot of the depth contours for from 0 to 1500 m by 100 meter intervals for the whole of the North Atlantic. Original data is from General Bathymetric Chart of the Oceans. All that was needed was some lines of code in R to convert point values to contour and then generate the kml-file. And adding on top of that the ICES statistical rectangles as a layer, coloured in red.

Zooming in or out in Google Earth environment is a breeze. An example on the macroscale, showing the Porupine bank is shown below.

Wandering if one should not use the Google Earth as the base platform to publish all fisheries and management information. Both the global as well as the local.


  1. Howdy! I can clearly notice the fact that you deeply understand what you are speaking about over here. Do you own a degree or an education which is linked with the theme of the entry? Thanks a lot in advance for your answer.

  2. Fascinating post! Any chance you could send me the kml file of the depth contours? Via dropbox perhaps? It would be fun to try to reproduce this.

    Mark Huisjes,

  3. Very interesting, I also want to reproduce this and incorporate similar content to google earth.

    Can I maybe receive the r scripts

    Kind regards Maarten de Jong

  4. i think the R-code was something like this:

    require(rgdal) # check
    options(max.contour.segments = 100000)


    xs <- seq(from= -45,to= 0,by=sec30)
    xs <- xs[-length(xs)]+sec15 # get the mid point
    ys <- seq(from=30,to=60,by=sec30)
    ys <- ys[-length(ys)]+sec15 # get the mid point
    xss = rep(xs,times=length(ys)) # expand
    yss = rep(ys,each=length(xs))
    spp = SpatialPoints(data.frame(x=xss,y=yss))
    gt = points2grid(spp) # Create grid topology

    z = get.var.ncdf(d,"z")
    sgdf = SpatialGridDataFrame(gt, data=as.data.frame(z),CRS("+proj=longlat +datum=WGS84"))

    dLevels <- c(-50,seq(-100,-1400,by = -100)) # could not get further down
    aDepth = contourLines(as.image.SpatialGridDataFrame(sgdf),
    aDepth = ContourLines2SLDF(aDepth)
    proj4string(aDepth)=CRS("+proj=longlat +datum=WGS84")

    writeOGR(aDepth, "aDepth.kml",c(as.character(dLevels)),driver="KML")