CN103759714A - 一种三线阵卫星影像区域网平差方法 - Google Patents

Abstract

本发明涉及一种三线阵卫星影像区域网平差方法，包括以下步骤：1）读取卫星影像，并将各类点的坐标信息存成相应的点文件；2）读取各文件，得到卫星影像的RPC模型以及控制点和检查点的地面坐标、影像坐标和连接点的像点坐标；3）确定地面点坐标和像方变换模型的初值；4）逐点构建误差方程并进行法化；5）对法方程利用谱修正法迭代估计求解改正参数和地面点坐标；6）计算此次平差迭代后所能达到的物方精度；7）计算此次平差迭代后所能达到的像方精度，通过像方残差并利用信息扩散模型确定下一次迭代像点坐标观测值的权值P；8）当满足定向参数中的平移参数均小于阈值时时，平差迭代结束；当不满足预设条件时继续迭代，直至迭代结束；9）当整个平差迭代结束时，输出最平差物方精度和像方精度的精度报告。

Description

一种三线阵卫星影像区域网平差方法

技术领域

本发明涉及摄影测量技术领域，特别是关于一种三线阵卫星影像区域网平差方法。

背景技术

高分辨率推扫式传感器可以获取高分辨率遥感影像，当前高分辨率卫星遥感影像基本都是基于这种特殊的成像方式来获取。进入21世纪以来，越来越多高分辨率商业卫星的成功发射，为广大用户提供了多种数据源的选择。但是由于各个卫星传感器设计的不同，使得获取的卫星影像的实用性受到了较大的限制。因此，为了建立影像坐标与地面坐标之间严密一一对应关系，需要建立影像的成像几何模型。目前常用的成像几何模型主要分为两种：严密成像几何模型和通用成像几何模型。

由于严密成像几何模型完整地反映了从像方坐标系到物方坐标系之间一系列严密的转换，一直是摄影测量相关领域专家学者的首选。然而，一些星载传感器系统虽然实现了商业化，但所涉及的卫星平台和有效载荷等关键技术参数仍然是航天大国的技术机密，高分辨率传感器的核心信息和卫星轨道参数并未公开，如IKONOS等卫星。因此，在缺少这些重要模型参数时，是无法利用严密成像几何模型进行卫星影像的几何处理；同时，传感器成像方式的多样化也对摄影测量软件提出了模型通用化的要求，使得通用成像几何模型得到越来越多的关注与应用。这样在软件框架设计时就不必为现有的或将来可能会出现的各种类型的传感器一一确定其严格几何模型的形式，而是可以统一采用通用成像几何模型进行处理，从而大大降低了程序设计的复杂性，更易于软件升级和维护，特别在同时处理多源卫星遥感影像数据时，通用成像几何模型更能显示出其优势。

常见的通用成像几何模型有以下几种：多项式、直接线性变换、仿射变换、有理多项式模型（RPC）等。RPC模型是卫星遥感影像通用成像几何模型的一种，它适用于各类传感器包括最新的航空和航天传感器。建立RPC模型并不要求了解传感器的实际特性和成像过程，是一种能高精度拟合严密成像模型的形式简单的广义成像模型。基于RPC模型具有模拟精度高，通用性好，应用方便，计算量小等等优点，但是，该模型同时存在一大缺点，即参数没有严格的几何意义，在进行平差处理时，无法建立起对应具有物理意义的严密模型。

发明内容

针对上述问题，本发明的目的是提供一种能够提高定位精度，而且能够在较少控制点甚至无控制点的情况下进行三线阵卫星影像区域网平差的方法。

