Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Split polygon parts of a single SpatialPolygons Object

In R, I have single SpatialPolygons object (i.e. multi-polygons) containing several hundred polygons. I would like to split this SpatialPolygons object into a list of Polygons (i.e. holes should remain attached to the parent polygon).

Any idea how to do this?

EDITED:

Using the following example provided in the sp package:

# simple example, from vignette("sp"):
Sr1 = Polygon(cbind(c(2,4,4,1,2),c(2,3,5,4,2)))
Sr2 = Polygon(cbind(c(5,4,2,5),c(2,3,2,2)))
Sr3 = Polygon(cbind(c(4,4,5,10,4),c(5,3,2,5,5)))
Sr4 = Polygon(cbind(c(5,6,6,5,5),c(4,4,3,3,4)), hole = TRUE)

Srs1 = Polygons(list(Sr1), "s1")
Srs2 = Polygons(list(Sr2), "s2")
Srs3 = Polygons(list(Sr3, Sr4), "s3/4")
SpP = SpatialPolygons(list(Srs1,Srs2,Srs3), 1:3)

Then running out = lapply(SpP@polygons, slot, "Polygons"). I get a list of three Polygons (i.e. Srs1, Srs2, Srs3).

However, the case I am trying to solve is a bit different from this example. The SpatialPolygons object I am trying to split is the result of a geometric union done with the gUnaryUnion function (in the RGEOS package). If I apply out <- lapply(merged.polygons@polygons, slot, "Polygons"), I get a unique list of Polygon objects (n.b. not a list of Polygons objects). In other words, each polygon is separated from its hole(s).

Running topol <- sapply(unlist(out), function(x) x@hole)

I get:

> length(topol)
[1] 4996


> sum(topol, na.rm=TRUE)
[1] 469

According to the RGEOS v0.3-2 manual (http://cran.r-project.org/web/packages/rgeos/rgeos.pdf):

In order for rgeos to function properly it is necessary that all holes within a given POLYGON or MULTIPOLYGON geometry must belong to a specific polygon. The SpatialPolygons class implementation does not currently include this information. To work around this limitation rgeos uses an additional comment attribute on the Polygons class that indicates which hole belongs to which polygon. Under the current implementation this comment is a text string of numbers separated by spaces where the order of the numbers corresponds to the order of the Polygon objects in the Polygons slot of the Polygons object. A 0 implies the Polygon object is a polygon, a non-zero number implies that the Polygon object is a hole with the value indicating the index of the Polygon that “owns” the hole.

So the createSPComment() function in RGEOS is likely to be a workaround to reaggregate Polygons and holes.

like image 881
jatobat Avatar asked Oct 11 '13 10:10

jatobat


People also ask

How do you split a polygon?

To split a polygon, use the Cut Polygons tool, then draw a line across the polygon. The cut operation updates the shape of the existing feature and creates one or more new features. If there is no domain assigned to a field, the attribute values are copied from the original feature to the new feature.

What is the split polygon tool?

Used to split a selected polygon into two polygons. If the polygons have properties, then the new polygons will have those same properties. Any annotation for the initial polygon will be deleted.


1 Answers

As I understand it, the OP wants to convert a SpatialPolygons object into a list of Polygons, preserving holes if present.

The SpP object created by the OP contains three polygons, the third of which has an associated hole.

You can use lapply to cycle through each polygon in SpP, returning a list of SpatialPolygons. The difference between a Polygons and SpatialPolygons object is the addition of plot order information. Since each resulting SpatialPolygons is of length = 1, however, this information is superfluous.

n_poly <- length(SpP)

out <- lapply(1:n_poly, function(i) SpP[i, ])

lapply(out, class)

> lapply(out, class)
   [[1]]
   [1] "SpatialPolygons"
   attr(,"package")
   [1] "sp"

   [[2]]
   [1] "SpatialPolygons"
   attr(,"package")
   [1] "sp"

   [[3]]
   [1] "SpatialPolygons"
   attr(,"package")
   [1] "sp"

plot(out[[3]]) # Hole preserved

If a list of Polygons is needed, simply pull the appropriate slot from the SpatialPolygons object:

out <- lapply(1:n_poly, function(i) SpP[i, ]@polygons[[1]])

lapply(out, class)

> lapply(out, class)
[[1]]
[1] "Polygons"
attr(,"package")
[1] "sp"

[[2]]
[1] "Polygons"
attr(,"package")
[1] "sp"

[[3]]
[1] "Polygons"
attr(,"package")
[1] "sp"
like image 122
adamdsmith Avatar answered Sep 19 '22 11:09

adamdsmith