1.GDAL库
GDAL（Geospatial Data Abstraction Library）是一个开源的栅格与矢量地理空间数据转换库，在功能上类似于arcgis系列的arcpy包。GDAL是C++编写的库，但它提供了非常完善的Python绑定，可以使用python调用。
1.1GDAL库安装
Windows系统下GDAL库的安装：建议使用conda install -c conda-forge gdal进行安装。不使用pip是因为，GDAL本质上是一个C++库，而不是python库，直接用pip下载会让PyPi下载C++的源码包，而不是.py的文件，需要额外的多余操作。使用-c conda-forge在conda-forge社区频道下载的包版本更新，更稳定。
1.2GDAL的安装检验
from osgeo import gdal, ogr, osrprint("GDAL 安装验证")
print(f"GDAL 版本: {gdal.__version__}")# 获取 GDAL 数据目录
data_dir = gdal.GetConfigOption('GDAL_DATA')
print(f"GDAL 数据目录: {data_dir}")# 获取驱动数量
driver_count = gdal.GetDriverCount()
print(f"GDAL 驱动数量: {driver_count}")# OGR 驱动数量
ogr_driver_count = ogr.GetDriverCount()
print(f"OGR 驱动数量: {ogr_driver_count}")print("模块 (gdal, ogr, osr) 导入成功")
2栅格数据
2.1读取栅格数据文件
from osgeo import gdal # 主要处理栅格数据gdal.UseExceptions() # 开启异常模式（gdal4之后默认）：出错时抛出异常
# gdal.DontUseExceptions() 不开启异常模式if __name__=="__main__":# 读取影像文件raster_path = r""ds = gdal.Open(raster_path)# 关闭影像文件# GDAL 的 Python 绑定是通过 SWIG 生成的，其资源释放机制依赖于引用计数，而不是 Python 的 delds = None
2.2获取栅格数据文件基本信息
from osgeo import gdal # 主要处理栅格数据# ===================================
gdal.UseExceptions() # 开启异常模式（gdal4之后默认）：出错时抛出异常
# gdal.DontUseExceptions() 不开启异常模式# ---------- 检查栅格----------
def check_raster(ds):print("\n【栅格数据检查】")if ds is None:print("ds为None")else:print(f"文件路径:{ds.GetDescription()}")print(f"栅格尺寸: {ds.RasterXSize} 列 × {ds.RasterYSize} 行")print(f"波段数量: {ds.RasterCount}")# 获取投影信息proj = ds.GetProjection()print(f"投影信息: {proj[:200]}..." if len(proj) > 200 else f"投影信息: {proj}")# 获取地理变换参数gt = ds.GetGeoTransform()if gt:print(f"左上角X坐标: {gt[0]:.6f}") # 左上角像元的X坐标print(f"像元宽度: {gt[1]:.8f}") # 像元宽度（X 方向分辨率）print(f"行旋转参数：{gt[2]:6f}") # 行旋转参数（X 方向随行的变化）print(f"左上角Y坐标: {gt[3]:.6f}") # 左上角像元的 Y 坐标print(f"列旋转参数：{gt[4]:6f}") # 列旋转参数（Y 方向随列的变化）print(f"像元高度: {gt[5]:.8f}") # 像元高度（Y 方向分辨率），通常为负值表示Y轴向下# 计算并打印坐标范围xmin = gt[0]xmax = gt[0] + ds.RasterXSize * gt[1]ymax = gt[3] # 左上角Y是最大值ymin = gt[3] + ds.RasterYSize * gt[5] # gt[5]是负数print(f"坐标范围: X[{xmin:.6f}, {xmax:.6f}], Y[{ymin:.6f}, {ymax:.6f}]")if __name__ == "__main__":# 栅格数据文件读取raster_path = r""ds = gdal.Open(raster_path)# 基础信息读取check_raster(ds)# GDAL 的 Python 绑定是通过 SWIG 生成的，其资源释放机制依赖于引用计数，而不是 Python 的 delds = None # 关闭文件
2.3转换栅格数据文件坐标系（重投影）
if __name__ == "__main__":# 栅格数据文件读取raster_path = r""output_raster = r""target_srs = "EPSG:4326"gdal.Warp(output_raster, # 输出文件路径raster_path, # 输入文件路径dstSRS=target_srs, # 目标坐标系format="GTiff", # 输出格式resampleAlg=gdal.GRA_Bilinear # 重采样算法)
3.矢量数据
3.1读取矢量数据文件
from osgeo import gdal # 主要处理栅格数据
from osgeo import ogr # 主要处理矢量数据
gdal.UseExceptions() #开启异常模式if __name__=="__main__":# ========== 路径 ==========vector_path = r""output_vector = r""ds = ogr.Open(vector_path)check_voctor(ds)ds = None # 关闭文件
3.2输出主要信息
from osgeo import gdal # 主要处理栅格数据
from osgeo import ogr # 主要处理矢量数据
gdal.UseExceptions() #开启异常模式def check_voctor(ds):print("\n【矢量数据检查】")if ds is None:print("ds为None")else:layer = ds.GetLayer(0)print(f"文件路径: {ds.GetDescription()}")print(f"图层名称: {layer.GetName()}")print(f"要素数量: {layer.GetFeatureCount()} 个站点")# 获取投影信息spatial_ref = layer.GetSpatialRef()if spatial_ref:print(f"投影信息: {spatial_ref.GetName()}")else:print("投影信息: 无（可能是地理坐标系）")# 获取范围extent = layer.GetExtent()print(f"坐标范围: X[{extent[0]:.6f}, {extent[1]:.6f}], Y[{extent[2]:.6f}, {extent[3]:.6f}]")# 获取字段信息layer_defn = layer.GetLayerDefn()print("字段列表:")for i in range(layer_defn.GetFieldCount()):field_defn = layer_defn.GetFieldDefn(i)print(f" [{i}] {field_defn.GetName()}")if __name__=="__main__":# ========== 路径 ==========vector_path = r""output_vector = r""ds = ogr.Open(vector_path)check_voctor(ds)ds = None # 关闭文件
3.3坐标系转换
target_srs = "EPSG:32649"gdal.VectorTranslate(output_vector,vector_path,format="ESRI Shapefile",layerCreationOptions=["ENCODING=UTF-8"], # 避免中文乱码dstSRS=target_srs # 目标坐标系)