gdal.Warp可以对栅格数据依据输入的矢量进行裁剪，还可以进行栅格重采样，重投影等操作，关键在于options的参数设置。Options的参数设置可以参见网址：https://gdal.org/python/，其中在左侧可以选择Warp和WarpOptions两个模块查看参数。
1. Warp参数解释
warp是扭曲或者变形的意思，该函数是对一个数据或者多个数据集进行变形，用法如下：
Warp(destNameOrDestDS, srcDSOrSrcDSTab, **kwargs)
参数包括：
- destNameOrDestDS --- 输出数据集名称或对象
- srcDSOrSrcDSTab --- 数据集对象或文件名的数组，或数据集对象或文件名的数组
关键字参数是：
- options --- 返回gdal.WarpOptions()，字符串或字符串数组
- gdal.WarpOptions()的其他关键词参数
如果选项是作为gdal.WarpOptions()对象提供的，则会忽略其他关键字。
2. WarpOptions参数解释（可以略过，直接看代码实例）
warp工具实现的关键在于参数options的设置， WarpOptions的功能即是创建一个可以传递给gdal.Warp()的WarpOptions()对象。下面分析WarpOptions参数，其参数较多，如下所示：
WarpOptions(options=None, format=None, outputBounds=None, outputBoundsSRS=None, xRes=None, yRes=None, targetAlignedPixels=False, width=0, height=0, srcSRS=None, dstSRS=None, coordinateOperation=None, srcAlpha=False, dstAlpha=False, warpOptions=None, errorThreshold=None, warpMemoryLimit=None, creationOptions=None, outputType=gdalconst.GDT_Unknown, workingType=gdalconst.GDT_Unknown, resampleAlg=None, srcNodata=None, dstNodata=None, multithread=False, tps=False, rpc=False, geoloc=False, polynomialOrder=None, transformerOptions=None, cutlineDSName=None, cutlineLayer=None, cutlineWhere=None, cutlineSQL=None, cutlineBlend=None, cropToCutline=False, copyMetadata=True, metadataConflictValue=None, setColorInterpretation=False, overviewLevel='AUTO', callback=None, callback_data=None)
挨个对参数进行下解释。
- options --- 可以是一个字符串数组，一个字符串或让空和填充从其他关键字。
- format --- 输出格式(“GTiff”等)
- outputBounds --- 在输出的空间参照系统中的输出边界(minX, minY, maxX, maxY)
- outputBoundsSRS --- 表示输出边界的SRS以防它们没有在dstSRS中表示，其中SRS是spatial reference system的缩写，即空间参照系统
- xRes, yRes --- 目标SRS的空间分辨率
- targetAlignedPixels --- 是否强制输出边界为输出分辨率的倍数
- width --- 输出光栅的像素宽度
- height --- 输出光栅的像素高度
- srcSRS --- 输入的 SRS
- dstSRS --- 输出的 SRS
- coordinateOperation -- 坐标操作是 PROJ string 或者是 WKT string
- srcAlpha --- 是否强制将输入数据集的最后一个频带视为alpha频带
- dstAlpha --- 是否强制创建输出波段
- outputType --- 输出类型 (gdalconst.GDT_Byte, etc...)
- workingType --- 操作类型 (gdalconst.GDT_Byte, etc...)
- warpOptions --- warping options 列表
- errorThreshold --- 近似转换器的误差阈值 (in pixels)
- warpMemoryLimit --- 工作缓冲区的大小，以MB为单位
- resampleAlg --- 重采样模式
- creationOptions --- creation options 列表
- srcNodata --- source nodata value(s)
- dstNodata --- output nodata value(s)
- multithread --- 是否多线程计算和I/O操作
- tps --- 是否使用薄板样条GCP转换器
- rpc --- 是否使用RPC转换器
- geoloc --- 是否使用GeoLocation阵列转换器
- polynomialOrder --- 多项式GCP插值阶
- transformerOptions --- transformer options 列表
- cutlineDSName --- 切割线的数据集的名称
- cutlineLayer --- 切割线的图层名称
- cutlineWhere --- 切割线的 WHERE 从句
- cutlineSQL --- 切割线的 SQL 语句
- cutlineBlend --- 以像素为单位的切割线混合距离
- cropToCutline --- 是否使用切割线范围作为输出边界
- copyMetadata --- 是否拷贝源元数据
- metadataConflictValue --- 元数据数据冲突值
- setColorInterpretation --- 是否强制输入波段的颜色解释到输出波段
- overviewLevel --- 指定必须使用的源文件的概述级别
- callback --- 回调函数
- callback_data --- 用于回调的用户数据
3. 栅格数据矢量裁剪（批量）
在矢量裁剪中，主要用到的是cutline相关的参数，包括切割后的数据集（cutlineDSName），切割图层（cutlineLayer ），是否使用切割线范围作为输出边界（cropToCutline）等，下面展示实例代码。
from osgeo import gdal
ds = gdal.Warp(tif_out, tif_input, format='GTiff', cutlineDSName=shp_boundary, cropToCutline=True, dstNodata=0)
4. 遥感影像重投影
对坐标系进行重投影也很简单，主要是坐标系的设置，如果是非标准的自定义的坐标系，可以参见方法ESRI中WKT格式的空间投影坐标转换为GDAL中的Proj4格式转化为proj4的格式。具体实例代码如下：
from osgeo import gdal
ds = gdal.Warp(tif_out, tif_input, srcSRS='EPSG:4326', dstSRS=srs_proj4) # srs_proj4是自定义的坐标系
5. 遥感影像重采样（修改分辨率）
栅格数据或者说遥感影像的重采样也比较简单，重采样方法的确定在gdal中可以通过ReprojectImage()函数查看，参见https://gdal.org/python/。
采用Warp函数的话，主要修改的resampleAlg参数，默认的最近邻重采样方法，查看源码发现其他可选参数如下。具体重采样方法的差异不再赘述。
if resampleAlg is not None:
mapMethodToString = {
gdalconst.GRA_NearestNeighbour: 'near',
gdalconst.GRA_Bilinear: 'bilinear',
gdalconst.GRA_Cubic: 'cubic',
gdalconst.GRA_CubicSpline: 'cubicspline',
gdalconst.GRA_Lanczos: 'lanczos',
gdalconst.GRA_Average: 'average',
gdalconst.GRA_RMS: 'rms',
gdalconst.GRA_Mode: 'mode',
gdalconst.GRA_Max: 'max',
gdalconst.GRA_Min: 'min',
gdalconst.GRA_Med: 'med',
gdalconst.GRA_Q1: 'q1',
gdalconst.GRA_Q3: 'q3',
gdalconst.GRA_Sum: 'sum',
}
代码如下：
from osgeo import gdal
ds = gdal.Warp(tif_out, tif_input, resampleAlg=gdalconst.GRA_NearestNeighbour, xRes=500, yRes=500)
至此，栅格数据（如遥感影像）的裁剪、重投影、重采样都可以用Warp函数实现，而且可以在同一个函数下实现所有的操作，只需要把各个options放在一起即可。
初次使用gdal完毕，撒花~
网友评论