利用GDAL库实现像素坐标与地理坐标的转换

在地理信息系统中，经常需要将像素坐标转换为地理坐标，或者将地理坐标转换为像素坐标。这对于遥感图像处理、地理信息系统开发以及空间数据可视化等任务都是至关重要的。GDAL（Geospatial Data Abstraction Library）是一个开源的栅格和矢量地理空间数据转换库，它提供了丰富的功能来处理这类转换。

下面，我们将展示如何使用GDAL库来实现像素坐标与地理坐标之间的转换，并将这些操作封装成一个Python函数，以便在项目中重复使用。

代码

首先，确保已经安装了GDAL库。

然后，我们可以编写一个Python函数来实现坐标转换的功能。

根据地理转换参数和像素坐标计算地理坐标

from osgeo import gdal  

def pixel_to_geo(col, row, gt):  
    """  
    根据地理转换参数和像素坐标计算地理坐标。

    参数:  
    col (int): 像素列坐标。
    row (int): 像素行坐标。
    gt (tuple): 地理转换参数，包含(a, b, c, d, e, f)。

    返回:  
    tuple: 包含(x_geo, y_geo)的地理坐标。
    """  
    a, b, c, d, e, f = gt  
    x_geo = a * col + b * row + a * 0.5 + b * 0.5 + c  
    y_geo = d * col + e * row + d * 0.5 + e * 0.5 + f  
    return x_geo, y_geo

根据地理坐标计算像素坐标（即从地理坐标转换到像素坐标，也称为逆地理转换）

def geo_to_pixel(geotransform, geo_x, geo_y):  
    """  
    根据地理转换参数和地理坐标计算像素坐标。
      
    参数:  
    geotransform (tuple): GDAL地理转换参数，通常是一个包含6个元素的元组。
    geo_x (float): 地理X坐标（经度）。
    geo_y (float): 地理Y坐标（纬度）。
      
    返回:  
    tuple: 包含像素坐标（X, Y）的元组。
    """  
    pixel_x = int((geo_x - geotransform[0]) / geotransform[1])  
    pixel_y = int((geo_y - geotransform[3]) / geotransform[5])  
    return pixel_x, pixel_y

在上面的代码中，我们定义了两个函数pixel_to_geo和geo_to_pixel，分别用于将像素坐标转换为地理坐标和将地理坐标转换为像素坐标。这两个函数都接受地理转换参数gt作为输入，这是通过调用gdal.Open方法读取GeoTIFF文件后，调用GetGeoTransform方法获得的。

注意，在geo_to_pixel函数中，我们使用了pixel_offset参数来考虑像素中心的位置。在遥感应用中，通常将像素的中心点作为地理坐标的参考点，而不是像素的左上角角点。因此，在反向转换时，我们通常需要考虑这个偏移量。

测试

在envi打开某景Landsat8影像，进行代码测试。这景Landsat8影像的投影为utm。

影像的测试点信息：

行列号为5469.4409, 4638.8315

X,Y为374068.2275， 2694650.0537

下图是我们的代码计算的地理坐标。

下图是我们的代码反算的屏幕坐标。

计算结果与envi的结果保持一致！

完成测试！

完！