如何使用 sf 和 R 将线串分割成线串块？

Question

问题已得到解答，但我仍然无法应用它。 我想使用 sf 包将线串分割成给定数量的块。

地图上显示为蓝色的道路是一个名为“line”的线串，我从中采样了5个点

line <- st_cast(line, "LINESTRING")
line <- st_transform(line, 3857)
pts <- st_line_sample(line,n = 5)

为了绘制它们，我将它们转换回纬度/经度，然后我将这些点捕捉到线上，以确保它们确实在线上：

pts <- st_cast(pts, "POINT")
pts <- st_snap(pts, line, tol=1e-9)
pts <- st_transform(pts,4326)
line <- st_transform(line,4326)

我确保它们看起来不错

leaflet() %>% addTiles() %>% addPolylines(data=line) %>% addCircleMarkers(data=pts)

然后我把线分开

seg <- st_split(x = line, y = pts)
seg2 <- st_collection_extract(seg, "LINESTRING")
print(length(seg))
print(length(seg2))

无论我做什么，长度都保持为 1。包网页上给出的示例有效

https://r-spatial.github.io/lwgeom/reference/st_split.html

但在这种情况下它失败了。

这是这条线的坐标：

> st_coordinates(line)
           X        Y L1
  [1,] -73.47588 45.57844  1
  [2,] -73.47538 45.57834  1
  [3,] -73.47487 45.57824  1
  [4,] -73.47451 45.57817  1
  [5,] -73.47374 45.57799  1
  [6,] -73.47288 45.57779  1
  [7,] -73.47240 45.57767  1
  [8,] -73.47209 45.57760  1
  [9,] -73.47176 45.57752  1
 [10,] -73.47140 45.57744  1
 [11,] -73.47099 45.57734  1
 [12,] -73.47056 45.57725  1
 [13,] -73.47022 45.57716  1
 [14,] -73.46961 45.57703  1
 [15,] -73.46928 45.57696  1
 [16,] -73.46897 45.57689  1
 [17,] -73.46873 45.57684  1
 [18,] -73.46857 45.57681  1
 [19,] -73.46840 45.57678  1
 [20,] -73.46819 45.57674  1
 [21,] -73.46788 45.57668  1
 [22,] -73.46748 45.57660  1
 [23,] -73.46716 45.57654  1
 [24,] -73.46678 45.57647  1
 [25,] -73.46651 45.57642  1
 [26,] -73.46621 45.57636  1
 [27,] -73.46587 45.57630  1
 [28,] -73.46556 45.57623  1
 [29,] -73.46527 45.57617  1
 [30,] -73.46498 45.57611  1
 [31,] -73.46424 45.57598  1
 [32,] -73.46377 45.57588  1
 [33,] -73.46233 45.57557  1
 [34,] -73.46164 45.57541  1
 [35,] -73.46122 45.57531  1
 [36,] -73.45898 45.57478  1
 [37,] -73.45831 45.57463  1
 [38,] -73.45798 45.57456  1
 [39,] -73.45766 45.57449  1
 [40,] -73.45705 45.57436  1
 [41,] -73.45498 45.57391  1
 [42,] -73.45440 45.57379  1
 [43,] -73.45070 45.57300  1
 [44,] -73.44766 45.57235  1
 [45,] -73.44441 45.57166  1
 [46,] -73.44156 45.57105  1
 [47,] -73.43243 45.56910  1
 [48,] -73.43218 45.56905  1
 [49,] -73.43198 45.56901  1
 [50,] -73.43175 45.56896  1
 [51,] -73.43155 45.56892  1
 [52,] -73.43134 45.56889  1
 [53,] -73.43113 45.56885  1
 [54,] -73.43095 45.56882  1
 [55,] -73.43078 45.56880  1
 [56,] -73.43053 45.56877  1
 [57,] -73.43032 45.56874  1
 [58,] -73.43013 45.56871  1
 [59,] -73.42994 45.56869  1
 [60,] -73.42973 45.56867  1
 [61,] -73.42956 45.56866  1
 [62,] -73.42937 45.56864  1
 [63,] -73.42920 45.56863  1
 [64,] -73.42901 45.56862  1
 [65,] -73.42883 45.56861  1
 [66,] -73.42861 45.56860  1
 [67,] -73.42839 45.56859  1
 [68,] -73.42815 45.56858  1
 [69,] -73.42789 45.56858  1
 [70,] -73.42764 45.56858  1
 [71,] -73.42736 45.56858  1
 [72,] -73.42710 45.56859  1
 [73,] -73.42685 45.56860  1
 [74,] -73.42659 45.56861  1
 [75,] -73.42634 45.56862  1
 [76,] -73.42605 45.56864  1
 [77,] -73.42577 45.56865  1
 [78,] -73.42225 45.56885  1
 [79,] -73.42144 45.56890  1
 [80,] -73.41945 45.56900  1
 [81,] -73.41828 45.56907  1
 [82,] -73.41735 45.56912  1
 [83,] -73.41703 45.56914  1
 [84,] -73.41661 45.56916  1
 [85,] -73.41598 45.56920  1
 [86,] -73.41535 45.56923  1
 [87,] -73.41287 45.56937  1
 [88,] -73.41079 45.56948  1
 [89,] -73.41015 45.56952  1
 [90,] -73.40708 45.56969  1
 [91,] -73.40558 45.56977  1
 [92,] -73.40510 45.56980  1
 [93,] -73.40366 45.56988  1
 [94,] -73.40205 45.56996  1
 [95,] -73.40055 45.57005  1
 [96,] -73.39993 45.57008  1
 [97,] -73.39943 45.57011  1
 [98,] -73.39913 45.57013  1
 [99,] -73.39886 45.57015  1
