每个网格单元区域的 r 栅格统计
r raster statistics per grid cell zone
我想计算格线内区域中的网格单元统计信息(见图表)。
获取显示的网格单元格的 cellstats 并不困难,但是如何为所有网格单元格完成此操作?
library(raster)
filename <- system.file("external/test.grd", package="raster")
r=raster(filename)
plot(r)
e=extent(180000,181000,330000,331000)
plot(e,add=T)
grid()
x=crop(r,e)
cellStats(x,mean)
你的可重现示例对我不起作用,因为我不确定为什么使用 r[r>500]=1
和 r[r<=500]=NA
而你的原始情节不同?!另外,我不确定 rx
是什么?!一种方法如下:
library(raster)
filename <- system.file("external/test.grd", package="raster")
r=raster(filename)
# r[r>500]=1
# r[r<=500]=NA
plot(r)
xmin <- seq(178000, 181000, 1000)
xmax <- seq(179000, 182000, 1000)
ymin <- seq(329000, 333000, 1000)
ymax <- seq(330000, 334000, 1000)
zonal.mtx <- matrix(nrow=length(xmin), ncol=length(ymin))
library(spatialEco)
for (i in 1:length(xmin)){
for (j in 1:length(ymin)){
e=extent(xmin[i],xmax[i],ymin[j],ymax[j])
plot(e,add=T)
grid()
p <- as(e, 'SpatialPolygons')
crs(p) <- crs(r)
p$ID=j
p <- SpatialPolygonsDataFrame(p,p@data)
zonal.mtx[i,j] <- zonal.stats(x=p, y=r, stat=sum, trace=TRUE, plot=TRUE)
#print(zonal.mtx[i,j])
}
}
注意: 这会在计算区域统计数据时绘制数据的每个部分。我相信您已经知道您还可以定义自己的函数以在 zonal.stats
而不是 sum
、mean
等
中使用
我想计算格线内区域中的网格单元统计信息(见图表)。
获取显示的网格单元格的 cellstats 并不困难,但是如何为所有网格单元格完成此操作?
library(raster)
filename <- system.file("external/test.grd", package="raster")
r=raster(filename)
plot(r)
e=extent(180000,181000,330000,331000)
plot(e,add=T)
grid()
x=crop(r,e)
cellStats(x,mean)
你的可重现示例对我不起作用,因为我不确定为什么使用 r[r>500]=1
和 r[r<=500]=NA
而你的原始情节不同?!另外,我不确定 rx
是什么?!一种方法如下:
library(raster)
filename <- system.file("external/test.grd", package="raster")
r=raster(filename)
# r[r>500]=1
# r[r<=500]=NA
plot(r)
xmin <- seq(178000, 181000, 1000)
xmax <- seq(179000, 182000, 1000)
ymin <- seq(329000, 333000, 1000)
ymax <- seq(330000, 334000, 1000)
zonal.mtx <- matrix(nrow=length(xmin), ncol=length(ymin))
library(spatialEco)
for (i in 1:length(xmin)){
for (j in 1:length(ymin)){
e=extent(xmin[i],xmax[i],ymin[j],ymax[j])
plot(e,add=T)
grid()
p <- as(e, 'SpatialPolygons')
crs(p) <- crs(r)
p$ID=j
p <- SpatialPolygonsDataFrame(p,p@data)
zonal.mtx[i,j] <- zonal.stats(x=p, y=r, stat=sum, trace=TRUE, plot=TRUE)
#print(zonal.mtx[i,j])
}
}
注意: 这会在计算区域统计数据时绘制数据的每个部分。我相信您已经知道您还可以定义自己的函数以在 zonal.stats
而不是 sum
、mean
等