RPC | 高分三号雷达数据RPC校正代码
RPC | 高分三号雷达数据RPC校正代码
ytkzRPC（rational polynomial coefficients），有理多项式系数。RPC是传感器严格几何模型的拟合形式。是根据卫星平台载荷测量的平台运行轨迹参数、姿态参数、传感器安装参数及传感器内部几何参数构建的像-地关系几何模型，一共有90个参数。
相关理论不再展示，感兴趣去网上搜，现在是基于高分三号进行RPC校正，代码如下：
from osgeo import gdal
import re
import os
#这是一个读取RPC文件的函数
def read_rpb (rpbfile):
with open(rpbfile, 'r') as f:
buff = f.read()
# 关键字1,2(修改引号间的内容)
# 命名参考的网址：http://geotiff.maptools.org/rpc_prop.html
ERR_BIAS1 = 'errBias = ' # 错误-偏差。图像中所有点的RMS偏差误差（单位：米/水平轴）（-1.0，如果未知）
ERR_BIAS2 = ';'
ERR_RAND1 = 'errRand = ' # 错误-随机。图像中每个点的每个水平轴的RMS随机误差，以米为单位（如果未知，则为-1.0）
ERR_RAND2 = ';'
LINE_OFF1 = 'lineOffset = ' # 线偏移
LINE_OFF2 = ';'
SAMP_OFF1 = 'sampOffset = ' # 样品偏移
SAMP_OFF2 = ';'
LAT_OFF1 = 'latOffset = ' # 大地纬度偏移
LAT_OFF2 = ';'
LONG_OFF1 = 'longOffset = ' # 大地经度偏移
LONG_OFF2 = ';'
HEIGHT_OFF1 = 'heightOffset = ' # 大地高程偏移
HEIGHT_OFF2 = ';'
LINE_SCALE1 = 'lineScale = ' # 线刻度
LINE_SCALE2 = ';'
SAMP_SCALE1 = 'sampScale = ' # 样本量表
SAMP_SCALE2 = ';'
LAT_SCALE1 = 'latScale = ' # 大地纬度刻度
LAT_SCALE2 = ';'
LONG_SCALE1 = 'longScale = ' # 大地经度标度
LONG_SCALE2 = ';'
HEIGHT_SCALE1 = 'heightScale = ' # 大地高度标尺
HEIGHT_SCALE2 = ';'
LINE_NUM_COEFF1 = 'lineNumCoef = (' # 线分子系数。r n方程分子中的多项式的二十个系数。
LINE_NUM_COEFF2 = ');'
LINE_DEN_COEFF1 = 'lineDenCoef = (' # 线分母系数。r n方程的分母中的多项式的二十个系数。
LINE_DEN_COEFF2 = ');'
SAMP_NUM_COEFF1 = 'sampNumCoef = (' # 样本分子系数。c n方程分子中的多项式的二十个系数。
SAMP_NUM_COEFF2 = ');'
SAMP_DEN_COEFF1 = 'sampDenCoef = (' # 样本分母系数。c n方程的分母中的多项式的二十个系数。
SAMP_DEN_COEFF2 = ');'
# 正则化提取数值
pat_ERR_BIAS = re.compile(ERR_BIAS1 + '(.*?)' + ERR_BIAS2, re.S)
result_ERR_BIAS = pat_ERR_BIAS.findall(buff)
ERR_BIAS = result_ERR_BIAS[0]
pat_ERR_RAND = re.compile(ERR_RAND1 + '(.*?)' + ERR_RAND2, re.S)
result_ERR_RAND = pat_ERR_RAND.findall(buff)
ERR_RAND = result_ERR_RAND[0]
pat_LINE_OFF = re.compile(LINE_OFF1 + '(.*?)' + LINE_OFF2, re.S)
result_LINE_OFF = pat_LINE_OFF.findall(buff)
LINE_OFF = result_LINE_OFF[0]
pat_SAMP_OFF = re.compile(SAMP_OFF1 + '(.*?)' + SAMP_OFF2, re.S)
result_SAMP_OFF = pat_SAMP_OFF.findall(buff)
SAMP_OFF = result_SAMP_OFF[0]
pat_LAT_OFF = re.compile(LAT_OFF1 + '(.*?)' + LAT_OFF2, re.S)
result_LAT_OFF = pat_LAT_OFF.findall(buff)
LAT_OFF = result_LAT_OFF[0]
pat_LONG_OFF = re.compile(LONG_OFF1 + '(.*?)' + LONG_OFF2, re.S)
result_LONG_OFF = pat_LONG_OFF.findall(buff)
LONG_OFF = result_LONG_OFF[0]
pat_HEIGHT_OFF = re.compile(HEIGHT_OFF1 + '(.*?)' + HEIGHT_OFF2, re.S)
result_HEIGHT_OFF = pat_HEIGHT_OFF.findall(buff)
HEIGHT_OFF = result_HEIGHT_OFF[0]
pat_LINE_SCALE = re.compile(LINE_SCALE1 + '(.*?)' + LINE_SCALE2, re.S)
result_LINE_SCALE = pat_LINE_SCALE.findall(buff)
LINE_SCALE = result_LINE_SCALE[0]
pat_SAMP_SCALE = re.compile(SAMP_SCALE1 + '(.*?)' + SAMP_SCALE2, re.S)
result_SAMP_SCALE = pat_SAMP_SCALE.findall(buff)
SAMP_SCALE = result_SAMP_SCALE[0]
pat_LAT_SCALE = re.compile(LAT_SCALE1 + '(.*?)' + LAT_SCALE2, re.S)
result_LAT_SCALE = pat_LAT_SCALE.findall(buff)
LAT_SCALE = result_LAT_SCALE[0]
pat_LONG_SCALE = re.compile(LONG_SCALE1 + '(.*?)' + LONG_SCALE2, re.S)
result_LONG_SCALE = pat_LONG_SCALE.findall(buff)
LONG_SCALE = result_LONG_SCALE[0]
pat_HEIGHT_SCALE = re.compile(HEIGHT_SCALE1 + '(.*?)' + HEIGHT_SCALE2, re.S)
result_HEIGHT_SCALE = pat_HEIGHT_SCALE.findall(buff)
HEIGHT_SCALE = result_HEIGHT_SCALE[0]
pat_LINE_NUM_COEFF = re.compile(LINE_NUM_COEFF1 + '(.*?)' + LINE_NUM_COEFF2, re.S)
result_LINE_NUM_COEFF = pat_LINE_NUM_COEFF.findall(buff)
LINE_NUM_COEFF = result_LINE_NUM_COEFF[0]
LINE_NUM_COEFF3 = LINE_NUM_COEFF[0]
LINE_NUM_COEFF3 = LINE_NUM_COEFF3.strip('()')
LINE_NUM_COEFF3 = LINE_NUM_COEFF3.replace('\n', '')
LINE_NUM_COEFF3 = LINE_NUM_COEFF3.replace('\t', '')
LINE_NUM_COEFF3 = LINE_NUM_COEFF3.replace(',', ' ')
pat_LINE_DEN_COEFF = re.compile(LINE_DEN_COEFF1 + '(.*?)' + LINE_DEN_COEFF2, re.S)
result_LINE_DEN_COEFF = pat_LINE_DEN_COEFF.findall(buff)
LINE_DEN_COEFF = result_LINE_DEN_COEFF[0]
LINE_DEN_COEFF3 = LINE_DEN_COEFF[0]
LINE_DEN_COEFF3 = LINE_DEN_COEFF3.strip('()')
LINE_DEN_COEFF3 = LINE_DEN_COEFF3.replace('\n', '')
LINE_DEN_COEFF3 = LINE_DEN_COEFF3.replace('\t', '')
LINE_DEN_COEFF3 = LINE_DEN_COEFF3.replace(',', ' ')
pat_SAMP_NUM_COEFF = re.compile(SAMP_NUM_COEFF1 + '(.*?)' + SAMP_NUM_COEFF2, re.S)
result_SAMP_NUM_COEFF = pat_SAMP_NUM_COEFF.findall(buff)
SAMP_NUM_COEFF = result_SAMP_NUM_COEFF[0]
SAMP_NUM_COEFF3 = SAMP_NUM_COEFF[0]
SAMP_NUM_COEFF3 = SAMP_NUM_COEFF3.strip('()')
SAMP_NUM_COEFF3 = SAMP_NUM_COEFF3.replace('\n', '')
SAMP_NUM_COEFF3 = SAMP_NUM_COEFF3.replace('\t', '')
SAMP_NUM_COEFF3 = SAMP_NUM_COEFF3.replace(',', ' ')
pat_SAMP_DEN_COEFF = re.compile(SAMP_DEN_COEFF1 + '(.*?)' + SAMP_DEN_COEFF2, re.S)
result_SAMP_DEN_COEFF = pat_SAMP_DEN_COEFF.findall(buff)
SAMP_DEN_COEFF = result_SAMP_DEN_COEFF[0]
SAMP_DEN_COEFF3 = SAMP_DEN_COEFF[0]
SAMP_DEN_COEFF3 = SAMP_DEN_COEFF3.strip('()')
SAMP_DEN_COEFF3 = SAMP_DEN_COEFF3.replace('\n', '')
SAMP_DEN_COEFF3 = SAMP_DEN_COEFF3.replace('\t', '')
SAMP_DEN_COEFF3 = SAMP_DEN_COEFF3.replace(',', ' ')
rpc = ['ERR_BIAS='+ERR_BIAS, 'ERR_RAND='+ERR_RAND, 'LINE_OFF='+LINE_OFF, 'SAMP_OFF='+SAMP_OFF, 'LAT_OFF='+LAT_OFF, 'LONG_OFF='+LONG_OFF, 'HEIGHT_OFF='+HEIGHT_OFF, 'LINE_SCALE='+LINE_SCALE, 'SAMP_SCALE='+SAMP_SCALE, 'LAT_SCALE='+LAT_SCALE, 'LONG_SCALE='+LONG_SCALE, 'HEIGHT_SCALE='+HEIGHT_SCALE,'LINE_NUM_COEFF='+LINE_NUM_COEFF3,'LINE_DEN_COEFF='+LINE_DEN_COEFF3,'SAMP_NUM_COEFF='+SAMP_NUM_COEFF3,'SAMP_DEN_COEFF='+SAMP_DEN_COEFF3]
#rpc = ['ERR BIAS=' + ERR_BIAS, 'ERR RAND=' + ERR_RAND, 'LINE OFF=' + LINE_OFF, 'SAMP OFF=' + SAMP_OFF,'LAT OFF=' + LAT_OFF, 'LONG OFF=' + LONG_OFF, 'HEIGHT OFF=' + HEIGHT_OFF, 'LINE SCALE=' + LINE_SCALE,'SAMP SCALE=' + SAMP_SCALE, 'LAT SCALE=' + LAT_SCALE, 'LONG SCALE=' + LONG_SCALE,'HEIGHT SCALE=' + HEIGHT_SCALE, 'LINE NUM COEFF=' + LINE_NUM_COEFF3, 'LINE DEN COEFF=' + LINE_DEN_COEFF3,'SAMP NUM COEFF=' + SAMP_NUM_COEFF3, 'SAMP DEN COEFF=' + SAMP_DEN_COEFF3]
return rpc
if __name__ == "__main__":
## 参数设置
demfile = r'C:\ENVI53\data\GMTED2010.jp2' # 使用ENVI的DEM文件
filepath = 'E:\data\GF3_MYN_SL_021996_E113.7_N22.0_20201013_L1A_HH_L10005135856' # 解压后的文件夹
rpbfile = filepath+'\GF3_MYN_SL_021996_E113.7_N22.0_20201013_L1A_HH_L10005135856.rpc' #对应的RPC文件
infile = filepath+'\GF3_MYN_SL_021996_E113.7_N22.0_20201013_L1A_HH_L10005135856.tiff'
outfile= os.path.splitext(infile)[0]+'_out.tiff'
os.chdir (filepath)
# strip函数是剔除文本中的'\n',split函数是分割规则
rpc = []
rpc=read_rpb(rpbfile)
# dataset = rasterio.open(file)
dataset = gdal.Open(infile)
dataset.SetMetadata(rpc, 'RPC')
dst_ds = gdal.Warp('outfile.tif', dataset, dstSRS='EPSG:4326', xRes=0.0002, yRes=0.0002,
rpc=True, transformerOptions=[r'RPC_DEM=' + demfile])
测试1
原始数据：
rpc校正后数据:
小结
雷达影像的相机坐标与光学影像相机坐标是不一样的。
光学卫星相机坐标是符合正常的拍摄逻辑，但是雷达卫星的相机坐标是上下颠倒、左右颠倒。通过无GCP点的PRC校正，可以初步将地理坐标赋予给雷达影像。
代码中的参数需要根据不同电脑重新设置，其他无需改动即可进行RPC校正。
## 参数设置
demfile = r'C:\ENVI53\data\GMTED2010.jp2' # 使用ENVI的DEM文件
filepath = 'E:\data\GF3_MYN_SL_021996_E113.7_N22.0_20201013_L1A_HH_L10005135856' # 解压后的文件夹
rpbfile = filepath+'\GF3_MYN_SL_021996_E113.7_N22.0_20201013_L1A_HH_L10005135856.rpc' #对应的RPC文件
infile = filepath+'\GF3_MYN_SL_021996_E113.7_N22.0_20201013_L1A_HH_L10005135856.tiff'
outfile= os.path.splitext(infile)[0]+'_out.tiff'
