问题描述
我有多条多段线需要分割成最大 X 英尺,同时也不能短于 Y 英尺。请参阅下面的折线示例。
df<-new("SpatialLinesDataFrame",data = structure(list(Count = "87",length = 13443.0406016608,Lat_TC = 2.5),class = "data.frame",row.names = 1L),lines = list(new("Lines",Lines = list(new("Line",coords = structure(c(1279653.38073649,1275804.03323203,476893.543169454,489773.677656633),.Dim = c(2L,2L)))),ID = "1")),bBox = structure(c(1275804.03323203,1279653.38073649,2L),.Dimnames = list(c("x","y"),c("min","max"))),proj4string = new("CRS",projargs = "+proj=tmerc +lat_0=31 +lon_0=-104.333333333333 +k=0.999909091 +x_0=165000 +y_0=0 +datum=NAD83 +units=us-ft +no_defs"))
#Segment polyline function
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])
y <- coords[i,2])
coordsOut <- rbind(coordsOut,y))
break
}
coordsOut <- rbind(coordsOut,c(coords[i + 1,1],coords[i + 1,2]))
}
}
return(coordsOut)
}
CreateSegments <- function(coords,length = 10000,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,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,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)
}
MergeLast <- function(lst) {
l <- length(lst)
lst[[l - 1]] <- rbind(lst[[l - 1]],lst[[l]])
lst <- lst[1:(l - 1)]
return(lst)
}
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))
}
我想每 6000 英尺分割一条多段线,但是如果最后一条线段短于 6000,则与前一条线合并并将合并线分成两半。
split_poly=st_as_sf(SegmentSpatialLines(x,length=6000,merge.last = T))
这会创建两条线段,我想将第二条线段一分为二,因为它大于 6000 英尺。最终结果将是三个线段,长度分别为 6000、3721 和 3721。
解决方法
使用下面的循环解决
datalist= list()
for (i in 1:nrow(split_poly)) {
if(st_drop_geometry(split_poly[i,"length"])<6000) next
if(st_drop_geometry(split_poly[i,"length"])>6000){
split_poly2=st_as_sf(SegmentSpatialLines(as_Spatial(split_poly[i,]),n.parts=3,merge.last = F))
split_poly2=split_poly2%>%mutate(length=st_length(xx),id=i)
datalist[[i]] <- xx # add it to your list
}
}
final_df=do.call(rbind,datalist)