[R-sig-Geo] maptools/SpatialPolygons2PolySet produces ascending instead of descending POS numbers for polygons with holes

Roger Bivand Roger.Bivand at nhh.no
Wed Aug 6 15:02:03 CEST 2014


On Wed, 6 Aug 2014, Daniel Rodolphe Schlaepfer wrote:

> Thank you for the suggestion. Unfortunately, this does not work for my 
> application: I plot the polygon on top of other another plot where I 
> need the holes of the polygon to be transparent. For instance, I don’t 
> want the African countries in the following examples to be covered by 
> white:
>
> library(maps)
> map()
> plot(polySP, density=10, angle=45, pbg=“white")

Obviously not, these are base graphics, and you have to consider the 
plotting order:

plot(polySP, density=10, angle=45, pbg="white", axes=TRUE)
map(add=TRUE)

>
> However, this was not the point that I try to ask about. My point is 
> that maptools/SpatialPolygons2PolySet does not translate into correct 
> PolySet objects as illustrated in my previous example.

Since you no longer need PolySet objects, why? Maybe the coercion method 
can be fixed, but it is not a priority - patch welcomed.

Roger

>
> Thanks,
> Daniel
>
> On Aug 6, 2014, at 12:52 PM, Roger Bivand <Roger.Bivand at nhh.no> wrote:
>
>> On Wed, 6 Aug 2014, Daniel Rodolphe Schlaepfer wrote:
>>
>>> Hello all,
>>>
>>> I try to convert a sp SpatialPolygons object with a hole to a PBS PolySet object - my ultimate goal is to plot the polygon with hatching considering the hole.
>>
>> The problem is not that you need a PolySet representation, but that you need to set the polygon background explicitly to a value other than "transparent" using the pbg= argument:
>>
>> plot(polySP, density=10, angle=45, pbg="white")
>>
>> or adjust par("bg") to suit. When hatching is used, polypath is not used, so automatic handling of holes in the plot method is not available.
>>
>> Roger
>>
>>>
>>> My problem is that maptools/SpatialPolygons2PolySet only produces PolySet objects with increasing POS even for polygons with holes; however, the documentation of PolySet indicates that “We adopt the convention that POS goes from 1 to n along an outer boundary, but from n to 1 along an inner boundary, regardless of rotational direction.”
>>>
>>>
>>> #Create simple doughnut-shaped polygon
>>> library(sp)
>>> library(maptools)
>>>
>>> coords1 <- matrix(c(108, -54, -108, -54, -108, 54, 108, 54, 108, -54), ncol=2, byrow=TRUE)
>>> coords2 <- matrix(c(36, -18, -36, -18, -36, 18, 36, 18, 36, -18), ncol=2, byrow=TRUE)
>>>
>>> polySP <- SpatialPolygons(list(Polygons(list(Polygon(coords1, hole=FALSE), Polygon(coords2, hole=TRUE)), ID=1)), proj4string=CRS("+proj=longlat +datum=WGS84"))
>>>
>>> #Convert sp-SpatialPolygons to PBS-PolySet with maptools function
>>> polyPBS <- SpatialPolygons2PolySet(polySP)
>>>
>>> #-> POS for SID == 2 are increasing and thus do not reflect PBS standards for a polygon with a hole
>>>  PID SID POS    X   Y
>>> 1    1   1   1  108 -54
>>> 2    1   1   2 -108 -54
>>> 3    1   1   3 -108  54
>>> 4    1   1   4  108  54
>>> 5    1   1   5  108 -54
>>> 6    1   2   1   36 -18
>>> 7    1   2   2   36  18
>>> 8    1   2   3  -36  18
>>> 9    1   2   4  -36 -18
>>> 10   1   2   5   36 -18
>>>
>>> I believe that it would require only a small change to maptools/SpatialPolygons2PolySet to produce PolySet objects that meet the PBS standards also for polygons with holes, i.e., replace the line
>>> 	POS <- c(POS, 1:k)
>>> with
>>> 	POS <- if(slot(srs[[j]], "hole")) c(POS, k:1) else c(POS, 1:k)
>>> inside the loops: for (i in 1:n) … for (j in 1:m) …
>>>
>>>
>>>
>>> Sincerely,
>>> Daniel Schlaepfer
>>>
>>>
>>> My session infos:
>>> sessionInfo()
>>> 	R version 3.1.1 (2014-07-10)
>>> 	Platform: x86_64-apple-darwin13.3.0 (64-bit)
>>>
>>> 	locale:
>>> 	[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
>>>
>>> 	attached base packages:
>>> 	[1] stats     graphics  grDevices
>>> 	[4] utils     datasets  methods
>>> 	[7] base
>>>
>>> 	other attached packages:
>>> 	[1] PBSmapping_2.67.60 maptools_0.8-30
>>> 	[3] sp_1.0-15
>>>
>>> 	loaded via a namespace (and not attached):
>>> 	[1] foreign_0.8-61  grid_3.1.1
>>> 	[3] lattice_0.20-29
>>>
>>>
>>> -------------------------------------------------------
>>> Daniel Schlaepfer, PhD
>>> University of Wyoming
>>> Laramie, WY 82071
>>>
>>> _______________________________________________
>>> 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
>
>

-- 
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