[100,] -73.39862 45.57017  1
[101,] -73.39828 45.57020  1
[102,] -73.39802 45.57022  1
[103,] -73.39779 45.57025  1
[104,] -73.39755 45.57027  1
[105,] -73.39734 45.57030  1
[106,] -73.39711 45.57032  1
[107,] -73.39681 45.57036  1
[108,] -73.39648 45.57040  1
[109,] -73.39609 45.57046  1
[110,] -73.38743 45.57163  1
[111,] -73.38551 45.57189  1
[112,] -73.38328 45.57219  1
[113,] -73.38275 45.57228  1
[114,] -73.38220 45.57236  1
[115,] -73.38090 45.57253  1
[116,] -73.37786 45.57295  1
[117,] -73.37743 45.57301  1
[118,] -73.37270 45.57364  1
[119,] -73.37035 45.57396  1
[120,] -73.36803 45.57428  1
[121,] -73.36627 45.57452  1
[122,] -73.36414 45.57480  1
[123,] -73.36278 45.57499  1
[124,] -73.36221 45.57505  1
[125,] -73.36151 45.57514  1
[126,] -73.36100 45.57521  1
[127,] -73.35975 45.57538  1
[128,] -73.35831 45.57558  1
[129,] -73.35577 45.57592  1
[130,] -73.35447 45.57609  1
[131,] -73.35367 45.57620  1
[132,] -73.35283 45.57632  1
[133,] -73.35195 45.57644  1
[134,] -73.35133 45.57652  1
[135,] -73.35059 45.57662  1
[136,] -73.34979 45.57673  1
[137,] -73.34637 45.57719  1
[138,] -73.34551 45.57730  1
[139,] -73.34409 45.57750  1
[140,] -73.34146 45.57785  1
[141,] -73.33933 45.57814  1
[142,] -73.33868 45.57823  1

[143,] -73.33827 45.57829  1
[144,] -73.33789 45.57834  1
[145,] -73.33755 45.57839  1
[146,] -73.33727 45.57843  1
[147,] -73.33692 45.57848  1
[148,] -73.33661 45.57853  1
[149,] -73.33628 45.57858  1
[150,] -73.33584 45.57865  1
[151,] -73.33530 45.57874  1

Answer 1

我也无法使用积分来让

st_split

工作。但是，如果将样本点转换为无限小的线段，效果会很好：

library(sf)
library(lwgeom)

line <- st_linestring(coords) |>
  st_sfc(crs = "WGS84")

blades <- line |>
  st_sample(5) |>
  st_cast("POINT") |>
  st_coordinates() |>
  asplit(1) |>
  lapply(function(x) rbind(x + c(0, -1e-6), x + c(0, 1e-6))) |>
  st_multilinestring() |>
  st_sfc(crs = "WGS84")

segments <- st_split(line, blades) |> 
  st_collection_extract("LINESTRING")

我们可以看到

segments

包含 6 个线串，正如预期的那样：

segments
#> Geometry set for 6 features 
#> Geometry type: LINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: -73.47588 ymin: 45.56858 xmax: -73.3353 ymax: 45.57874
#> Geodetic CRS:  WGS 84
#> First 5 geometries:
#> LINESTRING (-73.47588 45.57844, -73.47538 45.57...
#> LINESTRING (-73.44202 45.57115, -73.44156 45.57...
#> LINESTRING (-73.41764 45.5691, -73.41735 45.569...
#> LINESTRING (-73.40452 45.56983, -73.40366 45.56...
#> LINESTRING (-73.40304 45.56991, -73.40205 45.56...

它们沿着原来的线条完美契合

plot(line)
for(i in seq(length(segments))) plot(segments[i], col = i, add = TRUE, lwd = 3)

^{创建于 2024 年 10 月 30 日，使用 reprex v2.1.0}