【Python&RS】基于GDAL给无人机图片定义坐标系

        前段时间有过一个想法，就是如果可以给无人机拍摄的图片定义坐标系，再使用GADL库里的镶嵌拼接函数，是不是就可以实现快速拼接影像。虽然结果不是正射影像，但效率比无人机厂家的软件提高了很多很多，主要还是看用途。

        有了这个想法后就要行动起来，定义一个坐标系，一般需要两个点。一是需要投影坐标系的参数，二就是仿射地理变换参数。投影坐标系非常容易得到，投影坐标系的参数GDAL有内置的一部分，当然也可以读取已有数据的坐标系。但仿射地理变换参数就比较麻烦了，我印象中定义仿射地理变换参数时只能通过图片左上角的点就是设置，不像ENVI可以任意选点。ENVI方便很多，有兴趣可以查看【RS】基于ENVI给图片/影像定义坐标系。所以我们只要能解决仿射地理变换参数就可以解决图片定义坐标系的问题。

        大致思路是：通过之前计算图片某点的投影坐标的方法，得到图片左上角的投影坐标和偏转角度，以此建立仿射地理变换参数，再通过这些创建新的数据资源将原来的数据写入。

一、获取图片中心点经纬度

        这一步之前的博文已经说过无数次了，这里就不解释。如果你只是娱乐一下，并不需要定义真实的坐标系，那么可以直接跳到最后一步，仿射地理变换参数直接设置0，1即可。

def Get_LatLon(path_image):
    """
    :param path_image: 输入图片路径
    :return: 返回图片中心点经纬度
    """
    # 获取图片的经纬度信息
    f = open(path_image, 'rb')
    contents = exifread.process_file(f)
    longitude = contents["GPS GPSLongitude"].values
    longitude_f = longitude[0].num/longitude[0].den + (longitude[1].num/longitude[1].den/60) + (longitude[2].num/longitude[2].den/3600)
    latitude = contents["GPS GPSLatitude"].values
    latitude_f = latitude[0].num/latitude[0].den + (latitude[1].num/latitude[1].den/60) + (latitude[2].num/latitude[2].den/3600)
    # print("经度:", longitude_f)  # contents['GPS GPSLatitudeRef']
    # print("纬度:", latitude_f)  # contents['GPS GPSLongitudeRef']
    f.close()
    return latitude_f, longitude_f

二、转换投影

        无人机拍摄的影像一般都是GPS经纬度坐标（WGS84），我这里要将他定义为UTM/WGS84 51N，所以需要将它的地理坐标系先转换成投影坐标系。

def LonLat_Meter(lat, lon):
    """
    :param lat: 图片中心点纬度
    :param lon: 图片中心点经度
    :return: 返回UTM投影坐标
    """
    source = osr.SpatialReference()
    # 初始化osr.SpatialReference对象以形成一个合法的坐标系统
    source.ImportFromEPSG(4326)
    # 向对象中写入Source_EPSG坐标系统
    target = osr.SpatialReference()
    target.ImportFromEPSG(32651)
    # 这里是用内置的EPSG对应的坐标系作为转换参数
    transform = osr.CoordinateTransformation(source, target)
    point = ogr.CreateGeometryFromWkt("POINT (%s %s)" % (lat, lon))
    # 报错的话，将经纬度翻过来
    point.Transform(transform)
    # print(point.GetX(), point.GetY())
    return point.GetX(), point.GetY()

三、获取图片偏转角

        无人机航拍时不可能一直与正北方向一致，所以还需要获取其成像时的偏转角。这一步是为了计算图片左上角的仿射地理变换参数。我这里使用的是大疆的无人机，其他厂商的无人机成像参数需要自己查阅。

def Get_Image_Yaw_angle(file_path):
    """
    :param file_path: 输入图片路径
    :return: 图片的偏航角
    """
    # 获取图片偏航角
    b = b