为实现上述目的，本发明采取以下技术方案：一种三线阵卫星影像区域网平差方法，其包括以下步骤：1）读取卫星CCD阵列影像，选择相邻卫星影像重叠区域的同名点作为影像间的连接点，根据点之记或者控制影像在卫星影像上刺出控制点和检查点，并将连接点、控制点和检查点的坐标信息存成相应的点文件；2）读取卫星CCD阵列影像的rcp参数文件、连接点文件、控制点文件和检查点文件，得到卫星CCD阵列影像的RPC模型以及控制点和检查点的地面坐标、影像坐标和连接点的像点坐标；3）确定地面点坐标和像方变换模型的初值；4）利用量测得到像点坐标和对应的地面点坐标，针对连接点和控制点分别对定向参数和地面点坐标求偏导并逐点构建误差方程并将构建的误差方程进行法化；5）对法方程利用谱修正法迭代估计求解改正参数和地面点坐标；6）每完成一次平差迭代时，计算出检查点对应的地面点坐标，并同时计算此次平差迭代后所能达到的物方精度；7）每完成一次平差迭代时，通过计算的检查点对应的地面点坐标反投影计算出像点坐标，并计算此次平差迭代后所能达到的像方精度，通过像方残差并利用信息扩散模型确定下一次迭代像点坐标观测值的权值P；8）当满足定向参数中的平移参数均小于阈值时，整个平差迭代结束；当不满足预设条件时，返回步骤5）继续迭代计算，直到满足迭代收敛条件；如果迭代次数达到预设迭代次数，仍然不能收敛，那么迭代结束；9）当整个平差迭代结束时，输出最后一次迭代时所计算得到的平差物方精度和像方精度的精度报告，根据精度报告对平差效果进行评价。

步骤4）利用量测得到像点坐标和对应的地面点坐标，针对连接点和控制点分别对定向参数和地面点坐标求偏导并逐点构建误差方程并将构建的误差方程进行法化，具体过程为：在模型的基础上建立仿射变换模型：

$\left.\begin{array}{c}{F}_y=\mathrm{\Δy}=a_1+a_2\·\mathrm{sample}+a_3\·\mathrm{line}\\ F_x=\mathrm{\Δx}=b_1+b_2\·\mathrm{sample}+b_3\·\mathrm{line}\end{array}\right\}$

式中，Δy和Δx为控制点在影像坐标系中的量测坐标与真实坐标的差值，即改正数；a₁,a₂,a₃和b₁,b₂,b₃是影像的定向参数，line和sample是控制点在影像坐标系中的行列号；建立误差方程：

$v_x=\left(\begin{array}{c}\frac{\∂F_x}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_x}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_x}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_x}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_x}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_x}{\∂b_3}\·\mathrm{\Δ}{b}_3\\ +\frac{\∂F_x}{\∂D_{\mathrm{lat}}}\·\mathrm{\Δ}{D}_{\mathrm{lat}}+\frac{\∂F_x}{\∂D_{\mathrm{lon}}}\·\mathrm{\Δ}{D}_{\mathrm{lon}}+\frac{\∂F_x}{\∂D_{\mathrm{hei}}}\·\mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)+F_{x0}$

$v_y=\left(\begin{array}{c}\frac{\∂F_y}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_y}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_y}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_y}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_y}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_y}{\∂b_3}\·\mathrm{\Δ}{b}_3\\ +\frac{\∂F_y}{\∂D_{\mathrm{lat}}}\·\mathrm{\Δ}{D}_{\mathrm{lat}}+\frac{\∂F_y}{\∂D_{\mathrm{lon}}}\·\mathrm{\Δ}{D}_{\mathrm{lon}}+\frac{\∂F_y}{\∂D_{\mathrm{hei}}}\·\mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)+F_{y0}$

式中，ΔD_lat、ΔD_lon、ΔD_hei为地面点坐标改正数；Δa₁，Δa₂，Δa₃，Δb₁，Δb₂，Δb₃为影像定向参数改正数；v_x,v_y为像点坐标改正数；F_x0,F_y0为像点坐标近似值与像点坐标观测值之差；

（i=1、2、3）为误差方程对定向参数所求的偏导数；

为误差方程对地面点坐标所求的偏导数；

记为：

V＝Bt+AX-l  P

同样可以对每个控制点列如下线性方程：

$v_x=(\frac{\∂F_x}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_x}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_x}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_x}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_x}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_x}{\∂b_3}\·\mathrm{\Δ}{b}_3)+F_{x0}$

$v_y=(\frac{\∂F_y}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_y}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_y}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_y}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_y}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_y}{\∂b_3}\·\mathrm{\Δ}{b}_3)+F_{y0}$

记为：

V＝Bt-l  P

其中：

