python - 将 GeoTIFF 图像划分为规则网格并计算平均值

Question

如何将 GeoTIFF(单波段)图像(例如 1000 x 1000 像素)子集化为规则的正方形网格(例如 50 x 50 像素；总共 400 个网格对象)；并计算每个对象的平均纬度、平均经度和平均 DN 值？任何想法将不胜感激。

输入数据只是图像，例如从以下位置下载: http://eoimages.gsfc.nasa.gov/images/imagerecords/57000/57752/land_shallow_topo_2048.tif

这可以使用 gdal 在 python 中打开，如下所示:

import gdal
geotiff = gdal.Open ('land_shallow_topo_2048.tif')
colum_numbers,row_numbers,band_numbers=geotiff.RasterXSize,geotiff.RasterYSize,geotiff.RasterCount
print (colum_numbers,row_numbers,band_numbers)
2048 1024 3

Answer 1

有一个 Numpy 整形技巧可以做到这一点。您可以将原始 2D 数组重新整形为 4D 数组，其中将落在单个 50*50 网格中的所有单元格都放入两个唯一的维度中。在这些轴上调用任何函数都会得到聚合结果。

让我们制作一个示例二维数组:

n = 1000
grid_size = 20

grid = np.arange(n*n).reshape(n,n)

然后计算聚合因子和两个维度的网格数:

factor = n / grid_size
yblocks, xblocks = np.array(grid.shape) / factor

然后您可以将原始 grid 数组重新整形为 4 个维度，并在第二个和第四个维度上应用 mean。

grid_small = grid.reshape(yblocks, factor, xblocks, factor).mean(axis=3).mean(axis=1)

您可以通过自己切片一些网格并对其应用均值来测试它:

assert(grid[:factor,:factor].mean() == grid_small[0,0])
assert(grid[-factor:,-factor:].mean() == grid_small[-1,-1])

可视化的结果和差异: