【问题标题】:How to split network in equal piece line segments in R如何在R中将网络分割成相等的线段
【发布时间】:2021-02-25 03:07:59
【问题描述】:

我从 x 和 y 坐标创建一个 SpatialLines。我想知道有没有办法在网络上选择 5、10 或 m 个彼此距离相等的点,至少 (2:m-1) 个点的距离相等。

我在想也许我可以计算每个点到前一个点的距离并获得累积长度并使用 seq(., ., length.out=m) 会给我相等的距离,但我不能得到 x 和 y坐标。

xy <- cbind(xco=c(172868.6, 172891, 172926.8, 172949.2, 172985, 173007.3, 173029.7, 173065.5, 173087.9, 173123.7, 173146, 173181.9, 173204.2, 173226.6, 173262.4, 173320.7, 173343, 173365.4, 173401.3, 173423.7, 173459.6, 173482, 173504.3, 173541.2, 173563.5, 173601.2, 173623.4, 173644.9, 173684.7, 173706.7, 173747.3, 173769.3, 173811, 173832.8, 173875.2, 173896.9, 173918.5, 173962, 173983.6, 174027.1, 174048.6, 174092.1, 174113.7, 174157.2, 174178.8, 174222.2, 174243.7, 174287.1, 174308.7, 174330.3, 174373.6, 174395.2, 174438.9, 174460.4, 174504.2, 174525.7, 174569.4, 174590.9, 174634.5, 174656, 174700.4, 174722, 174742.4, 174790.4, 174860.9, 174914.5, 174932.6, 174988.4, 175005.4, 175063.6, 175079.2, 175138.5, 175200.4, 175213.2, 175276, 175340.7, 175405.8, 175415.9, 175481, 175546.1, 175611.2, 175621.3, 175686.4, 175751.5, 175761.6, 175826.3, 175893.1, 175960.8, 176029.4, 176098.8, 176168.1, 176237.6, 176307.1, 176375.5, 176442.9, 176509.2, 176517.7, 176584.2, 176648.6, 176712.4),
yco=c(3376152, 3376130, 3376096, 3376074, 3376040, 3376019, 3375997, 3375963, 3375941, 3375907, 3375886, 3375851, 3375830, 3375808, 3375774, 3375718, 3375697, 3375676, 3375641, 3375620, 3375586, 3375564, 3375543, 3375509, 3375489, 3375454, 3375434, 3375415, 3375381, 3375362, 3375328, 3375310, 3375276, 3375259, 3375226, 3375209, 3375192, 3375159, 3375142, 3375109, 3375093, 3375060, 3375043, 3375010, 3374994, 3374960, 3374944, 3374911, 3374894, 3374878, 3374844, 3374828, 3374795, 3374778, 3374745, 3374729, 3374696, 3374680, 3374646, 3374630, 3374597, 3374581, 3374567, 3374536, 3374494, 3374464, 3374455, 3374428, 3374420, 3374395, 3374389, 3374366, 3374346, 3374341, 3374324, 3374307, 3374292, 3374289, 3374274, 3374258, 3374243, 3374240, 3374225, 3374209, 3374207, 3374192, 3374181, 3374172, 3374166, 3374162, 3374161, 3374163, 3374167, 3374175, 3374184, 3374197, 3374198, 3374213, 3374230, 3374248))

library(sp)
myLine190 <- Line(xy)
myL190 <- Lines(list(myLine190), ID = 1)
Spt1 <- SpatialLines(list(myL190))
Y1 <- coordinates (Spt1)  # the same as xy above

累计距离

DD1 = sqrt(rowSums(do.call("rbind", lapply(1:(dim(Y1)[1]-1), function(i){(Y1[i,])-(Y1[i+1,])}))^2))
cd1 = cumsum(c(0,DD1))
grid = seq(min(cd1), max(cd1), length.out = 20)

【问题讨论】:

  • 如果您提供Y1 作为矩阵,对您的帮助会容易得多。例如。 Y1 &lt;- cbind(x=c( ), y = c( ))
  • 我已将它们添加为Y1$xY1$y
  • 您也可以在第一个Y1 上使用此代码Y1&lt;- read.table(text = 'copy-paste data from the post', header = TRUE) 它会起作用
  • 谢谢,但这仍然不起作用。我应该能够将您的代码复制并粘贴到 R 中并运行它。当然,如果我真的想的话,我可以想办法做到这一点,但我不是那个寻求帮助的人。
  • Y134&lt;- read.table(text = ' xco yco 1 172868.6 3376152 2 172891.0 3376130 3 172926.8 3376096 4 172949.2 3376074 5 172985.0 3376040 6 173007.3 3376019 7 173029.7 3375997 8 173065.5 3375963 9 173087.9 3375941 10 173123.7 3375907 11 173146.0 3375886', header = TRUE) 当我以这种方式复制时,它对我有用。我可以通过电子邮件将代码发送给您吗?或问题。是的,我确实需要帮助。感谢您的宝贵时间。

标签: r raster spatial sf


【解决方案1】:

示例数据

