0

我想计算每个栅格单元内线串的总长度。我有一个繁琐的解决方法，但我认为有一种更有效的方法可以直接使用该rasterize函数执行此操作。下面的 MWE 显示了我想要什么以及我迄今为止的工作解决方案。

library(sf)
library(raster)

## Base Raster
base_rast <- raster( matrix(runif(100, max=10),10,10))

## LineStrings
ls <- st_linestring(rbind(c(0,0),c(1,1),c(2,1)))
mls <- st_multilinestring(list(rbind(c(2,2),c(1,3)), rbind(c(0,0),c(1,1),c(2,1))))
sfc <- st_sf( geometry=st_sfc(ls,mls) )
sfc$ID <- 1:length(sfc)

## Would Like To Do Something Like
new_fun <- function(j){ sum(sapply(j, st_length)) }
i_length <- rasterize(sfc, base_rast, fun=new_fun)

## Cumbersome Workaround
base_polygon <- rasterToPolygons(base_rast)
base_polygon <- st_as_sf(base_polygon)
base_polygon$baseID <- 1:nrow(base_polygon)
i_poly <- st_intersection(i_poly, sfc)
i_poly_list <- split(i_poly, i_poly$baseID)
i_length <- sapply(i_poly_list, function(j) sum(st_length(j)))
i_mat <- cbind(baseID=as.numeric(names(i_poly_list)), i_length)
lengths <- merge(base_polygon, i_mat, by='baseID', all=T)
i_df <- data.frame(
    st_coordinates(st_centroid(lengths))[,c('X','Y')],
    river_length=lengths$i_length)
coordinates(i_df) <- ~X+Y
gridded(i_df) <- TRUE
i_length <- raster(i_df)
4

0 回答 0