$B=\left(\begin{array}{cccccc}\frac{\∂F_x}{\∂a_1}& \frac{\∂F_x}{\∂a_2}& \frac{\∂F_x}{\∂a_3}& \frac{\∂F_x}{\∂b_1}& \frac{\∂F_x}{\∂b_2}& \frac{\∂F_x}{\∂b_3}\\ \frac{\∂F_y}{\∂a_1}& \frac{\∂F_y}{\∂a_2}& \frac{\∂F_y}{\∂a_3}& \frac{\∂F_y}{\∂b_1}& \frac{\∂F_y}{\∂b_2}& \frac{\∂F_y}{\∂b_3}\end{array}\right),t={\left(\begin{array}{cccccc}\mathrm{\Δ}{a}_1& \mathrm{\Δ}{a}_2& \mathrm{\Δ}{a}_3& \mathrm{\Δ}{b}_1& \mathrm{\Δ}{b}_2& \mathrm{\Δ}{b}_3\end{array}\right)}^T$

$A=\left(\begin{array}{ccc}\frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{lat}}}& \frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{lon}}}& \frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{hei}}}\\ \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{lat}}}& \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{lon}}}& \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{hei}}}\end{array}\right),X={\left(\begin{array}{ccc}\mathrm{\Δ}{D}_{\mathrm{lat}}& {\mathrm{\ΔD}}_{\mathrm{lon}}& \mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)}^T$

$l=\left(\begin{array}{c}-F_{x0}\\ -F_{y0}\end{array}\right),V=\left(\begin{array}{c}{v}_x\\ v_y\end{array}\right)$

根据最小二乘平差构建法方程：

$\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)\left(\begin{array}{c}t\\ X\end{array}\right)=\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right).$

步骤5）对法方程利用谱修正法迭代估计求解改正参数和地面点坐标，具体过程为：将法方程两边同时加上 $\left(\begin{array}{c}t\\ X\end{array}\right),$ 得到：

$[\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)+E]\left(\begin{array}{c}t\\ X\end{array}\right)=\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right)+\left(\begin{array}{c}t\\ X\end{array}\right)$

式中，E为和 $\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)$ 阶数相同的单位矩阵，由于式子两边都含有未知参数 $\left(\begin{array}{c}t\\ X\end{array}\right),$ 采用迭代的方法求解，其迭代公式为：

${\left(\begin{array}{c}t\\ X\end{array}\right)}^{\left(k\right)}={[\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)+E]}^{-1}+[\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right)+{\left(\begin{array}{c}t\\ X\end{array}\right)}^{(k-1)}]$

式中，k为迭代次数。

本发明由于采取以上技术方案，其具有以下优点：1、本发明采用的是通用传感器模型作为区域网平差的处理模型，形式简单，使用方便。因此，本发明可以有效地对卫星相关载荷、姿态、轨道、卫星相机等参数进行保密，同时让缺少专业知识的操作者可以方便的进行数据处理。2、本发明采用的是三线阵的区域网平差模式，较之传统的卫星影像平差方法而言，平差模型之间的立体约束能力更强，可以获得更好的平面和高程精度。本发明可以广泛应用于摄影测量的平差计算中。

附图说明

以下结合附图来对本发明进行详细的描绘。然而应当理解，附图的提供仅为了更好地理解本发明，它们不应该理解成对本发明的限制。

图1是本发明的方法流程示意图。

具体实施方式

下面结合附图和实施例对本发明进行详细的描述。

如图1所示，本发明的三线阵卫星影像区域网平差方法针对CCD线阵推扫式光学卫星影像，利用稀（缺）控制点进行区域网平差，在求解像方改正模型的同时，求解出所有待定点（检查点）的地面点坐标，其包括以下步骤：

1）读取卫星CCD阵列影像，利用匹配技术或者人工采集等方式选择相邻卫星影像重叠区域的同名点作为影像间的连接点，根据点之记或者控制影像在卫星影像上刺出控制点和检查点，并将连接点、控制点和检查点的坐标信息存成相应的点文件；其中，控制点布设的原则应当尽量能够分布在测区的四角。

2）读取卫星CCD阵列影像的rcp参数文件、连接点文件、控制点文件和检查点文件，得到卫星CCD阵列影像的RPC模型以及控制点和检查点的地面坐标、影像坐标和连接点的像点坐标。

