[R-sig-Geo] Clipping problem

[Edzer Pebesma]

And that is how it is documented: the help obtained by

methods?overlay

tells you about the method used:

     x = "SpatialPoints", y = "SpatialPolygons" returns a numeric
          vector of length equal to the number of points; the number is
          the index (number) of the polygon of ‘y’ in which a point
          falls; NA denotes the point does not fall in a polygon; if a
          point falls in multiple polygons, the last polygon is
          recorded.

so

plot(grd[!is.na(clip)])

should give you the pixels inside the polygon(s).

On 10/28/2010 10:23 PM, Peter Larson wrote:
> Hello all,
> 
> I have a grid of points and a spatial polygon shape file. I want to
> clip the grid so that I only have the point which lie within the
> shapefile.
> 
> I tried the following code:
> 
> ## create a grid onto which we will interpolate:
> ## first get the range in data
> x.range <- c(37,50)
> y.range <- c(29,38)
> 
> ## now expand to a grid with 500 meter spacing:
> grd <- expand.grid(x=seq(from=x.range[1], to=x.range[2], by=.01),
> y=seq(from=y.range[1], to=y.range[2], by=.01) )
> 
> ## convert to SpatialPixel class
> coordinates(grd) <- ~ x+y
> gridded(grd) <- TRUE
> 
> ## Read in Shapefile
> ccShapet=readShapePoly('IRQ_adm0') # here cbg00barncnty is an ArcGIS shapefile
> 
> ## Clip
> clip <- overlay(grd, ccShapet)
> 
> #Plot
> plot(clip)
> 
> BUT it does not work. I just get a vector of number and NA's, but not
> the clipped grid.
> 
> Any ideas?
> 
> Thanks
> 
> Pete
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

-- 
Edzer Pebesma
Institute for Geoinformatics (ifgi), University of Münster
Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251
8333081, Fax: +49 251 8339763  http://ifgi.uni-muenster.de
http://www.52north.org/geostatistics      e.pebesma at wwu.de