CN102519433A - 一种利用rpc反演星载线阵传感器几何定标参数方法 - Google Patents

Abstract

本发明涉及星载CCD传感器几何定标参数反演方法，该方法根据RPC反求像元在卫星本体坐标系中的视线方向角和几何定标参数，包括以下步骤：(1)内插影像中间行摄影时刻卫星状态矢量和姿态值；(2)设定高程值，根据RFM计算中间行所有像点对应的地面点精确坐标；(3)计算像元在地球空间中的视线向量；(4)计算像元在本体坐标系中的视线向量，生成像元视线角文件；(5)计算镜头倾角、焦距和像元在像平面坐标系中的坐标，生成几何定标数据文件。本发明可根据星载线阵影像提供的RPC、星历和姿态参数反演传感器几何定标参数，用于影像的严密几何定位，在遥感影像摄影测量领域具有较高的应用价值。

Description

一种利用RPC反演星载线阵传感器几何定标参数方法

技术领域

本发明属于摄影测量领域，特别涉及星载线阵传感器几何定标参数提取和星载线阵影像严密几何处理技术。 

背景技术

目前，遥感影像的几何处理主要有严密模型和RFM模型。严格传感器模型是星载线阵遥感影像进行摄影测量的基础公式，是以传感器的内方位和外方位元素为定向参数的共线方程基础上建立的。在知道传感器几何定标参数或像元在传感器坐标系中的坐标的情况下，通过下述方式实现影像的对地定位。 

(1)传感器坐标系下像元视线方向 

像空间坐标系(传感器坐标系)以传感器投影中心为原点，x轴与飞行方向一致，z轴由像主点指向镜头中心，以右手规则建立y轴。设此坐标系中像点的坐标为： 

${\stackrel{\→}{u}}_c^{\′}={\left[\begin{array}{ccc}x-x_0& y-y_0& -f\end{array}\right]}^T=\left[\begin{array}{ccc}{x}^{\′}& y^{\′}& z^{\′}\end{array}\right]$

其中(x₀，y₀)为像主点坐标。故在像空间坐标系中，视线向量的单位矢量为： 

${\stackrel{\→}{u}}_s=\left[\begin{array}{ccc}{x}^{\′}& y^{\′}& z^{\′}\end{array}\right]/\sqrt{x^{\′2}+y^{\′2}+z^{\′2}}$

(2)本体坐标系中的视线方向 

本体坐标系的原点在卫星的质心，X、Y、Z轴为卫星的三主惯量轴。传感器坐标系的三轴与载体坐标系的三轴通常是平行的。由于安装误差，或根据需要人为设计而存在不平行，但本体坐标系到传感器坐标系总可以通过三个偏置角的旋转来实现两者的平行，与姿态角的定义一致，设此三偏置角依次为 

则传感器系到载体系坐标的转换为： 

因此，本体坐标系中像元的视线方向为： 

${\stackrel{\→}{u}}_b={R_c}^b{\stackrel{\→}{u}}_c$

(3)轨道坐标系中像元视线方向 

在惯性地心坐标系中，轨道坐标系的Z轴方向由地心指向卫星的位置，X轴与卫星的飞行方向一致，Y轴按右手规则确定。设本体坐标系相对于轨道坐标系的三姿态角为ω， κ，则本体坐标系到轨道坐标系的转换矩阵为： 

轨道系中像元的视线方向 

${\stackrel{\→}{u}}_o={R_b}^o{\stackrel{\→}{u}}_b$

(4)地心坐标系中像元视线方向 

上式中，R_c ^b的计算与前文类似，R_O ^E为地固轨道坐标系到地固地心直角坐标系的转换，可由地固地心坐标系的卫星状态矢量数据计算。设地固地心直角坐标系中摄影时刻卫星轨道位置和速度向量分别为 

地固轨道坐标系三轴的数学定义表达式： 

