[R-sig-Geo] Holes in polygons disappear after performing checkPolygonsHoles from maptools
Roger Bivand
Roger.Bivand at nhh.no
Fri Apr 3 20:42:56 CEST 2015
On Fri, 3 Apr 2015, Grigory Alexandrovich wrote:
>
>
> Am 03.04.2015 um 13:16 schrieb 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|
>> >
>> >
> My session infos:
>
> sessionInfo()
> R version 3.1.1 (2014-07-10)
> Platform: x86_64-w64-mingw32/x64 (64-bit)
>
> ...
> attached base packages:
> [1] stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] rgdal_0.9-1 maptools_0.8-34 sp_1.0-17
>
> loaded via a namespace (and not attached):
> [1] foreign_0.8-61 grid_3.1.1 lattice_0.20-29 tools_3.1.1
>
>
> rgeos::version_GEOS()
> [1] "3.4.2-CAPI-1.8.2 r3921"
>
> rgdal::getGDALVersionInfo()
> [1] "GDAL 1.11.0, released 2014/04/16"
>
> rgdal::getPROJ4VersionInfo()
> [1] "Rel. 4.8.0, 6 March 2012, [PJ_VERSION: 480]"
>
>
> But I also obtained the same result on another windows machine with R 3.1.2
> (dont have session info for it)
How often do you run update.packages() ? This can be done any time in the
3.1 release cycle, as packages are compatible with R releases at that
level unless otherwise stated. rgdal_0.9-1 was built linking against
sp_1.0-15, as we don't lock-step packages depending on sp with its
releases. So some of the code for creating sp objects being used in rgeos
and rgdal at the C level is possibly out of date. Could you please try
update.packages() with regard to the last problem which I cannot
reproduce?
Roger
>
> Grigory
>
>
>
--
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
More information about the R-sig-Geo
mailing list