Calculating azimuth of nearest features

Question

I'm calculating nearest points between two sf objects (point). I need to get the nearest hub's name, the distance and the azimut direction. At this moment I got the name and the distance but I can't figure it out how to get the direction between ben's points and each nearest hub.

library(sf)

#create the point layers
x_hub <- c(483049,482685)
y_hub <- c(4691026,4686387)
x_origin <- c(489597,489935,479537,476802,486356)
y_origin <- c(4696493,4695278,4692577,4684002,4690383)
p1 <- st_point(c(x_hub[1],y_hub[1]))
p2 <- st_point(c(x_hub[2],y_hub[2]))
p3 <- st_point(c(x_origin[1],y_origin[1]))
p4 <- st_point(c(x_origin[2],y_origin[2]))
p5 <- st_point(c(x_origin[3],y_origin[3]))
p6 <- st_point(c(x_origin[4],y_origin[4]))
p7 <- st_point(c(x_origin[5],y_origin[5]))
hub <- st_sfc(p1,p2, crs = 25831)
atr_hub <- data.frame(id= c(1,2))
hub <- st_sf(atr_hub, geometry = hub)
origin <- st_sfc(p3,p4,p5,p6,p7, crs=25831)
atr_origin <- data.frame(id= c(1:5))
origin <- st_sf(atr_origin, geometry = origin)


#find the hub name and save in origin
hubname <- st_nearest_feature(origin,hub)
origin$hubname <- hubname

#distance matrix between origin and hub
dist_matrix <- st_distance(origin,hub, by_element=FALSE)

#get the minimum value in each row and save in origin
dist <- 1
for (i in 1:nrow(dist_matrix)){
  print(i)
  dist[i] <- min(dist_matrix[i,])
}
origin$dist <- dist

Thank you.

user_1339626 · Answer

Please see the modified script below to get the (planar) azimuth to nearest point, using package nngeo. The development version needs to be installed from GitHub for this to work:

remotes::install_github("michaeldorman/nngeo")

The script:

library(sf)
library(nngeo)

#create the point layers
x_hub <- c(483049,482685)
y_hub <- c(4691026,4686387)
x_origin <- c(489597,489935,479537,476802,486356)
y_origin <- c(4696493,4695278,4692577,4684002,4690383)
p1 <- st_point(c(x_hub[1],y_hub[1]))
p2 <- st_point(c(x_hub[2],y_hub[2]))
p3 <- st_point(c(x_origin[1],y_origin[1]))
p4 <- st_point(c(x_origin[2],y_origin[2]))
p5 <- st_point(c(x_origin[3],y_origin[3]))
p6 <- st_point(c(x_origin[4],y_origin[4]))
p7 <- st_point(c(x_origin[5],y_origin[5]))
hub <- st_sfc(p1,p2, crs = 25831)
atr_hub <- data.frame(id= c(1,2))
hub <- st_sf(atr_hub, geometry = hub)
origin <- st_sfc(p3,p4,p5,p6,p7, crs=25831)
atr_origin <- data.frame(id= c(1:5))
origin <- st_sf(atr_origin, geometry = origin)

#find the hub name and save in origin
hubname <- st_nearest_feature(origin,hub)
origin$hubname <- hubname

##############

# Find distances and azimuths
origin$dist = st_distance(origin, hub[hubname, ], by_element = TRUE)
origin$az = st_azimuth(origin, hub[hubname, ])

# Plot
plot(st_union(st_geometry(origin), st_geometry(hub)), col = NA)
plot(st_geometry(origin), col = "grey", add = TRUE)
plot(st_geometry(hub), col = "red", add = TRUE)
plot(st_connect(origin, hub), col = "grey", add = TRUE)
text(st_coordinates(origin), as.character(round(origin$az)))

And here is the plot produced by the last part of the script, showing the direction towards the nearest point along with the azimuth:

Calculating azimuth of nearest features

1 Answers