面积计算和像素数的差异

Question

我有一个被河流网络以二进制形式掩盖的栅格图层。在计算编号时。 R (freq(raster)) 和使用 r.report 的 QGIS 中的像素数我在两者中发现了相同的数字。但是在以平方公里为单位计算面积时，我发现由 (tapply(area(raster), raster[], sum)) 和 QGIS 计算的面积 R 存在差异。但是，我遇到的主要问题是为什么面积计算与像素数不符？栅格的分辨率为 30 秒（约 1km*1km），因此像素数必须约等于以平方公里为单位的面积。栅格具有地理坐标系 OGC:CRS84 - WGS 84 (CRS84) - Geographic，格式为 .grd。我还将它投影到 QGIS 的 UTM，这略微增加了面积，但差别不大。 我还在下面发布来自 R 和 QGIS 的报告，如果您也想查看栅格，请点击下面的链接。我想要区域中的值，所以我真的不知道是否应该将像素数转换为以平方公里为单位的区域。在这种情况下，它应该相等或使用 R 或 QGIS 的答案之一。

在 QGIS 区域的 r.report 中，以平方公里为单位与像素： 0 值：222、520 与 290,767 1 个值：81,653 与 106,934

在 Rstudio 区域（以平方公里为单位）与像素： 0 值：222,068.53 与 290,767 1 个值：81,484.18 与 106,934

https://drive.google.com/drive/folders/1pBba0ejIc4t9ayKl36nyIo3vbgD4BKCT?usp=sharing

Answer 1

使用 R raster 包，经度/纬度栅格单元的面积以 m2 为单位计算。另请注意，在赤道，一个 30 秒的单元大约是 860 平方米，并且这个面积随着朝向两极而减小。

给定单元格数量n，以及您的栅格lat的平均纬度，该区域应该是大约

n * 860 * cos(lat*pi/180)

在raster 中显示正确（内存安全）方法的示例

#example data
library(raster)
r <- raster()
set.seed(67)
values(r) <- 0
r[sample(ncell(r), 1000)] <- 1

#select the cells that are not 0
rr <- reclassify(r, cbind(0,NA))

a <- area(rr)
aa <- mask(a, rr)

cellStats(aa, "sum")

---- 使用您的数据

r <- raster("raster.grd")
rr <- reclassify(r, cbind(0, NA))
a <- area(r)
aa <- mask(a, rr)
cellStats(aa, "sum")
#81484.18

即约8.1万平方公里