$\left\{\begin{array}{c}{\stackrel{\→}{Z}}_O=[{\left(Z_O\right)}_X,{\left(Z_O\right)}_Y,{\left(Z_O\right)}_Z]=\stackrel{\→}{P}\left(t\right)/\left|\right|\stackrel{\→}{P}\left(t\right)\left|\right|\\ {\stackrel{\→}{Y}}_O=[{\left(Y_O\right)}_X,{\left(Y_O\right)}_Y,{\left(Y_O\right)}_Z]=({\stackrel{\→}{Z}}_O\×{\stackrel{\→}{V}}_O)/\left|\right|({\stackrel{\→}{Z}}_O\×{\stackrel{\→}{V}}_O)\left|\right|\\ {\stackrel{\→}{X}}_O=[{\left(X_O\right)}_X,{\left(X_O\right)}_Y,{\left(X_O\right)}_Z]={\stackrel{\→}{Y}}_O\×{\stackrel{\→}{Z}}_O\end{array}\right.$

地固轨道坐标系到地固地心直角坐标系的转换矩阵为： 

${R_O}^E=\left[\begin{array}{ccc}{\left(X_O\right)}_X& {\left(X_O\right)}_Y& {\left(X_O\right)}_Z\\ {\left(Y_O\right)}_X& {\left(Y_O\right)}_Y& {\left(Y_O\right)}_Z\\ {\left(Z_O\right)}_X& {\left(Z_O\right)}_Y& {\left(Z_O\right)}_Z\end{array}\right],$

地心坐标系中，像元的视线向量可表示为： 

${\stackrel{\→}{u}}_E={R_o}^E{\stackrel{\→}{u}}_o$

卫星上有GPS接收机时，记录的数据为GPS天线所在的位置。 

(5)像点坐标对应地面点坐标的计算 

得到地心坐标系中的视线方向，就可以计算视线与椭球体上地面的交点。 

在知道地面高程值h的情况下，地球椭球模型 

$\frac{X^2+Y^2}{{(a+h)}^2}+\frac{Z^2}{{(b+h)}^2}=1$

地心坐标系中，像元的视线向量： 

$\stackrel{\→}{r}=\stackrel{\→}{R}+d{\stackrel{\→}{u}}_E$

假设交点M的地心坐标为(X，Y，Z)，代入上式有： 

$\left\{\begin{array}{c}X=X_P+\mathrm{\μ}\×{\left(u_3\right)}_X\\ Y=Y_P+\mathrm{\μ}\×{\left(u_3\right)}_Y\\ Z=Z_P+\mathrm{\μ}\×{\left(u_3\right)}_Z\end{array}\right.$

联合地球椭球方程和像元的视线向量方程，即可求解地面点在地心直角坐标系中的坐标(X，Y，Z)，实现像点坐标的对地几何定位。 

从定位过程可以看出，像元在传感器坐标系中的精确坐标或者像元的视线向量，是影像进行严格几何定位的基础。 

有些卫星影像厂商将传感器参数当作机密，没有向社会公开，但在影像的辅助文件中提供了影像的有理多项式模型系数(RPC)，根据RPC，可计算像点对应的地面点坐标，实现像点坐标的对地定位。 

根据像点坐标和RPC参数，给地物点高程H值时，计算相应的物方空间坐标， 

其算法主要为RFM的正解公式 

$\left.\begin{array}{c}{x}_n={\mathrm{RFM}}_X(P,L,H)=\frac{f_1(P_n,L_n,H_n)}{f_2(P_n,L_n,H_n)}\\ y_n={\mathrm{RFM}}_Y(P,L,H)=\frac{f_3(P_n,L_n,H_n)}{f_4(P_n,L_n,H_n)}\end{array}\right\}$

式中： 

x_n＝(x-x₀)/x_s，y_n＝(y-y₀)/y_s

$P_n=\frac{P-P_0}{P_S},$ $L_n=\frac{L-L_0}{L_S},$ $H_n=\frac{H-H_0}{H_S}$

f_i(P_n，L_n，H_n)＝c_i1+c_i2L_n+c_i3P_n+c_i4H_n+c_i5P_nL_n+c_i6L_nH_n+c_i1H_nP_n+c_i8L_n ²+c_i9P_n ²+c_i10H_n ²+c_i11P_nL_nH_n+c_i12L_n ³+c_i13P_n ²L_n+c_i14H_n ²L_n+c_i15L_n ²P_n+c_i16P_n ³+c_i17Z_n ²P_n+c_i18L_n ²H_n+c_i19P_n ²H_n+c_i20H_n ³

(P，L)为像点的像素坐标，(P_n，L_n)为归一化像素坐标，(P，L，H)为地面点坐标，(P_n，L_n，H_n)为归一化地面坐标，(x₀，y₀)(x_s，y_s)(P₀，L₀，H₀)(P_s，L_s，H_s)为归一化参数。f_i(P_n，L_n，H_n)为一般多项式，且多项式中每个坐标分量(P_n，L_n，H_n)的幂最大不超过3，每一坐标分量幂的总和也不超过3，提供RPC参数的影像，各多项式系数c_ij(i＝1，2，3，4；j＝1，2，...20)和归一化参数包含在RPC参数文件中。 

将RFM模型表达成： 

$\left\{\begin{array}{c}{x}_n{f}_2-f_1=0\\ y_n{f}_4-f_3=0\end{array}\right.$

建立误差方程： 

$\left\{\begin{array}{c}{V}_1=x_n{f}_2-f_1\\ V_2=y_n{f}_4-f_3\end{array}\right.$

当正则化像点坐标和地面点的大地高给定后，求解正则化地面点坐标[X_n，Y_n]，可将其作为未知数，并线性化： 

$\left\{\begin{array}{c}{V}_1=f_{1X}\mathrm{\δ}{P}_n+f_{1Y}\mathrm{\δ}{L}_n-L_1\\ V_2=f_{2X}\mathrm{\δ}{P}_n+f_{2Y}\mathrm{\δ}{L}_n-L_2\end{array}\right.$

通过最小二乘迭代求解，计算到精确的正则化像点地面坐标后，计算相应的真实地面点坐标。 

P＝P_nP_S+P₀，L＝L_nL_S+L₀，H＝H_nH_S+H₀

将地理坐标(P，L，H)根据坐标转换公式转换为地心直角坐标系中的坐标(X，Y，Z)，或高斯投影坐标系中的坐标(x，y，h)。 

RFM是严密模型的拟合，可以实现遥感影像的对地定位，但一些宽覆盖高分辨率CCD传感器，由多条线阵拼接成一条长CCD线阵，由于拼接和制造原因，像元的排列常是不规则的，复杂的结构会影响拼接后影像RFM模型的精度。有些卫星影像提供了各个拼接线阵影像 块的RPC来提高定位精度，但一景影像被分成了多个影像块独立处理影像了数据处理的效率。一景完整的影像进行统一处理，并不降低定位的精度，是本发明的能解决的问题，同时实现有RPC但无几何定标参数的传感器影像的严密几何模型处理。 

发明内容

1、一种星载线阵传感器像元视线向量重建和传感器几何定标参数提取方法，其特征是：根据RPC参数和星历参数重建本体坐标系中像元的视线向量和传感器几何定标数据，即采用以下步骤：(1)读入并提取影像辅助文件中的RPC、星历和姿态参数，根据星历参数内插影像中间行摄影时刻的卫星位置、速度和姿态数据；(2)设定地面高程，利用RPC依次计算中间行所有像点对应的地面点在地心直角坐标系中的坐标；(4)根据中间行所有地面点坐标和传感器位置依次计算每个像元在地心直角坐标系中的视线向量；(5)计算像元在本体坐标系中的视线向量，生成像元沿轨和垂轨方向的视线角并保存在文件中；(6)计算镜头倾角、焦距和所有像元在传感器像平面坐标系中的坐标，生成几何定标数据文件。 

进一步地，优选的方法是，步骤(1)中，中间行摄影时刻t的卫星位置和速度根据的星历数据利用拉格朗日插值函数插值得到，即： 

$\stackrel{\→}{P}\left(t\right)=\underset{j=1}{\overset{m}{\mathrm{\Σ}}}\left[\underset{i=1}{\overset{m}{\mathrm{\Π}}}\frac{\left({t-t}_j\right)}{(t_j-t_i)}\right]{\stackrel{\→}{P}}_j,i\≠j$

$\stackrel{\→}{V}\left(t\right)=\underset{j=1}{\overset{m}{\mathrm{\Σ}}}\left[\underset{i=1}{\overset{m}{\mathrm{\Π}}}\frac{(t-t_j)}{(t_j-t_i)}\right]{\stackrel{\→}{V}}_j,i\≠j$

式中t₁、t₂、t₃...为星历时刻， 为相应时刻的卫星位置， 

为相应时刻的速度，t为中间行影像摄影时刻、 

为t时刻的卫星位置(X_s，Y_s，Z_s)和速度(V_Xs，V_Ys，V_Zs)矢量。 

3、进一步地，优选的方法是，步骤(1)中，中间行摄影时刻t的卫星姿态角按线性内插的算法计算得到，即： 

$\mathrm{\α}\left(t\right)=\mathrm{\α}\left(t_i\right)+(\mathrm{\α}\left(t_{i+1}\right)-\mathrm{\α}\left(t_i\right))\×\frac{t-t_i}{t_{i+1}-t_i}$

其中：i代表离插值时间最近且小于插值时间的对应姿态序号(t_i≤t≤t_i+1)；α(t)代表t时刻的pitch(ω)角、 角和yaw(κ)角。 

进一步地，优选的方法是，步骤(2)中，地面点的高程H设定为RPC中高程正则化参数H₀，根据影像中间行每一列像点的正则化影像坐标(x_n，y_n)和给定高程的正则化值0，按公式 

$\left.\begin{array}{c}{x}_n={\mathrm{RFM}}_X(P,L,H)=\frac{f_1(P_n,L_n,H_n)}{f_2(P_n,L_n,H_n)}\\ y_n={\mathrm{RFM}}_Y(P,L,H)=\frac{f_3(P_n,L_n,H_n)}{f_4(P_n,L_n,H_n)}\end{array}\right\}$

通过解方程组计算正则化的地理坐标(P_n，L_n)，并根据RPC的正则化参数将(P_n，L_n，0)转化为地理坐标的正常值(P，L，H)，进而将正常值转化为地心直角坐标系坐标值(X，Y，Z)； 

进一步地，优选的方法是，步骤(3)中，所有像元在地心直角坐标系中的视线向量按公式[X-X_s Y-Y Z-Z_s]^T计算； 

进一步地，优选的方法是，步骤(4)中，，像元在本体坐标系中的视线向量计算公式为： 

$\left[\begin{array}{c}x\\ y\\ z\end{array}\right]=\mathrm{\λ}{R_O}^b{R_E}^O\left[\begin{array}{c}X-X_s\\ Y-Y_s\\ Z-Z_s\end{array}\right]$

式中，R_E ^O轨道坐标系到地固地心直角坐标系的转换矩阵，R_O ^b为轨道坐标系到本体坐标系的转换矩阵，即： 

${R_O}^E=\left[\begin{array}{ccc}{\left(X_O\right)}_X& {\left(X_O\right)}_Y& {\left(X_O\right)}_Z\\ {\left(Y_O\right)}_X& {\left(Y_O\right)}_Y& {\left(Y_O\right)}_Z\\ {\left(Z_O\right)}_X& {\left(Z_O\right)}_Y& {\left(Z_O\right)}_Z\end{array}\right],$

ω、 

κ为步骤(2)中内插获得的姿态角，R_O ^E矩阵元素根据步骤(2)中内插获得摄影时刻卫星轨道位置和速度向量为 

计算： 

$\left\{\begin{array}{c}{\stackrel{\→}{Z}}_O=[{\left(Z_O\right)}_X,{\left(Z_O\right)}_Y,{\left(Z_O\right)}_Z]=\stackrel{\→}{P}\left(t\right)/\left|\right|\stackrel{\→}{P}\left(t\right)\left|\right|\\ {\stackrel{\→}{Y}}_O=[{\left(Y_O\right)}_X,{\left(Y_O\right)}_Y,{\left(Y_O\right)}_Z]=({\stackrel{\→}{Z}}_O\×{\stackrel{\→}{V}}_O)/\left|\right|({\stackrel{\→}{Z}}_O\×{\stackrel{\→}{V}}_O)\left|\right|\\ {\stackrel{\→}{X}}_O=[{\left(X_O\right)}_X,{\left(X_O\right)}_Y,{\left(X_O\right)}_Z]={\stackrel{\→}{Y}}_O\×{\stackrel{\→}{Z}}_O\end{array}\right.$

进一步地，优选的方法是，步骤(5)中，，像元的视线角计算通过坐标向量[x，y，z]^T乘以系数1/z，将z值归化为1，即： 

[x′y′1]^T＝(1/z)[x y z]^T

并根据公式： 

ψ_Y＝a tan(x′) 

ψ_X＝-atan(y′) 

计算各像元在本体坐标系中沿轨和垂轨方向的视线向量角ψ_Y，ψ_X，并记录在文件中。 

进一步地，优选的方法是，步骤(6)中，，计算传感器几何定标数据时所有像元在传感器坐标系中的视线向量[x″y″z″]为： 

$\left[\begin{array}{c}{x}^{\′\′}\\ y^{\′\′}\\ z^{\′\′}\end{array}\right]=-{R_b}^c\left[\begin{array}{c}{x}^{\′}\\ y^{\′}\\ 1\end{array}\right]$

R_b ^c为本体坐标系到传感器坐标系的转换矩阵，即： 

R_b ^c＝R^T _Y(θ) 

θ为传感器镜头的前视或后视的偏置角，偏置角的值按传感器提供的理论值代入，在没有提供的情况下，对单镜头线阵传感器，θ取值0，对多镜头多线阵传感器，值按下式计算： 

$\mathrm{\θ}=-a\tan \frac{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x^{\′}}_i}{n}$

将所有像元的向量元素乘以系数g，获得像元在传感器坐标系中的坐标[x”’，y”’，f]，即：[x″′y″′f]＝g·[x″y″z″]＝[g·x″g·y″g·z″]系数g根据第一个和最后一个像点的坐标(x″₁，y″₁)、(x″_n，y″_n)、像元大小μ＝6.5×10^-5按公式 

$g=\mathrm{n\μ}/\sqrt{{({x^{\′\′}}_n-{x^{\′\′}}_1)}^2+{({y^{\′\′}}_n-{y^{\′\′}}_1)}^2}$ 计算。 

记录镜头倾角θ、镜头焦距f、以及每个像元的传感器坐标系中的坐标[x″′y″′]，生成摄影测量几何定标文件。 

本发明基于RPC参数和星历参数实现像元视线向量重建和等效传感器几何定标参数提取方法，根据计算获取的像元视线向量参数或传感器前后侧视角θ、像元坐标可正常应用于星载CCD传感器遥感影像的严密几何处理中，解决了一些高分辨率星载CCD影像缺乏内定向参数而不能进行严密定位的缺陷。对于一些多条线阵拼接成的长线阵CCD相机(如ALOS/PRISM传感器)，提供了每条CCD影像块的RPC参数和拼接后影像的RPC参数，多个RPC参数的综合应用会提高影像定位的进度，但提高了复杂性。这种情况下本发明通过多组RPC参数提取一套传感器定向参数，不仅能够实现严密定位，还能实现与多组RPC参数进行定位一致的精度，简化了定位，RPC与严格模型转换精度能优于千分之一像素。 

附图说明

下面结合附图对本发明进行进一步详细的描述，以使得本发明的上述优点更加明确。 

图1是本发明所述的基于RPC的星载CCD传感器内定向参数重建方法的流程图。 

图2(a)为ALOS/PRISM正下视传感器各条CCD获取的影像块；图2(b)为拼接后的影像；图2(c)为根据本发明重建后的像元在传感器坐标系像平面上的排列示意图，分别按影像辅助文件提供的各线阵影像块的RPC数据(CCD1-CCD4)和拼接后的影像RPC数据(CCD1-4)计算得到；图2(d)中为澳大利亚墨尔本大学以大量的高精度地面控制点对传感器进行检校得到的像元排列示意图。 

具体实施方式

下面结合附图1对本发明作进一步详细的描述。 

该种星载CCD传感器像可定标参数反演方法，依靠RPC参数和星历参数，计算传感 器每个像元在卫星本体坐标系中的视线方向角，或像元在传感器坐标系中的坐标，具体实施可通过计算机完成。 

步骤1，根据影像辅助文件说明提取该景影像的RPC参数、星历数据和姿态数据，内插出中间行影像摄影时刻卫星的空间位置、速度和姿态。 

一般来说，我们是按照下述方法来进行插值的： 

中间行摄影时刻t的卫星位置和速度根据的星历数据利用拉格朗日插值函数插值得到，即： 

$\stackrel{\→}{P}\left(t\right)=\underset{j=1}{\overset{m}{\mathrm{\Σ}}}\left[\underset{i=1}{\overset{m}{\mathrm{\Π}}}\frac{\left({t-t}_j\right)}{(t_j-t_i)}\right]{\stackrel{\→}{P}}_j,i\≠j$

$\stackrel{\→}{V}\left(t\right)=\underset{j=1}{\overset{m}{\mathrm{\Σ}}}\left[\underset{i=1}{\overset{m}{\mathrm{\Π}}}\frac{(t-t_j)}{(t_j-t_i)}\right]{\stackrel{\→}{V}}_j,i\≠j$

式中t₁、t₂、t₃...为星历时刻， 

为相应时刻的卫星位置， 

为相应时刻的速度，t为中间行影像摄影时刻、 

为t时刻的卫星位置(X_s，Y_s，Z_s)和速度(V_Xs，V_Ys，V_Zs)。 

中间行摄影时刻t的卫星姿态角按线性内插的算法计算得到，即： 

$\mathrm{\α}\left(t\right)=\mathrm{\α}\left(t_i\right)+(\mathrm{\α}\left(t_{i+1}\right)-\mathrm{\α}\left(t_i\right))\×\frac{t-t_i}{t_{i+1}-t_i}$

其中：i代表离插值时间最近且小于插值时间的对应姿态序号(t_i≤t≤t_i+1)；α(t)代表t时刻的pitch(ω)角、 

角和yaw(κ)角。 

步骤2，根据像点坐标和RPC参数，给定物方坐标H值为RPC参数中正则化参数H₀，根据影像中间行每一列像点的影像坐标(x，y)和给定的地面高程H，按下述步骤计算相应像点在地心直角坐标系中的坐标。 

2.1)根据公式 

x_n＝(x-x₀)/x_s

y_n＝(y-y₀)/y_s

H_n＝0 

计算正则化像点坐标和高程坐标。 

2.2)根据正则化像点坐标和高程坐标求解正则化地理坐标(P_n，L_n) 

其算法主要为RFM的正解公式 

$\left.\begin{array}{c}{x}_n={\mathrm{RFM}}_X(P,L,H)=\frac{f_1(P_n,L_n,H_n)}{f_2(P_n,L_n,H_n)}\\ y_n={\mathrm{RFM}}_Y(P,L,H)=\frac{f_3(P_n,L_n,H_n)}{f_4(P_n,L_n,H_n)}\end{array}\right\}$

