Python处理遥感影像：GDAL库使用教程

要处理遥感影像，需掌握GDAL的几个关键技巧：1. 读取基本信息，包括尺寸、波段数、地理变换和投影信息；2. 读写波段数据，使用NumPy数组进行计算并保存结果；3. 裁剪与重采样操作，通过gdal.Warp实现区域裁剪和分辨率调整；4. 注意数据类型、NoData值处理、内存管理和资源释放。这些步骤构成了Python中利用GDAL处理遥感影像的核心流程。

处理遥感影像在地理信息、环境监测、城市规划等领域非常常见，而Python凭借其强大的生态和易用性，成为很多人的首选语言。如果你打算用Python来处理遥感数据，GDAL（Geospatial Data Abstraction Library） 是一个绕不开的工具。

这篇文章会从实用角度出发，介绍几个你最可能需要掌握的 GDAL 使用技巧，帮助你快速上手处理遥感影像。

读取遥感影像的基本信息

在开始处理之前，了解影像的基本信息是必须的，比如波段数、分辨率、坐标范围、投影信息等。

使用 GDAL 打开文件后，可以通过 Dataset 对象获取这些内容：

from osgeo import gdal

dataset = gdal.Open("your_image.tif")
print("行数：", dataset.RasterYSize)
print("列数：", dataset.RasterXSize)
print("波段数：", dataset.RasterCount)

geotransform = dataset.GetGeoTransform()
projection = dataset.GetProjection()

print("地理变换参数：", geotransform)
print("投影信息：", projection)

• RasterYSize 和 RasterXSize 是图像尺寸。
• GetGeoTransform() 返回的是仿射变换参数，包括左上角坐标和像元大小。
• 投影信息一般是以 WKT 格式返回的字符串。

小提示：如果读取失败，可能是路径不对或者文件格式不支持，记得检查一下扩展名是否正确。

读取与写入波段数据

遥感影像通常包含多个波段（如红绿蓝、近红外等），我们经常需要单独操作某个波段，比如提取 NDVI 或进行分类。

读取第一个波段的数据可以这样写：

band = dataset.GetRasterBand(1)
array = band.ReadAsArray()

• ReadAsArray() 返回的是 NumPy 数组，可以直接参与计算。
• 如果你想把处理后的数组保存为新文件，可以用下面的方法创建并写入：

driver = gdal.GetDriverByName("GTiff")
out_dataset = driver.Create("output.tif", array.shape[1], array.shape[0], 1, band.DataType)
out_dataset.SetGeoTransform(geotransform)
out_dataset.SetProjection(projection)
out_band = out_dataset.GetRasterBand(1)
out_band.WriteArray(array)
out_band.FlushCache()

• 注意输出数据类型要和原始一致，否则可能会出错。
• 写完之后别忘了释放资源：del out_dataset。

常见预处理操作：裁剪与重采样

实际工作中，我们经常需要对影像进行裁剪或调整分辨率，这两个功能 GDAL 都原生支持。

裁剪影像

可以使用 gdal.Warp() 函数配合裁剪范围（bounding box）完成：

options = gdal.WarpOptions(outputBounds=[minX, minY, maxX, maxY])
gdal.Warp("clipped.tif", dataset, options=options)

• 其中 minX, maxX 等是你希望保留区域的边界坐标。
• 这个方法也适用于格式转换，只需要修改输出文件名的后缀即可。

重采样影像

改变影像分辨率也很简单，只需指定 xRes 和 yRes：

resampled = gdal.Warp("resampled.tif", dataset, xRes=30, yRes=30)

• 上面的例子将影像重采样到 30 米分辨率。
• 默认使用最近邻插值，也可以通过 resampleAlg="bilinear" 指定其他方式。

一些容易忽略但重要的细节

• 数据类型：不同传感器输出的数据类型可能不一样，比如有的是 uint16，有的是 float32，在做归一化或计算前最好先确认。

• NoData 值：有些影像有无效像素，处理时要记得设置或跳过它们：

  nodata = band.GetNoDataValue()
  if nodata is not None:
      array[array == nodata] = 0  # 示例替换为0

• 内存管理：大影像一次性加载可能吃不消，建议分块读取或使用 ReadAsArray(xoff, yoff, xsize, ysize) 来控制内存。

• 关闭数据集：处理完成后记得释放资源，避免锁住文件或内存泄漏。

基本上就这些。掌握了这几个关键步骤，你就能用 Python + GDAL 完成大多数遥感影像的基础处理了。