“Hole” argument ignored in call to Polygon function

Question

I am working with a dataset consisting of a list of coordinates for a set of polygons. Some of these polygons are holes. However, the dataset does not follow the convention of listing the coordinates of holes in counter-clockwise direction. Therefore, I have to check for each polygon whether they are a hole, and use the result as argument to the Polygon function. Here is a simple dummy dataset:

library(sp)

# Coordinates for two dummy polygons (squares), one of which is fully contained within the other

x1 <- c(5,9,9,5,5)
y1 <- c(5,5,1,1,5)
x2 <- c(6,8,8,6,6)
y2 <- c(4,4,2,2,4)

tables <- list()    # Will contain the coordinates of each object, each in a data.frame
tables[[1]] <- data.frame(x1,y1)
tables[[2]] <- data.frame(x2,y2)

nobjs <- 2  # Number of objects

Next, I implemented an ugly but functional way of determining whether a polygon is a hole (i.e. fully contained within another polygon), by checking whether all of the nodes of a certain polygon are contained in any other polygon of the dataset:

holes <- logical(nobjs)
# Check for polygons that are holes
for(i in 1:nobjs){
    vectTF <- vector("logical",length(tables[[i]]$V1))
    for(j in 1:length(tables[[i]][,1])){
        vectTFj <- vector("logical",length(tables))
        for(k in 1:length(tables)){
            if(k==i){
                vectTFj[k] <- FALSE
            } else {
            vectTFj[k] <- point.in.polygon(tables[[i]][j,1],tables[[i]][j,2],tables[[k]][,1],tables[[k]][,2])
            }
        }
        vectTF[j] <- as.logical(sum(vectTFj))
    }
    holes[i] <- as.logical(sum(vectTF))
}

In this example, the vector holes correctly tells me that my second object is contained within another polygon, while the first one is not:

holes

[1] FALSE TRUE

Next, I try to build a Polygons object using the vector holes as the hole argument to the Polygon function:

objs <- list()
id <- integer(nobjs)
for(i in 1:nobjs){
    objs[[i]] <- Polygons(list(Polygon(tables[[i]][,1:2],hole=holes[i])),as.character(i))
    id[i] <- i
}

Now, when I check the structure of the new Polygons object, I see that none of my polygons are marked as being holes:

str(objs)
List of 2
$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. ..@ Polygons :List of 1
.. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. ..@ labpt : num [1:2] 7 3
.. .. .. .. ..@ area : num 16
.. .. .. .. ..@ hole : logi FALSE
.. .. .. .. ..@ ringDir: int 1
.. .. .. .. ..@ coords : num [1:5, 1:2] 5 9 9 5 5 5 5 1 1 5
.. .. .. .. .. ..- attr(, "dimnames")=List of 2
.. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. ..$ : chr [1:2] "x1" "y1"
.. ..@ plotOrder: int 1
.. ..@ labpt : num [1:2] 7 3
.. ..@ ID : chr "1"
.. ..@ area : num 16
$ :Formal class 'Polygons' [package "sp"] with 5 slots
.. ..@ Polygons :List of 1
.. .. ..$ :Formal class 'Polygon' [package "sp"] with 5 slots
.. .. .. .. ..@ labpt : num [1:2] 7 3
.. .. .. .. ..@ area : num 4
.. .. .. .. ..@ hole : logi FALSE
.. .. .. .. ..@ ringDir: int 1
.. .. .. .. ..@ coords : num [1:5, 1:2] 6 8 8 6 6 4 4 2 2 4
.. .. .. .. .. ..- attr(, "dimnames")=List of 2
.. .. .. .. .. .. ..$ : NULL
.. .. .. .. .. .. ..$ : chr [1:2] "x2" "y2"
.. ..@ plotOrder: int 1
.. ..@ labpt : num [1:2] 7 3
.. ..@ ID : chr "2"
.. ..@ area : num 4

Can anyone tell me what I'm doing wrong?

Answer 1

This is because a hole should always belong to a non-hole - i.e. you cannot have a Polygons where the only ring is a hole.

Note that the hole argument is not ignored by Polygon(), it is overridden and the coordinate order reversed by ?Polygons() - "In Polygons, if all of the member Polygon objects are holes, the largest by area will be converted to island status."

Please note that "rgeos::createSPComment" is for this purpose, it's discussed also in ?Polygons.

I think you are meaning to do this:

```R

id <- integer(nobjs)
objs <- list()
for(i in 1:nobjs){
    objs[[i]] <- Polygon(tables[[i]][,1:2],hole=holes[i])
    id[i] <- i
}

Polygons(objs, "1")

```

Now it's a hole, in the sense that it's a non-parent part in a multi-part object.

Note that there's no "nesting" that makes a hole a child of another part, it's just a flat list and only one is not allowed be a hole, arbitrarily set as the largest ring - that means you can give it an island and call it a hole and the interpretation is nonsensical. The validation of this is left to other tools in rgeos, and presumably on export via rgdal and so on - see rgeos::createSPComment.

Note that assigning "hole" status here reverses the coordinate order, if it's determined to be clockwise. The rules are "clockwise is an island, anti- a hole" and "assigning it different to the coordinate order results in reverse-winding".

Once it goes into Polygons that overrides any assignment for the largest in the list (the authors of sp have noted that it's in need of updating given the OGC sfr standards).