[R-sig-Geo] Problems with image2Grid

Michael Sumner mdsumner at gmail.com
Thu Sep 24 12:01:30 CEST 2009


Hello, Roger's reply to this made me realize I replied only to
fernando, so I wanted to correct my first message and add some more.

This reminds me that image2Grid needs check for irregular x or y
values - it assumes the relationship of the first 2 in each define
them all, so it should error in this case.

I'll try to put together an example that doe the resampling in R with
overlay, or its underlying functions, to avoid external calls to GDAL.

The y values in the example data are irregular and so this is not supported by
Spatial[Grid/Pixels]DataFrame - although this is supported by image()
for xyz lists of this type.

range(diff(met$y))
[1] 0.09527307 0.13665998

You can convert this to points like this:

library(sp)
x <- data.frame(expand.grid(x = met$x, y = met$y), z = as.vector(met$z))
coordinates(x) <- ~x+y

If you know the projection and there is a regular transform that works
- this is common in NetCDF files and environmental models where an
underlying Mercator grid is use in lat/long form - only Y is
"irregular". Having a wild guess:

library(rgdal)
 proj4string(x) <- CRS("+proj=longlat +ellps=WGS84")
x2 <- spTransform(x, CRS("+proj=merc +ellps=WGS84"))
gridded(x2) <- TRUE

image(x2)

That works to create a regular grid and you could transform companion
datasets in the same way, but may not be useful depending on your
application. If you need to warp the grid to regular raster you can do
that externally - with the GDAL command line utilities, which I have
in my path. E.g.

writeGDAL(x2, "file2convert.tif")
system("gdalwarp file2convert.tif LL.tif -t_srs \"+proj=longlat +ellps=WGS84\"")

ll.x <- readGDAL("LL.tif")
summary(ll.x)
Object of class SpatialGridDataFrame
Coordinates:
       min       max
x -77.35186 -64.60699
y -58.38184 -41.03831
Is projected: FALSE
proj4string : [+proj=longlat +ellps=WGS84 +no_defs]
Number of points: 2
Grid attributes:
 cellcentre.offset  cellsize cells.dim
x         -77.28616 0.1313904        97
y         -58.31615 0.1313904       132
Data attributes:
   Min.  1st Qu.   Median     Mean  3rd Qu.     Max.
0.003633 0.438200 0.438200 0.426500 0.438200 0.983700

image(ll.x)
contour(met, add = TRUE)

HTH

Regards, Mike.



More information about the R-sig-Geo mailing list