我认为解决此问题的一种方法是获取加纳道路的形状文件并执行地理空间操作以找到最近的道路。从那里您应该能够使用 Google 的 API 来获取行驶距离
此解决方案中的步骤是
- 下载道路形状文件
- 在我们的出发地和目的地之间找到一条“直线”
- 找出与这条直线相交的所有道路
- 在 Google 路线查询中使用最近的十字路口作为新起点
在本例中,我使用的 shapefile 来自 here
library(sf) ## geospatial operations
library(googleway) ## plotting and google API
注意:
googleway 是我的包,你需要一个 Google API 密钥才能使用它
## Setting API keys I'll be using, one for the maps and one for directions
set_key("my_map_key"), api = "map")
set_key("my_other_api_key")
## Ghana roads
ghanaRoads <- sf::st_read("~/Downloads/gha_roads_dcw/GHA_rds_1m_dcw.shp")
## origin piont
df <- data.frame(lat = 7.9254948, lon = -0.6283887)
## destination point
dest <- data.frame(lat = 5.591622, lon = -0.187677)
google_map() %>%
add_markers(data = df) %>%
add_markers(data = dest) %>%
add_polylines(ghanaRoads)
您在此处使用的确切方法可能会有所不同。但在这个例子中,我使用两个坐标之间的一条线来合理猜测行进方向,以及起点应该在哪里。
## convert origin into an 'sf' object
sf_origin <- sf::st_sf(geometry = sf::st_sfc(sf::st_point(x = c(-0.6283887, 7.9254948))))
## create a line between the origin and destination
m <- matrix(c(-0.6283887, 7.9254948, -0.187677, 5.591622), ncol = 2, byrow = T)
sf_line <- sf::st_sf(geometry = sf::st_sfc(sf::st_linestring(x = m)))
## The coordinate reference system needs to match between the two for
## spatial operations
sf::st_crs(sf_line) <- sf::st_crs(ghanaRoads)
sf::st_crs(sf_origin) <- sf::st_crs(ghanaRoads)
## find all the intersecting points
sf_intersections <- sf::st_intersection(ghanaRoads, sf_line)
google_map() %>%
add_markers(data = df) %>%
add_markers(data = dest) %>%
add_polylines(data = ghanaRoads) %>%
add_markers(data = sf_intersections)
我们可以使用离原点最近的相交点作为路线查询的起点。
注意
找到最近的sf 方法是使用sf::st_distance,但这取决于您的系统上是否安装了lwgeom,但我在安装它时遇到问题,所以我不得不使用不同的方法
我正在使用我编写的 for this answer 的 data.table 函数来计算从每个点到原点的半正弦距离。然后我选择距离最小的那个。
library(data.table)
coords <- matrix(unlist(sf_intersections$geometry), ncol = 2, byrow = T)
## Taking a fucntion I wrote for this answer
## https://stackoverflow.com/a/42014364/5977215
dt.haversine <- function(lat_from, lon_from, lat_to, lon_to, r = 6378137){
radians <- pi/180
lat_to <- lat_to * radians
lat_from <- lat_from * radians
lon_to <- lon_to * radians
lon_from <- lon_from * radians
dLat <- (lat_to - lat_from)
dLon <- (lon_to - lon_from)
a <- (sin(dLat/2)^2) + (cos(lat_from) * cos(lat_to)) * (sin(dLon/2)^2)
return(2 * atan2(sqrt(a), sqrt(1 - a)) * r)
}
dt <- as.data.table(coords)
dt[, `:=`(origin_lon = df$lon, origin_lat = df$lat)]
dt[, distance := dt.haversine(origin_lat, origin_lon, V1, V2)]
## min distance
sf_nearest <- dt[order(-distance)][1, .(lon = V1, lat = V2)]
sf_nearest <- sf::st_point(c(sf_nearest$lon, sf_nearest$lat))
sf_nearest <- sf::st_sf(geometry = sf::st_sfc(sf_nearest))
sf_nearest$colour <- "green"
google_map() %>%
add_markers(data = df) %>%
add_markers(data = dest) %>%
add_markers(data = sf_nearest, colour = "colour")
我们可以在路线查询中使用这个绿色标记
orig <- sf_nearest$geometry[[1]]
orig <- as.numeric(orig)
df_orig <- data.frame(lat = orig[2], lon = orig[1])
google_map() %>%
add_markers(df_orig)
res <- google_directions(origin = df_orig,
destination = dest)
## all the api results are now stored in the 'res' object.
direction_legs(res)$distance
# text value
# 1 397 km 396829
## you can look at the route through the polygoin
df_route <- data.frame(polyline = direction_polyline(res))
google_map() %>%
add_markers(data = df_orig) %>%
add_markers(data = dest) %>%
add_polylines(data = df_route, polyline = "polyline")
dt[order(-distance)] 给出了从原点到新原点的距离,
dt[order(-distance)][1, distance]
# [1] 1329904
这是以米为单位。假设平均步行速度为 4kph,您可以将其添加到总时间中
替代方法
根据 cmets 的要求,另一种寻找最近道路的方法是在原点周围画一个缓冲区并找到任何相交的道路
sf_buffer <- sf::st_buffer(sf_origin, dist = 0.5)
sf::st_crs(sf_buffer) <- sf::st_crs(ghanaRoads)
google_map() %>%
add_polylines(ghanaRoads) %>%
add_polygons(sf_buffer)
然后您可以找到与此缓冲区相交的所有线
sf_intersection <- sf::st_intersection(sf_buffer, ghanaRoads)
google_map() %>%
add_markers(data = df) %>%
add_polylines(sf_intersection)
您可以在距离计算中使用这个新的sf_intersection 对象。