留言板

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

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

利用三视匹配元进行多视影像批处理重建

卢俊 张保明 郭海涛 张宏伟

卢俊, 张保明, 郭海涛, 张宏伟. 利用三视匹配元进行多视影像批处理重建[J]. 武汉大学学报 ● 信息科学版, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
引用本文: 卢俊, 张保明, 郭海涛, 张宏伟. 利用三视匹配元进行多视影像批处理重建[J]. 武汉大学学报 ● 信息科学版, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
LU Jun, ZHANG Baoming, GUO Haitao, ZHANG Hongwei. A Batch Reconstruction Algorithm of Multi-view Images Using Image Triplets[J]. Geomatics and Information Science of Wuhan University, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
Citation: LU Jun, ZHANG Baoming, GUO Haitao, ZHANG Hongwei. A Batch Reconstruction Algorithm of Multi-view Images Using Image Triplets[J]. Geomatics and Information Science of Wuhan University, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672

利用三视匹配元进行多视影像批处理重建

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

国家973计划 2012CB720000

详细信息
    作者简介:

    卢俊, 博士生, 讲师, 研究方向为数字摄影测量和计算机视觉。ljhb45@126.com

  • 中图分类号: P237.3;TP751

A Batch Reconstruction Algorithm of Multi-view Images Using Image Triplets

Funds: 

The National 973 Program of China 2012CB720000

