[R-sig-Geo] Union/Overlay single spatialPolygonsDataFrame containing overlapping polygonsto create polygons with all unique combinations of attributes

[nevil amos]

I am attempting to combine overlapping polygons in a
single spatialPolygonsDataFrame so that the output is a set of polygons
such that where there is no overlap one polygon, of the non-intersecting
area will be formed, and retain the attributes of the original polygon.
 where polygons intersect  new polygons, of the intersecting area only will
be formed, one with the attributes of each of the original polygons
intersecting in that area.

The example includes a single overlap, the real data may contain many tens
or hundreds of intersecting polygons, from among tens of thousands of
polygons.

the step that I cannot work out is "myUnion"  it is aiming to achieve the
equivalent of a self- union in ESRI geoprocessing.

I realize that I can achieve a similar result with raster but would prefer
to stay with polygons for this process.

many thanks in advance

Nevil Amos

Trivial example:
library(rgeos)
library(maptools)
library (tidyr)
library(plyr)
rm(list = ls())
par(mfrow=c(2,2))
cols=c("red", "black", "blue", "yellow")
cols<-adjustcolor(cols,alpha.f=.5)

p2<-
SpatialPolygons(list(Polygons(list(Polygon(cbind(c(1,4,4,1,1),c(1,1,4,4,1)),hole
= F)), 1 ),

 Polygons(list(Polygon(cbind(c(3,5,5,3,3),c(3,3,5,5,3)),hole = F)), 2)))

ID<-c(1,2)
YEAR<-c(1999,2002)
Type<-c(1,2)
mydf2<-data.frame(ID,YEAR,Type)
mySPDF2<-SpatialPolygonsDataFrame(p2,mydf2)
mySPDF2$Area<-gArea(mySPDF2, byid=TRUE)
print("mySPDF2 at data")
print(mySPDF2 at data)
plot(mySPDF2,axes=T,col=cols)
mtext("mySPDF2",cex=2)
mtext(paste("colors =","red", "yellow", "blue", "green,"), outer = TRUE,
cex = 1.5)
#this is the step that does not give the desired result. It has the same
two polygons ( areas 9 and 4) as the input, no intersection and clipping of
the polygons occurs   I am trying to acheive the equivalent of a "self
union in ESRI ARCGIS.  The result I would like be four polygons: two "L
shaped plygons formed from the non- intersecting parts of the input
polygons (Areas 8 and 3 respecively), and two square polygons ( coincident
both of Area = )  I manually create my desired output in the next section
to illustrate what I would like to achieve
myUnion<-unionSpatialPolygons(mySPDF2,IDs=mySPDF2$ID)
print("myUnion Areas:")
print (gArea(myUnion, byid=TRUE))
plot(myUnion,col=cols,axes=T)
mtext("myUnion",cex=2)
#Create example of desired output of "union" operation
p4<-SpatialPolygons(list(Polygons(list(Polygon(cbind(c(1,4,4,3,3,1,1),c(1,1,3,3,4,4,1)),hole
= F)), 1),

 Polygons(list(Polygon(cbind(c(4,4,5,5,3,3,4),c(4,3,3,5,5,4,4)),hole = F)),
2),

 Polygons(list(Polygon(cbind(c(3,4,4,3,3),c(3,3,4,4,3)),hole = F)), 3),

 Polygons(list(Polygon(cbind(c(3,4,4,3,3),c(3,3,4,4,3)),hole = F)), 4)))
ID<-c(1,2,3,4)
YEAR<-c(1999,2002,1999,2002)
Type<-c(1,2,1,2)
mydf4<- data.frame(ID,YEAR,Type)
mySPDF4<-SpatialPolygonsDataFrame(p4,mydf4)


centroids<-coordinates(mySPDF4)
mySPDF4$X_Centroid<-centroids[,1]
mySPDF4$Y_Centroid<-centroids[,2]
mySPDF4$MergeID<-paste(mySPDF4$X_Centroid,mySPDF4$Y_Centroid)
mySPDF4$MergeID<-unclass(as.factor(mySPDF4$MergeID))
mySPDF4$Area<-gArea(mySPDF4, byid=TRUE)
print("mySPDF4 at data")
print(mySPDF4 at data)
plot(mySPDF4, col=cols,axes=T)
mtext("mySPDF4",cex=2)
# the ultimate aim is to get a Flat ( ie non-overlapping polygons)  with
attributes that can then be joined back wide form table of to the YEAR and
Type of mySPDF4

#Flatten on "MergeID"
myUnion2<-unionSpatialPolygons(mySPDF4,IDs=mySPDF4$MergeID)
myDF<-as.data.frame(mySPDF4)
myDF<-arrange(myDF,MergeID,desc(YEAR))

#add sequnce number for subsequent conversion to wide form
for( row in 1:nrow(myDF)){
  thisVal = myDF$MergeID[row]
  if (row ==1) Seq = 1 else {if (thisVal==lastVal) Seq =Seq+1 else Seq =1}
  myDF$Sequence[row]<-Seq
  lastVal=thisVal
}
#make YEAR and Type into wide form
yearDF<-spread(myDF[,c("MergeID","Sequence",'YEAR')],Sequence,YEAR)
names(yearDF)<-c("ID","YEAR1","YEAR2")
typeDF<-spread(myDF[,c("MergeID","Sequence",'Type')],Sequence,Type)
names(typeDF)<-c("ID","Type1","Type2")
flatDF<-join(yearDF,typeDF)
rownames(flatDF)<-flatDF[,"ID"]
#append data frame to myUnion2 spatialpolygons
myUnion2<-SpatialPolygonsDataFrame(myUnion2,flatDF)
plot(myUnion2,col=cols,axes=T)
mtext("myUnion2",cex=2)
myUnion2$Area<-gArea(myUnion2, byid=TRUE)
print("myUnion2 at data")
print(myUnion2 at data)

	[[alternative HTML version deleted]]