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

[Roger Bivand]

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

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no