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

Add s2_convex_hull(), add na.rm argument to s2_convex_hull_agg() #163

Merged
merged 6 commits into from
Jan 25, 2022

Conversation

paleolimbot
Copy link
Collaborator

This PR adds a feature-wise convex-hull generator based on the implementation by @spiry34 in #151. It also adds na.rm, since this is how the other _agg() functions are implemented. This PR completes #150.

@spiry34 you're the most familiar with this code...is there anything I've missed or should have considered? You mentioned that the convex hull for lines is different and I haven't checked this in GEOS. I should probably check this in BigQuery to make sure that this implementation is consistent with that (since that's how I've tested the other functions).

library(s2)

s2_convex_hull(
  c(
    "GEOMETRYCOLLECTION(POINT(3.6 43.2), POINT (0 0), POINT(3.61 43.21))", 
    NA
  )
)
#> <s2_geography[2]>
#> [1] <POLYGON ((0 0, 3.61 43.21, 3.6 43.2, 0 0))>
#> [2] <NA>

s2_convex_hull_agg(c("POINT (0 0)", NA), na.rm = FALSE)
#> <s2_geography[1]>
#> [1] <NA>
s2_convex_hull_agg(
  c("GEOMETRYCOLLECTION(POINT(3.6 43.2), POINT (0 0), POINT(3.61 43.21))", NA),
  na.rm = TRUE
)
#> <s2_geography[1]>
#> [1] <POLYGON ((0 0, 3.61 43.21, 3.6 43.2, 0 0))>

Created on 2022-01-16 by the reprex package (v2.0.1)

@codecov-commenter
Copy link

codecov-commenter commented Jan 16, 2022

Codecov Report

Merging #163 (7569004) into main (90fa128) will decrease coverage by 0.02%.
The diff coverage is 94.87%.

Impacted file tree graph

@@            Coverage Diff             @@
##             main     #163      +/-   ##
==========================================
- Coverage   94.78%   94.75%   -0.03%     
==========================================
  Files          42       42              
  Lines        3200     3220      +20     
==========================================
+ Hits         3033     3051      +18     
- Misses        167      169       +2     
Impacted Files Coverage Δ
src/geography.h 32.00% <0.00%> (-2.79%) ⬇️
R/s2-transformers.R 100.00% <100.00%> (ø)
src/geography-collection.h 73.60% <100.00%> (+0.42%) ⬆️
src/s2-transformers.cpp 96.70% <100.00%> (+0.14%) ⬆️

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 90fa128...7569004. Read the comment docs.

@paleolimbot
Copy link
Collaborator Author

It's a good point about the convex hull of a segment/point. I don't know enough about the application here to know how important this is...it certainly could be included. I can see how it's difficult/verbose to replicate (my attempt is below):

library(s2)

fix <- function(result, original, tolerance = 1e-14) {
  simple <- s2_rebuild(
    result,
    options = s2_options(
      snap_radius = tolerance
    )
  )
  
  simple_boundary <- s2_rebuild(
    s2_boundary(result),
    options = s2_options(
      snap_radius = tolerance,
      edge_type = "undirected"
    )
  )
  
  simple_empty <- s2_is_empty(simple)
  simple_boundary_empty <- s2_is_empty(simple_boundary)
  simple[simple_boundary_empty] <- original[simple_boundary_empty]
  simple[simple_empty & !simple_boundary_empty] <- 
    original[simple_empty & !simple_boundary_empty]
  
  simple
}

geos::geos_convex_hull("MULTIPOINT (3.6 43.2, 0 0, 3.61 43.21)")
#> <geos_geometry[1]>
#> [1] <POLYGON ((0 0, 3.6 43.2, 3.61 43.21, 0 0))>
s2_convex_hull("MULTIPOINT (3.6 43.2, 0 0, 3.61 43.21)")
#> <s2_geography[1]>
#> [1] <POLYGON ((0 0, 3.61 43.21, 3.6 43.2, 0 0))>
fix(
  s2_convex_hull("MULTIPOINT (3.6 43.2, 0 0, 3.61 43.21)"), 
  "MULTIPOINT (3.6 43.2, 0 0, 3.61 43.21)"
)
#> <s2_geography[1]>
#> [1] <POLYGON ((0 0, 3.61 43.21, 3.6 43.2, 0 0))>


geos::geos_convex_hull("POINT (0 0)")
#> <geos_geometry[1]>
#> [1] <POINT (0 0)>
s2_convex_hull("POINT (0 0)")
#> <s2_geography[1]>
#> [1] <POLYGON ((0 0, -5.72949748e-14 3.03663366e-16, -3.03663366e-16 -5.72949748e-14, 0 0))>
fix(s2_convex_hull("POINT (0 0)"), "POINT (0 0)")
#> <s2_geography[1]>
#> [1] <POINT (0 0)>


geos::geos_convex_hull("LINESTRING (0 0, 3.6 43.2)")
#> <geos_geometry[1]>
#> [1] <LINESTRING (0 0, 3.6 43.2)>
s2_convex_hull("LINESTRING (0 0, 3.6 43.2)")
#> <s2_geography[1]>
#> [1] <POLYGON ((1.51774324 21.6094428, 3.6 43.2, 0 0, 1.51774324 21.6094428))>
fix(s2_convex_hull("LINESTRING (0 0, 3.6 43.2)"), "LINESTRING (0 0, 3.6 43.2)")
#> <s2_geography[1]>
#> [1] <LINESTRING (0 0, 3.6 43.2)>

Created on 2022-01-16 by the reprex package (v2.0.1)

@paleolimbot
Copy link
Collaborator Author

BigQuery does this:

SELECT ST_CONVEXHULL(ST_GEOGFROMTEXT("POINT (0 1)"))
POINT(-5.73037018596824e-14 1)
SELECT ST_CONVEXHULL(ST_GEOGFROMTEXT("LINESTRING (0 0, 3.6 43.2)"))
LINESTRING(3.6 43.2, 1.51774323905594 21.6094428417663, 0 0)

Merge commit '90fa1280a77c53cb81d71ae51642c2504c7b670e'

Conflicts:
	man/s2_boundary.Rd
@spiry34
Copy link
Contributor

spiry34 commented Jan 25, 2022

I really appreciate your work, and returning a polygon is definitely the most stable option for many usages.
I will still try to add a dedicated function able to return either a point/line/polygon depending on s2loop contents but I didn't succeed yet (I was just working one more bit on this). This will allow to fit standards. Anyway this will be just an add-on and will not modify your current implementation.
Glad to see soon this new version on the CRAN :-)
Many thanks for your contribution, all the best,
Sylvain

@paleolimbot
Copy link
Collaborator Author

I'm about to merge this one, but if you do get a working implementation feel free to PR it in!

@paleolimbot paleolimbot merged commit 0c74988 into main Jan 25, 2022
@paleolimbot paleolimbot deleted the convex-hull branch January 26, 2022 00:33
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.

3 participants