pts <- cbind(xco=c(172868.6, 172891, 172926.8, 172949.2, 172985, 173007.3, 173029.7, 173065.5, 173087.9, 173123.7, 173146, 173181.9, 173204.2, 173226.6, 173262.4, 173320.7, 173343, 173365.4, 173401.3, 173423.7, 173459.6, 173482, 173504.3, 173541.2, 173563.5, 173601.2, 173623.4, 173644.9, 173684.7, 173706.7, 173747.3, 173769.3, 173811, 173832.8, 173875.2, 173896.9, 173918.5, 173962, 173983.6, 174027.1, 174048.6, 174092.1, 174113.7, 174157.2, 174178.8, 174222.2, 174243.7, 174287.1, 174308.7, 174330.3, 174373.6, 174395.2, 174438.9, 174460.4, 174504.2, 174525.7, 174569.4, 174590.9, 174634.5, 174656, 174700.4, 174722, 174742.4, 174790.4, 174860.9, 174914.5, 174932.6, 174988.4, 175005.4, 175063.6, 175079.2, 175138.5, 175200.4, 175213.2, 175276, 175340.7, 175405.8, 175415.9, 175481, 175546.1, 175611.2, 175621.3, 175686.4, 175751.5, 175761.6, 175826.3, 175893.1, 175960.8, 176029.4, 176098.8, 176168.1, 176237.6, 176307.1, 176375.5, 176442.9, 176509.2, 176517.7, 176584.2, 176648.6, 176712.4),
yco=c(3376152, 3376130, 3376096, 3376074, 3376040, 3376019, 3375997, 3375963, 3375941, 3375907, 3375886, 3375851, 3375830, 3375808, 3375774, 3375718, 3375697, 3375676, 3375641, 3375620, 3375586, 3375564, 3375543, 3375509, 3375489, 3375454, 3375434, 3375415, 3375381, 3375362, 3375328, 3375310, 3375276, 3375259, 3375226, 3375209, 3375192, 3375159, 3375142, 3375109, 3375093, 3375060, 3375043, 3375010, 3374994, 3374960, 3374944, 3374911, 3374894, 3374878, 3374844, 3374828, 3374795, 3374778, 3374745, 3374729, 3374696, 3374680, 3374646, 3374630, 3374597, 3374581, 3374567, 3374536, 3374494, 3374464, 3374455, 3374428, 3374420, 3374395, 3374389, 3374366, 3374346, 3374341, 3374324, 3374307, 3374292, 3374289, 3374274, 3374258, 3374243, 3374240, 3374225, 3374209, 3374207, 3374192, 3374181, 3374172, 3374166, 3374162, 3374161, 3374163, 3374167, 3374175, 3374184, 3374197, 3374198, 3374213, 3374230, 3374248))

这是一个可以改进的解决方案

# interval
x <- 500

d <- raster::pointDistance(pts[-nrow(pts),], pts[-1,], lonlat=FALSE)
cd <- cumsum(d)
n <- trunc(sum(d) / x)
s <- (1:n) * x

xy <- matrix(ncol=2, nrow=length(s))  
for (i in 1:length(s)) {
    j <- which(cd >= s[i])[1]
    pd <- raster::pointDistance(pts[j,], pts[j+1,], lonlat=FALSE)
    r <- (s[i] - cd[j-1]) / pd
    xy[i,] <- c(((1 - r) * pts[j,1] + r * pts[j+1,1]), 
                ((1 - r) * pts[j,2] + r * pts[j+1,2]))
}

看看

plot(pts, type="l")
points(xy)

这可以通过在循环的最后一行中使用 geosphere::destPoint 来针对 lon/lat 数据进行调整。

spsample("regular")findInterval 可以让你接近。

sp <- raster::spLines(xy)
s <- sp::spsample(sp, 5, "regular")
plot(sp)
points(s)

但这会选择节点,而不是节点之间的点。

【讨论】:

  • sp &lt;- spLines(xy) 此代码对我不起作用。它来自哪个图书馆?
  • 对此感到抱歉。我添加了raster:: 它来自哪里
【解决方案2】:

这是一个可能的解决方案。

df 坐标是你的 Y1 df。

dist.matrix <- as.matrix(dist(coordinates, method = "euclidean"))

创建一个距离矩阵并将其拆分为对角线。详情请见here

d <- row(dist.matrix) - col(dist.matrix)
diagonals <- split(dist.matrix, d)

这是对角线下方的第一个对角线

dist.points <- diagonals$`1`

我们的 cum sum 向量。

dist.cum <- cumsum(0, dist.points)

我们的距离段向量。

dist.segments <- seq(min(dist.cum), max(dist.cum), length.out = 20)

这个布尔向量列表告诉您可以到达 dist.points 对角线中的哪个位置,使得 它在一个段内。因此,您也知道您在原始数据框中的哪一点。因为 entry (2, 1) 告诉您从第二点到第一点的距离。也就是说,这是关键。

segments <- lapply(dist.segments,
                   FUN = function(dist.segment)
                   {
                     return(dist.cum <= dist.segment)
                   })


从那里你只需要编写一个返回点和下一个点的函数(通过使用布尔向量的segments列表)并计算你需要获得精确cumsum的额外距离。并且因为知道下一个点,所以可以轻松计算 (x, y) 坐标。

HTH!

【讨论】:

  • 这几乎就是我所做的。好吧,我无法从距离到 x 和 y 坐标,这就是我发布问题的第一个地方。
猜你喜欢
  • 2023-03-23
  • 1970-01-01
  • 2016-11-20
  • 1970-01-01
  • 1970-01-01
  • 2017-05-25
  • 1970-01-01
  • 1970-01-01
  • 2020-11-22
相关资源
最近更新 更多