development of alpha hull algorithm

alphadev.R

nc <- sf::st_read(system.file("shape/nc.shp", package = "sf"))

geo <- sf::st_centroid(nc) |> 
  sf::st_geometry()

tris <- geo |> 
  geos::geos_make_collection() |> 
  geos::geos_delaunay_triangles()

# this is what i would've preferred (and now do)
x <- geos::geos_unnest(tris)

# Identify the circumscribed circle
circumscribe_circle <- function(geom) {
  tri_xy <- sf::st_coordinates(sf::st_as_sf(geom))[1:3,1:2]
  circum_circ <- tripack::circumcircle(tri_xy[,1], tri_xy[,2])
  circum_circ
}


# apply the above function to all triangles
all_circs <- do.call("rbind", lapply(x, circumscribe_circle)) 

# create new wk circles object
circs <- wk::new_wk_crc(
  x = list(x = unlist(all_circs$x),
           y = unlist(all_circs$y),
           r = unlist(all_circs$r)),
  crs = attr(x, "crs")
  )

# pull the radii from the circles
radii <- unclass(circs)[["r"]]

# sort 
radii_sorted <- sort(radii, TRUE)

# calculate alpha
alpha_og <- 1 / radii_sorted[1]

# plot tesselations
plot(x[radii < alpha_og])

for alphahull

nc <- sf::st_read(system.file("shape/nc.shp", package = "sf"))

geo <- sf::st_centroid(nc) |> 
  sfheaders::sf_to_df() |> dplyr::select(x, y)


a <- alphahull::ahull(geo[,1], geo[,2], 1)
plot(a)

Right, alphahull is a good choice. However, I suspect it assumes planar coordinates:

geo <- sf::st_centroid(nc) |> sf::st_transform("EPSG:32119") |>
  sfheaders::sf_to_df() |> dplyr::select(x, y)
a <- alphahull::ahull(geo[,1], geo[,2], 1000000)
plot(a, asp=1)

where changing the tension on alpha gives varying effects. Please also note that tripack should be replaced by interp or deldir as it has an encumbered license.

@rsbivand & @mdsumner thank you! I've found that using sf::st_jitter()enables alphahull to work. This is after either using unique() on a matrix, or distinct on a data frame.

Roger please ... if course it does, it's not rocket science

use a projection, we're just getting basic shit to work the actual reality is immaterial at this level, where it's pure geometry

of course jitter would work, you'll get unique vertices which is exactly what the message says is required

don't get sucked into thinking the longlat thing has anything to do with it, or any of the other specious sf bs

Let's not get worked up over triangles and circles. :) I've tried both long lat and projected coords, both work after jittering.

I'm sorry Jos, you know it's just numbers ... I'm not going to participate in any nonsense about it being "non-planar" (when that, while real, is irrelevant)

I don't have any patience for that idiocy and that's why I'll just dip out 🙏