[R-sig-Geo] extracting information from a raster file

[Roger Bivand]

On Sun, 2 Nov 2008, Harry Kim wrote:

> Dear R-sig-geo users,
>
>    I have a question about extracting information from a raster file
> after I read it in using rgdal package, and I would be infinitely
> grateful if someone could help me.

The overlay() method is what you need:

library(sp)
data(meuse.grid)
data(meuse)
coordinates(meuse.grid) <- c("x", "y")
coordinates(meuse) <- c("x", "y")
gridded(meuse.grid) <- TRUE
gatp <- overlay(meuse.grid, as(meuse, "SpatialPoints"))
plot(meuse.grid[gatp,])
plot(meuse, pch=4, col="red", add=TRUE)

overlay() gives the indices you are looking for, so you can either index 
the whole object to get a similar object back, or index the column in your 
object:

xxx$band1[gatp]

Hope this helps,

Roger

>    The summary of raster file looks like this:
>
> #summary of raster file
>> summary(xxx)
> Object of class SpatialGridDataFrame
> Coordinates:
>   min max
> x -180 180
> y  -90  90
> Is projected: FALSE
> proj4string :
> [+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs +towgs84=0,0,0]
> Number of points: 2
> Grid attributes:
>  cellcentre.offset cellsize cells.dim
> x           -179.10      1.8       200
> y            -89.55      0.9       200
> Data attributes:
>     Min.   1st Qu.    Median      Mean   3rd Qu.      Max.      NA's
>    0.000     0.000     0.000    47.010     3.264 13540.000 26563.000
>
> I would like to write a code that extracts the value of associated
> cell given a specific location in latitude and longitude. I've noticed
> that xxx at data contains the values as 40000 x 1 matrix. I've wrote the
> following code assuming that the values are stored as in C programing
> (column wise seperation):
>
>
> #initiate values
> x_inc=xxx at grid@cellsize[1]
> y_inc=xxx at grid@cellsize[2]
>
> x_start=xxx at grid@cellcentre.offset[1]
> y_start=xxx at grid@cellcentre.offset[2]
>
> x_dim=xxx at grid@cells.dim[1]
> y_dim=xxx at grid@cells.dim[2]
>
> #given coordinate
> x=26
> y=45
>
> #find the index
> index=floor((x-x_start)/x_inc)*y_dim + floor( (y-y_start)/y_inc )
> xxx at data[index,1]
>
> The result does not seem to match what I should get. Could anybody
> explain to me how the data is stored in SpatialGridDataFrame? Your
> help would be much much appreciated.
>
> Thank you and have a good night
> Harry
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no