3）确定地面点坐标和像方变换模型的初值，即：所有地面点的坐标可以通过直接前方交会，得到地面点坐标(lat,lon,h)作为地面点初值，像方变换模型初值可以设定为 $\left(\begin{array}{ccc}0& 1& 0\\ 0& 0& 1\end{array}\right);$ 其中，lat表示纬度，lon表示经度，h表示高程。

4）利用量测得到像点坐标和对应的地面点坐标，针对连接点和控制点分别对定向参数和地面点坐标求偏导并逐点构建误差方程并将构建的误差方程进行法化，具体过程为：

在上述RPC模型的基础上建立仿射变换模型：

$\left.\begin{array}{c}{F}_y=\mathrm{\Δy}=a_1+a_2\·\mathrm{sample}+a_3\·\mathrm{line}\\ F_x=\mathrm{\Δx}=b_1+b_2\·\mathrm{sample}+b_3\·\mathrm{line}\end{array}\right\}$

式中，Δy和Δx为控制点在影像坐标系中的量测坐标与真实坐标的差值，即改正数；a₁,a₂,a₃和b₁,b₂,b₃是影像的定向参数，line和sample是控制点在影像坐标系中的行列号。

在此基础上，建立误差方程：

$v_x=\left(\begin{array}{c}\frac{\∂F_x}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_x}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_x}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_x}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_x}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_x}{\∂b_3}\·\mathrm{\Δ}{b}_3\\ +\frac{\∂F_x}{\∂D_{\mathrm{lat}}}\·\mathrm{\Δ}{D}_{\mathrm{lat}}+\frac{\∂F_x}{\∂D_{\mathrm{lon}}}\·\mathrm{\Δ}{D}_{\mathrm{lon}}+\frac{\∂F_x}{\∂D_{\mathrm{hei}}}\·\mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)+F_{x0}$

$v_y=\left(\begin{array}{c}\frac{\∂F_y}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_y}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_y}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_y}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_y}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_y}{\∂b_3}\·\mathrm{\Δ}{b}_3\\ +\frac{\∂F_y}{\∂D_{\mathrm{lat}}}\·\mathrm{\Δ}{D}_{\mathrm{lat}}+\frac{\∂F_y}{\∂D_{\mathrm{lon}}}\·\mathrm{\Δ}{D}_{\mathrm{lon}}+\frac{\∂F_y}{\∂D_{\mathrm{hei}}}\·\mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)+F_{y0}$

式中，ΔD_lat、ΔD_lon、ΔD_hei为地面点坐标改正数；Δa₁，Δa₂，Δa₃，Δb₁，Δb₂，Δb₃为影像定向参数改正数；v_x,v_y为像点坐标改正数；F_x0,F_y0为像点坐标近似值与像点坐标观测值之差；

（i=1、2、3）为误差方程对定向参数所求的偏导数；为误差方程对地面点坐标所求的偏导数。

记为：

V＝Bt+AX-l  P

同样可以对每个控制点列如下线性方程:

$v_x=(\frac{\∂F_x}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_x}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_x}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_x}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_x}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_x}{\∂b_3}\·\mathrm{\Δ}{b}_3)+F_{x0}$

$v_y=(\frac{\∂F_y}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_y}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_y}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_y}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_y}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_y}{\∂b_3}\·\mathrm{\Δ}{b}_3)+F_{y0}$

记为：

V＝Bt-l  P

其中：

$B=\left(\begin{array}{cccccc}\frac{\∂F_x}{\∂a_1}& \frac{\∂F_x}{\∂a_2}& \frac{\∂F_x}{\∂a_3}& \frac{\∂F_x}{\∂b_1}& \frac{\∂F_x}{\∂b_2}& \frac{\∂F_x}{\∂b_3}\\ \frac{\∂F_y}{\∂a_1}& \frac{\∂F_y}{\∂a_2}& \frac{\∂F_y}{\∂a_3}& \frac{\∂F_y}{\∂b_1}& \frac{\∂F_y}{\∂b_2}& \frac{\∂F_y}{\∂b_3}\end{array}\right),t={\left(\begin{array}{cccccc}\mathrm{\Δ}{a}_1& \mathrm{\Δ}{a}_2& \mathrm{\Δ}{a}_3& \mathrm{\Δ}{b}_1& \mathrm{\Δ}{b}_2& \mathrm{\Δ}{b}_3\end{array}\right)}^T$

