[R-sig-Geo] gUnaryUnion Not Dissolving Correctly

[Michael Sumner]

Thanks for all this detail Roger, is there a way to "re-build" a spatial
object so that the given scale setting is applied? Are there any general
rounding or "orthogonalize" functions in the Spatial suite?

Cheers, Mike.

On Wed, 4 Nov 2015 at 18:16 Roger Bivand <Roger.Bivand at nhh.no> wrote:

> On Tue, 3 Nov 2015, Matt Strimas-Mackey wrote:
>
> > I'm working with a regular hexagonal grid stored as SPDF. At some
> > point I subset this SPDF, then want to combine all adjacent hexagons
> > together so that each contiguous set of hexagons is a single polygon.
> > I'm doing this last step using gUnaryUnion (or gUnionCascaded, not
> > clear what the different is actually). The problem is that in some
> > cases boundaries between clearly adjacent polygons are not dissolved.
> >
> > Example:
> >
> > ## Create three adjacent hexagons
> > library(sp)
> > library(rgeos)
> > p1 <- Polygon(cbind(
> >      c(1276503.26781119, 1281876.11747031, 1287248.96712942,
> >        1287248.96712942, 1281876.11747031, 1276503.26781119,
> 1276503.26781119),
> >      c(204391.40834643, 207493.42454344, 204391.40834643,
> 198187.37595242,
> >        195085.35975541, 198187.37595242, 204391.40834643)))
> > p2 <- Polygon(cbind(
> >      c(1287248.96712943, 1292621.81678854, 1297994.66644766,
> >        1297994.66644766, 1292621.81678854, 1287248.96712943,
> 1287248.96712943),
> >      c(204391.40834643, 207493.42454344, 204391.40834643,
> 198187.37595242,
> >        195085.35975541, 198187.37595242, 204391.40834643)))
> > p3 <- Polygon(cbind(
> >      c(1281876.11747031, 1287248.96712943, 1292621.81678854,
> >        1292621.81678854, 1287248.96712943, 1281876.11747031,
> 1281876.11747031),
> >      c(213697.45693745, 216799.47313446, 213697.45693745,
> 207493.42454344,
> >        204391.40834643, 207493.42454344, 213697.45693745)))
> > spoly <- SpatialPolygons(list(Polygons(list(p1, p2, p3), 's1')))
> > plot(gUnaryUnion(spoly))
>
> No, this is just numerical fuzz:
>
> plot(spoly)
> plot(gUnaryUnion(spoly), border="green", lty=3, lwd=2, add=TRUE)
> oS <- getScale()
> # default 1e+8
> setScale(1e+4)
> plot(gUnaryUnion(spoly), border="orange", lwd=2, add=TRUE)
> setScale(oS)
>
> JTS, GEOS, and consequently rgeos by default shift all coordinates to an
> integer grid after multiplying by a scale factor (finding integer matches
> is much easier than real matches). If the scaling is too detailed (in some
> cases), the operations do not give the expected outcomes.
>
> There is work in progress in GEOS and JTS to provide other scaling options
> and models, and to permit iteration over scaling values until a "clean"
> result is obtained (for some meanings of clean).
>
> gUnionCascaded() was the only possible function for GEOS < 3.3.0, from
> GEOS 3.3.0 gUnaryUnion() is available and the prefered and more efficient
> route. This is explained on the help page.
>
> Hope this clarifies,
>
> Roger
>
> >
> > Note that p2 and p3 are dissolved together, but p1 is separate. The
> > shared edge of p1 and p2 is:
> > p1:
> > [2,] 1281876 207493.4
> > [3,] 1287249 204391.4
> > p2:
> > [5,] 1287249 204391.4
> > [6,] 1281876 207493.4
> >
> > So, exactly the same apart from the order. I originally thought this
> > difference in order might be the problem, but this doesn't seem to be
> > an issue with in this example, where order is also flipped:
> > sss <- rasterToPolygons(raster(nrow=2, ncol=2, xmn=0, xmx=2, ymn=0,
> > ymx=2, vals=1:4))
> > lapply(sss at polygons, function(x) slot(x, 'Polygons')[[1]]@coords)
> > plot(sss)
> > plot(gUnaryUnion(sss))
> >
> > Session Info:
> > R version 3.2.2 (2015-08-14)
> > Platform: x86_64-apple-darwin13.4.0 (64-bit)
> > Running under: OS X 10.10.5 (Yosemite)
> >
> > Message when rgeos is loaded:
> >
> > rgeos version: 0.3-14, (SVN revision 511)
> > GEOS runtime version: 3.5.0-CAPI-1.9.0 r0
> > Linking to sp version: 1.2-0
> > Polygon checking: TRUE
> >
> > Any help on how to get these polygons to dissolve is appreciated.
> >
> > M
> >
> > _______________________________________________
> > 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
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

	[[alternative HTML version deleted]]