[R-sig-Geo] intersecting two polygons

[Jeff Laake]

I ran into the same problem with spatstat and then discovered the gpclib 
package.  I wrote a simple function below that converts an owin poly to 
a gpc class (see below).  The snippet of code below shows how to use 
intersect and get.pts from gpclib and then to construct a poly for owin 
that uses the intersection points to create a poly dataframe for owin in 
spatstat.  I've only used this with intersection of rectangles with 
other shapes and have not tested it broadly but I don't know of any 
reason it will not work with 2 general polygons. I mentioned this in 
correspondence with Adrian Baddeley but didn't send him this code yet. 
This worked for me but there may be a better solution.  I'm a newcomer 
to this forum. --jeff

      gpc.area=owin.gpc.poly(study.area)
      b=as(data.frame(x=x,y=y),"gpc.poly")
      inside.poly=get.pts(intersect(b,gpc.area))[[1]]
      xdf= data.frame(x=rev(inside.poly$x),y=rev(inside.poly$y))
      poly=xdf[!duplicated(xdf),]

owin.gpc.poly=function(window)
##################################################################################
# Converts an owin class composed of a single polygon to a gpc.poly
#
# Arguments:  window  - an owin class
#
# Value    :  gpc.poly from first polygon in owin
#
# Jeff Laake
# 18 April 2008
##################################################################################
{
if(is.null(window$bdry))
  return(as(cbind(c(window$xrange,rev(window$xrange)),
             c(rep(window$yrange[1],2),rep(window$yrange[2],2))),
             "gpc.poly"))
else
  return(as(cbind(window$bdry[[1]]$x,window$bdry[[1]]$y),"gpc.poly"))
}