留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

一种求解单像空间后方交会的监督学习方法

李佳田 贾成林 牛一如 阿晓荟 高鹏 晏玲

李佳田, 贾成林, 牛一如, 阿晓荟, 高鹏, 晏玲. 一种求解单像空间后方交会的监督学习方法[J]. 武汉大学学报 ● 信息科学版, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
引用本文: 李佳田, 贾成林, 牛一如, 阿晓荟, 高鹏, 晏玲. 一种求解单像空间后方交会的监督学习方法[J]. 武汉大学学报 ● 信息科学版, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
LI Jiatian, JIA Chenglin, NIU Yiru, A Xiaohui, GAO Peng, YAN Ling. A Supervised Learning Method for Solving Space Resection of Single Image[J]. Geomatics and Information Science of Wuhan University, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
Citation: LI Jiatian, JIA Chenglin, NIU Yiru, A Xiaohui, GAO Peng, YAN Ling. A Supervised Learning Method for Solving Space Resection of Single Image[J]. Geomatics and Information Science of Wuhan University, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292

一种求解单像空间后方交会的监督学习方法

doi: 10.13203/j.whugis20170292
基金项目: 

国家自然科学基金 41561082

国家自然科学基金 41161061

详细信息
    作者简介:

    李佳田, 博士, 教授, 主要研究方向为数值最优化方法与机器场景理解。ljtwcx@163.com

    通讯作者: 贾成林, 硕士生。1774356658@qq.com
  • 中图分类号: P231

A Supervised Learning Method for Solving Space Resection of Single Image

Funds: 

The National Natural Science Foundation of China 41561082

The National Natural Science Foundation of China 41161061

