在地理數據處理中，矢量裁剪柵格是一個非常重要的操作，它可以幫助我們提取感興趣的區域并獲得更精確的分析結果。其重要性包括：

1. 區域限定：地球科學研究通常需要關注特定的地理區域。通過矢量裁剪柵格，我們可以將柵格數據限定在特定的研究區域內，去除不相關的數據，從而減少計算量和提高分析效率。
2. 空間一致性：柵格數據通常具有均勻的空間分辨率，但我們可能只對特定區域感興趣。通過矢量裁剪柵格，我們可以將柵格數據限定在感興趣的區域內，確保分析結果在這個區域內具有更高的精度和一致性。
3. 數據匹配：在地球科學研究中，我們通常需要將不同來源的數據進行比較和分析。通過矢量裁剪柵格，我們可以將不同數據源的柵格數據裁剪為相同的區域，使它們具有相同的空間范圍和分辨率，從而方便進行數據匹配和比較。
4. 空間分析：地球科學研究中常常需要進行各種空間分析，如土地利用變化、水文模擬、地質建模等。通過矢量裁剪柵格，我們可以將柵格數據限定在特定區域內，以便更準確地進行空間分析，并獲得與實際情況更符合的結果。

簡單來說，制作數據的人要盡可能的考慮到多個用戶的需求，所以往往數據的范圍比較大，比如全國范圍、全球范圍，但是使用數據的人往往有自己聚焦的研究區，所以其 “研究區” 更明確。所以，在進行科研或者實際工作時，我們一般會裁剪出自己的目標區域數據，減少數據量，也方便分析和挖掘數據。

🌰?假設你獲取到了全中國多年的氣溫柵格數據，而你的研究區是北京市，老板讓你求北京市的平均溫度序列，后續還要用北京市的數據進行其他分析。這時候我們第一步自然是先把我們需要的數據裁剪出來，既方便處理，也方便出圖。

裁剪柵格數據一般有兩種情況：

1. 依據范圍裁剪：有大概的范圍，比如有研究區的經緯度范圍，需要根據相應的經緯度范圍來進行裁剪
2. 基于矢量數據裁剪：有研究區的矢量數據，那么我們可以使用矢量數據來進行裁剪

可能用到的包

import rasterio as rio
from rasterio.plot import show
from rasterio.windows import Window
from rasterio.transform import from_origin
import geopandas as gpd
from rasterio.mask import mask

讀取柵格數據并簡單可視化：

dataset = rio.open('/home/mw/input/GeoPy17233/GeoPy1/temp_grid_84/tem2015_84.tif')
show(dataset)

我們讀取的數據是整個中國范圍的，可以通過下面的代碼查看數據的范圍：

dataset.bounds

BoundingBox(left=65.45935893761649, bottom=13.120243022441713, right=137.45935893761649, top=56.01024302244171)

依據經緯度范圍裁剪

假設我們的研究區是湖南省，借助搜索引擎，發現其經緯度范圍為東經 108.786106 ~ 114.256514，北緯 24.643089 ~ 30.1287。

我們用下面的代碼對上面的數據進行裁剪并且保存：

# 定義經緯度范圍,湖南省范圍
min_lon = 108
max_lon = 115
min_lat = 24
max_lat = 31# 打開柵格數據
with rio.open('/home/mw/input/GeoPy17233/GeoPy1/temp_grid_84/tem2015_84.tif') as src:# 獲取柵格數據的元數據和變換信息profile = src.profiletransform = src.transform# 計算裁剪后的柵格數據的行列范圍col_start, row_start = ~transform * \(min_lon, max_lat)  # ~transform表示坐標轉換的操作符號col_end, row_end = ~transform * (max_lon, min_lat)# 計算裁剪后的柵格數據的寬度和高度width = int(col_end - col_start)height = int(row_end - row_start)# 創建裁剪后的柵格數據數組cropped_data = src.read(window=Window(col_start, row_start, width, height))# 更新元數據profile.update(width=width, height=height, transform=from_origin(min_lon, max_lat, profile['transform'][0], -profile['transform'][4]))# 創建新的柵格數據文件并寫入裁剪后的數據with rio.open('/home/mw/temp/tem2015_84_hunan.tif', 'w', **profile) as dst:dst.write(cropped_data)# 可視化裁剪后的數據
dataset = rio.open('/home/mw/temp/tem2015_84_hunan.tif')
show(dataset)

從上面結果可以看出我們已經成功裁剪出了我們需要范圍的數據。

# 定義經緯度范圍,湖南省范圍
min_lon = 108
max_lon = 115
min_lat = 24
max_lat = 31# 打開柵格數據
with rio.open('/home/mw/input/GeoPy17233/GeoPy1/temp_grid_84/tem2015_84.tif') as src:# 獲取柵格數據的元數據和變換信息profile = src.profiletransform = src.transform# 計算裁剪后的柵格數據的行列范圍col_start, row_start = ~transform * \(min_lon, max_lat)  # ~transform表示坐標轉換的操作符號col_end, row_end = ~transform * (max_lon, min_lat)# 計算裁剪后的柵格數據的寬度和高度width = int(col_end - col_start)height = int(row_end - row_start)# 創建裁剪后的柵格數據數組cropped_data = src.read(window=Window(col_start, row_start, width, height))# 更新元數據profile.update(width=width,height=height,transform=from_origin(min_lon,max_lat,profile['transform'][0],profile['transform'][4]))  # 提示：仔細觀察這一行~# 創建新的柵格數據文件并寫入裁剪后的數據with rio.open('/home/mw/temp/tem2015_84_hunan_wrong.tif', 'w', **profile) as dst:dst.write(cropped_data)# 可視化裁剪后的數據
dataset = rio.open('/home/mw/temp/tem2015_84_hunan_wrong.tif')
show(dataset)

用矢量數據來裁剪

下面我們再來學習如何用矢量文件來進行裁剪。

以北京市為例：

# 打開柵格數據
with rio.open('/home/mw/input/GeoPy17233/GeoPy1/temp_grid_84/tem2015_84.tif') as src:# 打開矢量文件beijing = gpd.read_file('/home/mw/input/GeoPy17233/GeoPy1/beijing.shp')# 對柵格數據進行矢量裁剪clipped_data, clipped_transform = mask(src, beijing.geometry, crop=True)# 獲取裁剪后的元數據clipped_profile = src.profileclipped_profile.update(transform=clipped_transform,width=clipped_data.shape[2],height=clipped_data.shape[1])# 創建新的柵格數據文件并寫入裁剪后的數據with rio.open('/home/mw/temp/tem2015_84_beijing.tif', 'w', **clipped_profile) as dst:dst.write(clipped_data)# 可視化結果
dataset = rio.open('/home/mw/temp/tem2015_84_beijing.tif')
show(dataset)

?可以看出我們已經成功的用北京市的矢量裁剪出了我們想要的結果。