[R-sig-Geo] gIntersection with drop_lower_td=T drops entire intersection

Pius Korner pius.korner at bluewin.ch
Thu Oct 1 18:11:37 CEST 2015


Dear all

 

I intersect a SpatialPolygons-Object with two polygons 1a and 1b with a
polygon 2. Polygon 1a intersects with polygon 2 in one area but also only
touches polygon 2 at another position.

Here a schematic case for illustration:

 

poly1 <-
SpatialPolygons(list(Polygons(list(Polygon(coords=matrix(c(0,0,2,2,0,1,1,0),
ncol=2,byrow=F))),ID=c("a")),

 
Polygons(list(Polygon(coords=matrix(c(0,0,2,2,2,3,3,2),ncol=2,byrow=F))),ID=
c("b"))))

poly2 <-
SpatialPolygons(list(Polygons(list(Polygon(coords=matrix(c(0,0,2,2,1,1,1,3,3
,0,0,2),ncol=2,byrow=F))),ID=c("c"))))

plot(poly1,border="orange")

plot(poly2,border="blue",add=T,lty=2,density=8,angle=30,col="blue")

 

I'm interested in the area of overlap of polygons 1a and 1b with polygon 2.

gIntersection(poly1,poly2,byid=T)

does not work because of the different subgeometries that result from the
intersection of 1a with 2; they cannot be collected.

gIntersection(poly1,poly2,byid=T,drop_lower_td=T)

drops the entire intersection between 1a with 2, hence I get no value for
1a, even though it does overlap with polygon 2.

 

Is there a way so that only the Point-subgeometries (with area=0) are
dropped, but the Polygon-subgeometry is retained and returned for the
intersection of 1a and 2?

 

my temporary workaround was to very slightly buffer poly1 (width=-0.0001),
to prevent the point-subgeometry.

I could create a loop and intersect each polygon in poly1 separately, but
using byid=T would be nicer (faster..)

 

many thanks

Pius


	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list