计算最近要素的方位角

Question

我正在计算两个点的最近距离，需要获取最近枢纽名称、距离和方位角方向。目前我已获得了名称和距离，但我不知道如何获取本点与每个最近枢纽之间的方向。

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

谢谢。

- Josep Pueyo

---

你好，你的示例不可复制，因为只有你可以访问你的.shp文件。你能否修改你的帖子，提供一个可复制的示例？ - bretauv

我编辑了一些我的数据点。现在它是可重复的。感谢您的建议。 - Josep Pueyo

Answer 1

请查看以下修改后的脚本，使用nngeo包获取到最近点的（平面）方位角。需要从GitHub安装开发版本才能正常工作：

remotes::install_github("michaeldorman/nngeo")

脚本：

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)))

这是脚本最后部分生成的图表，显示了朝向最近点的方向以及方位角：

- Michael Dorman

---

非常感谢你。还有让我的代码更高效和优雅。我是一个自学者，我知道我写了很多无用的代码行。一步一步来。 - Josep Pueyo