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

本公开提供了一种基于不完全分解预处理的共轭梯度法光束法平差方法，包括：步骤S1，构建光束法平差的法方程的统一形式，采用多类未知数交替趋近解法的多类法方程；步骤S2，采用对角线矩阵块作为预处理矩阵进行不完全Cholesky分解预处理，降低法方程系数阵的条件数；步骤S3，给定共轭梯度法光束法平差的迭代初值，依次沿着相互共轭的方向进行相关目标函数解的搜索，并将目标函数改造为能够直接求得原法方程的解X的值的形式；步骤S4，给定一小值ε作为迭代收敛条件进行迭代计算，直至满足收敛条件迭代终止，输出方程的最优解。本公开能够有效地提高区域网平差的稳定性和计算效率，更好地适应大区域影像产品生产和三维重建等工作。

【技术实现步骤摘要】
基于不完全分解预处理的共轭梯度法光束法平差方法
本专利技术涉及摄影测量与遥感
，尤其涉及一种基于不完全分解预处理的共轭梯度法光束法平差方法。
技术介绍
基于共线条件方程式的光束法平差解法，是一种把控制点的像点坐标、待定点的像点坐标以及其他内业、外业量测数据的一部分或全部均视作观测值，以整体且同时解求它们的最或是值和待定点空间坐标的解算方法。在航天航空摄影测量、低空摄影测量以及近景摄影测量等区域网平差工作中，一方面由于摄影方式呈现多样化的特点，另一方面也常常需要引入多种类型的附加参数(如定位和定姿设备的安装矩阵、自检校模型参数和附加几何约束条件等)，法方程系数阵的稀疏带状结构遭到了破坏，难以确定带宽且为加边矩阵。这给法方程求解快速算法的设计带来了显著的困难。针对这一问题，国内外学者提出了多种解决方案。如McGloneC.(1996)考虑各种几何约束条件，构建统一模式的带附加参数的条件平差模型，区域网平差求解时采用整体广义逆解法；WangX.和ClarkeT.A.(2001)基于分组最小二乘估计方式，对像方参数和物方坐标进行相对独立估计的分组平差；BartoliA.(2002)采用通常的代价函数，将不同类型的模型统一到同一框架下进行准线性优化；冯其强,等(2008)根据工业摄影测量对运算速度的需求，提出逐点修正未知近似值的自检较点松弛迭代解法；NiK.,等(2007)将全局的非线性代价函数转化为各自坐标系下的子问题，尽可能利用外存空间和实现并行处理等。这些方法都在一定程度上兼顾了解算精度和运算效率的要求，但模型的设计和解算方法相对复杂，且不容易通过程序实现。共轭梯度法是一种在可预测模式中不需要非零输入的，处理大规模稀疏系统的非常有用的迭代估计方法。该方法由最速下降法演化而来，最早由HestenesM.R.和StiefelE.(1952)提出，沿着梯度所构建的相互共轭的方向计算，从而最快速度得到目标函数的最小值。但由于对舍入误差的敏感性问题，算法的应用受到了极大的限制。近些年来，该问题通过与各种预处理算法相结合而得到解决。由于算法简便、存储需求小等优势，共轭梯度法已经成为求解大规模、稀疏线性方程组最常用的一类方法，并且开始引入到测量平差领域。因而，利用共轭梯度法在正交基下的不变性进行基于共轭梯度法的光束法平差算法设计，并采用基于Cholesky分解的预处理方法减低法方程系数阵的条件数，能够有效地提高区域网平差的稳定性和计算效率，更好地适应大区域影像产品生产和三维重建等工作。公开内容(一)要解决的技术问题本公开提供了一种基于不完全分解预处理的共轭梯度法光束法平差方法，以至少部分解决以上所提出的技术问题。(二)技术方案根据本公开的一个方面，提供了一种基于不完全分解预处理的共轭梯度法光束法平差方法，包括：步骤S1，构建光束法平差的法方程的统一形式，并构建采用多类未知数交替趋近解法的多类法方程；步骤S2，采用对角线矩阵块构造与光束法平差法方程系数阵近似的系数阵形式，将其作为预处理矩阵进行不完全Cholesky分解预处理，构造光束法平差新的法方程结构，降低法方程系数阵的条件数；步骤S3，给定共轭梯度法光束法平差的迭代初值，依次沿着相互共轭的方向进行相关目标函数解的搜索，并将目标函数改造为能够直接求得原法方程的解X的值的形式；步骤S4，给定一小值ε作为迭代收敛条件进行迭代计算，直至满足收敛条件迭代终止，输出基于不完全分解预处理的共轭梯度法光束法平差法方程的最优解。在本公开一些实施例中，所述步骤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类法方程的常数项矩阵。在本公开一些实施例中，所述n类未知数包括：影像外方位元素的改正数、物方点空间坐标的改正数，以及各非摄影测量附件参数的改正数。在本公开一些实施例中，所述步骤S2包括：子步骤S201，采用对角线矩阵块构造与光束法平差法方程系数阵近似的系数阵形式，子步骤S103中多类未知数交替趋近解法，则交替一步的过程用如下方程组作为近似表示：其中，M为子步骤S103中n类未知数的法方程式的系数矩阵的对角阵；子步骤S202，将子步骤S201中的对角线矩阵块构造的系数阵作为共轭梯度法光束法平差的预处理矩阵，进行不完全Cholesky分解预处理，得到：M＝LLT其中，L为一个下三角矩阵；子步骤S203，对系数阵M不完全Cholesky分解后，将子步骤S102中的法方程转换为如下形式：[L-1N(L-1)T]LTX＝L-1W其中，法方程的未知数X替换为Y＝LTX；法方程的系数阵N替换为F＝L-1N(L-1)T，仍满足正定对称矩阵条件；法方程的常数项W替换为G＝L-1W，将上式简化为如下形式：FY＝G经过变换，将子步骤S102中共轭梯度法光束法平差原法方程的求解问题转化为对新法方程的求解问题，其中，系数阵F的条件数远低于原法方程系数阵N的条件数。在本公开一些实施例中，所述步骤S3包括：子步骤S301，给定子步骤S102中法方程的初值为X0，将预处理后法方程的初值替换为Y0＝LTX0，设定共轭梯度法光束法平差法的迭代第一步为k＝0；子步骤S302：在n维向量空间中，共轭梯度法光束法平差原法方程依次沿着相互共轭的方向X1,X2,…Xn，进行相关目标函数解的搜索；当搜索至第k步时的解为Xk＝α1X1+α2X2+...+αkXk，其中，相关目标函数的解满足Yk＝LTXk；子步骤S303：为使共轭梯度法光束法平差能够直接求得原法方程的解X的值，将目标函数中的Xk替换为(L-1)TYk。在本公开一些实施例中，所述步骤S302包括：对于新的法方程，引入残差变量Rk：Rk＝G-FYk＝L-1[W-N(L-1)TLTXk]＝L-1rk其中，Yk＝LTXk，rk是针对原法方程引入的一个残差变量，满足：rk＝W-N(L-1)TLTXk；再计算残差变量Rk的初值：R0＝L-1[W-N(L-1)TLTX0]＝L-1r0。在本公开一些实施例中，所述步骤S303包括：引入一个变量Zk：Zk＝(L-1)TL-1rk＝M-1rk解算方程MZ0＝R0得到变量Zk的初值Z0，并令X1＝Z0。在本公开一些实施例中，所述步骤S4包括：子步骤S401：给定一小值，设定迭代收敛条件对平差是否继续迭代进行判断，若满足收敛条件则转至子步骤S404，若不满足则转至子步骤S402至子步骤S403；子步骤S402：计算迭代至该步时与目标函数取值相关的参数α、方程组的搜索结果X、残差变量R本文档来自技高网...

