[geos-devel] rgeos interface to R classes (sp, others)

Martin Davis mbdavis at refractions.net
Mon Feb 1 12:48:57 EST 2010


Roger,

The technique of using buffer(0) to union polygons is now deprecated in 
favour of using Unary Union (Geometry.union() - not sure what the exact 
GEOS signature is).  Unary Union is usually faster and more robust than 
the previous technique.  You might want to check this out.

Martin

Roger Bivand wrote:
> Hi,
>
> I've felt that I've been making reasonable progress with interfacing 
> GEOS geometries and methods for the R language (cran.r-project.org), 
> in a draft contributed package rgeos:
>
> https://r-forge.r-project.org/projects/rgeos/
>
> However, I encountered a problem that I do not understand, and would 
> be very grateful if I could be put back on track. The specific problem 
> is that an R function uses the C API to dissolve polygon borders for 
> adjacent polygons sharing a category given as its second argument.
>
> I've used the buffer technique from the JTS documentation, and all was 
> well until I tried to dissolve borders between touching squares when 
> the coordinate measures were small (square side 0.1). When the squares 
> are 100 units, all is well, and GEOMTouches is TRUE. But for 0.1, 
> GEOMTouches is FALSE, and no dissolving takes place.
>
> This can be reproduced (I'm on RHEL5, x86_64) by installing R, 
> contributed packages sp and maptools from CRAN, GEOS (3.1.1 or 3.2.0), 
> and installing rgeos from R-Forge. Start R, say
>
> library(rgeos) example(unionSpatialPolygonsGEOS)
>
> and look for undissolved borders in:
>
> image(grd, axes=TRUE)
> plot(spol1, add=TRUE)
>
> but dissolved in:
>
> image(grdx, axes=TRUE)
> plot(spol1x, add=TRUE)
>
> I started on an alternative implementation using GEOMTouches in 
> unionSpatialPolygonsGEOS(..., buffer=FALSE) output to console, where 
> one sees in the example output:
>
>> spol1F <- unionSpatialPolygonsGEOS(as(spol, "SpatialPolygons"),
> +   as.character(spol$xvs), buffer=FALSE)
> # 4 squares within (0,0), (0.2,0.2) NE, NW, SW share category, SE doesn't
> npls: 3, nnpls: 3
> type[0] Polygon
> i: 0, j: 1, touches: 0
> i: 0, j: 2, touches: 0
> type[1] Polygon
> i: 1, j: 2, touches: 0
> type[2] Polygon
> out of function
> npls: 1, nnpls: 1
> type[0] Polygon
> out of function
>> spol1xF <- unionSpatialPolygonsGEOS(as(spolx, "SpatialPolygons"),
> +   as.character(spolx$xvs), buffer=FALSE)
> # 4 squares within (0,0), (200,200), same categories
> npls: 3, nnpls: 3
> type[0] Polygon
> i: 0, j: 1, touches: 1
> i: 0, j: 2, touches: 1
> type[1] Polygon
> i: 1, j: 2, touches: 1
> type[2] Polygon
> out of function
> npls: 1, nnpls: 1
> type[0] Polygon
> out of function
>
> The C code is in:
>
> https://r-forge.r-project.org/plugins/scmsvn/viewcvs.php/pkg/src/rgeos_sp.c?rev=44&root=rgeos&view=log 
>
>
> around line 466 rgeos_SpatialPolygonsUnion(), calling 
> rgeos_plspairUnion() - desperate test framework, or 
> rgeos_plsbufUnion() which I had thought worked, but which clearly 
> doesn't dissolve small squares.
>
> Very grateful for any help, this is a show-stopper, and I had hoped to 
> release rgeos before February (small chance now!)
>
> Roger
>

-- 
Martin Davis
Senior Technical Architect
Refractions Research, Inc.
(250) 383-3022



More information about the geos-devel mailing list