从道路形状文件计算道路密度栅格

Question

我希望使用 R 将美国东半部的道路形状文件（包括每条道路的长度列）转换为 1x1 公里道路密度的栅格。

我在 Arcmap 中找不到直接的方法（线密度与单元中心的半径一起工作，而不仅仅是单元格）。

Answer 1

这是一个从栅格单元创建多边形的解决方案（改编自我的回答 here）。您可能需要对数据集的子集进行此操作，然后合并。

示例数据

library(terra)
v <- vect(system.file("ex/lux.shp", package="terra"))
roads <- as.lines(v)
rs <- rast(v)

解决方案

values(rs) <- 1:ncell(rs)
names(rs) <- "rast"    
rsp <- as.polygons(rs)

rp <- intersect(roads, rsp)

rp$length <- perim(rp) / 1000 #km
x <- tapply(rp$length, rp$rast, sum)

r <- rast(rs)
r[as.integer(names(x))] <- as.vector(x)

plot(r)
lines(roads)