式中： 

f_i(P_n，L_n，H_n)＝c_i1+c_i2L_n+c_i3P_n+c_i4H_n+c_i5P_nL_n+c_i6L_nH_n+c_i1H_nP_n+c_i8L_n ²+c_i9P_n ²+c_i10H_n ²+c_i11P_nL_nH_n+c_i12L_n ³+c_i13P_n ²L_n+c_i14H_n ²L_n+c_i15L_n ²P_n+c_i16P_n ³+c_i17Z_n ²P_n+c_i18L_n ²H_n+c_i19P_n ²H_n+c_i20H_n ³C_ij为RPC参数。 

将RFM模型表达成： 

$\left\{\begin{array}{c}{x}_n{f}_2-f_1=0\\ y_n{f}_4-f_3=0\end{array}\right.$

建立误差方程： 

$\left\{\begin{array}{c}{V}_1=x_n{f}_2-f_1\\ V_2=y_n{f}_4-f_3\end{array}\right.$

当正则化像点坐标和地面点的大地高给定后，求解正则化地面点坐标x_n，y_n，可将其作为未知数，并线性化： 

$\left\{\begin{array}{c}{V}_1=f_{1X}\mathrm{\δ}{P}_n+f_{1Y}\mathrm{\δ}{L}_n-L_1\\ V_2=f_{2X}\mathrm{\δ}{P}_n+f_{2Y}\mathrm{\δ}{L}_n-L_2\end{array}\right.$

通过解方程组，得到像点的正则化地面坐标(P_n，L_n). 

步骤2.3)计算相应的真实地面点坐标(P，L，H)。 

P＝P_nP_S+P₀，L＝L_nL_S+L₀，H＝H_nH_S+H₀

步骤2.4)将地理坐标(P，L，H)根据坐标转换公式转换为地心直角坐标系中的坐标(X，Y，Z)。 

步骤3)计算像元在物方空间坐标系中的视线向量为[X-X_s Y-Ys Z-Zs]^T； 

