gIntersection with drop_lower_td=T drops entire intersection

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]]