(no subject)
On Tue, 5 Jan 2010, laure velez wrote:
Dear list member, Although the answer to this might be straightforward I'm struggling with this issue for weeks now : I have a shapefile (only one polygon) and a raster (.tif). I would like to know the value of each pixel on the raster which overlap with the polygon in the shapefile. I had various tries but none did work : library(rgdal) pol = readOGR(?shape.shp?, layer = ?shape?) grid = readGDAL(?raster.tif?) overlay(grid,pol) With this I obtain ?NA? for all pixels Could you help me ?
Please try: image(grid, axes=TRUE) plot(pol, add=TRUE) If the polygon doesn't show up in the correct place, do examine summary() output for both objects. I suspect that the coordinate values of the objects - seem in the bounding boxes in the summary output - do not overlap, probably because the coordinate reference systems differ. Roger
Thanks in advance, Laure Please find attached the shape and the raster Vous ?tes ?tudiants ? B?n?ficiez de Windows 7 ? 35? ! D?couvrez notre offre !
_________________________________________________________________ [[elided Hotmail spam]] [[alternative HTML version deleted]]
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