步骤4)计算像元在传感器坐标系中的视线向量，按下述步骤实现 

4.1)计算轨道坐标系到本体坐标系的转换矩阵R_o ^b

上式中，ω、 

κ角为三姿态角。 

4.2)计算地心直角坐标系到轨道系的转换矩阵，R_E ^b该矩阵由摄影时刻的卫星状态矢量计算： 

${R_O}^E=\left[\begin{array}{ccc}{\left(X_O\right)}_X& {\left(X_O\right)}_Y& {\left(X_O\right)}_Z\\ {\left(Y_O\right)}_X& {\left(Y_O\right)}_Y& {\left(Y_O\right)}_Z\\ {\left(Z_O\right)}_X& {\left(Z_O\right)}_Y& {\left(Z_O\right)}_Z\end{array}\right]$

4.3)计算像元在本体坐标系中的视线向量计算公式为： 

$\left[\begin{array}{c}x\\ y\\ z\end{array}\right]=\mathrm{\λ}{R_O}^b{R_E}^O\left[\begin{array}{c}{X}_n-X_s\\ Y_n-Y_s\\ Z_n-Z_s\end{array}\right]$

4.4)对[x，y，z]^T同时乘以-1/z，将z值归化为-1，即： 

[x′y′1]^T＝(1/z)[x y z]^T

