基于不完全分解预处理的共轭梯度法光束法平差方法与流程

技术特征：

1.一种基于不完全分解预处理的共轭梯度法光束法平差方法，包括：

步骤S1，构建光束法平差的法方程的统一形式，并构建采用多类未知数交替趋近解法的多类法方程；

步骤S2，采用对角线矩阵块构造与光束法平差法方程系数阵近似的系数阵形式，将其作为预处理矩阵进行不完全Cholesky分解预处理，构造光束法平差新的法方程结构，降低法方程系数阵的条件数；

步骤S3，给定共轭梯度法光束法平差的迭代初值，依次沿着相互共轭的方向进行相关目标函数解的搜索，并将目标函数改造为能够直接求得原法方程的解X的值的形式；

步骤S4，给定一小值ε作为迭代收敛条件进行迭代计算，直至满足收敛条件迭代终止，输出基于不完全分解预处理的共轭梯度法光束法平差法方程的最优解。

2.根据权利要求1所述的共轭梯度法光束法平差方法，所述步骤S1包括：

子步骤S101，获取区域内具有重叠度的多景光学影像和对应的成像外方位参数，包括：与成像姿态相关的3个外方位角元素：俯仰角，翻滚角和偏航角，与成像位置相关的3个外方位线元素：像主点相对于摄影测量坐标系在X，Y和Z方向的坐标；

子步骤S102，通过多景光学影像和对应的成像参数，利用共线方程和最小二乘平差方法，构建光束法平差法方程的统一形式：

NX＝W

其中，X为法方程的未知数矩阵；N为法方程的系数矩阵，满足n阶对阵正定矩阵；W为法方程的常数项矩阵；

子步骤S103，构造采用多类未知数交替趋近解法的n类未知数的法方程式：

其中，n≥2，所述n类未知数包括：影像外方位元素的改正数、物方点空间坐标的改正数，Xn为第n类法方程的未知数矩阵；Nnn为第n类法方程的系数矩阵，满足n阶对阵正定矩阵；Wn为第n类法方程的常数项矩阵。

3.根据权利要求2所述的共轭梯度法光束法平差方法，所述n类未知数包括：影像外方位元素的改正数、物方点空间坐标的改正数，以及各非摄影测量附件参数的改正数。

4.根据权利要求1所述的共轭梯度法光束法平差方法，所述步骤S2包括：

子步骤S201，采用对角线矩阵块构造与光束法平差法方程系数阵近似的系数阵形式，子步骤S103中多类未知数交替趋近解法，则交替一步的过程用如下方程组作为近似表示：

其中，M为子步骤S103中n类未知数的法方程式的系数矩阵的对角阵；

子步骤S202，将子步骤S201中的对角线矩阵块构造的系数阵作为共轭梯度法光束法平差的预处理矩阵，进行不完全Cholesky分解预处理，得到：

M＝LL^T

其中，L为一个下三角矩阵；

子步骤S203，对系数阵M不完全Cholesky分解后，将子步骤S102中的法方程转换为如下形式：

[L^-1N(L^-1)^T]L^TX＝L^-1W

其中，法方程的未知数X替换为Y＝L^TX；法方程的系数阵N替换为F＝L^-1N(L^-1)^T，仍满足正定对称矩阵条件；法方程的常数项W替换为G＝L^-1W，将上式简化为如下形式：

FY＝G

经过变换，将子步骤S102中共轭梯度法光束法平差原法方程的求解问题转化为对新法方程的求解问题，其中，系数阵F的条件数远低于原法方程系数阵N的条件数。

5.根据权利要求1所述的共轭梯度法光束法平差方法，所述步骤S3包括：

子步骤S301，给定子步骤S102中法方程的初值为X⁰，将预处理后法方程的初值替换为Y⁰＝L^TX⁰，设定共轭梯度法光束法平差法的迭代第一步为k＝0；

子步骤S302：在n维向量空间中，共轭梯度法光束法平差原法方程依次沿着相互共轭的方向X1,X2,…Xn，进行相关目标函数解的搜索；当搜索至第k步时的解为X^k＝α1X1+α2X2+...+αkXk，其中，相关目标函数的解满足Y^k＝L^TX^k；

子步骤S303：为使共轭梯度法光束法平差能够直接求得原法方程的解X的值，将目标函数中的X^k替换为(L^-1)^TY^k。

6.根据权利要求5所述的共轭梯度法光束法平差方法，所述步骤S302包括：

对于新的法方程，引入残差变量R^k：

R^k＝G-FY^k＝L^-1[W-N(L^-1)^TL^TX^k]＝L^-1r^k

其中，Y^k＝L^TX^k，r^k是针对原法方程引入的一个残差变量，满足：r^k＝W-N(L^-1)^TL^TX^k；

再计算残差变量R^k的初值：

R⁰＝L^-1[W-N(L^-1)^TL^TX⁰]＝L^-1r⁰。

7.根据权利要求6所述的共轭梯度法光束法平差方法，所述步骤S303包括：

引入一个变量Z^k：

Z^k＝(L^-1)^TL^-1r^k＝M^-1r^k

解算方程MZ⁰＝R⁰得到变量Z^k的初值Z⁰，并令X1＝Z⁰。

8.根据权利要求7所述的共轭梯度法光束法平差方法，所述步骤S4包括：

子步骤S401：给定一小值，设定迭代收敛条件对平差是否继续迭代进行判断，若满足收敛条件则转至子步骤S404，若不满足则转至子步骤S402至子步骤S403；

子步骤S402：计算迭代至该步时与目标函数取值相关的参数α、方程组的搜索结果X、残差变量R和参数Z的取值；

子步骤S403：计算迭代至该步时与目标函数取值相关的参数β以及下一步的迭代方向；

子步骤S404：迭代终止，输出目标函数取得最小值时X^k的值，作为共轭梯度法光束法平差法方程的最优解。

9.根据权利要求8所述的共轭梯度法光束法平差方法，所述子步骤S401包括：

给定一小值ε，进行判断，若||R^k||≥ε则执行子步骤S402至子步骤S403，并进行迭代；否则迭代中止，执行子步骤S404。

10.根据权利要求8所述的共轭梯度法光束法平差方法，所述子步骤S402包括：

令k＝k+1，在第k步沿着方向Xk进行搜索，Xk同时满足与全部的X1,X2,…,Xk-1均共轭，计算方程组解的搜索结果：

X^k＝X^k-1+αkXk

其中，与目标函数取值相关的参数αk由以下公式计算得出：

再计算变量Z^k：

Z^k＝(L^-1)^TL^-1R^k＝M^-1R^k，

其中，残差变量R^k＝R^k-1-αkNXk；

所述子步骤S403包括：

计算第k+1步的搜索方向Xk+1：

Xk+1＝Z^k+βkXk

其中，参数βk满足：