如何测量 r 中 sf 多线串对象的两个特征之间插入的 sf 线串部分的长度?

How to measure length of a section sf linestring inserted between two features from sf multilinestring object in r?

我有一个线串格式的数据集和一个海洋测深多线串(几条带有深度信息的线)。我的目标是在每个深度范围(两个多线串之间)测量这些线串的长度部分,并获得每个间隔的总长度。

预期输出示例

       interval        length 
1  -2250 and -2500    5200.56 [m]
2  -2500 and -2750    xxxxxxx [m]
3  -2750 and -3000    xxxxxxx [m]
4  -3000 and -3250    xxxxxxx [m]

但是,当我使用 st_intersection() 时,我可以获得这段线的深度,但是除了关于深度的信息是例如 -1500 米而不是间隔,(例如: -1500 m 到 -1750 m) 我无法测量长度,因为 st_intersection() 只计算交点。有没有办法在 r 中做这个?

数据集示例

library(sf)
library(dplyr)


#creating dataset
id <- c("A","A", "B","B","C","C")
lat <- c(-25.31157, -25.42952, -25.4253, -25.19177, -25.18697, -25.12748)
long <- c(-41.39523, -39.99665, -41.00311, -41.29756, -41.30314, -39.37707)

df <- dplyr::tibble(id = as.factor(id), lat, long)

#convert sf
df.sf = sf::st_as_sf(df,
                      coords = c("long","lat"), 
                      crs    = 4326)

#creating linestrings
lines <- df.sf %>% 
  dplyr::group_by(id) %>%
  dplyr::summarise(do_union = FALSE) %>%
  sf::st_cast("LINESTRING")

######attempt at something similar:

#intersect
intersection <- sf::st_intersection (bathy, lines) %>% sf::st_as_sf() %>%
  dplyr::mutate(length=sf::st_length(.)) %>% 
  sf::st_drop_geometry()

#sum length per depth 
output <- intersection %>% 
  group_by(depth) %>% 
  summarise(length = sum(length)) %>%
  arrange(desc(length)) 

> output
   depth   length[m]
1  -2250     0 [m]
2  -2500     0 [m]
3  -2750     0 [m]
4  -3000     0 [m]

不幸的是,我无法在此处重新创建多线串对象的子集,因此我放入 github 以下载 shapefile 格式的子集(如有必要),只需单击“代码”,然后单击“下载”邮编”enter link description here。对象是这样的:

> bathy
Simple feature collection with 15690 features and 1 field
Geometry type: MULTILINESTRING
Dimension:     XY
Bounding box:  xmin: -53.37185 ymin: -35.77762 xmax: -25 ymax: 7.081626
Geodetic CRS:  WGS 84
First 10 features:
   PROFUNDIDA                       geometry
1         -50 MULTILINESTRING ((-52.4267 ...
2         -50 MULTILINESTRING ((-52.77632...
3         -75 MULTILINESTRING ((-51.04274...
4         -75 MULTILINESTRING ((-52.38656...
5        -100 MULTILINESTRING ((-51.07005...
6        -100 MULTILINESTRING ((-52.18633...
7        -200 MULTILINESTRING ((-51.97665...
8        -300 MULTILINESTRING ((-51.95862...
9        -400 MULTILINESTRING ((-51.94465...
10       -500 MULTILINESTRING ((-51.93161...

这里是一个粗略的答案,可能会有帮助。如果您的数据较小,则提供 reprex 会更容易。它有助于只选择绝对最小量的数据来重现您的问题。

在这种情况下,您只需要一个路径线串,以及一对测深线串。我建议对您的数据进行大量子集化,然后使用 datapasta 包或 reprex 制作一个可以轻松发布的示例。

假设 intersection 只是 1 LINESTRING 与深海数据的交集,这可能有效:

# st_intersects returns both POINT and MULTIPOINT geometries.
# It returns MULTIPOINT when the line x intersects the feature y 
# multiple times 
# We want just POINT geometries
# this is an unfortunately complex:
intersection %>% 
  group_by(geomtype = st_geometry_type(geometry)) %>%
  group_modify(~ st_cast(.x, "POINT")) %>% 
  ungroup() %>% 
  select(-geomtype) %>% 
  st_as_sf() ->
  intersection_points

# Once we have one POINT for each intesection, 
# along with the associated information 
# (Depth in this example) we can calculate the distance
# between adjacent points using dplyr::lead and sf::st_distance
intersection_points %>% 
  mutate(geometry2 = lead(geometry)) %>%
  mutate(length = st_distance(geometry, geometry2, by_element = TRUE)) ->
intersection_points

# Add Depth ranges:
intersection_points %>% 
   mutate(interval = paste(as.character(PROFUNDIDA),
                           "and",
                           as.character(lead(PROFUNDIDA))))