4.5)根据公式 

ψ_Y＝a tan(x′) 

ψ_X＝-atan(y′) 

计算像元在本体坐标系中的视线向量。 

步骤5，按下述步骤计算传感器几何定标数据： 

步骤5.1)将所有像元视线向量转换到传感器坐标系中，转换公式为： 

$\left[\begin{array}{c}{x}^{\′\′}\\ y^{\′\′}\\ z^{\′\′}\end{array}\right]=-{R_b}^c\left[\begin{array}{c}{x}^{\′}\\ y^{\′}\\ 1\end{array}\right]$

R_b ^c为本体坐标系到传感器坐标系的转换矩阵，即 

R_b ^c＝R^T _Y(θ) 

θ为传感器镜头的前视或后视的偏置角，偏置角的值按提供的值代入，如果没有提供，对单镜头线阵传感器，θ取值0，对应多镜头多线阵传感器，按下式计算： 

$\mathrm{\θ}=-a\tan \frac{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x^{\′}}_i}{n}$

步骤5.2)按下式计算传感器真实空间归化系数g： 

$g=\mathrm{n\μ}/\sqrt{{({x^{\′\′}}_n-{x^{\′\′}}_1)}^2+{({y^{\′\′}}_n-{y^{\′\′}}_1)}^2}$

