r-spatial / s2

Spherical Geometry Operators Using the S2 Geometry Library
https://r-spatial.github.io/s2/
71 stars 9 forks source link

invalid polygon #222

Open edzer opened 1 year ago

edzer commented 1 year ago

I was hoping that this polygon:

x

would be made valid by:

library(s2)
x = s2_geog_from_text("POLYGON ((0 0, 3 0, 3 2, 1 1, 2 1, 0 2, 0 0))", check = FALSE)
o = s2_options(split_crossing_edges = TRUE)
s2_rebuild(x, o) |> s2_is_valid()
# [1] FALSE

Zooming in, the resulting geometry

library(sf)
# Linking to GEOS 3.11.1, GDAL 3.6.2, PROJ 9.1.1; sf_use_s2() is TRUE
s2_rebuild(x, o) |> st_as_sfc() -> y
y[[1]]
# POLYGON ((0 0, 3 0, 3 2, 1.5 1.250243, 0 2, 0 0), (2 1, 1 1, 1.5 1.250243, 2 1))

seems valid (inner ring CW, touches outer ring in a point), and

st_crs(y) = NA # so that GEOS is_valid is called:
st_is_valid(y)
# [1] TRUE

is considered valid by GEOS. Does s2geometry have a different concept of valid, or does the package handle this?

edzer commented 1 year ago

Interestingly,

library(s2)
x = s2_geog_from_text("POLYGON ((0 0, 3 0, 3 2, 1 1, 2 1, 0 2, 0 0))", check = FALSE)
o = s2_options(split_crossing_edges = TRUE, snap = s2_snap_precision(0))
s2_rebuild(x, o) |> s2_is_valid()
# [1] TRUE
# Warning messages:
# 1: In cpp_s2_rebuild(as_s2_geography(x), options) :
#   NAs introduced by coercion to integer range
# 2: In cpp_s2_rebuild(as_s2_geography(x), options) :
#   NAs introduced by coercion to integer range

but the warnings here seem to indicate that something went wrong, as x didn't change:

library(sf)
# Linking to GEOS 3.11.1, GDAL 3.6.2, PROJ 9.1.1; sf_use_s2() is TRUE
s2_rebuild(x, o) |> st_as_sfc() -> y
# Warning messages:
# 1: In cpp_s2_rebuild(as_s2_geography(x), options) :
#   NAs introduced by coercion to integer range
# 2: In cpp_s2_rebuild(as_s2_geography(x), options) :
#   NAs introduced by coercion to integer range
y[[1]]
# POLYGON ((0 0, 3 0, 3 2, 2 1, 1 1, 0 2, 0 0))
paleolimbot commented 1 year ago

First, the warning is coming from s2_snap_precision(0) (i.e., log10(0)). Did you mean s2_snap_distance(0)?

Maybe a unary union is needed to fix this one?

library(s2)

x = s2_geog_from_text("POLYGON ((0 0, 3 0, 3 2, 1 1, 2 1, 0 2, 0 0))", check = FALSE)
o = s2_options(split_crossing_edges = TRUE)
x |> 
  s2_is_valid()
#> [1] FALSE

x |> 
  s2_rebuild(o) |> 
  s2_is_valid()
#> [1] FALSE

x |> 
s2_rebuild(o) |> 
  s2_union() |> 
  s2_is_valid()
#> [1] TRUE

Created on 2023-03-30 with reprex v2.0.2