Rasterio簡介

有沒有覺得用GDAL的Python繫結書寫的程式碼很不Pythonic，強迫症的你可能有些忍受不了。不過，沒關係，MapBox旗下的開源庫Rasterio幫我們解決了這個痛點。

Rasterio是基於GDAL庫二次封裝的更加符合Python風格的主要用於空間柵格資料處理的Python庫。

Rasterio中柵格資料模型基本和GDAL類似，需要注意的是：

在Rasterio 1.0以後，對於GeoTransform的表示棄用了GDAL風格的放射變換，而使用了Python放射變換的第三方庫affine庫的風格。

對於放射變換

affine.Affine(a,b,c,d,e,f)

GDAL中對應的引數順序是：(c,a,f,e)

採用新的放射變換模型的好處是，如果你需要計算某個行列號的地理座標，直接使用行列號跟給放射變換物件相乘即可，完全符合數學上矩陣乘法的操作，更加直觀和方便。

柵格資料讀取程式碼示例

下面的示例程式中演示瞭如何讀取一個GeoTIFF檔案並獲取相關資訊，需要注意的是：

1、rasterio使用rasterio.open()函式開啟一個柵格檔案

2、rasterio使用read()函式可以將資料集轉為numpy.ndarray，該函式如果不帶引數，將把資料的所有波段做轉換（第一維是波段數），如果指定波段，則只取得指定波段對應的資料（波段索引從1開始）

3、資料的很多元資訊都是以資料集的屬性進行表示的

import rasterio

with rasterio.open('example.tif') as ds:
 print('該柵格資料的基本資料集資訊（這些資訊都是以資料集屬性的形式表示的）：')
 print(f'資料格式：{ds.driver}')
 print(f'波段數目：{ds.count}')
 print(f'影像寬度：{ds.width}')
 print(f'影像高度：{ds.height}')
 print(f'地理範圍：{ds.bounds}')
 print(f'反射變換引數（六引數模型）：\n {ds.transform}')
 print(f'投影定義：{ds.crs}')
 # 獲取第一個波段資料，跟GDAL一樣索引從1開始
 # 直接獲得numpy.ndarray型別的二維陣列表示，如果read()函式不加引數，則得到所有波段（第一個維度是波段）
 band1 = ds.read(1)
 print(f'第一波段的最大值：{band1.max()}')
 print(f'第一波段的最小值：{band1.min()}')
 print(f'第一波段的平均值：{band1.mean()}')
 # 根據地理座標得到行列號
 x,y = (ds.bounds.left + 300,ds.bounds.top - 300) # 距離左上角東300米，南300米的投影座標
 row,col = ds.index(x,y) # 對應的行列號
 print(f'(投影座標{x},{y})對應的行列號是({row},{col})')
 # 根據行列號得到地理座標
 x,y = ds.xy(row,col) # 中心點的座標
 print(f'行列號({row},{col})對應的中心投影座標是({x},{y})')
 # 那麼如何得到對應點左上角的資訊
 x,y = (row,col) * ds.transform
 print(f'行列號({row},{col})對應的左上角投影座標是({x},{y})')

輸出如下：

該柵格資料的基本資料集資訊（這些資訊都是以資料集屬性的形式表示的）：
資料格式：GTiff
波段數目：3
影像寬度：4800
影像高度：4800
地理範圍：BoundingBox(left=725385.0,bottom=2648415.0,right=869385.0,top=2792415.0)
反射變換引數（六引數模型）：
 | 30.00,0.00,725385.00|
| 0.00,-30.00,2792415.00|
| 0.00,1.00|
投影定義：CRS({'init': 'epsg:32649'})
第一波段的最大值：5459
第一波段的最小值：-313
第一波段的平均值：489.80300625
(投影座標725685.0,2792115.0)對應的行列號是(10,10)
行列號(10,10)對應的中心投影座標是(725700.0,2792100.0)
行列號(10,10)對應的左上角投影座標是(725685.0,2792115.0)

以上這篇使用Rasterio讀取柵格資料的例項講解就是小編分享給大家的全部內容了，希望能給大家一個參考，也希望大家多多支援我們。