(x″₁，y″₁)、(x″_n，y″_n)为第一个和最后一个像元的视线向量值，μ为传感器像元大小，一般为μ＝6.5×10^-5。 

步骤5.3)计算传感器坐标系中的像元坐标，用系数g乘以步骤5.1中得到的坐标，即g·[x″y″z″]＝[g·x″g·y″g·z″]。其中，[g·x″g·y″]为像元在像平面坐标系中的坐标，f＝g·z″即为传感器的焦距。 

对所有像元执行上述操作，得到每个像元的传感器空间坐标系中的坐标，生成摄影测量几何定标文件。 

提取的星载CCD传感器内定向参数，能正常应用于相应传感器影像的严密几何处理中，在有些情况下可以获得不低于RFM模型的几何处理精度。 

以日本ALOS卫星上的PRISM传感器为例，其为大幅宽高分辨率三镜头三线阵传感器，每镜头的焦平面上安装了多条线阵拼接而成的一组长线阵，影像由各条线阵独立成像拼接而成。图2(a)为ALOS/PRISM正下视传感器各条CCD获取的影像块；图2(b)为拼接后的影像；图2(c)为根据本发明重建后的像元在传感器坐标系像平面上的排列示意图，分别按影像辅助 文件提供的各线阵影像块的RPC数据(CCD1-CCD4)和拼接后的影像RPC数据(CCD1-4)计算得到；图2(d)中为澳大利亚墨尔本大学以大量的高精度地面控制点对传感器进行检校得到的像元排列示意图。从图2(c)图形上看，利用本发明方法反演获取的CCD1-CCD4曲线与图1(d)是一致的。 

