如何獲取geotif中單元格的座標？

Question

我有一個地理信息tif。使用gdal，我可以將光柵文件轉換爲數組（numpy）。

如何獲取該數組中一個條目的座標？

Answer 1

使用仿射變換矩陣，它將像素座標映射到世界座標。例如，使用affine軟件包。 （還有其他的方法做同樣的，使用簡單的數學。）

from affine import Affine 
fname = '/path/to/raster.tif' 

這裏有兩種方式來獲得仿射變換矩陣，T0。例如，使用GDAL/Python的：

from osgeo import gdal 
ds = gdal.Open(path, gdal.GA_ReadOnly) 
T0 = Affine.from_gdal(*ds.GetGeoTransform()) 
ds = None # close 

例如，使用rasterio：

import rasterio 
with rasterio.open(fname, 'r') as r: 
    T0 = r.affine 

13759 GDAL（T0）用於變換陣列的約定是參考像素的角落。您可能要代替參考像素中心，所以它需要50％翻譯：

T1 = T0 * Affine.translation(0.5, 0.5) 

而且現在從像素座標轉換爲世界座標，乘以矩陣，這是可以做到的座標一個簡單的函數：

rc2xy = lambda r, c: (c, r) * T1 

現在，得到一個光柵的座標在第一行，第二列（索引[0, 1]）：

print(rc2xy(0, 1)) 

而且，請注意，如果您需要從世界座標中獲取像素座標，則可以使用倒立仿射變換矩陣~T0。