Expose edge tessellator as a wk filter

[paleolimbot]

The S2EdgeTessellator is a C++ class that can be applied to import and export to approximate implicit cartesian edges (on import) or implicit geodesic edges (on export). The phrasing of these as a wk filter written in C requiring a C API for the edge tessellator seems a bit obscure...I promise it will pay off! Even without any additional infrastructure, one could use these to expose the tessellator to sf:

library(s2)
library(sf)
#> Warning: package 'sf' was built under R version 4.0.5
#> Linking to GEOS 3.8.1, GDAL 3.2.0, PROJ 7.2.0

tessellate_geod <- function(x, distance, radius = s2_earth_radius_meters()) {
  coords_tes <- wk::wk_handle(
    x,
    s2_unprojection_filter(
      s2_projection_filter(
        wk::sfc_writer(),
        tessellate_tol = distance / radius
      )
    )
  )
  
  coords_tes %>% st_set_crs(st_crs(x))
}

tessellate_cartesian <- function(x, distance, radius = s2_earth_radius_meters()) {
  coords_tes <- wk::wk_handle(
    x,
    s2_unprojection_filter(
      s2_projection_filter(wk::sfc_writer()),
      tessellate_tol = distance / radius
    )
  )
  
  coords_tes %>% st_set_crs(st_crs(x))
}

plot(tessellate_geod(sf::st_as_sfc("LINESTRING (0 0, 0 45, -60 45)"), 100))

plot(tessellate_cartesian(sf::st_as_sfc("LINESTRING (0 0, 0 45, -60 45)"), 100))

^{Created on 2021-05-17 by the reprex package (v0.3.0)}

[codecov-commenter]

Codecov Report

Merging #115 (9086f5e) into master (8dbd305) will increase coverage by 0.48%.
The diff coverage is 99.21%.

❗ Current head 9086f5e differs from pull request most recent head c56d300. Consider uploading reports for the commit c56d300 to get more accurate results

@@            Coverage Diff             @@
##           master     #115      +/-   ##
==========================================
+ Coverage   93.55%   94.04%   +0.48%     
==========================================
  Files          37       40       +3     
  Lines        2731     2987     +256     
==========================================
+ Hits         2555     2809     +254     
- Misses        176      178       +2     

Impacted Files	Coverage Δ
src/wk-c-utils.c	98.82% <98.82%> (ø)
R/wk-utils.R	100.00% <100.00%> (ø)
src/s2-c-api.cpp	100.00% <100.00%> (ø)

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 8dbd305...c56d300. Read the comment docs.

[edzer]

This addresses the invalid Sudan geometry in rnaturalearth, and importing of GeoJSON in general, right?

[paleolimbot]

It doesn't happen to fix Sudan in this case, but it's one tool in the toolbox (along with rebuild and unary union). It does the same thing as st_segmentize() or lwgeom::st_geod_segmentize() but potentially faster and skipping the addition of unnecessary points (e.g., if the edge is already a geodesic edge).

sudan <- rnaturalearth::ne_countries(country = "Sudan", returnclass = "sf")
s2_is_valid(sudan)
s2_is_valid(tessellate_cartesian(sudan, 0.01))

[edzer]

For some reason I can't get naturalearth's Sudan boundaries to be valid, I tried:

library(sf)
# Linking to GEOS 3.11.1, GDAL 3.6.2, PROJ 9.1.1; sf_use_s2() is TRUE
library(s2)
sudan <- rnaturalearth::ne_countries(country = "Sudan", returnclass = "sf")
st_as_s2(sudan, planar = TRUE, rebuild = TRUE) |> st_as_sfc() |> st_make_valid() |> st_is_valid()
# [1] FALSE


st_geometry(sudan) |> st_as_binary() |> s2_geog_from_wkb(planar = TRUE, check = FALSE) |> s2_rebuild() |> s2_is_valid()
# [1] FALSE

o = s2_options(snap = s2_snap_precision(1e7)) # also tried 1e5, 1e6, 1e8, 1e9
st_geometry(sudan) |> st_as_binary() |> s2_geog_from_wkb(planar = TRUE, check = FALSE) |> s2_rebuild(options = o) |> s2_is_valid()
# [1] FALSE

What am I missing?

> sessionInfo()
R version 4.2.2 Patched (2022-11-10 r83330)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 22.04.1 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.20.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_GB.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_GB.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] s2_1.1.2  sf_1.0-10

loaded via a namespace (and not attached):
 [1] rnaturalearth_0.3.2 Rcpp_1.0.10         magrittr_2.0.3     
 [4] units_0.8-1         tidyselect_1.2.0    lattice_0.20-45    
 [7] R6_2.5.1            rlang_1.0.6         fansi_1.0.4        
[10] httr_1.4.4          dplyr_1.0.10        wk_0.7.1.9000      
[13] tools_4.2.2         grid_4.2.2          KernSmooth_2.23-20 
[16] utf8_1.2.3          cli_3.6.0           e1071_1.7-12       
[19] DBI_1.1.3           class_7.3-20        assertthat_0.2.1   
[22] tibble_3.1.8        lifecycle_1.0.3     vctrs_0.5.2        
[25] glue_1.6.2          sp_1.6-0            proxy_0.4-27       
[28] compiler_4.2.2      pillar_1.8.1        generics_0.1.3     
[31] classInt_0.4-8      jsonlite_1.8.4      pkgconfig_2.0.3    

[edzer]

It doesn't happen to fix Sudan

Ah, missed that - it's not supposed to fix it.