Claims (8)

1.一种星载线阵(CCD)传感器几何定标参数反演方法，其特征是：根据RPC参数和星历参数重建本体坐标系中像元的视线向量，反演传感器几何定标参数，即采用以下步骤：(1)读入并提取影像辅助文件中的RPC、星历和姿态参数，根据星历和姿态参数内插影像中间行摄影时刻的卫星位置、速度和姿态值；(2)设定地面高程，利用RFM模型和RPC参数依次计算中间行所有像点对应的地面点在地心直角坐标系中的坐标；(3)根据中间行所有地面点坐标和传感器位置依次计算每个像元在地心直角坐标系中的视线向量；(4)计算像元在本体坐标系中的视线向量，生成像元沿轨和垂轨方向的视线角并保存在几何定标文件中；(5)计算镜头倾角、焦距和所有像元在传感器像平面坐标系中的坐标，生成几何定标数据并保存在几何定标文件中。

2.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，所述步骤(1)中间行摄影时刻t的卫星位置和速度

根据的星历数据利用拉格朗日插值函数插值得到，t时刻的姿态角pitch(ω)、

和yaw(κ)角按线性内插的算法计算得到。

3.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，所述步骤(2)中，地面点的高程H设定为RPC中高程正则化参数H₀，依次按下式计算影像中间行所有不同列像点的正则化影像坐标(x_n，y_n)和正则化高程值H₀＝0，利用RFM模型

$\left.\begin{array}{c}{x}_n={\mathrm{RFM}}_X(P,L,H)=\frac{f_1(P_n,L_n,H_n)}{f_2(P_n,L_n,H_n)}\\ y_n={\mathrm{RFM}}_Y(P,L,H)=\frac{f_3(P_n,L_n,H_n)}{f_4(P_n,L_n,H_n)}\end{array}\right\}$

通过解方程组计算正则化的地理坐标(P_n，L_n)，并根据RPC的正则化参数将正则化坐标(P_n，L_n，0)转化为地理坐标的正常值(P，L，H)，进而将正常值转化为地心直角坐标系坐标值(X，Y，Z)。

4.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，所述步骤(3)中，所有像元在地心直角坐标系中的视线向量依次按公式[X-X_s Y-Y Z-Z_s]^T计算。