$A=\left(\begin{array}{ccc}\frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{lat}}}& \frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{lon}}}& \frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{hei}}}\\ \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{lat}}}& \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{lon}}}& \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{hei}}}\end{array}\right),X={\left(\begin{array}{ccc}\mathrm{\Δ}{D}_{\mathrm{lat}}& {\mathrm{\ΔD}}_{\mathrm{lon}}& \mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)}^T$

$l=\left(\begin{array}{c}-F_{x0}\\ -F_{y0}\end{array}\right),V=\left(\begin{array}{c}{v}_x\\ v_y\end{array}\right)$

根据最小二乘平差构建法方程：

$\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)\left(\begin{array}{c}t\\ X\end{array}\right)=\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right).$

5）对法方程利用谱修正法迭代估计求解改正参数t和地面点坐标X，具体过程为：将法方程两边同时加上 $\left(\begin{array}{c}t\\ X\end{array}\right),$ 得到：

$[\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)+E]\left(\begin{array}{c}t\\ X\end{array}\right)=\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right)+\left(\begin{array}{c}t\\ X\end{array}\right)$

式中，E为和 $\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)$ 阶数相同的单位矩阵，由于式子两边都含有未知参数 $\left(\begin{array}{c}t\\ X\end{array}\right),$ 所以只能采用迭代的方法求解，其迭代公式为：

${\left(\begin{array}{c}t\\ X\end{array}\right)}^{\left(k\right)}={[\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)+E]}^{-1}+[\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right)+{\left(\begin{array}{c}t\\ X\end{array}\right)}^{(k-1)}]$

式中，k为迭代次数。

6）每完成一次平差迭代时，计算出检查点对应的地面点坐标，并同时计算出检查点的地面点坐标和已知的检查点地面坐标之差，即为检查点精度，也就是此次平差迭代后所能达到的物方精度。

7）每完成一次平差迭代时，通过计算的检查点对应的地面点坐标反投影计算出像点坐标，并计算像方残差，也就是此次平差迭代后所能达到的像方精度，通过像方残差并利用信息扩散模型确定下一次迭代像点坐标观测值的权值P。

8）当满足定向参数中的平移参数a₁，b₁均小于阈值时（本实施例为0.1个像元pixel，但是不限于此）时，整个平差迭代结束；当不满足预设条件时，返回步骤5）继续迭代计算，直到满足迭代收敛条件，如果迭代次数达到预设迭代次数（本实施例设定的迭代次数为20次，但是不限于此），仍然不能收敛，那么迭代结束。

9）当平差迭代结束时，输出最后一次迭代时所计算得到的平差物方精度和像方精度的精度报告，根据精度报告可以对平差效果进行评价，如果通过精度报告得知平差精度较高则说明平差效果较好，平差参数可以用于后续的生产和测绘中，如果通过报告得知平差精度较低，可以从精度报告中查找原因，进而避免对后续的应用产生影响。

上述各实施例仅用于说明本发明，其中方法的各实施步骤等都是可以有所变化的，凡是在本发明技术方案的基础上进行的等同变换和改进，均不应排除在本发明的保护范围之外。

Claims (3)

1.一种三线阵卫星影像区域网平差方法，其包括以下步骤：

1）读取卫星CCD阵列影像，选择相邻卫星影像重叠区域的同名点作为影像间的连接点，根据点之记或者控制影像在卫星影像上刺出控制点和检查点，并将连接点、控制点和检查点的坐标信息存成相应的点文件；

2）读取卫星CCD阵列影像的rcp参数文件、连接点文件、控制点文件和检查点文件，得到卫星CCD阵列影像的RPC模型以及控制点和检查点的地面坐标、影像坐标和连接点的像点坐标；

3）确定地面点坐标和像方变换模型的初值；

4）利用量测得到像点坐标和对应的地面点坐标，针对连接点和控制点分别对定向参数和地面点坐标求偏导并逐点构建误差方程并将构建的误差方程进行法化；

5）对法方程利用谱修正法迭代估计求解改正参数和地面点坐标；