More Information
    Author Bio:

    LI Jiatian, PhD, professor, specializes in numerical optimization and scene understanding for robot. E-mail:ljtwcx@163.com

    Corresponding author: JIA Chenglin, postgraduate. E-mail: 1774356658@qq.com
  • 摘要: 单像空间后方交会可描述为非线性最小二乘问题,不可导、法方程系数矩阵病态以及陷入局部极值是造成其数值过程不收敛的主要原因。不同地区的控制点空间分布不具相似性,若把同一地区同一组控制点之下数张已知外方位元素的像片看作一个样本集,则在给定每个外方位元素初值的前提下,可通过监督学习方法求取外方位元素的整体下降方向;而对于单像空间后方交会中因前述原因不收敛的情况,则可采用整体下降方向近似解算。以此为出发点,提出一种单像空间后方交会求解的监督学习方法,主要过程是:①训练阶段,利用监督学习过程,对同一测区内不同姿态像片所组成的样本集进行整体外方位元素的求解,得到该测区外方位元素的整体下降方向集合;②测试阶段,对该测区的任意像片,给定外方位元素的初值,直接采用训练阶段得到的整体下降方向集合进行外方位元素的迭代求解。对比试验表明,该方法在数值过程收敛性与初值依赖性上均表现出较强的优势,并能克服欧拉角法因法方程系数矩阵病态而无法收敛的情况。
  • 图  1  监督学习方法示意图

    Figure  1.  Schematic Diagram of Supervised Learning Method

    图  2  牛顿法和监督学习方法对比图

    Figure  2.  Comparison of Newton's Method and Supervised Learning Method

    图  3  理想垂直摄影示意图

    Figure  3.  Schematic Diagram of Ideal Vertical Photography

    图  4  多重交向摄影试验示意图

    Figure  4.  Schematic Diagram of Multi-intersection Photography

    图  5  计算值与真值对比曲线图

    Figure  5.  Comparison Curve Between Calculated Value and True Value

    图  6  3种方法计算结果对比曲线图

    Figure  6.  Comparison of the Results of Three Methods

    图  7  3种方法迭代次数对比图

    Figure  7.  Comparison of Iterations of Three Methods

    表  1  理想垂直摄影外方位元素计算结果

    Table  1.   Result of Calculation of Ideal Vertical Photography

    元素 项目 真值 欧拉角法 监督学习方法
    线元素/m XS 0.000 000 不收敛 0.001 565
    YS 0.000 000 -0.002 479
    ZS 150.000 000 150.003 025
    角元素/rad φ 0.000 000 不收敛 0.000 009
    ω 0.000 000 -0.000 002
    κ 0.000 000 0.000 012
    下载: 导出CSV

    表  2  地面控制点坐标

    Table  2.   Coordinates of GCPs

    点号 X/m Y/m Z/m
    1 10 10 11
    2 160 40 24
    3 320 40 34
    4 10 160 11
    5(A) 160 160 0
    6 320 160 28
    7 40 320 14
    8 160 320 29
    9 320 320 27
    下载: 导出CSV
  • [1] 张祖勋, 张剑清.数字摄影测量学[M].武汉:武汉大学出版社, 1996

    Zhang Zuxun, Zhang Jianqing. Digital Photogrammetry[M]. Wuhan:Wuhan University Press, 1996
    [2] 王之卓.摄影测量原理[M].武汉:武汉大学出版社, 2007

    Wang Zhizhuo. The Principles of Photogrammetry[M]. Wuhan:Wuhan University Press, 2007
    [3] 袁修孝.航空摄影测量影像定向的若干探讨[J].地球科学进展, 2007, 22(8):828-834 doi:  10.3321/j.issn:1001-8166.2007.08.008

    Yuan Xiuxiao. Some Investigations of Image Orientation in Aerial Photogrammetry[J]. Advances in Earth Science, 2007, 22(8):828-834 doi:  10.3321/j.issn:1001-8166.2007.08.008
    [4] 唐利民.非线性最小二乘的不适定性及算法研究[D].长沙: 中南大学, 2011 http://cdmd.cnki.com.cn/Article/CDMD-10533-1011177506.htm

    Tang Limin. Research on the Ill-posed and Solving Methods of Nonlinear Least Squares Problem[D]. Changsha: Central South University, 2011 http://cdmd.cnki.com.cn/Article/CDMD-10533-1011177506.htm
    [5] Haralick B M, Lee C N, Ottenberg K, et al. Review and Analysis of Solutions of the Three Point Perspective Pose Estimation Problem[J]. International Journal of Computer Vision, 1994, 13(3):331-356 doi:  10.1007/BF02028352
    [6] 冯文灏.近景摄影测量[M].武汉:武汉大学出版社, 2002

    Feng Wenhao. Close Range Photogrammetry[M]. Wuhan:Wuhan University Press, 2002
    [7] 付仲良, 周凡, 俞志强.综合多种特征的后方交会法[J].测绘学报, 2014, 43(8):827-834 http://www.cnki.com.cn/Article/CJFDTOTAL-CHXB201408009.htm

    Fu Zhongliang, Zhou Fan, Yu Zhiqiang. A Space Resection Synthesized the Multiple Features[J]. Acta Geodaetica et Cartographica Sinica, 2014, 43(8):827-834 http://www.cnki.com.cn/Article/CJFDTOTAL-CHXB201408009.htm
    [8] 李斐, 郝卫峰, 王文睿, 等.非线性病态问题解算的扰动分析[J].测绘学报, 2011, 40(1):5-9 http://d.old.wanfangdata.com.cn/Periodical/chxb201101002

    Li Fei, Hao Weifeng, Wang Wenrui, et al. The Perturbation Analysis of Nonlinear Ill-conditioned Solution[J]. Acta Geodaetica et Cartographica Sinica, 2011, 40(1):5-9 http://d.old.wanfangdata.com.cn/Periodical/chxb201101002
    [9] 郭忠磊, 翟京生, 邓凯亮.基于LM方法的大姿态角影像空间后方交会[J].武汉大学学报·信息科学版, 2016, 41(4):565-568 http://ch.whu.edu.cn/CN/abstract/abstract5432.shtml

    Guo Zhonglei, Zhai Jingsheng, Deng Kailiang. Space Resection of Big Attitude Image Based on Levenberg-Marquardt Method[J]. Geomatics and Information Science of Wuhan University, 2016, 41(4):565-568 http://ch.whu.edu.cn/CN/abstract/abstract5432.shtml
    [10] 何敬, 李永树.基于特征点和最优路径的无人机影像拼接方法[J].遥感技术与应用, 2012, 27(2):173-219 http://d.old.wanfangdata.com.cn/Periodical/ygjsyyy201202003

    He Jing, Li Yongshu. A Mosaic Method in Unmanned Aerial Vehicle Images Based on Feature Points and Optimal Path[J]. Remote Sensing Technology and Application, 2012, 27(2):173-219 http://d.old.wanfangdata.com.cn/Periodical/ygjsyyy201202003
    [11] Ji Q, Costa M S, Haralick R M, et al. A Robust Linear Least-Squares Estimation of Camera Exterior Orientation Using Multiple Geometric Features[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2000, 55(2):75-93 doi:  10.1016/S0924-2716(00)00009-5
    [12] 官云兰, 程效军, 周世健, 等.基于单位四元数的空间后方交会解算[J].测绘学报, 2008, 37(1):30-35 doi:  10.3321/j.issn:1001-1595.2008.01.006

    Guan Yunlan, Cheng Xiaojun, Zhou Shijian, et al. A Solution to Space Resection Based on Unit Quaternion[J]. Acta Geodaetica et Cartographica Sinica, 2008, 37(1):30-35 doi:  10.3321/j.issn:1001-1595.2008.01.006
    [13] Lepetit V, Moreno-Noguer F, Fua P. EPnP:An Accurate O(n)Solution to the PnP Problem[J]. International Journal of Computer Vision, 2009, 81(2):155-166 doi:  10.1007/s11263-008-0152-6
    [14] Penatesanchez A, Andradecetto J, Morenonoguer F. Exhaustive Linearization for Robust Camera Pose and Focal Length Estimation[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2013, 35(10):2387-2400 doi:  10.1109/TPAMI.2013.36
    [15] 龚辉, 姜挺, 江刚武, 等.一种基于四元数的空间后方交会全局收敛算法[J].测绘学报, 2011, 40(5):639-645 http://www.cnki.com.cn/Article/CJFDTOTAL-CHXB201105020.htm

    Gong Hui, Jiang Ting, Jiang Gangwu, et al. A Globally Convergent Algorithm of Space Resection Based on Quaternion[J]. Acta Geodaetica et Cartographica Sinica, 2011, 40(5):639-645 http://www.cnki.com.cn/Article/CJFDTOTAL-CHXB201105020.htm
    [16] 江刚武, 姜挺, 王勇, 等.基于单位四元数的无初值依赖空间后方交会[J].测绘学报, 2007, 36(2):169-175 doi:  10.3321/j.issn:1001-1595.2007.02.010

    Jiang Gangwu, Jiang Ting, Wang Yong, et al. Space Resection Indepentent of Initial Value Based on Unit Quaternions[J]. Acta Geodaetica et Cartographica Sinica, 2007, 36(2):169-175 doi:  10.3321/j.issn:1001-1595.2007.02.010
    [17] 张永军, 吴磊, 林立文, 等.摄影测量中病态问题的条件数指标分析[J].武汉大学学报·信息科学版, 2010, 35(3):308-312 http://ch.whu.edu.cn/CN/abstract/abstract882.shtml

    Zhang Yongjun, Wu Lei, Lin Liwen, et al. Condition Numbers for Evaluation of Ill-posed Problems in Photogrammetry[J]. Geomatics and Information Science of Wuhan University, 2010, 35(3):308-312 http://ch.whu.edu.cn/CN/abstract/abstract882.shtml
    [18] 王海亮, 向茂生, 李杨, 等.两种单像空间后方交会模型分析[J].光学技术, 2010, 36(1):14-19 http://d.old.wanfangdata.com.cn/Periodical/gxjs201001004

    Wang Hailiang, Xiang Maosheng, Li Yang, et al. The Analysis of Two Single Image Space Resection Models[J]. Optical Technique, 2010, 36(1):14-19 http://d.old.wanfangdata.com.cn/Periodical/gxjs201001004
    [19] 王振杰.测量中不适定问题的正则化解法[M].北京:科学出版社, 2006

    Wang Zhenjie. Regularization of Ill-posed Problems in Surverying[M]. Beijing:Science Press, 2006
    [20] 葛旭明, 伍吉仓.病态总体最小二乘问题的广义正则化[J].测绘学报, 2012, 41(3):372-377 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201203011

    Ge Xuming, Wu Jicang. Generalized Regularization to Ill-posed Total Least Squares Problem[J]. Acta Geodaetica et Cartographica Sinica, 2012, 41(3):372-377 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=chxb201203011
    [21] 欧吉坤.测量平差中不适定问题的统一表达和选权拟合法[J].测绘学报, 2004, 33(4):283-288 doi:  10.3321/j.issn:1001-1595.2004.04.001

    Ou Jikun. Uniform Expression of Solutions of Illposed Problems in Surveying Adjustment and the Fitting Method by Selection of the Parameter Weights[J]. Acta Geodaetica et Cartographica Sinica, 2004, 33(4):283-288 doi:  10.3321/j.issn:1001-1595.2004.04.001
    [22] 蒋涛, 李建成, 王正涛, 等.航空重力向下延拓病态问题求解[J].测绘学报, 2011, 40(6):684-689 http://cdmd.cnki.com.cn/Article/CDMD-10732-1014421889.htm

    Jiang Tao, Li Jiancheng, Wang Zhengtao, et al. Solution of Ill-posed Problem in Downward Continuation of Airborne Gravity[J]. Acta Geodaetica et Cartographica Sinica, 2011, 40(6):684-689 http://cdmd.cnki.com.cn/Article/CDMD-10732-1014421889.htm
    [23] 闫利, 聂前, 赵展.利用单位四元数描述线阵CCD影像的空间后方交会[J].武汉大学学报·信息科学版, 2010, 35(2):201-204 http://ch.whu.edu.cn/CN/abstract/abstract846.shtml

    Yan Li, Nie Qian, Zhao Zhan. Space Resection of Line Scanner CCD Image Based on the Description of Quaternions[J]. Geomatics and Information Science of Wuhan University, 2010, 35(2):201-204 http://ch.whu.edu.cn/CN/abstract/abstract846.shtml
    [24] 武汉大学测绘学院测量平差学科组.误差理论与测量平差基础[M]. 2版.武汉:武汉大学出版社, 2009

    The Subject Group of Surveying Adjustment of Institute of Surveying and Mapping of Wuhan University. Error Theory and Fundation of Surveying Adjustment[M]. 2nd ed. Wuhan:Wuhan University Press, 2009
    [25] 席少霖.非线性最优化方法[M].北京:高等教育出版社, 1992

    Xi Shaolin. The Method of Nonlinear Optimization[M]. Beijing:Higher Education Press, 1992
    [26] 张剑清, 潘励, 王树根.摄影测量学原理[M]. 2版.武汉:武汉大学出版社, 2009

    Zhang Jianqing, Pan Li, Wang Shugen. The Principles of Photogrammetry[M]. 2nd ed. Wuhan:Wuhan University Press, 2009
  • [1] 孙振, 曲国庆, 苏晓庆, 杜存鹏, 邓晓景.  测距定位方程参数估计的Frozen-Barycentre算法 . 武汉大学学报 ● 信息科学版, 2020, 45(9): 1478-1484. doi: 10.13203/j.whugis20180129
    [2] 姚宜斌, 黄书华, 张良, 胡羽丰, 李国平.  求解三维坐标转换参数的整体最小二乘新方法 ! . 武汉大学学报 ● 信息科学版, 2015, 40(7): 853-857. doi: 10.13203/j.whugis20130302
    [3] 姚宜斌, 黄书华, 陈家君.  求解自回归模型参数的整体最小二乘新方法 . 武汉大学学报 ● 信息科学版, 2014, 39(12): 1463-1466.
    [4] 姚宜斌, 黄书华, 孔建, 何军泉.  空间直线拟合的整体最小二乘算法 . 武汉大学学报 ● 信息科学版, 2014, 39(5): 571-574. doi: 10.13203/j.whugis20120104
    [5] 闫利, 胡晓斌.  利用GA求解卫星影像的空间后方交会 . 武汉大学学报 ● 信息科学版, 2013, 38(11): 1286-1289.
    [6] 向巍, 郭际明, 傅露.  基于垂直距离最小二乘拟合的双曲线沉降模型 . 武汉大学学报 ● 信息科学版, 2013, 38(5): 571-574.
    [7] 游为, 范东明, 傅淑娟.  同伦函数与填充函数相结合的非线性最小二乘平差模型 . 武汉大学学报 ● 信息科学版, 2010, 35(2): 185-188.
    [8] 邱卫宁, 齐公玉, 田丰瑞.  整体最小二乘求解线性模型的改进算法 . 武汉大学学报 ● 信息科学版, 2010, 35(6): 708-710.
    [9] 邓兴升, 花向红.  动态最小二乘支持向量机学习算法 . 武汉大学学报 ● 信息科学版, 2008, 33(11): 1122-1125.
    [10] 陈义, 陆珏, 郑波.  总体最小二乘方法在空间后方交会中的应用 . 武汉大学学报 ● 信息科学版, 2008, 33(12): 1271-1274.
    [11] 胡志刚, 花向红, 李昭, 韩红超.  基于同伦方法的非线性测量模型参数估计 . 武汉大学学报 ● 信息科学版, 2008, 33(9): 930-933.
    [12] 张双成, 杨元喜, 张勤, 高为广.  一种基于Bancroft算法的GPS动态抗差自适应滤波 . 武汉大学学报 ● 信息科学版, 2007, 32(4): 309-311.
    [13] 康志忠, 张祖勋, 阳凡林.  沿主光轴影像的单像空间后方交会方法研究 . 武汉大学学报 ● 信息科学版, 2007, 32(4): 293-296.
    [14] 宁伟, 陶华学, 卿熙宏.  广义非线性最小二乘测量参数平差的快速差分迭代解算 . 武汉大学学报 ● 信息科学版, 2005, 30(7): 617-620.
    [15] 张勤, 陶本藻.  基于同伦法的非线性最小二乘平差统一模型 . 武汉大学学报 ● 信息科学版, 2004, 29(8): 708-710.
    [16] 陶本藻, 张勤.  GPS非线性数据处理的同伦最小二乘模型 . 武汉大学学报 ● 信息科学版, 2003, 28(S1): 115-118.
    [17] 张祖勋, 张剑清, 胡翔云.  基于物方空间几何约束最小二乘匹配的建筑物半自动提取方法 . 武汉大学学报 ● 信息科学版, 2001, 26(4): 290-295.
    [18] 李桂苓, 万剑华, 陶华学.  用差商代替导数的非线性最小二乘估计 . 武汉大学学报 ● 信息科学版, 2001, 26(2): 118-121.
    [19] 白亿同.  非线性最小二乘平差迭代解法的收敛性 . 武汉大学学报 ● 信息科学版, 1991, 16(2): 92-95.
    [20] 白亿同.  从黎曼流形的观点看非线性最小二乘平差 . 武汉大学学报 ● 信息科学版, 1991, 16(4): 78-84.
  • 加载中
图(7) / 表(2)
计量
  • 文章访问数:  1018
  • HTML全文浏览量:  95
  • PDF下载量:  161
  • 被引次数: 0
出版历程
  • 收稿日期:  2018-05-16
  • 刊出日期:  2019-08-05

一种求解单像空间后方交会的监督学习方法

doi: 10.13203/j.whugis20170292
    基金项目:

    国家自然科学基金 41561082

    国家自然科学基金 41161061

    作者简介:

    李佳田, 博士, 教授, 主要研究方向为数值最优化方法与机器场景理解。ljtwcx@163.com

    通讯作者: 贾成林, 硕士生。1774356658@qq.com
  • 中图分类号: P231

摘要: 单像空间后方交会可描述为非线性最小二乘问题,不可导、法方程系数矩阵病态以及陷入局部极值是造成其数值过程不收敛的主要原因。不同地区的控制点空间分布不具相似性,若把同一地区同一组控制点之下数张已知外方位元素的像片看作一个样本集,则在给定每个外方位元素初值的前提下,可通过监督学习方法求取外方位元素的整体下降方向;而对于单像空间后方交会中因前述原因不收敛的情况,则可采用整体下降方向近似解算。以此为出发点,提出一种单像空间后方交会求解的监督学习方法,主要过程是:①训练阶段,利用监督学习过程,对同一测区内不同姿态像片所组成的样本集进行整体外方位元素的求解,得到该测区外方位元素的整体下降方向集合;②测试阶段,对该测区的任意像片,给定外方位元素的初值,直接采用训练阶段得到的整体下降方向集合进行外方位元素的迭代求解。对比试验表明,该方法在数值过程收敛性与初值依赖性上均表现出较强的优势,并能克服欧拉角法因法方程系数矩阵病态而无法收敛的情况。

English Abstract

李佳田, 贾成林, 牛一如, 阿晓荟, 高鹏, 晏玲. 一种求解单像空间后方交会的监督学习方法[J]. 武汉大学学报 ● 信息科学版, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
引用本文: 李佳田, 贾成林, 牛一如, 阿晓荟, 高鹏, 晏玲. 一种求解单像空间后方交会的监督学习方法[J]. 武汉大学学报 ● 信息科学版, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
LI Jiatian, JIA Chenglin, NIU Yiru, A Xiaohui, GAO Peng, YAN Ling. A Supervised Learning Method for Solving Space Resection of Single Image[J]. Geomatics and Information Science of Wuhan University, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
Citation: LI Jiatian, JIA Chenglin, NIU Yiru, A Xiaohui, GAO Peng, YAN Ling. A Supervised Learning Method for Solving Space Resection of Single Image[J]. Geomatics and Information Science of Wuhan University, 2019, 44(8): 1144-1152. doi: 10.13203/j.whugis20170292
  • 单像空间后方交会是摄影测量领域的基础性课题[1-2],它是指利用一定数量的地面控制点及像点求解像片外方位元素的过程。摄影测量中常借助快速、正确地恢复图像获取时的空间方位来确定地面目标的空间位置[3]。因此,对单像空间后方交会的解算进行研究具有典型的实际意义。空间后方交会是一个非线性模型,一般采用迭代方法进行解算[4],典型的解法包括欧拉角法[1-2]、角锥体法[5]及直接线性变换法[6],其中角锥体法对于外方位线元素初值依赖较强,直接线性变换法不适合高精度的测量任务[7]。欧拉角法将成像姿态用欧拉角描述,对共线方程线性化,并利用最小二乘法迭代求解[2]。在应用中,欧拉角法的限制主要表现为:①需要已知良好的外方位元素初值, 当初值与真值差距较大时,收敛速度慢,甚至不收敛;②对原始非线性模型用线性模型近似代替,有一定的模型误差,随着观测精度的提高,因线性近似所产生的模型误差与观测误差相当,有的甚至超过观测误差[8];③其迭代过程中需假定法方程系数矩阵满秩,否则将可能失效[9]。虽然外方位元素初值可从全球定位系统(Global Positioning System,GPS)、惯性导航系统和轨迹数据中获取,但无人机等平台自带GPS及姿态量角系统本身的精确性不高,并存在记录滞后的现象,通常不能直接使用[10]

    为减少初值依赖及线性近似的影响,学者们提出了许多有参考价值的方法,主要分为直接解法和迭代解法两类[11]。直接解法无需外方位元素初值,但解算精度与鲁棒性较低,文献[12]利用单位四元数分别解算外方位线元素与角元素。代表性的迭代解法如文献[13-14]以线性解为初值进行高斯-牛顿全局优化;文献[15-16]采用对偶四元数描述坐标系之间的平移和旋转。迭代解法虽具有较高精度和鲁棒性,但计算复杂度高且对参数初值依然具有依赖性。文献[17]指出,在控制点间高差较小时,虽空间后方交会模型为良态,但其法方程系数矩阵也会存在严重的病态问题;文献[18]提出当像平面与控制点平面平行时,法方程系数矩阵的病态性将大大加强。对于病态方程,通常采用正则化方法求解[19-22]。文献[23]采用谱修正迭代算法来改善方程的病态性,用以保持法方程在解算过程中的数值稳定性;文献[9]则利用Levenberg-Marquardt(LM)方法来避免法方程系数矩阵病态时造成的不收敛问题。

    为减少线性化影响及加快收敛速度,考虑用牛顿法代替欧拉角法中的最速下降法,即对共线方程二阶泰勒近似,用二次函数模型近似原始非线性模型。非线性迭代求解的本质是在给定初值的前提下,计算梯度矩阵或海森矩阵决定下降方向以达到函数极值点,此时,梯度矩阵和海森矩阵的病态性强弱将决定法方程计算时是否考虑系数矩阵病态性的影响。为此,作为一种方法上的有效补充,本文提出一种求解单像空间后方交会的监督学习方法,在摄站或摄影基线及其附近一定范围内模拟生成样本集,并用样本集的整体下降方向来代替牛顿法中的牛顿方向,从而避免梯度矩阵和海森矩阵的直接求解,以消除牛顿法中因梯度矩阵和海森矩阵病态性对法方程求解的影响。此外,将函数凹凸性隐藏在学习过程中,同时达到增强初值选取鲁棒性的目的。

    • 设$\hat \theta $为待估计参数,V为改正数,Y为观测值,则原始误差方程为:

      $$ \mathit{\boldsymbol{V}} = f\left( {\hat \theta } \right) - \mathit{\boldsymbol{Y}} $$ (1)

      式中,f为非线性函数,非线性最小二乘法就是在式(1)的条件下,求解以下目标函数[24]

      $$ \begin{array}{*{20}{c}} {\hat \theta = \arg \min {\mathit{\boldsymbol{V}}^{\rm{T}}}\mathit{\boldsymbol{PV}} = }\\ {\arg \min \left( {{{\left( {f\left( {\hat \theta } \right) - Y} \right)}^{\rm{T}}}\mathit{\boldsymbol{P}}\left( {f\left( {\hat \theta } \right) - Y} \right)} \right)} \end{array} $$ (2)

      式中,P为权矩阵,当同精度观测时为单位矩阵,即P=I。令函数F(θ)=VTV,则式(2)可化简为:

      $$ \begin{array}{*{20}{c}} {\hat \theta = \arg \min {\mathit{\boldsymbol{V}}^{\rm{T}}}\mathit{\boldsymbol{V}} = \arg \min F\left( {\hat \theta } \right) = }\\ {\arg \min \left( {{{\left( {f\left( {\hat \theta } \right) - Y} \right)}^{\rm{T}}}\left( {f\left( {\hat \theta } \right) - Y} \right)} \right)} \end{array} $$ (3)

      牛顿法是解决式(3)中非线性最小二乘的常用工具,诸多其他方法均以其为基础进行改进。它假设在最小值的邻域内,函数F(θ)可近似为一个二次函数,即对式(3)按二阶泰勒展开:

      $$ \begin{array}{*{20}{c}} {F\left( \theta \right) \approx F\left( {\hat \theta } \right) + \nabla _F^{\rm{T}}\left( {\hat \theta } \right)\left( {\theta - \hat \theta } \right) + }\\ {\frac{1}{2}{{\left( {\theta - \hat \theta } \right)}^{\rm{T}}}{H_F}\left( {\hat \theta } \right)\left( {\theta - \hat \theta } \right)} \end{array} $$ (4)

      式中,${\nabla _F}\left( {\hat \theta } \right)$与${\mathit{\boldsymbol{H}}_F}\left( {\hat \theta } \right)$是函数F在估计值${\hat \theta }$处的梯度矩阵和海森矩阵。式(4)对(θ-${\hat \theta }$)求偏导并令其为0,若海森矩阵正定,则可得到非线性最小二乘牛顿法迭代解算更新公式:

      $$ {\hat \theta _{k + 1}} - {\hat \theta _k} = - H_F^{ - 1}\left( {{{\hat \theta }_k}} \right){\nabla _F}\left( {{{\hat \theta }_k}} \right) $$ (5)

      牛顿法只是局部收敛[25],并且在每次迭代中要计算函数在当前估计值处的梯度矩阵和海森矩阵。

      结合式(3)、式(5),将式(5)中的${\nabla _F}\left( {{{\hat \theta }_k}} \right)$按链式法则展开:

      $$ {\hat \theta _{k + 1}} - {\hat \theta _k} = - 2H_F^{ - 1}\left( {{{\hat \theta }_k}} \right){\nabla _f}\left( {{{\hat \theta }_k}} \right)\left( {f\left( {{{\hat \theta }_k}} \right) - Y} \right) $$ (6)

      为提高计算精度,借鉴最小二乘法多次迭代收敛的思想,可得到监督学习方法解决非线性最小二乘问题第k次迭代的更新公式:

      $$ {\hat \theta _{k + 1}} - {\hat \theta _k} = {R_k}\left( {f\left( {{{\hat \theta }_k}} \right) - Y} \right) $$ (7)

      式中,${\mathit{\boldsymbol{R}}_k} = - 2\mathit{\boldsymbol{H}}_F^{ - 1}({\hat \theta _k}){\nabla _f}({\hat \theta _k})$, 称为第k次迭代的整体下降方向,其通过最小化m个样本构成的训练集中所有样本第k+1次迭代的预测值和真值之间的差得到,即:

      $$ \begin{array}{l} \min \sum\limits_{i = 1}^m {{{\left\| {\theta _ * ^i - \hat \theta _{k + 1}^i} \right\|}^2}} = \\ \min \sum\limits_{i = 1}^m {{{\left\| {\theta _ * ^i - \left( {\hat \theta _k^i + {R_k}\left( {f\left( {\hat \theta _k^i} \right) - Y_ * ^i} \right)} \right.} \right\|}^2}} \end{array} $$ (8)

      式中,$\theta _*^i$为第i个样本的待估计参数真值(训练阶段已知);$\hat \theta _k^i$为第i个样本第k次迭代待估计参数的估计值;$\mathit{\boldsymbol{Y}}_*^i = f\left( {\theta _*^i} \right)$。

      由式(7)、式(8)可知,监督学习方法的解算过程如图 1所示,具体步骤为:①训练阶段,给定由多个样本组成的训练集,其中每个样本均已知待估计参数的真值与初值,监督学习方法从训练数据中学习一个使得所有样本都由初值下降到真值的整体下降方向Rk,即最小化式(8)。为提高预测精度,采用多次迭代,得到Rk,将其代入式(7)更新后得到的${\hat \theta _{k + 1}}$作为新的初值,再利用新初值学习一个新的整体下降方向Rk+1,经k次迭代后,得到整体下降方向组成的集合{Rk}。②测试阶段,给定待估计参数初值,依次按训练阶段得到的整体下降方向集合{Rk},迭代式(7)进行待估计参数的优化。

      图  1  监督学习方法示意图

      Figure 1.  Schematic Diagram of Supervised Learning Method

      监督学习方法与牛顿法解决非线性最小二乘的区别如图 2所示。图 2中,牛顿法每次的迭代方向由函数在当前参数估计值处的梯度矩阵和海森矩阵共同决定。监督学习方法是用训练阶段学习得到的整体下降方向进行优化,即是沿着初值到真值的理想方向进行,如图 2(b)所示,第0个样本第k+1次迭代中,初值到真值的理想方向和优化更新方向分别用虚线箭头和实线箭头表示。由图 2可知,实线箭头沿着虚线箭头方向进行优化。

      图  2  牛顿法和监督学习方法对比图

      Figure 2.  Comparison of Newton's Method and Supervised Learning Method

      由上述可知,监督学习方法较之牛顿法有以下优点:

      1) 无需计算函数的梯度矩阵和海森矩阵,因此,对函数模型是否可导无要求。

      2) 式(5)成立的前提是函数模型在当前参数估计值处的海森矩阵正定,但对于其海森矩阵非正定的情况,牛顿法无法直接通过式(5)迭代求解。而由式(7)、式(8)可知,监督学习方法在求解整体下降方向Rk以及预测新估计值的过程中只用到了当前估计值的函数值,所以求解过程不会受到梯度矩阵和海森矩阵病态性的影响。因此,对于牛顿法中因法方程系数矩阵病态而无法解决的问题,监督学习方法可用从训练数据中学习到的整体下降方向进行待估计参数的优化。

      3) 无初值依赖性,由式(7)可知待估计参数当前值${\hat \theta _k}$与计算值${\hat \theta _{k + 1}}$由整体下降方向Rk决定。而式(8)中Rk的计算不依赖于函数的梯度矩阵与海森矩阵,即与函数的凹凸性质无关。因此,不会出现牛顿法中因初值选取不当而造成无法收敛的情况。

    • 当使用欧拉角描述摄站姿态时,设像片内方位元素已知,将待求解外方位元素采用列向量形式表示为$\mathit{\boldsymbol{p = }}{\left( {{X_S}{Y_S}{Z_S}\varphi \omega \kappa } \right)^{\rm{T}}}$,则共线方程表示为[26]

      $$ \left\{ \begin{array}{l} {h_x}\left( p \right) = x - {x_0} = - f\frac{{{a_1}\left( {X - {X_S}} \right) + {b_1}\left( {Y - {Y_S}} \right) + {c_1}\left( {Z - {Z_S}} \right)}}{{{a_3}\left( {X - {X_S}} \right) + {b_3}\left( {Y - {Y_S}} \right) + {c_3}\left( {Z - {Z_S}} \right)}}\\ {h_y}\left( p \right) = y - {y_0} = - f\frac{{{a_2}\left( {X - {X_S}} \right) + {b_2}\left( {Y - {Y_S}} \right) + {c_2}\left( {Z - {Z_S}} \right)}}{{{a_3}\left( {X - {X_S}} \right) + {b_3}\left( {Y - {Y_S}} \right) + {c_3}\left( {Z - {Z_S}} \right)}} \end{array} \right. $$ (9)

      式中,$\left( {{x_0},{y_0},f} \right)$为像片内方位元素;(x, y)为像点坐标;(X, Y, Z)为地面坐标;(XS, YS, ZS)为摄站在地面坐标系中的坐标;${a_i},{b_i},{c_i}\left( {i = 1,2,3} \right)$为变换前后两坐标轴相应夹角的余弦。

      令$\mathit{\boldsymbol{h}} = {({h_x}{h_y})^{\rm{T}}}$,$\mathit{\boldsymbol{l}} = {(x - {x_0}y - {y_0})^{\rm{T}}}$,则可得每个点的误差方程:

      $$ \mathit{\boldsymbol{v}} = \mathit{\boldsymbol{h}}\left( {\mathit{\boldsymbol{\hat p}}} \right) - \mathit{\boldsymbol{l}} $$ (10)

      为提高外方位元素的求解精度,常有多余观测。若有n个控制点,则可依据式(10)建立n组误差方程$\mathit{\boldsymbol{V = }}{\left( {{v_1}{v_2} \cdots {v_n}} \right)^{\rm{T}}}$,因此,单像空间后方交会问题可描述为如下最小二乘过程,令$F\left( {\mathit{\boldsymbol{\hat p}}} \right) = {\mathit{\boldsymbol{V}}^{\rm{T}}}\mathit{\boldsymbol{V}}$,则:

      $$ \begin{array}{l} \mathit{\boldsymbol{\hat p}} = \arg \min F\left( {\hat p} \right) = \\ \;\;\;\;\;\;\arg \min \left( {{{\left( {H\left( {\hat p} \right) - L} \right)}^{\rm{T}}}\left( {H\left( {\hat p} \right) - L} \right)} \right) \end{array} $$ (11)

      式中,$\mathit{\boldsymbol{H = }}{\left( {{\mathit{\boldsymbol{h}}_1}{\mathit{\boldsymbol{h}}_2} \cdots {\mathit{\boldsymbol{h}}_n}} \right)^{\rm{T}}}$; $\mathit{\boldsymbol{L = }}{\left( {{\mathit{\boldsymbol{l}}_1}{\mathit{\boldsymbol{l}}_2} \cdots {\mathit{\boldsymbol{l}}_n}} \right)^{\rm{T}}}$。依据1.1节非线性最小二乘的描述可知,若函数在当前外方位元素估计值${\mathit{\boldsymbol{\hat p}}_k}$处的海森矩阵${\mathit{\boldsymbol{H}}_F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)$正定,则可得单像空间后方交会的牛顿法迭代更新公式:

      $$ {\hat p_{k + 1}} - {\hat p_k} = - H_F^{ - 1}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right){\nabla _F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) $$ (12)

      式中,${\nabla _F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)$和${\mathit{\boldsymbol{H}}_F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) = \nabla \left( {{\nabla _F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)} \right)$是函数F在当前外方位元素估计值处的梯度矩阵和海森矩阵。将${{\nabla _F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)}$用链式法则展开:

      $$ \begin{array}{*{20}{c}} {{\nabla _F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) = 2{\nabla _V}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)\mathit{\boldsymbol{V}} = }\\ {2{\nabla _H}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)\left( {H\left( {\mathit{\boldsymbol{\hat p}}} \right) - L} \right)} \end{array} $$ (13)

      式中,

      $$ \begin{array}{l} {\nabla _H}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) = \left( {{\nabla _{{h_1}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right){\nabla _{{h_2}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) \cdots {\nabla _{{h_n}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)} \right) = \\ \;\;\;\;\;\;\left( {{\nabla _{{h_{1x}}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right){\nabla _{{h_{1y}}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) \cdots {\nabla _{{h_{nx}}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right){\nabla _{{h_{ny}}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)} \right) \end{array} $$ (14)

      其中,

      $$ {\nabla _{{h_{1x}}}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) = {\left( {\frac{{\partial {h_{1x}}\left( {\mathit{\boldsymbol{\hat p}}} \right)}}{{\partial {X_S}}}\frac{{\partial {h_{1x}}\left( {\mathit{\boldsymbol{\hat p}}} \right)}}{{\partial {Y_S}}}\frac{{\partial {h_{1x}}\left( {\mathit{\boldsymbol{\hat p}}} \right)}}{{\partial {Z_S}}}\frac{{\partial {h_{1x}}\left( {\mathit{\boldsymbol{\hat p}}} \right)}}{{\partial \varphi }}\frac{{\partial {h_{1x}}\left( {\mathit{\boldsymbol{\hat p}}} \right)}}{{\partial \omega }}\frac{{\partial {h_{1x}}\left( {\mathit{\boldsymbol{\hat p}}} \right)}}{{\partial \kappa }}} \right)^{\rm{T}}} $$ (15)

      将式(13)${{\nabla _F}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)}$的值代入式(12),并令${\mathit{\boldsymbol{R}}_k} = - 2\mathit{\boldsymbol{H}}_F^{ - 1}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right){\nabla _H}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right)$,可得监督学习方法第k次迭代下单像空间后方交会解算的更新公式:

      $$ {\mathit{\boldsymbol{\hat p}}_{k + 1}} - {\mathit{\boldsymbol{\hat p}}_k} = {\mathit{\boldsymbol{R}}_k}\left( {\mathit{\boldsymbol{H}}\left( {{{\mathit{\boldsymbol{\hat p}}}_k}} \right) - L} \right) $$ (16)

      由1.1节监督学习方法解算过程可知,式(16)中的Rk可通过最小化m个样本所构成的训练集中,所有样本第k+1次迭代的预测值与真值之间的差得到,即:

      $$ \begin{array}{l} \min \sum\limits_{i = 1}^m {{{\left\| {\mathit{\boldsymbol{p}}_ * ^i - \mathit{\boldsymbol{\hat p}}_{k + 1}^i} \right\|}^2}} = \\ \;\;\;\;\;\min \sum\limits_{i = 1}^m {{{\left\| {\mathit{\boldsymbol{p}}_ * ^i - \left( {\mathit{\boldsymbol{\hat p}}_k^i + {\mathit{\boldsymbol{R}}_k}\left( {H\left( {\hat p_k^i} \right) - L_ * ^i} \right)} \right)} \right\|}^2}} \end{array} $$ (17)

      式中,p*i表示第i个样本,即第i个外方位元素的真值;$\mathit{\boldsymbol{\hat p}}_k^i$为第i个样本第k次迭代外方位元素的估计值;$\mathit{\boldsymbol{L}}_*^i = H\left( {\mathit{\boldsymbol{p}}_*^i} \right)$表示控制点在外方位元素真值p*i下的像点坐标。

      由式(16)、式(17)可知,采用监督学习方法求解单像空间后方交会的两个条件为:

      1) 训练阶段需要同一测区内多组已知外方位元素的像片组成训练集对外方位元素的整体下降方向进行学习;

      2) 因不同测区地面控制点分布一般不具有相似性,为此,训练阶段与测试阶段为同一地区同一组控制点。

      在给定一测区内一组地面控制点及像片内方位元素的前提下,考虑采用由该测区多组模拟的外方位元素作为样本组成训练集进行外方位元素整体下降方向的学习。其具体计算过程如下:

      1) 样本数据准备阶段

      ① 原始数据,包括像片内方位元素、地面控制点在地面坐标系中的坐标MRnn表示地面控制点个数。

      ② 随机产生m个外方位元素${\mathit{\boldsymbol{P}}_*} = \left( {\mathit{\boldsymbol{p}}_*^1\mathit{\boldsymbol{p}}_*^2 \cdots \mathit{\boldsymbol{p}}_*^m} \right)$,P*Rm。其中p*mR6为${\left( {{X_S}{Y_S}{Z_S}\varphi \omega \kappa } \right)^{\rm{T}}}$,表示第m个外方位元素,并通过共线方程分别计算Mm个外方位元素下的像点坐标${\mathit{\boldsymbol{U}}_*} = \left( {\mathit{\boldsymbol{L}}_*^1\mathit{\boldsymbol{L}}_*^2 \cdots \mathit{\boldsymbol{L}}_*^m} \right)$,U*R2n×m,其中$\mathit{\boldsymbol{L}}_*^m = \mathit{\boldsymbol{H}}(\mathit{\boldsymbol{p}}_*^m) \in {\bf{R}^{2\mathit{\boldsymbol{n}}}}$为${({x_1}{y_1} \cdots {x_n}{y_n})^{\rm{T}}}$,表示外方位元素p*m下的像点坐标。

      2) 训练阶段

      ① 初始化外方位元素${\mathit{\boldsymbol{P}}_k} = (\mathit{\boldsymbol{p}}_k^1\mathit{\boldsymbol{p}}_k^2 \cdots \mathit{\boldsymbol{p}}_k^m) = 0$,此时迭代次数k=0,PkRm,其中pkmR6为${\left( {{X_S}{Y_S}{Z_S}\varphi \omega \kappa } \right)^{\rm{T}}}$,表示第m个外方位元素第k次迭代的计算值。

      ② 通过式(9)计算地面控制点在当前外方位元素Pk下的像点坐标${\mathit{\boldsymbol{U}}_k} = (\mathit{\boldsymbol{L}}_k^1\mathit{\boldsymbol{L}}_k^2 \cdots \mathit{\boldsymbol{L}}_k^m)$,UkR2n×m,其中$\mathit{\boldsymbol{L}}_k^m = \mathit{\boldsymbol{H}}(\mathit{\boldsymbol{p}}_k^m) \in {\bf{R}^{2n}}$为${({x_1}{y_1} \cdots {x_n}{y_n})^{\rm{T}}}$,表示外方位元素pkm下的像点坐标,计算U*-UkP*-Pk的值。通过解线性方程组${\mathit{\boldsymbol{P}}_*} - {\mathit{\boldsymbol{P}}_k} = {\mathit{\boldsymbol{R}}_k}\left( {{\mathit{\boldsymbol{U}}_*} - {\mathit{\boldsymbol{U}}_k}} \right)$得到整体下降方向Rk,其中RkR6×2n

      ③ 用学习得到的Rk,代入式(16)计算Pk+1k=k+1返回步骤②,直到满足指定的迭代次数或达到收敛条件为止。

      ④ 输出整体下降方向集合{Rk}。

      3) 测试阶段

      ① 待估计外方位元素像片的像点观测值${\mathit{\boldsymbol{L}}_*} = {({x_1}{y_1} \cdots {x_n}{y_n})^{\rm{T}}}$。

      ② 初始化外方位元素pk=0∈R6,此时迭代次数k=0。

      ③ 由式(9)计算地面控制点在当前外方位元素pk下的像点${\mathit{\boldsymbol{L}}_k} = {(x_1^ky_1^k \cdots x_n^ky_n^k)^{\rm{T}}}$及Lk-L*。代入式(16),用{Rk}中Rk计算pk+1,并将pk+1作为外方位元素新的估计值。

      ④ 迭代次数k=k+1,返回步骤③,直至{Rk}遍历完为止。

      ⑤ 输出外方位元素估计值。

    • 为验证本文监督学习方法解决单像空间后方交会的有效性,分别进行病态性试验与初值依赖性试验。

    • 假设存在理想情况,如图 3所示。已知摄像机焦距f为35 mm,像平面坐标系的o'-xy平面与世界坐标系O-XYZO-XY平面平行。航高为500 m,摄影中心S在世界坐标系的坐标为(0,0,150),单位为m,像片外方位角元素的真实值为(0,0,0),单位为rad。假设地面控制点的坐标分别为A(-50,50,0),B(-50,-50,0),C(50,-50,0),D(50,50,0),单位为m,则根据三角形相似性,可得对应像点坐标a(-3.5,3.5),b(-3.5,-3.5),c(3.5,-3.5),d(3.5,3.5),单位为mm。根据地面控制点及其对应像点坐标,采用欧拉角法计算像片的外方位元素。由于像平面与控制点平面平行,导致欧拉角法中的法方程系数矩阵病态,因此,无法直接通过法方程解算此时的外方位元素[17]

      图  3  理想垂直摄影示意图

      Figure 3.  Schematic Diagram of Ideal Vertical Photography

      以线元素均值为0、方差为0.1,以角元素均值为0、方差为0.05的高斯噪声,在外方位元素真实值(0,0,150,0,0,0)的基础上扰动产生多个外方位元素样本组成监督学习方法的训练集,训练整体下降方向并求解以上问题中的外方位元素,计算结果如表 1所示。

      表 1  理想垂直摄影外方位元素计算结果

      Table 1.  Result of Calculation of Ideal Vertical Photography

      元素 项目 真值 欧拉角法 监督学习方法
      线元素/m XS 0.000 000 不收敛 0.001 565
      YS 0.000 000 -0.002 479
      ZS 150.000 000 150.003 025
      角元素/rad φ 0.000 000 不收敛 0.000 009
      ω 0.000 000 -0.000 002
      κ 0.000 000 0.000 012

      表 1可知,当欧拉角法中因法方程系数矩阵病态而无法求解时,监督学习方法仍能准确解算,且能保持较高的计算精度。因为其通过训练集来学习获得整体下降方向,不受到梯度矩阵或海森矩阵病态的限制。

    • 给定地面控制点坐标,如表 2所示。

      表 2  地面控制点坐标

      Table 2.  Coordinates of GCPs

      点号 X/m Y/m Z/m
      1 10 10 11
      2 160 40 24
      3 320 40 34
      4 10 160 11
      5(A) 160 160 0
      6 320 160 28
      7 40 320 14
      8 160 320 29
      9 320 320 27

      试验方案如图 4所示。已知像片参数为:镜头焦距f为24 mm,相对航高H为114 m,x0=y0=0,以表 2中第5个控制点(即图 4A)为被摄地面的中心点,坐标为A(XA, YA, 0)。基于计算方便及像片上控制点全部可见两点考虑,设摄影基线的方程为:Y=YA, Z=H

      图  4  多重交向摄影试验示意图

      Figure 4.  Schematic Diagram of Multi-intersection Photography

      为保证定位精度最佳,选取φ∈[-45°, 45°]。则图 4BCD的外方位元素分别为B(XA YA H 0 0 0),C(XA - H YA H φmax 0 0),D(XA + H YA H φmin 0 0)。以B为基准,沿摄影基线以Δx=±13 m为步长、CD为终点产生如图 4所示的系列基准摄站,每个摄像机光轴对准地面中心点A。即B左边第i个摄站外方位元素为(XA-i×Δx YA H、$\arctan \frac{{i \times \Delta x}}{H}$ 0 0),右边第i个摄站外方位元素为(XA+i×Δx YA H-$\arctan \frac{{i \times \Delta x}}{H}$ 0 0)。再以每个基准摄站为基准,线元素增量范围为$\left[ { - \frac{{\Delta x}}{2},\frac{{\Delta x}}{2}} \right]$,角元素增量范围为[-5°, 5°],随机产生新的外方位元素。试验所有基准摄站共产生了1 000个外方位元素作为训练集进行训练。

      测试阶段,为使训练得到的整体下降方向集合能准确解算测试阶段的外方位元素,假设测试阶段的外方位线元素与角元素的最大范围不超过训练阶段的训练样本。具体为:以Δx=±8 m产生系列基准摄站,在每个基准摄站的基础上,以线元素增量范围为$\left[ { - \frac{{\Delta x}}{2},\frac{{\Delta x}}{2}} \right]$,角元素增量范围为[-3°, 3°], 随机产生新的外方位元素作为测试样本。试验共产生301个外方位元素。

      计算地面控制点在每个测试样本下的像点坐标,采用训练阶段得到的整体下降方向计算得到外方位元素的计算值,并与其随机产生的真实值对比,验证本文方法的正确性。考虑到计算效率,采用Eigen进行编程实现,训练过程花费时间为2 521 ms,平均计算每个测试样本的时间花费约为5 ms。为方便比较本文方法的准确性,将计算结果投影至二范数空间进行比较。为验证本文方法的初值依赖性,将训练阶段与测试阶段外方位元素的初值设为0,即$X_S^0 = Y_S^0 = Z_S^0 = 0$,${\varphi ^0} = {\omega ^0} = {\kappa ^0} = 0$。剔除计算错误的点后,所训练的整体下降方向准确解算测试样本数达283个,占全部测试样本的94%,计算结果如图 5所示。

      图  5  计算值与真值对比曲线图

      Figure 5.  Comparison Curve Between Calculated Value and True Value

      图 5可知,计算值曲线与真值曲线基本重合。经统计,测试集中摄站位置最大差值为0.016 m,摄影姿态最大差值为0.001°,满足测量精度要求[15],表明该方法的正确性,可用于求解单像空间后方交会。

      为验证本文方法的初值依赖性,分别采用本文方法、欧拉角法和谱修正方法解算上述的283个测试样本,并将计算结果投影至二范数空间进行对比。由于试验采用交向摄影,部分倾角较大,无法直接采用欧拉角法进行计算,因此,采用直接线性变换方法提供初值,3种方法的计算结果误差曲线如图 6所示,迭代次数如图 7所示。

      图  6  3种方法计算结果对比曲线图

      Figure 6.  Comparison of the Results of Three Methods

      图  7  3种方法迭代次数对比图

      Figure 7.  Comparison of Iterations of Three Methods

      图 6图 7可知,在计算精度上,欧拉角法、监督学习方法与谱修正方法对于外方位元素的求解精度相当。但欧拉角法依赖于迭代初值,会因迭代初值的选取不当而无法收敛。因此,在不提供初值的前提下,较难解算大姿态下像片的外方位元素。在迭代次数上,监督学习方法总体优于其他两种算法,具有较快的收敛速度。虽然监督学习方法的迭代更新公式是在牛顿法的基础上推导的,但较之牛顿法却无较强的初值依赖且能保持较快的收敛速度。其原因在于,监督学习方法是通过大量的样本,在已知外方位元素初值和真值的前提下,通过每次迭代中最小化所有样本当前外方位元素初值到真值之差,最终得到整体下降方向集合{Rk},并以此解算测试阶段待求像片的外方位元素。因此,只需测试样本与训练样本近似,且外方位元素的初值相同,则利用训练阶段得到的整体下降方向集合就可以准确解算待求的外方位元素。

      从测试数据解算结果可知,监督学习方法具有解算精度较高、收敛速度较快及初值依赖性弱3个特点,但在测试阶段要求待求像片的所有像点均可见,且要求在同一航线上对测区进行不同倾角的重复观测,需采用多重交向摄影的方式进行试验。

    • 本文提出了一种求解单像空间后方交会的监督学习方法,并在牛顿法的基础上推导了迭代更新公式。对比试验验证了监督学习方法克服了法方程系数矩阵病态、快速收敛及初值依赖性弱的缺点,与此同时,能够保持较高的解算精度与效率。由于监督学习方法是对摄影基线上系列基准摄站扰动产生训练样本以获得整体下降方向,因此,目前只能将其应用于摄影基线附近一定范围内的像片外方位元素的解算。文中试验采用的摄影基线比较简单,较为复杂的摄影基线下外方位元素的解算有待进一步研究。此外,采用交向摄影构成的测量网与常规航空摄影测量不同,并且线阵CCD(charge coupled devices)比面阵CCD的姿态运动方程复杂,能否将本文方法分别应用于光束法平差和线阵CCD像片的外方位元素求解是今后的研究方向。

参考文献 (26)

目录

    /

    返回文章
    返回