5.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，所述步骤(4)中，所有像元在本体坐标系中视线向量依次按下式计算：

$\left[\begin{array}{c}x\\ y\\ z\end{array}\right]=\mathrm{\λ}{R_O}^b{R_E}^O\left[\begin{array}{c}X-X_s\\ Y-Y_s\\ Z-Z_s\end{array}\right]$

式中，R_E ^O轨道坐标系到地固地心直角坐标系的转换矩阵，R_O ^b为轨道坐标系到本体坐标系的转换矩阵，即：

${R_O}^E=\left[\begin{array}{ccc}{\left(X_O\right)}_X& {\left(X_O\right)}_Y& {\left(X_O\right)}_Z\\ {\left(Y_O\right)}_X& {\left(Y_O\right)}_Y& {\left(Y_O\right)}_Z\\ {\left(Z_O\right)}_X& {\left(Z_O\right)}_Y& {\left(Z_O\right)}_Z\end{array}\right],$

ω、

κ为步骤(1)中内插获得的姿态角，R_O ^E矩阵元素根据步骤(1)中内插获得的摄影时刻卫星位置和速度向量为

按下式计算：

$\left\{\begin{array}{c}{\stackrel{\→}{Z}}_O=[{\left(Z_O\right)}_X,{\left(Z_O\right)}_Y,{\left(Z_O\right)}_Z]=\stackrel{\→}{P}\left(t\right)/\left|\right|\stackrel{\→}{P}\left(t\right)\left|\right|\\ {\stackrel{\→}{Y}}_O=[{\left(Y_O\right)}_X,{\left(Y_O\right)}_Y,{\left(Y_O\right)}_Z]=({\stackrel{\→}{Z}}_O\×{\stackrel{\→}{V}}_O)/\left|\right|({\stackrel{\→}{Z}}_O\×{\stackrel{\→}{V}}_O)\left|\right|\\ {\stackrel{\→}{X}}_O=[{\left(X_O\right)}_X,{\left(X_O\right)}_Y,{\left(X_O\right)}_Z]={\stackrel{\→}{Y}}_O\×{\stackrel{\→}{Z}}_O\end{array}\right.$

6.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，所述步骤(5)中，像元的视线角计算通过向量[x，y，z]^T乘以系数1/z，将z值归化为1，即：[x′y′1]^T＝(1/z)[x y z]^T

并根据公式：

ψ_Y＝atan(x′)

ψ_X＝-atan(y′)

计算得到各像元在本体坐标系中沿轨和垂轨方向的视线向量角ψ_Y，ψ_X，并记录在定标文件中。

7.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，所述步骤(6)中，计算传感器几何定标数据时所有像元在传感器坐标系中的坐标，首先按下式计算传感器坐标系中的视线向量[x″y″z″]：

$\left[\begin{array}{c}{x}^{\′\′}\\ y^{\′\′}\\ z^{\′\′}\end{array}\right]=-{R_b}^c\left[\begin{array}{c}{x}^{\′}\\ y^{\′}\\ 1\end{array}\right]$

R_b ^c为本体坐标系到传感器坐标系的转换矩阵，按下式计算：

R_b ^c＝R^T _Y(θ)

θ为传感器镜头的前视或后视的偏置角，偏置角的值按传感器提供的理论值代入，在没有提供的情况下，对单镜头线阵传感器，θ取值0，对多镜头多线阵传感器，值按下式计算：

$\mathrm{\θ}=-a\tan \frac{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{x^{\′}}_i}{n}$

将所有像元的视线向量元素乘以系数g，获得像元在传感器坐标系中的坐标[x”’，y”’，f]，即：[x″′y″′f]＝g·[x″y″z″]＝[g·x″g·y″g·z″]系数g根据第一个和最后一个像点的坐标(x″₁，y″₁)、(x″_n，y″_n)、像元大小μ，按公式 $g=\mathrm{n\μ}/\sqrt{{({x^{\′\′}}_n-{x^{\′\′}}_1)}^2+{({y^{\′\′}}_n-{y^{\′\′}}_1)}^2}$ 计算。

记录镜头倾角θ、镜头焦距f、以及每个像元在传感器坐标系中的坐标[x″′y″′]，生成摄影测量几何定标文件。

8.根据权利1关于星载CCD传感器几何定标参数反演方法，其特征在于，根据计算获取的像元视线向量参数或传感器前后侧视角θ、焦距f以及像元在传感器坐标系中的坐标[x″′y″′]，可正常应用于星载CCD传感器遥感影像的摄影测量严密几何处理中。

Images