6）每完成一次平差迭代时，计算出检查点对应的地面点坐标，并同时计算此次平差迭代后所能达到的物方精度；

7）每完成一次平差迭代时，通过计算的检查点对应的地面点坐标反投影计算出像点坐标，并计算此次平差迭代后所能达到的像方精度，通过像方残差并利用信息扩散模型确定下一次迭代像点坐标观测值的权值P；

8）当满足定向参数中的平移参数均小于阈值时，整个平差迭代结束；当不满足预设条件时，返回步骤5）继续迭代计算，直到满足迭代收敛条件；如果迭代次数达到预设迭代次数，仍然不能收敛，那么迭代结束；

9）当整个平差迭代结束时，输出最后一次迭代时所计算得到的平差物方精度和像方精度的精度报告，根据精度报告对平差效果进行评价。

2.如权利要求1所述的一种三线阵卫星影像区域网平差方法，其特征在于：步骤4）利用量测得到像点坐标和对应的地面点坐标，针对连接点和控制点分别对定向参数和地面点坐标求偏导并逐点构建误差方程并将构建的误差方程进行法化，具体过程为：

在模型的基础上建立仿射变换模型：

$\left.\begin{array}{c}{F}_y=\mathrm{\Δy}=a_1+a_2\·\mathrm{sample}+a_3\·\mathrm{line}\\ F_x=\mathrm{\Δx}=b_1+b_2\·\mathrm{sample}+b_3\·\mathrm{line}\end{array}\right\}$

式中，Δy和Δx为控制点在影像坐标系中的量测坐标与真实坐标的差值，即改正数；a₁,a₂,a₃和b₁,b₂,b₃是影像的定向参数，line和sample是控制点在影像坐标系中的行列号；

建立误差方程：

$v_x=\left(\begin{array}{c}\frac{\∂F_x}{\∂a_1}\·{\mathrm{\Δa}}_1+\frac{{\∂F}_x}{{\∂a}_2}\·{\mathrm{\Δa}}_2+\frac{\∂F_x}{\∂a_3}\·{\mathrm{\Δa}}_3+\frac{\∂F_x}{{\∂b}_1}\·{\mathrm{\Δb}}_1+\frac{{\∂F}_x}{\∂b_2}\·{\mathrm{\Δb}}_2+\frac{\∂F_x}{\∂b_3}\·\mathrm{\Δ}{b}_3\\ +\frac{\∂F_x}{\∂D_{\mathrm{lat}}}\·\mathrm{\Δ}{D}_{\mathrm{lat}}+\frac{\∂F_x}{\∂D_{\mathrm{lon}}}\·\mathrm{\Δ}{D}_{\mathrm{lon}}+\frac{\∂F_x}{\∂D_{\mathrm{hei}}}\·\mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)+F_{x0}$

$v_y=\left(\begin{array}{c}\frac{\∂F_y}{\∂a_1}\·{\mathrm{\Δa}}_1+\frac{{\∂F}_y}{{\∂a}_2}\·{\mathrm{\Δa}}_2+\frac{\∂F_y}{\∂a_3}\·{\mathrm{\Δa}}_3+\frac{\∂F_y}{{\∂b}_1}\·{\mathrm{\Δb}}_1+\frac{{\∂F}_y}{\∂b_2}\·{\mathrm{\Δb}}_2+\frac{\∂F_y}{\∂b_3}\·\mathrm{\Δ}{b}_3\\ +\frac{\∂F_y}{\∂D_{\mathrm{lat}}}\·\mathrm{\Δ}{D}_{\mathrm{lat}}+\frac{\∂F_y}{\∂D_{\mathrm{lon}}}\·\mathrm{\Δ}{D}_{\mathrm{lon}}+\frac{\∂F_y}{\∂D_{\mathrm{hei}}}\·\mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)+F_{y0}$

式中，ΔD_lat、ΔD_lon、ΔD_hei为地面点坐标改正数；Δa₁，Δa₂，Δa₃，Δb₁，Δb₂，Δb₃为影像定向参数改正数；v_x,v_y为像点坐标改正数；F_x0,F_y0为像点坐标近似值与像点坐标观测值之差；

（i=1、2、3）为误差方程对定向参数所求的偏导数；