【技术保护点】
1.一种基于不完全分解预处理的共轭梯度法光束法平差方法，包括：步骤S1，构建光束法平差的法方程的统一形式，并构建采用多类未知数交替趋近解法的多类法方程；步骤S2，采用对角线矩阵块构造与光束法平差法方程系数阵近似的系数阵形式，将其作为预处理矩阵进行不完全Cholesky分解预处理，构造光束法平差新的法方程结构，降低法方程系数阵的条件数；步骤S3，给定共轭梯度法光束法平差的迭代初值，依次沿着相互共轭的方向进行相关目标函数解的搜索，并将目标函数改造为能够直接求得原法方程的解X的值的形式；步骤S4，给定一小值ε作为迭代收敛条件进行迭代计算，直至满足收敛条件迭代终止，输出基于不完全分解预处理的共轭梯度法光束法平差法方程的最优解。

【技术特征摘要】
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＝LLT其中，L为一个下三角矩阵；子步骤S203，对系数阵M不完全Cholesky分解后，将子步骤S102中的法方程转换为如下形式：[L-1N(L-1)T]LTX＝L-1W其中，法方程的未知数X替换为Y＝LTX；法方程的系数阵N替换为F＝L-1N(L-1)T，仍满足正定对称矩阵条件；法方程的常数项W替换为G＝L-1W，将上式简化为如下形式：FY＝G经过变换，将子步骤S102中共轭梯度法光束法平差...

【专利技术属性】
技术研发人员：胡堃，张永军，雷斌，黄旭，
申请(专利权)人：中国科学院电子学研究所，
类型：发明
国别省市：北京,11