如何对面积变量的矢量进行栅格化？
我多次使用了同样的光栅化过程，效果相当好：
raster <- rasterize(vect(shapefile.shp), base_grid, "my_variable")其中光栅是光栅化的形状文件，shapefile.shp是原始矢量，base_grid是光栅骨架，"my“是要考虑的变量。对于与多边形面积无关的变量，这种方法是令人满意的，因为它使用平均计算来重新排列数据(例如:人口、生产、产量、温度、降水量)。
然而，现在我试图从矢量转换为栅格多边形，这些多边形具有可变的收获面积，这不是严格意义上的多边形的区域，而是可以被认为是与整个多边形面积成正比的区域。上述方法在考虑收获面积时会产生膨胀值(是对应向量的3-4倍)，这可能是因为多边形通常比网格单元大。因此，一个包含100个单元的大多边形被划分成10个网格单元，每个单元将给出100个单元，而我希望它们每个都有10个单元(因为它是一个区域)。
我认为我的方法是这样的：“在每个网格单元中，加权平均所有多边形的值与它们在网格单元中的存在成正比。”
但我要寻找的是：“对于每个网格单元中的每个多边形部分，计算网格单元内多边形的比例(wrt与总多边形面积)，并将网格单元内的所有值之和(因为它是一个区域单位)”。
任何帮助都是非常感谢的。
更新：
矢量数据的视图。光栅其实很多，因为我有好几年的时间：
Simple feature collection with 9382 features and 3 fields
Geometry type: MULTIPOLYGON
Dimension: XYZ
Bounding box: xmin: -67.38379 ymin: -41.03791 xmax: -53.63737 ymax: -21.99877
z_range: zmin: 0 zmax: 0
Geodetic CRS: WGS 84
First 10 features:
ADM2_REF anio my_variable geometry
2 Tres Arroyos 1978 180 MULTIPOLYGON Z (((-60.16947...
3 Tres Arroyos 1979 0 MULTIPOLYGON Z (((-60.16947...
4 Tres Arroyos 1988 1000 MULTIPOLYGON Z (((-60.16947...
5 Tres Arroyos 1989 1000 MULTIPOLYGON Z (((-60.16947...
6 Tres Arroyos 1990 3000 MULTIPOLYGON Z (((-60.16947...
7 Tres Arroyos 1991 1500 MULTIPOLYGON Z (((-60.16947...
8 Tres Arroyos 1992 2800 MULTIPOLYGON Z (((-60.16947...
9 Tres Arroyos 1993 2800 MULTIPOLYGON Z (((-60.16947...
10 Tres Arroyos 1994 2500 MULTIPOLYGON Z (((-60.16947...
11 Tres Arroyos 1995 1250 MULTIPOLYGON Z (((-60.16947...将上面的数据文件转换为光栅的步骤如下：
baserast <- rast(nrows=nrows, ncol=nrows,
extent= extent,
crs="+proj=longlat +datum=WGS84",
vals=NA)
rasters <- rast(lapply(1978:2019,
function(x)
rasterize(vect(shp.soy.yld.arg %>%
filter(anio==x)), baserast, "my variable")))链接一年的数据.gpkg (对所有年份来说都太大了)。
回答 2
Stack Overflow用户
发布于 2022-04-11 08:51:02
在这种情况下，你应该用栅格化密度而不是面积。
示例数据：
library(terra)
f <- system.file("ex/lux.shp", package="terra")
v <- vect(f)
# all areas have 100 ha of the crop
v$crop_area <- 100
r <- rast(v, res=.01) 计算密度
e <- expanse(v, unit="ha") # or "km" to avoid small numbers
v$density <- v$crop_area / e拉斯特尔
x <- rasterize(v, r, "density")返回地区
ra <- cellSize(r, unit="ha")
area <- x * ra检查数字是否合理(大面积/小单元的误差应最小)。每个多边形的期望值为100。
extract(area, v, sum, na.rm=TRUE, exact=TRUE) |> round()
# ID density
# [1,] 1 99
# [2,] 2 101
# [3,] 3 99
# [4,] 4 95
# [5,] 5 99
# [6,] 6 98
# [7,] 7 96
# [8,] 8 99
# [9,] 9 98
#[10,] 10 98
#[11,] 11 101
#[12,] 12 100下面，我将展示如何在多年的循环中完成此操作。
示例数据：
library(terra)
f <- system.file("ex/lux.shp", package="terra")
v <- vect(f)
# all areas have 100 ha of the crop
v$crop_area <- 100
v$year <- rep(c(1990,1991, 1992), each=4)
r <- rast(v, res=.01) 解决方案
ra <- cellSize(r, unit="ha")
e <- expanse(v, unit="ha")
v$density <- v$crop_area / e
years <- unique(v$year)
out <- list()
for (i in 1:length(years)) {
vv <- v[v$year == years[i], ]
x <- rasterize(vv, r, "density")
out[[i]] <- x * ra
}
out <- rast(out)
names(out) <- years
out
#class : SpatRaster
#dimensions : 73, 78, 3 (nrow, ncol, nlyr)
#resolution : 0.01, 0.01 (x, y)
#extent : 5.74414, 6.52414, 49.44781, 50.17781 (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 (EPSG:4326)
#sources : memory
# memory
# memory
#names : 1990, 1991, 1992
#min values : 0.2544571, 0.3028134, 0.3200223
#max values : 1.0492719, 0.6249076, 0.4335355 Stack Overflow用户
发布于 2022-04-11 06:17:48
如果我很好地理解你的问题(一个可重复的例子会很好)，你会希望所有的像素在栅格化多边形之和的收获值("my_variable“在您的代码)。
在这里，我创建了一个玩具示例，向您展示我的推理：
- 首先加载库
- 用一个示例创建玩具数据总数和收获面积
- 计算多边形覆盖的每个像素的分数。
- 将每个覆盖部分除以多边形的总面积，再乘以收获的面积。 图书馆(Sf)图书馆(栅格)图书馆(正式版) rast <-光栅:光栅(矩阵(rep(1,100)，ncol=10)，xmn=0，ymn=0，xmx=10，ymx=10) pol <- ymn=0 c(1，0，4，1) pol <- st_sf(data.frame(area = st_area(pol)，harvest=0.7，(geom=pol)光栅：：情节( cov_frac )光栅：：情节(波尔，波尔)cov_frac <- exactextractr::coverage_fraction(rast，pol)[1]光栅:情节(Cov_frac)光栅:情节(pol，add=T)结果<- cov_frac/st_area(Pol)*$收获费(值(结果))
正如您所看到的，在栅格化多边形中所有像素的和等于所获取的面积。
https://stackoverflow.com/questions/71825642
复制
腾讯云开发者
