使用一个shapefile (具有多个多边形)批量处理/提取一个栅格的原始数据？
如果我们有一个栅格，比方说一个国家的整数高程数据，和一个多边形形状文件，比方说在那个国家离散了300个河流流域，每个流域都有一个唯一的名称，我们如何才能最容易地为它们得到像这样的输出呢？
basinID, gridcellelev
a, 320
a, 321
a, 320
b, 17
b, 18
b, 19最繁琐的方法似乎是过滤/转换单个shapefile为300个shapefile，将栅格裁剪300次为300个uniqueID栅格，读回它们，为每个流域生成单独的表格，然后将它们全部rbind()在一起。
另一方面，理想的方式似乎是跳过文件生成，而不是保存xy数据，只使用一个栅格和一个shapefile创建相同的表-通过某种方式选择流域中的单元格，用basinID标记它们，创建一个表，丢失坐标，并不断迭代和追加该表，直到第300个流域。
我没有寻找任何统计数据，只是网格单元高程的原始数据列表，这将是一些标准裁剪方法的一部分。我相信ArcMap输出的标准栅格裁剪属性表不是原始数据，而是单元格的计数/频率。这也是可行的。
我不知道最低限度地再现一个栅格和一个多边形shapefile，所以如果有任何提示/库/函数/示例，我会非常感激。作为起点：
library(tidyverse)
library(raster)
library(rgdal)
library(sf)
elev_raster <- raster("spain_elev_meters.tif") #integer raster
basins <- readOGR("spainbasins.shp", "spainbasins") %>% st_as_sf() #unique basin ID column: `basinID` 除非建议使用其他起始格式。
任何建议都非常感谢！
回答 1
Stack Overflow用户
发布于 2019-07-19 22:31:07
如果你搜索“提取栅格shapefile”，我想你会找到答案的。
示例数据：
library(raster)
p <- shapefile(system.file("external/lux.shp", package="raster"))
r <- raster(p)
values(r) <- 1:ncell(r)p有12个多边形
p
#class : SpatialPolygonsDataFrame
#features : 12
#extent : 5.74414, 6.528252, 49.44781, 50.18162 (xmin, xmax, ymin, ymax)
#crs : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
#variables : 5
#names : ID_1, NAME_1, ID_2, NAME_2, AREA
#min values : 1, Diekirch, 1, Capellen, 76
#max values : 3, Luxembourg, 12, Wiltz, 312解决方案1.提取值并应用函数。例如mean
x <- extract(r, p, fun=mean, na.rm=TRUE)或者像这样
x <- extract(r, p)
v <- sapply(x, mean)每个多边形一个值
v
#[1] 14.00000 43.40000 49.00000 36.00000 29.50000 59.00000 91.00000 71.83333
#[9] 73.50000 87.16667 78.50000 59.57143解决方案2.获得您所要求的结构
x <- extract(r, p)
z <- do.call(rbind, lapply(1:length(x), function(i) cbind(i, x[[i]])))或者像这样
i <- sapply(x, length)
j <- rep(1:length(i), i)
z <- cbind(j, unlist(x))
colnames(z) = c("ID", "value")
head(z)
# ID value
#[1,] 1 4
#[2,] 1 5
#[3,] 1 13
#[4,] 1 14
#[5,] 1 15
#[6,] 1 22
tail(z)
# ID value
#[51,] 12 55
#[52,] 12 56
#[53,] 12 57
#[54,] 12 64
#[55,] 12 65
#[56,] 12 66稍后(参见讨论中的问题)
要获取频率，可以使用提取的值列表和table函数。
使用新的示例值来获得频率的一些变化
values(r) <- rep(1:4, 25)现在就做吧
f <- extract(r, p, fun=function(i,...) table(i)) 或者分两步执行：
x <- extract(r, p)
ff <- lapply(x, table)每个图层都有一个表
ff[1:2]
#[[1]]
#1 2 3 4
#3 2 1 1
#[[2]]
#1 2 3 4
#1 1 2 1 获得此结果的另一种方法是使用我上面所做的
x <- extract(r, p)
z <- do.call(rbind, lapply(1:length(x), function(i) cbind(i, x[[i]])))紧接着是
fff <- tapply(z[,2], z[,1], table)如果您想要进一步操作这些函数，请为此编写一个函数，并将其与lapply或for循环一起使用。
https://stackoverflow.com/questions/57107258
复制相似问题
腾讯云开发者
Copyright © 2013 - 2026 Tencent Cloud. All Rights Reserved. 腾讯云 版权所有
深圳市腾讯计算机系统有限公司 ICP备案/许可证号：粤B2-20090059 
粤公网安备44030502008569号
腾讯云计算（北京）有限责任公司 京ICP证150476号 | 京ICP备11018762号