More Information
    Author Bio:

    LU Jun, PhD candidate, lecturer, specializes in digital photogrammetry and computer vision. E-mail: ljhb45@126.com

  • 摘要: 无序多视影像的三维重建对噪声非常敏感,错误的匹配关系会影响重建的精度,甚至直接导致重建失败。提出了一种稳健的批处理重建算法,首先利用回路闭合约束剔除可能存在误匹配的三视匹配元,然后以三视匹配元中的三焦张量约束代替传统算法的核线约束来计算所有影像旋转矩阵和相机中心位置的全局最优解。重建过程中引入高效的并查集算法来提取多视匹配点,并利用迭代线性三角形算法计算空间点的三维坐标。实验结果表明,所提算法在重建效率和计算精度方面都能取得较好的结果。
  • 图  1  三视匹配元中的几何关系

    Figure  1.  Geometric Relationship in a Triplet

    图  2  基于并查集的多视匹配点提取

    Figure  2.  Multi-view Matching Points Extraction Based on UF Method

    图  3  回路闭合差统计

    Figure  3.  Closed Cycle Error Statistics

    表  1  批处理计算结果

    Table  1.   Calculation Results of Batch Method

    影像集 影像数量 三视匹配
    元数量
    回路闭合
    检测后
    多视匹配
    点数
    平差前δ 第一次平差 第二次平差
    迭代次数 δ 迭代次数 δ
    FountainP11 11 69 69 3 717 0.726 233 2 0.277 844 2 0.268 669
    HerzJesusP25 25 421 400 5 251 0.756 457 2 0.517 371 4 0.498 520
    CastleP30 30 491 362 6 573 1.544 170 2 0.450 979 3 0.321 487
    下载: 导出CSV

    表  2  计算时间

    Table  2.   Computing Time of Batch Method and Bundle

    影像集 t1 t2 t3 t4 t5 t6 t7 t Bundler
    FountainP11 0.07 15.23 0.02 0.06 1.11 3.77 5.43 25.69 37.35
    HerzJesusP25 0.12 29.09 0.19 0.27 1.85 9.41 14.79 55.72 176.81
    CastleP30 0.16 31.85 0.32 0.29 2.37 11.41 14.54 60.94 470.34
    下载: 导出CSV

    表  3  影像集精度统计

    Table  3.   PrecisionStatistics of 3 Image Sets

    影像集 ΔX/m ΔY/m ΔZ/m ΔC/m ΔRoll/(°) ΔPitch/(°) ΔYaw/(°) ΔR/(°)
    FountainP11 最大值 0.003 496 0.002 420 0.003 591 0.004 506 0.043 5 0.031 5 0.039 0 0.074 8
    均值 0.001 543 0.000 945 0.002 030 0.002 947 0.018 1 0.018 5 0.013 9 0.029 1
    均方根 0.001 832 0.001 230 0.002 267 0.003 164 0.021 9 0.020 1 0.017 6 0.039 0
    HerzJesusP25 最大值 0.008 147 0.007 807 0.010 924 0.012 594 0.306 2 0.041 9 0.337 5 0.095 6
    均值 0.002 805 0.003 363 0.004 790 0.007 199 0.060 4 0.014 9 0.081 2 0.039 1
    均方根 0.003 536 0.004 007 0.005 707 0.007 818 0.093 0 0.019 3 0.109 5 0.049 1
    CastleP30 最大值 0.051 092 0.047 150 0.028 073 0.067 435 0.384 0 0.101 0 0.357 1 0.119 0
    均值 0.013 625 0.014 809 0.012 217 0.026 935 0.063 7 0.035 8 0.047 4 0.052 1
    均方根 0.018 333 0.018 451 0.014 798 0.029 925 0.096 7 0.042 2 0.083 3 0.062 2
    下载: 导出CSV

    表  4  精度比较

    Table  4.   Precision Comparison

    FountainP11 HerzJesusP25 CastleP30
    ΔC/m ΔR/(°) ΔC/m ΔR/(°) ΔC/m ΔR/(°)
    文献[13] 0.007 2 0.112 0.030 8 0.110 0.300 0 0.279
    文献[18] 0.009 9 0.116 0.023 3 0.104 0.264 0 0.398
    文献[19] 0.131 7 - 0.253 8 - - -
    文献[4] 0.015 3 - 0.084 5 - - -
    文献[8] 0.014 0 0.195 0.064 0 0.188 0.235 0 0.480
    本文算法 0.003 2 0.039 0.007 8 0.049 0.029 9 0.062
    下载: 导出CSV
  • [1] Hartley R I, Trumpf J, Dai Y, et al. Rotation Averaging[J]. International Journal of Computer Vision, 2013, 103(3):267-305 doi:  10.1007/s11263-012-0601-0
    [2] Arie-Nachimson M, Kovalsky S Z, Kemelmacher-Shlizerman I, et al. Global Motion Estimation from Point Matches[C]. 3D Imaging, Modeling, Processing, Visualization and Transmission (3DIMPVT), Zurich, 2012
    [3] 郭复胜, 高伟.基于辅助信息的无人机图像批处理三维重建方法[J].自动化学报, 2013, 39(6):834-845 doi:  10.1016/S1874-1029(13)60056-7

    Guo Fusheng, Gao Wei. Batch Reconstruction from UAV Images with Prior Information[J]. Acta Automatica Sinica, 2013, 39(6):834-845 doi:  10.1016/S1874-1029(13)60056-7
    [4] Martinec D, Pajdla T. Robust Rotation and Translation Estimation in Multiview Reconstruction[C]. Computer Vision and Pattern Recognition, Minnesota, New York, 2007
    [5] Kahl F, Hartley R I. Multiple-view Geometry under the l∞-norm[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2008, 30(9):1603-1617 doi:  10.1109/TPAMI.2007.70824
    [6] 郑顺义, 张祖勋, 翟瑞芳.基于非量测相机的复杂物体三维重建[J].武汉大学学报·信息科学版, 2008, 33(5):446-449 http://ch.whu.edu.cn/CN/abstract/abstract1580.shtml

    Zheng Shunyi, Zhang Zuxun, Zhai Ruifang. 3D Reconstruction of Complex Objects Based on Non-metric Image[J]. Geomatics and Information Science of Wuhan University, 2008, 33(5):446-449 http://ch.whu.edu.cn/CN/abstract/abstract1580.shtml
    [7] Hartley R I, Zisserman A. Multiple View Geometry in Computer Vision[M]. New York:Cambridge University Press, 2003
    [8] Jiang Nianjuan, Cui Zhaopeng, Tan Ping. A Global Linear Method for Camera Pose Registration[C]. IEEE International Conference on Computer Vision, Sydney, 2013
    [9] Hartley R I, Sturm P. Triangulation[J]. Computer Vision and Image Understanding, 1997, 68(2):146-157 doi:  10.1006/cviu.1997.0547
    [10] Moulon P, Monasse P, Marlet R. Global Fusion of Relative Motions for Robust, Accurate and Scalable Structure from Motion[C]. IEEE International Conference on Computer Vision, Sydney, 2013
    [11] Fasano P, Forrest J. The COIN-OR Open Solver Interface:Technology Overview[EB/OL].https://projects.coin-or.org/Osi,2014
    [12] Galler B A, Fischer M J. An Improved Equivalence Algorithm[J]. Communications of the ACM, 1964, 7(5):301-303 doi:  10.1145/364099.364331
    [13] Snavely, Seitz S M, Szeliski R. Modeling The World from Internet Photo Collections[J]. International Journal of Computer Vision, 2008, 80(2):189-210 doi:  10.1007/s11263-007-0107-3
    [14] Lourakis M I A, Argyros A A. Sba:A Software Package for Generic Sparse Bundle Adjustment[J]. ACM Transactions on Mathematical Software, 2009, 36(1):1-30 https://dl.acm.org/purchase.cfm?id=1486527
    [15] Strecha C, van Hansen W, van Gool L J, et al. On Benchmarking Camera Calibration and Multi-View Stereo for High Resolution Imagery[C].Computer Vision and Pattern Recognition, Alaska, New York, 2008
    [16] Lowe D. Distinctive Image Features from Scale-Invariant Keypoints[J]. International journal of computer vision, 2004, 60(2):91-110 doi:  10.1023/B:VISI.0000029664.99615.94
    [17] Nistér D. An Efficient Solution to the Five-Point Relative Pose Problem[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2004, 26(6):756-770 doi:  10.1109/TPAMI.2004.17
    [18] Wu C. Towards Linear-Time Incremental Structure from Motion[C].3DTV-Conference, Aberdeen, 2013
    [19] Sinha S, Steedly D, Szeliski R. A Multi-stage Linear Approach to Structure from Motion[M]. Berlin:Springer, 2012
  • [1] 张红, 徐珊, 龚恩慧.  顾及实时路况的城市浪费性通勤测算 . 武汉大学学报 ● 信息科学版, 2021, 46(5): 650-658. doi: 10.13203/j.whugis20190363
    [2] 王玮琦, 游雄, 杨剑, 李钦.  一种改进匹配点对选取策略的ElasticFusion室内三维重建算法 . 武汉大学学报 ● 信息科学版, 2020, 45(9): 1469-1477. doi: 10.13203/j.whugis20180278
    [3] 姚国标, 邓喀中, 艾海滨, 杜全叶.  倾斜立体影像自动准稠密匹配与三维重建算法 . 武汉大学学报 ● 信息科学版, 2014, 39(7): 843-849.
    [4] 袁理, 陈庆虎, 廖海斌, 段柳云.  单视影像下的人脸快速三维重建 . 武汉大学学报 ● 信息科学版, 2012, 37(4): 487-491.
    [5] 孔云峰.  利用GIS与线性规划学校最优学区划分 . 武汉大学学报 ● 信息科学版, 2012, 37(5): 513-515.
    [6] 徐卫明, 暴景阳, 陆建华.  一种基于空时联合约束的相对定位方法 . 武汉大学学报 ● 信息科学版, 2009, 34(4): 495-498.
    [7] 张剑清, 孙明伟, 郑顺义, 季铮.  基于轮廓约束的摄影测量法元青花瓶数字三维重建 . 武汉大学学报 ● 信息科学版, 2009, 34(1): 7-10.
    [8] 杨小雄, 刘耀林, 王晓红, 段滔.  基于约束条件的元胞自动机土地利用规划布局模型 . 武汉大学学报 ● 信息科学版, 2007, 32(12): 1164-1167.
    [9] 郑肇葆.  基于Bayesian线性规划的影像纹理识别方法 . 武汉大学学报 ● 信息科学版, 2007, 32(3): 193-196.
    [10] 朱庆, 李逢春, 张叶廷.  一种改进的三维点集表面重建的区域生长算法 . 武汉大学学报 ● 信息科学版, 2006, 31(8): 667-670.
    [11] 李订芳, 章文, 牛艳庆.  基于粗糙集的缺失数据循环搜索重建算法 . 武汉大学学报 ● 信息科学版, 2005, 30(11): 1016-1019.
    [12] 邱卫宁.  具有稳健初值的选权迭代法 . 武汉大学学报 ● 信息科学版, 2003, 28(4): 452-454.
    [13] 耿红, 王泽民.  基于灰色线性规划的土地利用结构优化研究 . 武汉大学学报 ● 信息科学版, 2000, 25(2): 167-171.
    [14] 郑肇葆.  基于MRF参数的纹理识别的线性规划法 . 武汉大学学报 ● 信息科学版, 1996, 21(3): 228-231.
    [15] 王新生.  线性规划原理在城市道路控制点标高优化设计中的应用 . 武汉大学学报 ● 信息科学版, 1995, 20(3): 269-272.
    [16] 赵少荣.  边角网二类优化设计 . 武汉大学学报 ● 信息科学版, 1991, 16(1): 37-48.
    [17] 郑肇葆, 谭春健, 张润根.  按残差绝对值和最小原理的光束法区域网平差 . 武汉大学学报 ● 信息科学版, 1990, 15(2): 34-40.
    [18] 高扬.  高精度工程控制网的二类优化设计研究 . 武汉大学学报 ● 信息科学版, 1987, 12(1): 92-103.
    [19] 郑肇葆.  线性规划在摄影测量粗差检测中应用的尝试 . 武汉大学学报 ● 信息科学版, 1985, 10(3): 50-60.
    [20] 高贤君, 冉树浩, 张广斌, 杨元维.  基于多特征融合与对象边界联合约束网络的建筑物提取 . 武汉大学学报 ● 信息科学版, 0, 0(0): -. doi: 10.13203/j.whugis20210520
  • 加载中
图(3) / 表(4)
计量
  • 文章访问数:  1393
  • HTML全文浏览量:  90
  • PDF下载量:  403
  • 被引次数: 0
出版历程
  • 收稿日期:  2014-12-16
  • 刊出日期:  2017-01-05

利用三视匹配元进行多视影像批处理重建

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

    国家973计划 2012CB720000

    作者简介:

    卢俊, 博士生, 讲师, 研究方向为数字摄影测量和计算机视觉。ljhb45@126.com

  • 中图分类号: P237.3;TP751

摘要: 无序多视影像的三维重建对噪声非常敏感,错误的匹配关系会影响重建的精度,甚至直接导致重建失败。提出了一种稳健的批处理重建算法,首先利用回路闭合约束剔除可能存在误匹配的三视匹配元,然后以三视匹配元中的三焦张量约束代替传统算法的核线约束来计算所有影像旋转矩阵和相机中心位置的全局最优解。重建过程中引入高效的并查集算法来提取多视匹配点,并利用迭代线性三角形算法计算空间点的三维坐标。实验结果表明,所提算法在重建效率和计算精度方面都能取得较好的结果。

English Abstract

卢俊, 张保明, 郭海涛, 张宏伟. 利用三视匹配元进行多视影像批处理重建[J]. 武汉大学学报 ● 信息科学版, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
引用本文: 卢俊, 张保明, 郭海涛, 张宏伟. 利用三视匹配元进行多视影像批处理重建[J]. 武汉大学学报 ● 信息科学版, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
LU Jun, ZHANG Baoming, GUO Haitao, ZHANG Hongwei. A Batch Reconstruction Algorithm of Multi-view Images Using Image Triplets[J]. Geomatics and Information Science of Wuhan University, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
Citation: LU Jun, ZHANG Baoming, GUO Haitao, ZHANG Hongwei. A Batch Reconstruction Algorithm of Multi-view Images Using Image Triplets[J]. Geomatics and Information Science of Wuhan University, 2017, 42(1): 109-115. doi: 10.13203/j.whugis20140672
  • 多视三维重建的核心是在一个统一的全局坐标系中计算所有影像的绝对方位(包括绝对旋转和相机中心坐标),处理策略大体可以分为增量式算法、分步式算法和批处理算法(或全局算法)。增量式和分步式重建算法需要进行光束法平差的迭代运算,计算效率低下,稳定性不高,并且重建结果可能会产生漂移;而批处理重建算法不依赖于迭代式的优化构架,具有很高的计算效率和精度,近年来成为了三维重建领域的研究热点[1-3]

    由于相机旋转的计算精度和稳定性要高于相机位移,批处理重建算法通常将两者分开求解[1]。旋转平均是解算相机旋转矩阵的高效算法,能获得很高的精度,但受噪声影响较大。二阶锥规划(second order cone programming, SOCP)算法[3-5]将相机中心坐标的求解转化为l范数下的最小化问题,但对噪声非常敏感。文献[6]利用旋转平台获取的多视影像实现物体的三位重建,但算法依赖于先验信息。本文提出了一种稳健的批处理三维重建算法,首先利用回路闭合约束剔除可能存在误匹配的三视匹配元,然后计算每个三视匹配元中的一致性方位关系,最后以三视匹配元为基础,利用旋转平均算法计算旋转矩阵,通过最小化欧氏距离来计算相机中心位置的全局最优解。

    • 如果影像集S中两张影像ij之间存在足够多的匹配点,则认为这两张影像相关。S中的所有相关影像构成一个无向图G(VE),其中V={ < i > |i=1, 2, …, n}代表顶点集,E={ < i, j>|i, j=1, 2, …, n}代表边集。图G中的每一个顶点代表一张影像,每一条边代表两张影像间的匹配关系。如果影像ijk之间同时存在边e < i, j>、e < j, k>和e < i, k>,并且具有足够多的三视匹配点,则影像ijk构成三视匹配元τ < i, j, k>,影像集中的所有三视匹配元构成集合W= < i, j, k>|i, j, k=1, 2, …, n

    • G中的每一对相关影像ij,利用相对定向算法可以计算它们之间的相对旋转矩阵Rij和相对位移方向向量tij,假设相关影像ij对应的绝对旋转矩阵分别为RiRj(RiRj都为正交阵),则:

      $$ {\mathit{\boldsymbol{R}}_j} = {\mathit{\boldsymbol{R}}_{ij}}{\mathit{\boldsymbol{R}}_i}(\langle i, j\rangle \in E) $$ (1)
      $$ {\mathit{\boldsymbol{R}}_i}\mathit{\boldsymbol{R}}_i^T = \mathit{\boldsymbol{I}}, \left| {{\mathit{\boldsymbol{R}}_i}} \right| = 1(i \in V) $$ (2)

      Ri=[ri1 ri2  ri3],Rj=[rj1 rj2 rj3],则式(1)可表达为:

      $$ \mathit{\boldsymbol{r}}_j^k - {\mathit{\boldsymbol{R}}_{ij}}\mathit{\boldsymbol{r}}_i^k = {\left[ {\begin{array}{*{20}{c}} 0&0&0 \end{array}} \right]^{\rm{T}}} $$ (3)

      式中,k=1,2,3;rikrjk分别表示RiRj的第k列。令R1为单位阵,通过最小二乘方法可以计算任意一张影像i的近似一致性矩阵${\mathit{\boldsymbol{\hat R}}}$i,正交性约束式(2)可以利用奇异值分解(singular value decomposition,SVD)算法,将计算的矩阵进行投影变换,得到Frobenius范数下最接近的旋转矩阵[4]。首先对${\mathit{\boldsymbol{\hat R}}} $i进行SVD分解得到:

      $$ {{\mathit{\boldsymbol{\hat R}}}_i} = {\mathit{\boldsymbol{U}}_i}{\mathit{\boldsymbol{ \boldsymbol{\varSigma} }}_i}\mathit{\boldsymbol{V}}_{^i}^H $$ (4)

      式中, UiVi为3×3酉矩阵;Σi为对角阵,由${\mathit{\boldsymbol{\hat R}}}$i的奇异值组成。对${\mathit{\boldsymbol{\hat R}}} $i进行强制正交化可得:

      $$ {\mathit{\boldsymbol{R}}_i} = {\mathit{\boldsymbol{U}}_i}\mathit{\boldsymbol{V}}_{^i}^H $$ (5)
    • 首先利用§1.4节并查集算法搜索G中的所有三视匹配点数量不少于20的三视匹配元,构成集合W={〈i, j, k〉|i, j, k=1, 2,…,n},W中的每一个三视匹配元τi, j, k〉都构成一个闭合回路,理论上存在回路闭合约束:

      $$ {\mathit{\boldsymbol{R}}_{ki}}{\mathit{\boldsymbol{R}}_{jk}}{\mathit{\boldsymbol{R}}_{ij}} = \mathit{\boldsymbol{I}} $$ (6)

      但由于噪声的影响以及数值计算中的舍入误差,式(6)并不一定能严格成立。以RkiRjkRij和单位阵I之间的几何距离ΔRτ来表示三视匹配元的回路闭合差,即:

      $$ \Delta {\mathit{\boldsymbol{R}}^\tau } = {\rm{arccos}}(12({\rm{tr}}\left( {{\mathit{\boldsymbol{R}}_{ki}}{\mathit{\boldsymbol{R}}_{jk}}{\mathit{\boldsymbol{R}}_{ij}}} \right) - 1)) $$ (7)

      理论上ΔRτ的值为0,如果ΔRτ过大,则三视匹配元τi, j, k〉中可能存在误匹配,本文剔除W中的所有回路闭合差大于2°的三视匹配元。接下来计算每个三视匹配元τi, j, k〉中的局部一致性旋转矩阵和相机中心坐标。

      τ中建立一个局部坐标系,令其与影像k的相机坐标系一致,即Rkτ为单位阵,利用§1.1节中的算法可以求出RiτRjτ。假设三幅影像的投影中心在该局部坐标系下的坐标向量分别为CiτCjτCkτ(其中取Ckτ=[0  0  0]T),tij表示相关影像(ij)中的相对基线方向向量,Tijτ表示局部坐标系中基线方向向量,则存在关系式[7]

      $$ \mathit{\boldsymbol{T}}_{ij}^\tau = {\left( {\mathit{\boldsymbol{-R}}_j^\tau } \right)^{\rm{T}}}{\mathit{\boldsymbol{t}}_{ij}} $$ (8)

      通过相对定向可以计算tijtiktjk,然后利用式(8)可以算出TijτTikτTjkτ,理论上它们应该共面,但由于噪声的影响,通常不能满足共面条件,如图 1所示。

      图  1  三视匹配元中的几何关系

      Figure 1.  Geometric Relationship in a Triplet

      先假设Tijτ没有误差,线段ABTikτTjkτ的公垂线,并且和它们的交点分别为点A和点B,则Ckτ的最优解理论上应该位于线段AB的中点。利用文献[8]中最小化近似几何误差算法可得线性方程:

      $$ 2\mathit{\boldsymbol{C}}_k^\tau - \mathit{\boldsymbol{C}}_i^\tau - \mathit{\boldsymbol{C}}_j^\tau = \mathit{\boldsymbol{R}}_i^\tau \left( {{{\theta '}_j}} \right)\rho _{ij}^{ik}\left( {\mathit{\boldsymbol{C}}_j^\tau - \mathit{\boldsymbol{C}}_i^\tau } \right) + \mathit{\boldsymbol{R}}_j^\tau \left( { - {{\theta '}_j}} \right)\rho _{ij}^{jk}(\mathit{\boldsymbol{C}}_i^\tau - \mathit{\boldsymbol{C}}_j^\tau ) $$ (9)

      式中,θi为向量TikτTijτ的夹角; θj为向量TijτTjkτ的夹角; Riτ(φ)代表旋转矩阵Riτ绕轴Tijτ×Tikτ逆时针旋转φ角; ρijik代表基线‖Ciτ-Ckτ‖与基线‖Cjτ-Ckτ‖的比值。对于在三张影像上都可见的空间点X,假设在影像ij的基线长度为1时,通过影像(ij)重建得到其在影像i中的深度值为diij,在影像ik的基线长度为1时,通过影像(ij)重建得到其在影像i中的深度值为diik,则基线比为ρijjk=diij/diik。由于存在多个三视可见的空间点,基线比可在随机抽样一致(random sample consensus,RANSAC)引擎中通过计算平均值获得。

      同理, 假设TikτTjkτ没有误差时,可以得到:

      $$ 2\mathit{\boldsymbol{C}}_j^\tau-\mathit{\boldsymbol{C}}_i^\tau-\mathit{\boldsymbol{C}}_k^\tau = \mathit{\boldsymbol{R}}_i^\tau \left( {{{\theta '}_i}} \right)\rho _{ik}^{ij}\left( {\mathit{\boldsymbol{C}}_k^\tau-\mathit{\boldsymbol{C}}_i^\tau } \right) + \mathit{\boldsymbol{R}}_k^\tau \left( { - {{\theta '}_k}} \right)\rho _{ik}^{jk}(\mathit{\boldsymbol{C}}_i^\tau - \mathit{\boldsymbol{C}}_k^\tau ) $$ (10)
      $$ 2\mathit{\boldsymbol{C}}_i^\tau-\mathit{\boldsymbol{C}}_j^\tau-\mathit{\boldsymbol{C}}_k^\tau = \mathit{\boldsymbol{R}}_j^\tau \left( {{{\theta '}_j}} \right)\rho _{jk}^{ij}\left( {\mathit{\boldsymbol{C}}_k^\tau-\mathit{\boldsymbol{C}}_j^\tau } \right) + \mathit{\boldsymbol{R}}_k^\tau \left( { - {{\theta '}_k}} \right)\rho _{jk}^{jk}(\mathit{\boldsymbol{C}}_j^\tau - \mathit{\boldsymbol{C}}_k^\tau ) $$ (11)

      解算线性方程组式(9)、(10)、(11)可以得到每个三视匹配元中的局部一致性相机位置。

      通过三视匹配元τi, j, k〉中的局部一致性相机方位,本文利用三角形法[9]计算τ中所有三视匹配点{(xvu, yvu)}u∈{i, j, k}的空间坐标Xv,其中v代表点的序号,如果Xv的重投影误差大于5,则在τ中剔除该三视匹配点,如果三视匹配元中剩余的三视匹配点数量少于10,则在W中将该三视匹配元剔除。然后利用G中所有顶点和W中包含的边构建新的无向图G′(V′, E′),顶点集V′与G中顶点集一致,边集E′由W中所有三视匹配元包含的边构成。由于G′不一定能构成完整的联通关系,利用深度优先算法可以构建其中的最大连通图G(V, E),如果W中任意一个三视匹配元的任意一个顶点不属于V,则在W中剔除该三视匹配元。

    • G中的每一条边ei, j〉,寻找W中所有包含该边且具有最多空间点的三视匹配元τi, j, k〉,构建边和三视匹配元的对应列表U=〈e, τ〉。

      G中的每一张影像iV,首先利用§1.1算法计算其绝对旋转矩阵Ri,接下来利用G中所有匹配像对间的相对方位(Rij, tij)以及影像的绝对旋转矩阵Ri来计算每张影像的投影中心在全局坐标系中的坐标向量Ci。记影像(ij)投影中心间的绝对位移为Cij=Cj-CiTij表示全局坐标系中的基线方向向量,由于相对定向过程中基线方向向量只能确定到相差一个尺度因子,所以理论上TijCij成比例,即:

      $$ {\lambda _{ij}}{\mathit{\boldsymbol{T}}_{ij}} = {\mathit{\boldsymbol{C}}_{ij}} $$ (12)

      其中,λij为尺度系数。如果将§1.2节中局部坐标系改为全局坐标系,则式(8)可变换为Tij=-RjTtij,代入式(12)可得:

      $$ {\mathit{\boldsymbol{R}}_j}{\mathit{\boldsymbol{C}}_j}-{\mathit{\boldsymbol{R}}_j}{\mathit{\boldsymbol{C}}_i} + {\lambda _{ij}}{\mathit{\boldsymbol{t}}_{ij}} = 0 $$ (13)

      二阶锥规划算法的实质是将式(13)转换成最小化向量RjCi-RjCj和向量λijtij的角度误差,从而构建了一个拟凸的非线性系统,拟凸系统对噪声比较敏感,并且不能保证所有的尺度系数λij>0,即不一定能满足Cheirality约束[7]。为强制实现Cheirality约束,本文首先假设任意两张影像间的尺度系数不小于1,利用类似文献[10]中最小化欧氏距离的算法将式(13)转换为l范数下的线性优化问题:

      $$ \left| {{\mathit{\boldsymbol{R}}_j}{\mathit{\boldsymbol{C}}_j}-{\mathit{\boldsymbol{R}}_j}{\mathit{\boldsymbol{C}}_i} + {\lambda _{ij}}{\mathit{\boldsymbol{t}}_{ij}}} \right| \le \eta $$ (14)

      即在λij≥1条件下,对G中所有的相关影像,计算CiCjλij以及最小化参数η,令C1=[0  0  0]T,利用线性规划算法可以计算相机中心的全局最优解。本文利用三视匹配元中的约束关系来代替立体像对,通过对应列表U=〈e, τ〉,对任意一条边ei, j〉,可以获取与其对应的三视匹配元τ,则τ中的位移向量为:

      $$ \begin{align} & \mathit{\boldsymbol{C}}_{ij}^{\tau }=\mathit{\boldsymbol{C}}_{j}^{\tau }-\mathit{\boldsymbol{C}}_{i}^{\tau } \\ & \mathit{\boldsymbol{C}}_{jk}^{\tau }=\mathit{\boldsymbol{C}}_{k}^{\tau }-\mathit{\boldsymbol{C}}_{j}^{\tau }~ \\ & \mathit{\boldsymbol{C}}_{ik}^{\tau }=\mathit{\boldsymbol{C}}_{k}^{\tau }-\mathit{\boldsymbol{C}}_{i}^{\tau }~ \\ \end{align} $$ (15)

      τ中局部坐标系相对全局坐标系的旋转矩阵为Rτ,则这三个位移向量在全局坐标系下的坐标分别为RτTCijτRτTCjkτRτTCikτ,并且这些位移向量在τ中满足局部一致性,所以:

      $$ \begin{align} & {{\mathit{\boldsymbol{C}}}_{ij}}={{\lambda }^{\tau }}\mathit{\boldsymbol{R}}_{\tau }^{\rm{T}}\mathit{\boldsymbol{C}}_{ij}^{\tau } \\ & {{\mathit{\boldsymbol{C}}}_{jk}}={{\lambda }^{\tau }}\mathit{\boldsymbol{R}}_{\tau }^{\rm{T}}\mathit{\boldsymbol{C}}_{jk}^{\tau } \\ & {{\mathit{\boldsymbol{C}}}_{ik}}={{\lambda }^{\tau }}\mathit{\boldsymbol{R}}_{\tau }^{\rm{T}}\mathit{\boldsymbol{C}}_{ik}^{\tau } \\ \end{align} $$ (16)

      式中, λττ中的尺度系数。以三视匹配元中的局部一致性位移向量Cijτ和尺度系数λτ替代立体像对中的tijλij,于是可将式(14)改写为:

      $$ \left| {{\mathit{\boldsymbol{R}}_j}{\mathit{\boldsymbol{C}}_j}-{\mathit{\boldsymbol{R}}_j}{\mathit{\boldsymbol{C}}_i} + {\lambda ^\tau }\mathit{\boldsymbol{R}}_\tau ^{\rm{T}}\mathit{\boldsymbol{C}}_{ij}^\tau } \right| \le \eta $$ (17)

      假设:

      $$ \begin{align} & {{\mathit{\boldsymbol{R}}}_{j}}=\left( \begin{matrix} R_{j}^{00} & R_{j}^{01} & R_{j}^{02} \\ R_{j}^{10} & R_{j}^{11} & R_{j}^{12} \\ R_{j}^{20} & R_{j}^{21} & R_{j}^{22} \\ \end{matrix} \right), \\ & \mathit{\boldsymbol{R}}_{\tau }^{\rm{T}}\mathit{\boldsymbol{C}}_{ij}^{\tau }=\left( \begin{matrix} {{M}^{0}} \\ {{M}^{1}} \\ {{M}^{2}} \\ \end{matrix} \right), \\ & {{\mathit{\boldsymbol{C}}}_{j}}=\left( \begin{matrix} C_{j}^{x} \\ C_{j}^{y} \\ C_{j}^{z} \\ \end{matrix} \right), {{\mathit{\boldsymbol{C}}}_{i}}=\left( \begin{matrix} C_{i}^{x} \\ C_{i}^{y} \\ C_{i}^{z} \\ \end{matrix} \right) \\ \end{align} $$

      式(17)可以转化为标准的线性规划形式,其中目标函数为:

      $$ \text{minimize}\ \ \eta $$ (18)

      约束条件为:

      $$ \begin{align} & R_{j}^{00}C_{j}^{x}+R_{j}^{01}C_{j}^{y}+R_{j}^{02}C_{j}^{z}-R_{j}^{00}C_{i}^{x}-R_{j}^{01}C_{i}^{y}-R_{j}^{02}C_{i}^{z}+{{M}^{0}}{{\lambda }^{\tau }}-\eta \le 0 \\ & R_{j}^{10}C_{j}^{x}+R_{j}^{11}C_{j}^{y}+R_{j}^{12}C_{j}^{z}-R_{j}^{10}C_{i}^{x}-R_{j}^{11}C_{i}^{y}-R_{j}^{12}C_{i}^{z}+{{M}^{1}}{{\lambda }^{\tau }}-\eta \le 0 \\ & R_{j}^{20}C_{j}^{x}+R_{j}^{21}C_{j}^{y}+R_{j}^{22}C_{j}^{z}-R_{j}^{20}C_{i}^{x}-R_{j}^{21}C_{i}^{y}-R_{j}^{22}C_{i}^{z}+{{M}^{2}}{{\lambda }^{\tau }}-\eta \le 0 \\ & R_{j}^{00}C_{j}^{x}+R_{j}^{01}C_{j}^{y}+R_{j}^{02}C_{j}^{z}-R_{j}^{00}C_{i}^{x}-R_{j}^{01}C_{i}^{y}-R_{j}^{02}C_{i}^{z}+{{M}^{0}}{{\lambda }^{\tau }}+\eta \ge 0 \\ & R_{j}^{10}C_{j}^{x}+R_{j}^{11}C_{j}^{y}+R_{j}^{12}C_{j}^{z}-R_{j}^{10}C_{i}^{x}-R_{j}^{11}C_{i}^{y}-R_{j}^{12}C_{i}^{z}+{{M}^{1}}{{\lambda }^{\tau }}+\eta \ge 0 \\ & R_{j}^{20}C_{j}^{x}+R_{j}^{21}C_{j}^{y}+R_{j}^{22}C_{j}^{z}-R_{j}^{20}C_{i}^{x}-R_{j}^{21}C_{i}^{y}-R_{j}^{22}C_{i}^{z}+{{M}^{2}}{{\lambda }^{\tau }}+\eta \ge 0 \\ & \forall <i, j>\in E \\ \end{align} $$ (19)

      未知数界限条件为:

      $$ \begin{array}{l} \eta \ge 0\\ {\lambda ^\tau } \ge 1{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \forall \tau \in U\\ {\mathit{\boldsymbol{C}}_1} = {\left[{\begin{array}{*{20}{c}} 0&0&0 \end{array}} \right]^{\rm{T}}} \end{array} $$ (20)

      本文利用OSI线性规划算法[11]计算G中所有相机中心的全局最优解,由于三视匹配元比立体像对具有更多的内部约束和更稳健的空间几何关系,所以本文求解相机中心位置的算法更加稳健。

    • 为了恢复场景的三维结构,本文首先利用并查集(union find, UF)算法[12]提取影像集G中的多视匹配点。定义无向图P(F, U),F为顶点集,由G中所有影像的所有特征点组成,U为边集,每一对匹配点构成U中的一条边。

      首先令F为空集,然后依次查找无向图G中的每一条边e,将e对应的两幅影像包含的所有匹配点对加入F(F中不重复添加同一张影像的特征点),然后利用e中包含的特征点匹配关系,将F中对应的每一对特征点相连,构成U中的一条边,直到无向图G中的所有边完成查找。图 2是利用并查集算法提取多视匹配点的一个示意图,图 2中包含P1~P4共4张影像,而(P1P2)、(P1P3)、(P2P4)和(P3P4)分别构成匹配像对,依次将这些像对中各幅影像的匹配特征点对加入F,并连接匹配点对。当完成所有边的查找之后,得到的图P中每一条联通边对应的所有节点构成一个多视匹配点。联通边集对应G中的所有多视匹配点,如图 2中对应的多视匹配点为(1,4,9,13)、(2,5,10)、(3,6,12,15)、(7,14,11)和(8,16)。并查集算法是一种近似线性算法,复杂度为O((n)),其中α(n)为反Ackermann函数[12],增量法[13]提取多视匹配点的算法复杂度为O(nlgn),所以并查集算法具有更高的计算效率,并且理论上可以提取所有的多视匹配点。

      图  2  基于并查集的多视匹配点提取

      Figure 2.  Multi-view Matching Points Extraction Based on UF Method

      通过影像集中所有影像的全局一致性方位,本文利用三角形法计算每个多视匹配点的空间坐标。假设空间点Xm幅影像上可见,则每幅影像中都存在如下关系[7]

      $$ {\mathit{\boldsymbol{x}}_i} = {\mathit{\boldsymbol{P}}_i}\mathit{\boldsymbol{X}}, i = 1, \cdots m $$ (21)

      其中, xi=wi[ui  vi  1]TX=[Xx Xy  Xz  1]T分别为像点和空间点的齐次坐标向量;(ui, vi)为X在影像i上对应的像点;wi为尺度系数;Pi=KiRi[I|-Ci]为影像i对应的摄影机矩阵,Ki表示相机标定矩阵,I表示单位阵。若以Pik(k=0, 1, 2)表示Pi的第k行,则:

      $$ {w_i}{u_i} = \mathit{\boldsymbol{P}}_i^0\mathit{\boldsymbol{X}}, {w_i}{v_i} = \mathit{\boldsymbol{P}}_i^1\mathit{\boldsymbol{X}}, {w_i} = \mathit{\boldsymbol{P}}_i^2\mathit{\boldsymbol{X}} $$

      消去尺度系数wi可得:

      $$ {u_i}\mathit{\boldsymbol{P}}_i^2\mathit{X = }\mathit{\boldsymbol{P}}_i^0\mathit{\boldsymbol{X}}, {v_i}\mathit{\boldsymbol{P}}_i^2\mathit{\boldsymbol{X}} = \mathit{\boldsymbol{P}}_i^1\mathit{\boldsymbol{X}} $$ (22)

      以$\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} = }}{\left[{\begin{array}{*{20}{c}} {{X^x}}&{{X^y}}&{{X^z}} \end{array}} \right]^{\rm{T}}}$表示空间点的坐标向量,式(22)可以写成矩阵形式:$\mathit{\boldsymbol{A\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }} = l$,其中A为2 m×3的系数矩阵,利用最小二乘法可以解算${\mathit{\boldsymbol{\mathord{\buildrel{\lower3pt\hbox{$\scriptscriptstyle\frown$}} \over X} }}}$,但由于此时最小化的函数为uiPi2X-Pi0XviPi2X-Pi1X,而并非重投影误差,所以直接利用最小二乘不能获得最优解。本文利用文献[9]的迭代线性最小二乘三角形算法计算空间点的三维坐标,将式(22)改写为${u_i} = \frac{{\mathit{\boldsymbol{P}}_i^0\mathit{X}}}{{\mathit{\boldsymbol{P}}_i^2\mathit{X}}} $和${v_i} = \frac{{\mathit{\boldsymbol{P}}_i^1\mathit{X}}}{{\mathit{\boldsymbol{P}}_i^2\mathit{X}}} $,以式(22)的直接最小二乘解作为空间点坐标向量的初值X(0),通过如下的迭代过程可以获得更高的精度[9]

      $$ {u_i} = \frac{{\mathit{\boldsymbol{P}}_i^0{\mathit{\boldsymbol{X}}^{\left( {t + 1} \right)}}}}{{\mathit{\boldsymbol{P}}_i^2{\mathit{\boldsymbol{X}}^{\left( t \right)}}}}, {v_i} = \frac{{\mathit{\boldsymbol{P}}_i^1{\mathit{\boldsymbol{X}}^{\left( {t + 1} \right)}}}}{{\mathit{\boldsymbol{P}}_i^2{\mathit{\boldsymbol{X}}^{\left( t \right)}}}} $$ (23)

      式中, t代表迭代的次数。同时计算空间点的重投影误差:

      $$ \rho \left( \mathit{\boldsymbol{X}} \right) = \sqrt {\frac{{\sum\limits_{i = 1}^m {\left( {{{\left( {{u_i}-\frac{{\mathit{\boldsymbol{P}}_i^0\mathit{X}}}{{\mathit{\boldsymbol{P}}_i^2\mathit{X}}}} \right)}^2} + {{\left( {{v_i}-\frac{{\mathit{\boldsymbol{P}}_i^1\mathit{X}}}{{\mathit{\boldsymbol{P}}_i^2\mathit{X}}}} \right)}^2}} \right)} }}{m}} $$ (24)

      剔除所有重投影误差大于5的多视匹配点,然后利用光束法平差[14]对相机方位及空间点坐标同时进行优化得到最终的重建结果。

    • 本文利用文献[15]提供的Benchmark公开测试数据中的FountainP11、HerzJesuP25和CastleP30三组影像集进行了实验。这三组影像集中的每一幅影像都完成了严格的畸变差校正,同时,Benchmark数据集提供了每一幅影像的内方位元素及其在地面坐标系中的精确方位,包括相机中心坐标以及相片的旋转矩阵。实验所用的计算机配置为Intel (R) CORETMi7 2.1 GHz 16 G内存。

      实验中首先对影像集中所有可能的像对进行尺度不变特征转换(scale invariant feature transform, SIFT)特征匹配[16],如果获取的匹配点对数量不少于50,则认为这两张影像相关,所有相关影像构成无向图G; 然后在RANSAC引擎中利用Nistér五点算法[17]计算G中任意两张相关影像ij之间的相对旋转矩阵Rij和相对基线方向向量tij,如果内点数量少于30,则认为计算结果不可靠,在G中剔除边ei, j〉。然后对G中所有可能的三视匹配元进行回路闭合检测,图 3为所有三视匹配元回路闭合差的统计图。

      图  3  回路闭合差统计

      Figure 3.  Closed Cycle Error Statistics

      剔除影像集中所有回路闭合差大于2°的三视匹配元,然后利用本文算法计算影像的全局一致性旋转矩阵和相机中心位置,并利用迭代最小二乘三角形算法计算所有多视匹配点的三维坐标,实验中将迭代次数设置为4。剔除所有重投影误差大于5的多视匹配点之后利用光束法平差来对计算结果进行优化。本文实验中的平差过程分为两步进行:第一次平差时保持所有影像的旋转矩阵不变,对相机位置和空间点的三维坐标进行优化,第二次平差则同时对旋转矩阵、相机位置和空间点三维坐标进行优化。表 1记录了本文算法的计算结果。

      表 1  批处理计算结果

      Table 1.  Calculation Results of Batch Method

      影像集 影像数量 三视匹配
      元数量
      回路闭合
      检测后
      多视匹配
      点数
      平差前δ 第一次平差 第二次平差
      迭代次数 δ 迭代次数 δ
      FountainP11 11 69 69 3 717 0.726 233 2 0.277 844 2 0.268 669
      HerzJesusP25 25 421 400 5 251 0.756 457 2 0.517 371 4 0.498 520
      CastleP30 30 491 362 6 573 1.544 170 2 0.450 979 3 0.321 487

      表 1中, δ代表多视匹配点重投影误差的均方根。从表 1中可以看出,本文算法在平差前就可以获得较好的计算结果。将本文的批处理算法与增量法Bundler[13]的计算时间进行对比,对比结果见表 2

      表 2  计算时间

      Table 2.  Computing Time of Batch Method and Bundle

      影像集 t1 t2 t3 t4 t5 t6 t7 t Bundler
      FountainP11 0.07 15.23 0.02 0.06 1.11 3.77 5.43 25.69 37.35
      HerzJesusP25 0.12 29.09 0.19 0.27 1.85 9.41 14.79 55.72 176.81
      CastleP30 0.16 31.85 0.32 0.29 2.37 11.41 14.54 60.94 470.34

      表 2中,本文批处理算法的计算时间主要分为7个部分:回路闭合检测(t1),三视匹配元内的一致性相对方位计算(t2),绝对旋转矩阵估算(t3),相机中心坐标估算(t4),多视匹配点三维坐标估算(t5), 第一次光束法平差(t6),第二次光束法平差(t7),t代表总的计算时间。从表 2中可以看出,批处理算法的计算效率要远高于增量法,并且影像数量越多,效率优势体现得越明显。

      最后利用空间相似变换实现全局坐标系和地面坐标系之间的转换。实验中选取影像集中6张影像的相机中心作为控制点,以其余影像的相机中心作为检查点。

      表 3统计了3个影像集位置误差(ΔX,ΔY,ΔZ)、ΔC和方向误差(ΔRoll,ΔPitch,ΔYaw)、ΔR绝对值的最大值、均值和均方根。其中,位置误差只统计检查点的精度,由于空间相似变换过程中不涉及控制点对应影像的旋转矩阵,方向误差统计了所有影像的旋转矩阵的精度。

      表 3  影像集精度统计

      Table 3.  PrecisionStatistics of 3 Image Sets

      影像集 ΔX/m ΔY/m ΔZ/m ΔC/m ΔRoll/(°) ΔPitch/(°) ΔYaw/(°) ΔR/(°)
      FountainP11 最大值 0.003 496 0.002 420 0.003 591 0.004 506 0.043 5 0.031 5 0.039 0 0.074 8
      均值 0.001 543 0.000 945 0.002 030 0.002 947 0.018 1 0.018 5 0.013 9 0.029 1
      均方根 0.001 832 0.001 230 0.002 267 0.003 164 0.021 9 0.020 1 0.017 6 0.039 0
      HerzJesusP25 最大值 0.008 147 0.007 807 0.010 924 0.012 594 0.306 2 0.041 9 0.337 5 0.095 6
      均值 0.002 805 0.003 363 0.004 790 0.007 199 0.060 4 0.014 9 0.081 2 0.039 1
      均方根 0.003 536 0.004 007 0.005 707 0.007 818 0.093 0 0.019 3 0.109 5 0.049 1
      CastleP30 最大值 0.051 092 0.047 150 0.028 073 0.067 435 0.384 0 0.101 0 0.357 1 0.119 0
      均值 0.013 625 0.014 809 0.012 217 0.026 935 0.063 7 0.035 8 0.047 4 0.052 1
      均方根 0.018 333 0.018 451 0.014 798 0.029 925 0.096 7 0.042 2 0.083 3 0.062 2

      同时,将本文算法与部分经典的三维重建算法进行了比较,包括增量法(文献[13]、文献[18])、分步式算法(文献[19])、二阶锥规划批处理算法(文献[4])、线性批处理算法(文献[8])。表 4是对比结果,从表 4中可以看出,本文算法的位置精度和旋转精度较其他几种算法具有明显的优势。

      表 4  精度比较

      Table 4.  Precision Comparison

      FountainP11 HerzJesusP25 CastleP30
      ΔC/m ΔR/(°) ΔC/m ΔR/(°) ΔC/m ΔR/(°)
      文献[13] 0.007 2 0.112 0.030 8 0.110 0.300 0 0.279
      文献[18] 0.009 9 0.116 0.023 3 0.104 0.264 0 0.398
      文献[19] 0.131 7 - 0.253 8 - - -
      文献[4] 0.015 3 - 0.084 5 - - -
      文献[8] 0.014 0 0.195 0.064 0 0.188 0.235 0 0.480
      本文算法 0.003 2 0.039 0.007 8 0.049 0.029 9 0.062
    • 本文提出了一种基于三视匹配元的多视影像批处理三维重建算法,对Benchmark标准测试数据中的三组影像集进行了实验,并与部分其他重建算法进行了比较。实验结果表明,本文批处理算法可以避免迭代式的平差运算,计算效率较增量法有较大提高,由于多视影像中三焦张量比核线几何具有更强的约束关系,基于三视匹配元的重建算法精度更高,并且在三维重建过程中加入回路闭合检测可以有效避免错误的核线几何关系,提高计算结果的可靠性。

参考文献 (19)

目录

    /

    返回文章
    返回