对比结果

校正前：无经纬度信息

校正后：带经纬度，但对不齐，可能存在问题

完整代码

import numpy as np
from osgeo import gdal
import spectral

def read_hdr_file(hdr_file):
    """
    使用 spectral 读取 HDR 文件，并提取其中的元数据
    """
    hdr = spectral.io.envi.read_envi_header(hdr_file)
    return hdr

def extract_rpc_info(rpc_list):
    """
    格式化 RPC 信息
    """
    rpc_metadata = {
        'LINE_OFF': float(rpc_list[0]),
        'SAMP_OFF': float(rpc_list[1]),
        'LAT_OFF': float(rpc_list[2]),
        'LONG_OFF': float(rpc_list[3]),
        'HEIGHT_OFF': float(rpc_list[4]),
        'LINE_SCALE': float(rpc_list[5]),
        'SAMP_SCALE': float(rpc_list[6]),
        'LAT_SCALE': float(rpc_list[7]),
        'LONG_SCALE': float(rpc_list[8]),
        'HEIGHT_SCALE': float(rpc_list[9]),
        'LINE_NUM_COEFF': list(map(float, rpc_list[10:30])),
        'LINE_DEN_COEFF': list(map(float, rpc_list[30:50])),
        'SAMP_NUM_COEFF': list(map(float, rpc_list[50:70])),
        'SAMP_DEN_COEFF': list(map(float, rpc_list[70:90])),
    }
    return rpc_metadata

def convert_hdr_rpc_to_geotiff(raster_file, hdr_file, output_tiff):
    # 打开栅格数据集
    raster_dataset = gdal.Open(raster_file)

    # 使用 spectral 读取并提取 HDR 文件中的元数据
    metadata = read_hdr_file(hdr_file)
    rpc_list = metadata.get('rpc info', '')

    if not rpc_list:
        raise ValueError("HDR 文件中没有找到有效的 RPC 信息")

    # 提取并格式化 RPC 信息
    rpc_metadata = extract_rpc_info(rpc_list)

    # 设置 RPC 元数据
    for key, value in rpc_metadata.items():
        if isinstance(value, list):
            raster_dataset.SetMetadataItem(key, ' '.join(map(str, value)), "RPC")
        else:
            raster_dataset.SetMetadataItem(key, str(value), "RPC")

    # 创建 GeoTIFF 驱动
    driver = gdal.GetDriverByName("GTiff")

    # 使用 RPC 变换进行地理校正
    warp_options = gdal.WarpOptions(dstSRS='EPSG:4326', rpc=True, resampleAlg='bilinear')
    output_dataset = gdal.Warp(output_tiff, raster_dataset, options=warp_options)

    # 清理资源
    output_dataset = None
    raster_dataset = None

    print(f"导出完成：{output_tiff}")



# 输入栅格文件和 HDR 文件路径
raster_file = r'rad.dat'
hdr_file = r'_rad.hdr'
output_tiff = 'newtif.tif'

# 调用函数进行转换
convert_hdr_rpc_to_geotiff(raster_file, hdr_file, output_tiff)