python矢量裁剪栅格(按掩模提取)
保持好奇，持续热爱。
前言
矢量裁剪栅格，也就是将栅格中对应矢量区域的像素裁剪出来，形成一个新栅格的过程，对应Arcgis中空间分析工具里的按掩模提取功能。
过程
gdal方法实现
其中一种方法是利用gdal.Warp工具进行裁剪。
一般情况下，只需要用到几个参数，具体使用方式如下：
try:
import gdal
except:
from osgeo import gdal
def extract_by_shp_gdal(in_shp_path, in_raster_path, out_raster_path):
input_raster = gdal.Open(in_raster_path)
# 利用gdal.Warp进行裁剪
# https://gdal.org/api/python/osgeo.gdal.html#osgeo.gdal.Warp
result = gdal.Warp(
out_raster_path,
input_raster,
format = 'GTiff',
cutlineDSName = in_shp_path, # 用于裁剪的矢量
cropToCutline = True, # 是否使用cutlineDSName的extent作为输出的界线
dstNodata = 0 # 输出数据的nodata值
result.FlushCache()
del resultrasterio方法实现
gdal方法好处就是安装了gdal就行，但是偶尔可能会出现小bug。rasterio则比较稳健，非常好用。
import fiona
import rasterio
from rasterio.mask import mask
import gdal
except:
from osgeo import gdal
def extract_by_shp_rasterio(in_shp_path, in_raster_path, out_raster_path):
with fiona.open(in_shp_path, "r", encoding='utf-8') as shapefile:
# 获取所有要素feature的形状geometry
geoms = [feature["geometry"] for feature in shapefile]
with rasterio.open(in_raster_path) as src:
out_image, out_transform = mask(src, geoms, crop=True)
out_meta = src.meta.copy()
# 更新输出文件的元数据
out_meta.update({"driver": "GTiff",
"height": out_image.shape[1],