Skip to content

Instantly share code, notes, and snippets.

@tim-salabim

tim-salabim/nearest_neighbor_graph

Last active May 21, 2020 12:55
Show Gist options
  • Save tim-salabim/172efb527f0c8fa200c1e29209a11df0 to your computer and use it in GitHub Desktop.
Save tim-salabim/172efb527f0c8fa200c1e29209a11df0 to your computer and use it in GitHub Desktop.
Download ZIP
Raw
nearest_neighbor_graph
library(sf)
library(data.table)
library(mapview)
n = 1000
x = runif(n, 0, 2)
y = runif(n, 0, 1)
pts = st_as_sf(data.frame(id_pt = 1:n, x = x, y = y), coords = c("x", "y"))
pts_union = st_union(pts)
deltri = st_triangulate(pts_union, bOnlyEdges = TRUE)
deltri = st_as_sf(
data.frame(
id_ln = 1:(npts(deltri)/2),
geometry = st_cast(deltri, "LINESTRING")
)
)
pts_deltri_isec = as.data.table(st_join(deltri, pts))
lut = pts_deltri_isec[, (id_ln[which.min(sf:::CPL_length(geometry))]), by = "id_pt"]
shrtst = deltri[deltri$id_ln %in% lut$V1, ]
mapview(deltri, color = "#BEBEBE") +
mapview(shrtst, color = "red") +
mapview(pts, cex = 2, col.regions = "black")
@tim-salabim
Copy link
Author

shortest

@tim-salabim
Copy link
Author

constrain to certain distance:

library(sf)
library(data.table)
library(mapview)

n = 1000

x = runif(n, 0, 200)
y = runif(n, 0, 100)

pts = st_as_sf(data.frame(id_pt = 1:n, x = x, y = y), coords = c("x", "y"))
pts_union = st_union(pts)

deltri = st_triangulate(pts_union, bOnlyEdges = TRUE)

deltri = st_as_sf(
  data.frame(
    id_ln = 1:(npts(deltri)/2)
    , geometry = st_cast(deltri, "LINESTRING")
  )
)
deltri$lngth = st_length(deltri$geometry)
deltri = deltri[deltri$lngth <= 3, ]

mapview(deltri, color = "red", lwd = 2, legend = FALSE) +
  mapview(pts, cex = 2, col.regions = "black", legend = FALSE)

image

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment