Сегментация пространственных линий и сохранение атрибутов - PullRequest
0 голосов
/ 31 марта 2020

Мне нужно разбить наборы линий на сегменты и сохранить данные атрибутов исходных линий, а также пометить сегменты последовательными идентификаторами. Например, линия 1 в исходной коллекции разделена на 50-метровые сегменты и обозначена 1-1, 1-2, 1-3.

Мой подход состоит в том, чтобы использовать представленный здесь код, а затем присоединиться к исходным строкам, чтобы получить данные атрибутов.
http://rstudio-pubs-static.s3.amazonaws.com/10685_1f7266d60db7432486517a111c76ac8b.html

с использованием данные из исходного примера. Это то, что я сделал до сих пор

library(mvoutlier)
library(maptools)
data(kola.background)

#First, we need to define a function, which accepts two-column matrix of coordinates and lenghts from and to, which clip the segment. 
#Distances are always calculated from 0, therefore resulting segment will start at distance start from the beginning of the line and end at distance end from the beginning.
CreateSegment <- function(coords, from, to) {
  distance <- 0
  coordsOut <- c()
  biggerThanFrom <- F
  for (i in 1:(nrow(coords) - 1)) {
    d <- sqrt((coords[i, 1] - coords[i + 1, 1])^2 + (coords[i, 2] - coords[i + 
                                                                             1, 2])^2)
    distance <- distance + d
    if (!biggerThanFrom && (distance > from)) {
      w <- 1 - (distance - from)/d
      x <- coords[i, 1] + w * (coords[i + 1, 1] - coords[i, 1])
      y <- coords[i, 2] + w * (coords[i + 1, 2] - coords[i, 2])
      coordsOut <- rbind(coordsOut, c(x, y))
      biggerThanFrom <- T
    }
    if (biggerThanFrom) {
      if (distance > to) {
        w <- 1 - (distance - to)/d
        x <- coords[i, 1] + w * (coords[i + 1, 1] - coords[i, 1])
        y <- coords[i, 2] + w * (coords[i + 1, 2] - coords[i, 2])
        coordsOut <- rbind(coordsOut, c(x, y))
        break
      }
      coordsOut <- rbind(coordsOut, c(coords[i + 1, 1], coords[i + 1, 
                                                               2]))
    }
  }
  return(coordsOut)
}

#Now, lets say we want to cut the given line to several segments with given length. 
#Another option is to give number of desired parts (instead of their length) as an argument.
CreateSegments <- function(coords, length = 0, n.parts = 0) {
  stopifnot((length > 0 || n.parts > 0))
  # calculate total length line
  total_length <- 0
  for (i in 1:(nrow(coords) - 1)) {
    d <- sqrt((coords[i, 1] - coords[i + 1, 1])^2 + (coords[i, 2] - coords[i + 
                                                                             1, 2])^2)
    total_length <- total_length + d
  }

  # calculate stationing of segments
  if (length > 0) {
    stationing <- c(seq(from = 0, to = total_length, by = length), total_length)
  } else {
    stationing <- c(seq(from = 0, to = total_length, length.out = n.parts), 
                    total_length)
  }

  # calculate segments and store the in list
  newlines <- list()
  for (i in 1:(length(stationing) - 1)) {
    newlines[[i]] <- CreateSegment(coords, stationing[i], stationing[i + 
                                                                       1])
  }
  return(newlines)
}

#Since the actual length of line is rarely a multiple of given length, last segment is shorter.
#We can, however, very simply merge it with penultimate segment, if we wish.
MergeLast <- function(lst) {
  l <- length(lst)
  lst[[l - 1]] <- rbind(lst[[l - 1]], lst[[l]])
  lst <- lst[1:(l - 1)]
  return(lst)
}

#Attributes of lines, even in case they are given as SpatialLinesDataFrame, are not kept. 
#If desired length of segments is bigger than actual length of original line, the original line is returned.
SegmentSpatialLines <- function(sl, length = 0, n.parts = 0, merge.last = FALSE) {
  stopifnot((length > 0 || n.parts > 0))
  id <- 0
  newlines <- list()
  sl <- as(sl, "SpatialLines")
  for (lines in sl@lines) {
    for (line in lines@Lines) {
      crds <- line@coords
      # create segments
      segments <- CreateSegments(coords = crds, length, n.parts)
      if (merge.last && length(segments) > 1) {
        # in case there is only one segment, merging would result into error
        segments <- MergeLast(segments)
      }
      # transform segments to lineslist for SpatialLines object
      for (segment in segments) {
        newlines <- c(newlines, Lines(list(Line(unlist(segment))), ID = as.character(id)))
        id <- id + 1
      }
    }
  }
  return(SpatialLines(newlines))
}
#=========================================================================================
#Load the data (coast) from early example, plus data for this part of the example (coast2)
#========================================================
coast <- kola.background$coast[147:351, ]
coast2 <- kola.background$coast[353:416, ]
sl <- SpatialLines(list(Lines(list(Line(coords = coast)), ID = "1"), Lines(list(Line(coords = coast2)), ID = "2")))

#This is where my example diverges from orginal post. I need to turn sl into a spatiallinedataframe with routelength and Route_num
sldf<-data.frame(Distance.m = sapply(1:length(sl), function(i) gLength(sl[i, ])))
sldf$Route_num<-NA
sldf[1,]$Route_num<-"ONE"
sldf[2,]$Route_num<-"TWO"
sl.sldf<-SpatialLinesDataFrame(sl,data=sldf)
#========================================================
sl2 <- SegmentSpatialLines(sl, length = 1000, merge.last = TRUE)
df <- data.frame(Distance.m = sapply(1:length(sl2), function(i) gLength(sl2[i, ])))
rownames(df) <- sapply(1:length(sl2), function(i)sl2@lines[[i]]@ID)
#Create spatial lines dataframe.  This exports with the correct projection/datum.  Converting to SF didn't.  
sl2.sldf <- SpatialLinesDataFrame(sl2, data = df)
#========================================================
Join <- over(sl2.sldf, sl.sldf[,"Route_num"])

Однако, когда вы проверяете Join, вы можете видеть, что хотя sl2 (и впоследствии sl2.sldf) был создан из sl, там являются сегментами, которые не объединяются и, следовательно, Join $ Route_num равно "NA". Я считаю, что все сегменты должны правильно соединяться с исходной линией. # BANGHEADONWALL

Помогите, пожалуйста!

...