[R-sig-Geo] rgeos::gDifference

[Roger Bivand]

On Mon, 5 Jan 2015, Edzer Pebesma wrote:

> When running the following script:
>
> library(sp)
> gt = GridTopology(cellcentre.offset = c(0,0), cellsize = c(1,1),
> 	cells.dim = c(3,3))
> sgdf = SpatialGridDataFrame(SpatialGrid(gt),
>    data.frame(r = round(runif(9, max = 10), 1)))
>
> gt = GridTopology(cellcentre.offset = c(.25,.25), cellsize = c(1,1),
> 	cells.dim  = c(1,1))
> grd.red = SpatialGrid(gt)
>
> library(rgeos)
> sp1 = as(sgdf, "SpatialPolygonsDataFrame")
> sp2 = as(grd.red, "SpatialPolygons")
>
> plot(sp1, axes = TRUE)
> plot(sp2, border = 'red', add = TRUE)
>
> dif = gDifference(sp2, gUnionCascaded(sp1))
> dif = gDifference(sp2, sp1)
>
>
> the last command gives the following error:
>
> Error in RGEOSBinTopoFunc(spgeom1, spgeom2, byid, id, drop_not_poly,
> "rgeos_difference") :
>  TopologyException: no outgoing dirEdge found at -0.25 -0.25
> Calls: gDifference -> RGEOSBinTopoFunc -> .Call
> Execution halted
>
> What does this mean? Why does gUnionCascaded(sp1) not do this?
>

gUnionCascaded() is deprecated, use gUnaryUnion(). This isn't the problem, 
the argument order is probably the cause:

dif = gDifference(sp1, sp2, byid=c(TRUE, FALSE))
length(dif)

works, preserving the input geometries, as does:

dif = gDifference(gUnaryUnion(sp1), sp2)
length(dif)

if you want to aggregate.

This:

dif = gDifference(sp2, gUnionCascaded(sp1))

gives a NULL result, as sp2 is within sp1, as does:

dif = gDifference(sp2, gUnaryUnion(sp1))

but:

dif = gDifference(sp2, sp1, byid=c(FALSE, TRUE))

preserves the input geometries. If it was clear what you wanted as output 
from the operation, it would be easier to see what is going on.

Best wishes,

Roger

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