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

Robert J. Hijmans r.hijmans at gmail.com
Wed Aug 12 03:02:16 CEST 2015


I would think raster::union(polygons) should get you there.
Robert

On Tue, Aug 11, 2015 at 7:58 AM, Michael Sumner <mdsumner at gmail.com> wrote:
> On Mon, 10 Aug 2015 at 17:52 nevil amos <nevil.amos at gmail.com> wrote:
>
>> 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.
>>
>>
>
> Is raster::cover closer to what you want?
>
> library(raster)
>
> plot(mySPDF2)
>  plot(cover(mySPDF2[1,], mySPDF2[,2]), add = TRUE, col = rgb(c(0, 1), c(0,
> 0), c(1, 0), c(0.5)))
>  plot(mySPDF2)
>  plot(cover(mySPDF2[2,], mySPDF2[1,]), add = TRUE, col = rgb(c(0, 0), c(0,
> 1), c(1, 0), c(0.5)))
>
> You'll still need to take the non-overlapping piece from the complement
> overlay. I don't know if there's an easier way.
> See ?"raster-package" for a complete overview.
>
> HTH, Mike.
>
>
> 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]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo



More information about the R-sig-Geo mailing list