测量沿 R 中一条线的所有点的距离（sf 中的线串）答案

【问题标题】：Measuring the distance of all points along a line in R (linestring in sf)测量沿 R 中一条线的所有点的距离（sf 中的线串）
【发布时间】：2019-09-17 11:51:09
【问题描述】：

我有一条平滑线（海岸线的简化抽象），它是 linestring，我想沿着它以频繁的间隔测量线的长度。我可以创建平滑线，并测量它的长度：

library(raster)
library(sf)
library(tidyverse)
library(rnaturalearth)
library(smoothr)

# create bounding box for the line 
xmin=-80
xmax=-66
ymin=24
ymax=45
bbox <- extent(xmin, xmax, ymin, ymax)

# get coarse outline of the US 
usamap <- rnaturalearth::ne_countries(scale = "small", country = "united states of america", returnclass = "sf")[1] %>% 
  st_cast("MULTILINESTRING")

# crop to extent of the bbox and get rid of a border line that isn't coastal 
bbox2 <- st_set_crs(st_as_sf(as(raster::extent(-80, -74, 42, 45.5), "SpatialPolygons")), st_crs(usamap))
cropmap <- usamap %>% 
  st_crop(bbox) %>% 
  st_difference(bbox2)

# smooth the line 
smoothmap <- cropmap %>% 
  smoothr::smooth(method="ksmooth", smoothness=8)

# measure the line length
st_length(smoothmap) # I get 1855956m

我将把样本站点“捕捉”到沿这条线的点，我需要知道它们沿海岸线有多远。但是，我不知道如何沿着海岸线每隔一段时间测量线的长度。间隔大小并不是非常重要，只要它具有相对精细的分辨率，可能每 1 公里或每 0.01 度。

我想要生成的是一个数据框，其中 x,y 列包含沿线的点的纬度/经度，以及包含沿线距离（从原点到该点）的“长度”列。以下是我尝试过的一些事情：

迭代边界框。我尝试用越来越小的边界框裁剪线（使用纬度/经度的规则间隔），但是因为线在 -76、38 附近“向后”弯曲，一些裁剪框不包含我预期的完整线段。这种方法适用于该行的右上半部分，但不适用于左下半部分——它只返回零长度。

在实现平滑之前先裁剪范围，然后测量线条。 由于平滑功能如果只是测量原始线条的一部分，则不会产生相同的形状，因此这不会实际上沿着同一条线测量距离。

用st_coordinates获取linestring的坐标，剪掉一行（线上的一个点），并将剩余的坐标重铸为linestring。这种方法不产生一个linestring，而是一个点链（因为st_cast不知道如何再次连接它们），因此无法正常测量。

理想的做法是“编辑”smoothmap 的几何图形，一次删除一行坐标，重复测量直线，然后将终点坐标和直线长度写入数据帧。但是，我不确定是否可以编辑 sf 对象的坐标而不将其转换为数据框。

【问题讨论】：

标签： r spatial sf

【解决方案1】：

如果我理解你的问题，我认为你可以这样做：

x <- as(smoothmap, "Spatial")
g <- geom(x)
d <- pointDistance(g[-nrow(g), c("x", "y")], g[-1, c("x", "y")], lonlat=TRUE)
gg <- data.frame(g[, c('x','y')], seglength=c(d, 0))

gg$lengthfromhere <- rev(cumsum(rev(gg[,"seglength"])))

head(gg)

#          x        y seglength lengthfromhere
#1 -67.06494 45.00000 70850.765        1855956
#2 -67.74832 44.58805  2405.180        1785105
#3 -67.77221 44.57474  2490.175        1782700
#4 -67.79692 44.56095  2577.254        1780210
#5 -67.82248 44.54667  2666.340        1777633
#6 -67.84890 44.53189  2757.336        1774967

tail(gg)
#            x        y  seglength lengthfromhere
#539 -79.09383 33.34224   2580.481       111543.5
#540 -79.11531 33.32753   2542.557       108963.0
#541 -79.13648 33.31306   2512.564       106420.5
#542 -79.15739 33.29874   2479.949       103907.9
#543 -79.17802 33.28460 101427.939       101427.9
#544 -80.00000 32.68751      0.000            0.0

【讨论】：

确实可以运行，但不会产生准确的长度值。长度估计值沿着gg 的行上下波动（当它们应该随着行变短而稳步下降时），并且长度值远小于我从st_length() 得到的数量级。
它有节点之间的间隔，正如您要求的“以频繁的间隔测量线的长度”，但我现在添加了另一个变量“lengthfromhere”，将它们沿线求和。

【解决方案2】：

相信你需要sf::st_line_sample:

# Transform to metric
smoothmap_utm <- st_transform(smoothmap, 3857)

# Get samples at every kilometer
smoothmap_samples <- st_line_sample(smoothmap_utm, density = 1/1000)

# Transform back to a sf data.frame
smoothmaps_points <- map(smoothmap_samples, function(x) data.frame(geometry = st_geometry(x))) %>%
  map_df(as.data.frame) %>% st_sf() %>%
  st_cast("POINT") %>%
  st_set_crs(3857) %>%
  st_transform(4326)

mapview(smoothmaps_points) + mapview(smoothmap)

得到你想要的输出：

# Function to transform sf to lon,lat

sfc_as_cols <- function(x, names = c("lon","lat")) {
  stopifnot(inherits(x,"sf") && inherits(sf::st_geometry(x),"sfc_POINT"))
  ret <- sf::st_coordinates(x)
  ret <- tibble::as_tibble(ret)
  stopifnot(length(names) == ncol(ret))
  x <- x[ , !names(x) %in% names]
  ret <- setNames(ret,names)
  ui <- dplyr::bind_cols(x,ret)
  st_set_geometry(ui, NULL)
}

smoothmaps_points_xy <- sfc_as_cols(smoothmaps_points) %>%
  mutate(dist = cumsum(c(0, rep(1000, times = n() - 1))))

smoothmaps_points_xy

        lon      lat dist
1 -67.06836 44.99794    0
2 -67.07521 44.99383 1000
3 -67.08206 44.98972 2000
4 -67.08891 44.98560 3000
5 -67.09575 44.98149 4000
6 -67.10260 44.97737 5000

重要

但如果您的最终目标是获取路径中的点距离，我建议您查看rgeos::gProject。

【讨论】：