[R-sig-Geo] Holes in polygons disappear after performing checkPolygonsHoles from maptools

[Roger Bivand]

On Fri, 3 Apr 2015, Grigory Alexandrovich wrote:

> Am 02.04.2015 um 22:24 schrieb Roger Bivand:
>> On Thu, 2 Apr 2015, Roger Bivand wrote:
>>
>>> On Thu, 2 Apr 2015, Edzer Pebesma wrote:
>>>
>>>>
>>>>
>>>>  On 04/01/2015 08:03 PM, Grigory Alexandrovich wrote:
>>>>>  Hello all,
>>>>>>  I encountered the following unexpected result using
>>>> checkPolygonsHoles:
>>>>>>  # attach the worldmap as SpatialPolygonsDataFrame from the
>>>> package >  maptools
>>>>>  library(sp)
>>>>>  library(maptools)
>>>>>  data(wrld_simpl)
>>>>>>  # get a polygon with a hole
>>>>>  shape_with_hole <- wrld_simpl[5,]
>>>>>>  # plot it (hole is left white, surrounded by blue color)
>>>>>  plot(shape_with_hole, col = "blue")
>>>>>>  # perform checkPolygonsHoles
>>>>>  shape_with_hole at polygons <- lapply(shape_with_hole at polygons,
>>>>>  checkPolygonsHoles)
>>>>>>  # plot again, now holes aren't recognized as such
>>>>>  plot(shape_with_hole, col = "blue")
>>>>>>  # and even the original SpatialPolygonsDataFrame object is
>>>> changed !?
>>>>>  plot(wrld_simpl[5,], col = "blue")
>>>>>>  One irritating side effect here is that the original object
>>>> wrld_simpl
>>>>>  is also changed. This result looks to me like a bug, or have I
>>>> missed
>>>>>  something?
>>>>>
>>>>  It looks as if it is related to rgeos, as
>>>>
>>>>  avoidGEOS = FALSE
>>>>  library(sp)
>>>>  library(maptools)
>>>>  data(wrld_simpl)
>>>>  x = wrld_simpl
>>>>
>>>>  shape_with_hole <- wrld_simpl[5,]
>>>>
>>>>  if (avoidGEOS)
>>>>      gpclibPermit()
>>>>  shape_with_hole at polygons <- lapply(shape_with_hole at polygons,
>>>>  checkPolygonsHoles,
>>>>   avoidGEOS = avoidGEOS)
>>>>
>>>>  data(wrld_simpl)
>>>>  all.equal(x, wrld_simpl)
>>>
>>> No, but:
>>>
>>> library(sp)
>>> library(maptools)
>>> data(wrld_simpl)
>>> x = wrld_simpl
>>> shape_with_hole <- wrld_simpl[5,]
>>> pls <- slot(shape_with_hole, "polygons")
>>> pls1 <- lapply(pls, checkPolygonsHoles, avoidGEOS = FALSE)
>>> isTRUE(all.equal(pls, pls1))
>>> # TRUE
>>> gpclibPermit()
>>> pls1 <- lapply(pls, checkPolygonsHoles, avoidGEOS = TRUE)
>>> isTRUE(all.equal(pls, pls1))
>>> # FALSE
>>>
>>> where the Polygon objects are re-ordered and plotted out of order, so
>>> the deprecated gpclib code may be broken here. However, this isn't
>>> the problem,
>>
>> The gpclib code is OK, it re-orders the Polygon objects, but the plot
>> order is correct. The problem is as described below.
>>
>> Roger
>>
>>> but rather (perhaps) graphics::polypath. If it is used (by default it
>>> is), and with the default path fill mode:
>>>
>>> plot(x[5,], col="blue", usePolypath=TRUE, rule="winding")
>>>
>>> does not whiten the holes. However:
>>>
>>> plot(x[5,], col="blue", usePolypath=TRUE, rule="evenodd")
>>>
>>> does. Using:
>>>
>>> plot(x[5,], col="blue", pbg="white", usePolypath=FALSE)
>>>
>>> also works by overpainting the island with the holes in white.
>>>
>>> Please summarise the your SO thread.
>>>
>>> By the way, in R it is unusual for changes in one object to propagate
>>> in an unintended way to another object, and this is not happening
>>> here. The visual impression is being driven by polypath doing
>>> different things depending on its rule setting. Look at ?polypath for
>>> details.
>>>
>>> Well, this doesn't clarify, but it explains what is happeing, I hope ...
>>>
>>> Roger
>>>
>>>>
>>>>  returns a difference, but TRUE if avoidGEOS = TRUE.
>>>>
>>>
>
> |Thank you for your answers guys. But I am still not fully convinced
> that it is an issue of the plot function,
> because the object shape_with_hole edited with checkPolygonsHoles
> before, continues to behave in a strange way:
>
> # we check which polygons are flagged as holes.
> # The flags are still set properly,
> # although the plot() function didn't recognize them:
> sapply(shape_with_hole at polygons[[1]]@Polygons, slot, "hole")
>
> [1] FALSE  TRUE  TRUE  TRUE
>
> # load library rgdal for reprojection
> library(rgdal)
>
> # reproject with spTransform, just for testing
> shape_with_hole <- spTransform(shape_with_hole, CRS("+proj=longlat
> +ellps=WGS84 +datum=WGS84"))
>
> # after reprojection all flags are set to FALSE
> sapply(shape_with_hole at polygons[[1]]@Polygons, slot, "hole")
>
> [1] FALSE FALSE FALSE FALSE
>
>
> So it is still something wrong with shape_with_hole, or?

Please provide the output of sessionInfo(), rgeos::version_GEOS(), and 
rgdal::getGDALVersionInfo() with rgdal::getPROJ4VersionInfo(). I cannot 
reproduce the degradation of hole status following a no-op use of 
spTransform() with:

> sessionInfo()
R version 3.1.3 (2015-03-09)
Platform: x86_64-unknown-linux-gnu (64-bit)
Running under: Fedora 21 (Twenty One)
...
attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] rgdal_0.9-3     maptools_0.8-34 sp_1.0-17

loaded via a namespace (and not attached):
[1] foreign_0.8-63  grid_3.1.3      lattice_0.20-31 rgeos_0.3-8
[5] tools_3.1.3
> rgeos::version_GEOS()
[1] "3.4.2-CAPI-1.8.2 r3921"
> rgdal::getGDALVersionInfo()
[1] "GDAL 1.11.2, released 2015/02/10"
> rgdal::getPROJ4VersionInfo()
[1] "Rel. 4.9.1, 04 March 2015, [PJ_VERSION: 491]"

Roger


>
> Best
>
> Grigory|
>
>
> 	[[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
>

-- 
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 91 00
e-mail: Roger.Bivand at nhh.no