Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix st_difference on a single sfc so it works on entirely contained geometries correctly. #975

Merged
merged 2 commits into from Feb 18, 2019

Conversation

jmarshallnz
Copy link

In playing around with sf and some very simple geometries, I noticed that st_difference when run on a single sfc returned unexpected results in the case where earlier geometries are entirely contained in later geometries. e.g. in the existing test for st_difference with fully contained geometries, if we change the order of pl1 and pl2 (pl2 is entirely contained within pl1) then

library(sf)
pl1 = st_polygon(list(matrix(c(0, 0, 2, 0, 2, 2, 0, 2, 0, 0), byrow = TRUE, ncol=2)))
pl2 = st_polygon(list(matrix(c(0.5, 0.5, 1.5, 0.5, 1.5, 1.5, 0.5, 1.5, 0.5, 0.5), byrow = TRUE, ncol = 2)))
st_contains(pl1, pl2, sparse=FALSE)
out = st_difference(st_sfc(list(pl2, pl1)))
all.equal(out[[2]], pl1)

returns pl1 in it's second geometry, where we'd expect pl1 with a hole where pl2 is, as the usual st_difference gives the correct geometry with a hole:

out2 = st_difference(pl2, pl1)
all.equal(out2, out[[2]])

This is primarily due to GEOSOverlaps_r being called before we call GEOSDifference_r in CPL_nary_difference. Overlaps returns false in the case that one of the geometries are entirely contained (and in my understanding is symmetric, so that Overlaps(a,b) == Overlaps(b,a)).

The patch changes this to use GEOSIntersects_r which I'm guessing might be equivalent to Overlaps || Contains.

First commit adds to the existing test for entirely contained geometries, second fixes it :)

Jonathan Marshall added 2 commits February 15, 2019 09:56
…ns the first, which currently fails as st_difference calls st_overlaps which returns false for geometries that are entirely contained (regardless of order)
…es. Previously this code used GEOSOverlaps_r which returns false in the case that the geometries are entirely contained (regardless of which is the larger object). Instead, we need to use GEOSIntersects_r which returns true in this case, and in the case where they aren't entirely contained, yet intersect. An alternate might be to call both GEOSOverlaps and GEOSContains?
@edzer edzer merged commit 43836ba into r-spatial:master Feb 18, 2019
@edzer
Copy link
Member

edzer commented Feb 18, 2019

Thanks for the very clear report!

edzer added a commit that referenced this pull request Feb 18, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

2 participants