0

I have a file that defines a raster, although all the metadata was lost. I know the lat/long coordinates of the grid centroids, and that the grid is in a Mercator projection with a 1/12deg equatorial resolution. Is there any way to define the extent of the raster by just knowing the centroids and resolution? Currently I'm doing:

xmin = min(longitude)-1/24
xmax = max(longitude)+1/24
ymin = min(latitude)-1/24
ymax = max(latitude)+1/24
  • What software and version are you using? – juturna Jan 9 '15 at 16:20
  • I'm using various packages within R. – user13317 Jan 9 '15 at 17:11
  • What do you mean by the centroids? The center point of the corner cells? All cells? – mkennedy Jan 9 '15 at 17:22
  • @mkennedy Perhaps I should say the center point, and not centroid (???). I have the center points for all cells. – user13317 Jan 9 '15 at 19:15
2

The extent of the raster is equal to the extent of the cell centres, expanded by half the resolution.

Here's an example:

  1. Create a dummy raster with extent c(0, 1, 0, 1) and resolution c(0.1, 0.1):

    library(raster)
    r <- raster(res=0.1, xmn=0, xmx=1, ymn=0, ymx=1)
    
  2. Extract cell centres:

    p <- rasterToPoints(r)
    head(p)
    
    #         x    y
    # [1,] 0.05 0.95
    # [2,] 0.15 0.95
    # [3,] 0.25 0.95
    # [4,] 0.35 0.95
    # [5,] 0.45 0.95
    # [6,] 0.55 0.95
    
  3. Calculate the extent of the cell centres:

    e <- extent(p)
    
  4. Expand the extent by half the resolution. Note we don't need to specify 0.5 * res(r). The + method for extent objects yields a new extent that is res(r) units longer, overall, on each axis.

    e + res(r)
    
    # class       : Extent 
    # xmin        : 1.387779e-17 
    # xmax        : 1 
    # ymin        : -5.551115e-17 
    # ymax        : 1 
    

(There's minor floating point error in the recovered extent.)

However, I don't see any problem with your approach, either:

min(p[, 'x']) - xres(r)/2
# [1] 0

max(p[, 'x']) + xres(r)/2
# [1] 1

min(p[, 'y']) - xres(r)/2
# [1] -6.938894e-17

max(p[, 'y']) + yres(r)/2
# [1] 1
  • Nice answer! And thanks for confirming my initial thoughts! I'm not a GIS analyst, and I just wanted to make sure I wasn't making an unreasonable assumption. – user13317 Jan 12 '15 at 12:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.