为误差方程对地面点坐标所求的偏导数；

记为：

V＝Bt+AX-l  P

同样可以对每个控制点列如下线性方程：

$v_x=(\frac{\∂F_x}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_x}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_x}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_x}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_x}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_x}{\∂b_3}\·{\mathrm{\Δb}}_3)+F_{x0}$

$v_y=(\frac{\∂F_y}{\∂a_1}\·\mathrm{\Δ}{a}_1+\frac{\∂F_y}{\∂a_2}\·\mathrm{\Δ}{a}_2+\frac{\∂F_y}{\∂a_3}\·\mathrm{\Δ}{a}_3+\frac{\∂F_y}{\∂b_1}\·\mathrm{\Δ}{b}_1+\frac{\∂F_y}{\∂b_2}\·\mathrm{\Δ}{b}_2+\frac{\∂F_y}{\∂b_3}\·{\mathrm{\Δb}}_3)+F_{y0}$

记为：

V＝Bt-l  P

其中：

$B=\left(\begin{array}{cccccc}\frac{\∂F_x}{\∂a_1}& \frac{\∂F_x}{\∂a_2}& \frac{\∂F_x}{\∂a_3}& \frac{\∂F_x}{\∂b_1}& \frac{\∂F_x}{\∂b_2}& \frac{\∂F_x}{\∂b_3}\\ \frac{\∂F_y}{\∂a_1}& \frac{\∂F_y}{\∂a_2}& \frac{\∂F_y}{\∂a_3}& \frac{\∂F_y}{\∂b_1}& \frac{\∂F_y}{\∂b_2}& \frac{\∂F_y}{\∂b_3}\end{array}\right),t={\left(\begin{array}{cccccc}\mathrm{\Δ}{a}_1& \mathrm{\Δ}{a}_2& {\mathrm{\Δa}}_3& \mathrm{\Δ}{b}_1& \mathrm{\Δ}{b}_2& \mathrm{\Δ}{b}_3\end{array}\right)}^T$

$A=\left(\begin{array}{ccc}\frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{lat}}}& \frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{lon}}}& \frac{\∂F_x}{\mathrm{\Δ}{D}_{\mathrm{hei}}}\\ \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{lat}}}& \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{lon}}}& \frac{\∂F_y}{\mathrm{\Δ}{D}_{\mathrm{hei}}}\end{array}\right),X={\left(\begin{array}{ccc}\mathrm{\Δ}{D}_{\mathrm{lat}}& \mathrm{\Δ}{D}_{\mathrm{lon}}& \mathrm{\Δ}{D}_{\mathrm{hei}}\end{array}\right)}^T$

$l=\left(\begin{array}{c}-F_{x0}\\ -F_{y0}\end{array}\right),V=\left(\begin{array}{c}{v}_x\\ v_y\end{array}\right)$

根据最小二乘平差构建法方程：

$\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)\left(\begin{array}{c}t\\ X\end{array}\right)=\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right).$

3.如权利要求2所述的一种三线阵卫星影像区域网平差方法，其特征在于：步骤5）对法方程利用谱修正法迭代估计求解改正参数和地面点坐标，具体过程为：

将法方程两边同时加上 $\left(\begin{array}{c}t\\ X\end{array}\right),$ 得到：

$[\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)+E]\left(\begin{array}{c}t\\ X\end{array}\right)=\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right)+\left(\begin{array}{c}t\\ X\end{array}\right)$

式中，E为和 $\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)$ 阶数相同的单位矩阵，由于式子两边都含有未知参数 $\left(\begin{array}{c}t\\ X\end{array}\right),$ 采用迭代的方法求解，其迭代公式为：

${\left(\begin{array}{c}t\\ X\end{array}\right)}^{\left(k\right)}=\left[{\left(\begin{array}{cc}{A}^T\mathrm{PA}& A^T\mathrm{PB}\\ B^T\mathrm{PA}& B^T\mathrm{PB}\end{array}\right)+E]}^{-1}\right[\left(\begin{array}{c}{A}^T\mathrm{Pl}\\ B^T\mathrm{Pl}\end{array}\right)+{\left(\begin{array}{c}t\\ X\end{array}\right)}^{(k-1)}]$

式中，k为迭代次数。