论文 Structure-from-Motion Revisited 对ISFM改进的理解
这篇论文对应于增量式SFM开源软件COLMAP。对一般的ISFM做了5点改进,以进一步提高重建模型的完整性与准确性。现详细记录自己对每个改进的理解。
- Scene Graph Augmentation:引入几何验证策略来标记场景图中边的类型,以此提高初始化和三角化的鲁棒性;
- Next Best View Selection:选择下一个最佳视图时充分考虑观测到的地图点数目和地图点的均匀分布;
- Robust and Efficient Triangulation:一种计算成本更低,且鲁棒性更优的三角化方法;
- Bundle Adjustment:迭代BA,再三角化和异常值滤波策略能够减小误差累计,从而显著提高重建完整性和准确性;
- Redundant View Mining:从无序图像集集合挖掘相似视图,使之成组,从而减小BA计算量,提高BA准确性;
以上五点总结引用公众号:泡泡机器人SLAM 2017-10-21
Scene Graph Augmentation
基于多模型几何验证的场景图增强目的是为了使构建的correspondence graph更加稳健,从而提高后续初始化以及三角化的鲁棒性。
Scene Graph Augmentation主要在初始化时采用了该策略,在“two_view_geometry.cc”文件下的“EstimateCalibrated”函数中实现。该函数主要根据估计的两张影像之间的本质矩阵E、基础矩阵F、单应矩阵H及对应的内点数来确定求解相对位姿或运动关系的方式,即对应程序中确定变量“config”的设定,可选的“config”设定包括:UNDEFINED、DEGENERATE、CALIBRATED、UNCALIBRATED、PLANAR、PANORAMIC、PLANAR_OR_PANORAMIC、WATERMARK、MULTIPLE。利用前面求解的E、F、H矩阵及对应的内点,计算出几项指标值(后面会详细阐述),根据指标值满足的条件确定最终“config”的设定,从而确定在“two_view_geometry.cc”文件下的“EstimateRelativePose”函数中求解R、T矩阵的方式。比如经过计算,设定“config”的值为CALIBRATED,意味着将从本质矩阵中分解求得R、T矩阵。所以“EstimateCalibrated”函数包含6个参数,分别为存储两个相机内外参的变量,两张影像提取的特征点二维坐标,两张影像之间的匹配点对,以及对极几何的一些参数设定。
- 首先判断两张影像之间的匹配点对数量是否满足min_num_inliers阈值,不满足设置“config”的值为DEGENERATE;
- 接着提取两张影像之间的匹配点对,利用RANSAC算法估计E、F、H矩阵及对应的内点数,在此基础上如果有一种矩阵求解失败或者对应的内点数量不满足最少内点数量阈值要求,则设置“config”的值为DEGENERATE;
- 接着计算三个指标值:E_F_inlier_ratio、H_F_inlier_ratio、H_E_inlier_ratio,第一个为本质矩阵内点数比基础矩阵内点数,其它类似;
- 之后分五大类情况确定“config”的值,源码如下:
if (E_report.success && E_F_inlier_ratio > options.min_E_F_inlier_ratio &&
E_report.support.num_inliers >= options.min_num_inliers) {
// Calibrated configuration.
// Always use the model with maximum matches.
if (E_report.support.num_inliers >= F_report.support.num_inliers) {
num_inliers = E_report.support.num_inliers;
best_inlier_mask = &E_report.inlier_mask;
} else {
num_inliers = F_report.support.num_inliers;
best_inlier_mask = &F_report.inlier_mask;
}
if (H_E_inlier_ratio > options.max_H_inlier_ratio) {
config = PLANAR_OR_PANORAMIC;
if (H_report.support.num_inliers > num_inliers) {
num_inliers = H_report.support.num_inliers;
best_inlier_mask = &H_report.inlier_mask;
}
} else {
config = ConfigurationType::CALIBRATED;
}
} else if (F_report.success &&
F_report.support.num_inliers >= options.min_num_inliers) {
// Uncalibrated configuration.
num_inliers = F_report.support.num_inliers;
best_inlier_mask = &F_report.inlier_mask;
if (H_F_inlier_ratio > options.max_H_inlier_ratio) {
config = ConfigurationType::PLANAR_OR_PANORAMIC;
if (H_report.support.num_inliers > num_inliers) {
num_inliers = H_report.support.num_inliers;
best_inlier_mask = &H_report.inlier_mask;
}
} else {
config = ConfigurationType::UNCALIBRATED;
}
} else if (H_report.success &&
H_report.support.num_inliers >= options.min_num_inliers) {
num_inliers = H_report.support.num_inliers;
best_inlier_mask = &H_report.inlier_mask;
config = ConfigurationType::PLANAR_OR_PANORAMIC;
} else {
config = ConfigurationType::DEGENERATE;
return;
}
if (best_inlier_mask != nullptr) {
inlier_matches =
ExtractInlierMatches(matches, num_inliers, *best_inlier_mask);
if (options.detect_watermark &&
DetectWatermark(camera1, matched_points1, camera2, matched_points2,
num_inliers, *best_inlier_mask, options)) {
config = ConfigurationType::WATERMARK;
}
}
4.1.第一种情况,E矩阵估计成功,E_F_inlier_ratio满足阈值条件,E矩阵对应内点数满足阈值条件。这个大前提下,再比较E矩阵和F矩阵对应的内点数,哪个内点数大,就先暂时确定使用哪个矩阵用于后续的估计;同时判断H_E_inlier_ratio是否满足阈值条件,满足则说明场景更趋向于平面场景,从而设定“config”的值为PLANAR_OR_PANORAMIC,在此基础上则继续判断H矩阵对应的内点数量是否大于阈值,是则最终确定采用H矩阵估计R、T;如果H_E_inlier_ratio不满足阈值条件,则设置“config”的值为CALIBRATED,最终由前面确定的矩阵来估计R、T;
4.2.第二种情况,F矩阵估计成功,F矩阵对应内点数满足阈值条件。先暂时设定由F矩阵来估计R、T,之后判断H_F_inlier_ratio是否满足阈值,如果满足,设定“config”的值为PLANAR_OR_PANORAMIC,并进一步判断H矩阵对应的内点数是否大于之前F矩阵对应的内点数,大于就最终设定由H矩阵来估计R、T;如果H_F_inlier_ratio不满足阈值,则设置“config”的值为UNCALIBRATED;
4.3.第三种情况,H矩阵估计成功,H矩阵对应内点数满足阈值条件。设定“config”的值为PLANAR_OR_PANORAMIC,并且设定由H矩阵来估计R、T,内点数量为H矩阵对应的内点数量;
4.4.第四种情况,以上都不满足,设定“config”的值为DEGENERATE;
4.5.第五种情况,针对有水印的影像,设定“config”的值为WATERMARK;
5. 最后,提取内点对应的匹配点对,进行最终的R、T矩阵求解。
Next Best View Selection
选择最优的下一张影像进行register是为了最小化重建误差[1-2],选择了一张不合适的影像进行register,将会导致一连串的相机位姿错误求解以及三角化错误的空间点。或者说下一张需要register的影像质量直接影响位姿的估计以及三角化的精度和完整性。对于Internet photo collections来说,由于没有关于场景范围以及相机参数的先验信息,通常基于影像表面信息[1]、对应点数量、已有的增量式重建好的场景[2-3]来确定最优的下一张影像。[4]中本着最小化相机后方交会不确定的目标选择可以看到最多三角化点的影像;[2]中提出了一种uncertainty-driven的方式最小化重建误差。由于在求解确定最优的下一张影像位姿时采用的一般是PnP算法,[5]中实验探索了利用PnP算法求解影像位姿的精度与这张影像观测值的数量及这些观测值在影像分布情况之间的关系。主要结论总结来说就是2D-3D匹配点对数量越多,这些三维点在待求位姿的影像上的二维点分布越均匀,内外参数估计得就越准确[5][6]。
对于Internet photo collections来说,通常需要从很多候选影像中选出最优下一张影像。[2]中暴力协方差传播的方式由于需要对每一个候选影像进行计算和分析协方差,显然不可取。对于可以看到更多已有空间三维点并且这些可见三维点在待确定位姿的影像上的二维点分布越均匀,就给定一个更高的得分S[7]。论文的方法如下:
- 对于每一张候选影像构建一个L层的金字塔,将每一层金字塔影像划分为的格网,其中每一层中每一个格网设定一个权重:,(注意:论文中是,根据论文中给的图例,我认为应该是,下面说明);
- 统计每一张候选影像的分数。这个分数由金字塔影像所有层得分求和得到,每一层计算:如果一个格网中有一个空间三维点对应的二维点,S就增加一个对应的权重,否则不增加。这里以论文给出的图例来说明,构建金字塔的目的就是选取那些二维点分布最均匀的候选影像;
- 最后,根据这个S对候选影像进行排序,得分最高的就是最优的下一张影像。
这里有四张候选影像,每张影像构建三层金字塔,每层维度分别为2×2,4×4,8×8,对应的每层中每个格网的权重是2,4,8。
那么第一张影像得分计算如下:
首先第一层中,所有的可见三维空间点对应的二维点都在一个网格内,所以得分是1×2=2;
第二层中,所有的可见三维空间点对应的二维点分布在四个网格内,所以4×4=16;
第三层中,所有的可见三维空间点对应的二维点分布在六个网格内,所以6×8=48;
所以总分为:2+16+48=66.
同理对于第二张影像:
146=1×2+4×4+16×8;
对于第三张影像:
80=4×2+6×4+6×8;
对于第四张影像:
200=4×2+16×4+16×8。
本部分参考文献:
1 S. Chen, Y. F. Li, J. Zhang, and W. Wang. Active Sensor Planning for Multiview Vision Tasks. 2008.
2 S. Haner and A. Heyden. Covariance propagation and next best view planning for 3d reconstruction. ECCV, 2012.
3 N. Snavely, S. Seitz, and R. Szeliski. Photo tourism: exploring photo collections in 3d. ACM TOG, 2006.
4 N. Snavely. Scene reconstruction and visualization from internet photo collections. PhD thesis, 2008.
5 V. Lepetit, F. Moreno-Noguer, and P. Fua. EPnP: An accurate O(n) solution to the PnP problem. IJCV, 2009.
6 C. McGlone, E. Mikhail, and J. Bethel. Manual of photogrammetry. 1980.
7 A. Irschara, C. Zach, J.-M. Frahm, and H. Bischof. From structure-from-motion point clouds to fast location recognition. CVPR, 2009.
Robust and Efficient Triangulation
一般的图像特征匹配策略往往对场景相似度较高的影像对效果最优,但是这些影像对往往对应的相机中心连线,即基线较短,这将导致三角化失败。所以COLMAP的改进是利用传递性(transitivity)建立大基线影像之间的对应点关系,使得能够进行更准确的三角化。对于从带有噪声的观测值(有一定的误匹配点对)中进行多视角三角化,已有的方法包括[1-8],但是都无法处理噪声很高的情况。
论文分析了两点造成外点率高的原因,1)双目几何极线上错误的任意特征点匹配,比如四个特征点追踪了相同的距离,会导致3/4的外点率;2)不正确的相机位姿。Bundler取样多对匹配点对执行三角化,挑选三角化角足够大的解,再对所有匹配点对三角化;同时,所有的观测值需要通过cheirality constraint[9]。该方法对外点并不鲁棒,而且密集的三角化计算量大。本文的改进如下:
设输入的为feature track,,具有一定的外点率,其中𝑇𝑖包括归一化后的像平面坐标,对应的相机位姿,表示从世界坐标系到相机坐标系的变换,我们的目的是使满足下式表示的well-conditioned双目三角化的
越多越好,其中表示三角化方式,表示空间已经三角化的点;
而通过验证下面三个条件是否都成立来判断是否是一个well-conditioned双目三角化模型。
首先是足够大的交会角,按照下式计算:
接着是两个相机是否都具有正的深度,按照下式计算:
最后是该三维点在影像上的重投影误差小于设定的阈值,按照下式计算:
在sfm/incremental_mapper.cc的RegisterInitialImagePair函数中,我们可以看出如何实现更加鲁棒的三角化,源码如下:
Track track;
track.Reserve(2);
track.AddElement(TrackElement());
track.AddElement(TrackElement());
track.Element(0).image_id = image_id1;
track.Element(1).image_id = image_id2;
for (const auto& corr : corrs) {
const Eigen::Vector2d point1_N =
camera1.ImageToWorld(image1.Point2D(corr.point2D_idx1).XY());
const Eigen::Vector2d point2_N =
camera2.ImageToWorld(image2.Point2D(corr.point2D_idx2).XY());
const Eigen::Vector3d& xyz =
TriangulatePoint(proj_matrix1, proj_matrix2, point1_N, point2_N);
const double tri_angle =
CalculateTriangulationAngle(proj_center1, proj_center2, xyz);
if (tri_angle >= min_tri_angle_rad &&
HasPointPositiveDepth(proj_matrix1, xyz) &&
HasPointPositiveDepth(proj_matrix2, xyz)) {
track.Element(0).point2D_idx = corr.point2D_idx1;
track.Element(1).point2D_idx = corr.point2D_idx2;
reconstruction_->AddPoint3D(xyz, track);
}
}
可以看出首先调用TriangulatePoint函数进行三角化,之后判断交会角、相机深度来确定上一步三角化的每一个点是否能够加入到最终重建的三维点中。
本部分参考文献:
1 R. I. Hartley and P. Sturm. Triangulation. 1997.
2 F. Lu and R. Hartley. A fast optimal algorithm for L2 triangulation. ACCV, 2007.
3 C. Aholt, S. Agarwal, and R. Thomas. A QCQP Approach to Triangulation. ECCV, 2012.
4 R. Hartley and F. Schaffalitzky. L∞ minimization in geometric reconstruction problems. CVPR, 2004.
5 H. Li. A practical algorithm for L∞ triangulation with outliers. CVPR, 2007.
6 S. Agarwal, N. Snavely, and S. Seitz. Fast algorithms for L∞ problems in multiview geometry. CVPR, 2008.
7 C. Olsson, A. Eriksson, and R. Hartley. Outlier removal using duality. CVPR, 2010.
8 L. Kang, L. Wu, and Y.-H. Yang. Robust multi-view l2 triangulation via optimal inlier selection and 3d structure refinement. Pattern Recognition, 2014.
9 R. Hartley and A. Zisserman. Multiple view geometry in computer vision. 2003.
Bundle Adjustment
从COLMAP软件测试的多组数据来看,整体光束法平差的流程可以归纳为:
由上面的流程图可以看出整体平差策略是在每一张新增影像registration后,在包含该张新增影像在内的局部范围内进行光束法平差;在初始化以及再三角化后进行了全局光束法平差。虽然论文中提到“类似于VisualSFM,当重建模型规模增长一定比例后再进行全局光束法平差”,但实际软件操作中默认设置重建条件下,还是在每register一张新影像的最后进行了整体光束法平差。
除了以上整体策略外,还对参数的设定、外点(outliers)滤除、再三角化进行了说明。
参数的设定:
- 局部光束法平差中,使用Cauchy函数作为鲁棒的损失函数ρj;
- 对于含有几百个相机的重建场景,使用sparse direct solver(稀疏矩阵求解器);
- 对于更大的场景,使用PCG(共轭梯度法);
- 平差过程采用了Ceres Solver库[1],还可以选择任意数量图像间是否共享相机模型。 外点(outliers)滤除:
- 滤去重投影误差大的观测值[2][3];
- 对于每一个空间三维点,检查交会于该点的光束之间的交会角,是否满足well-conditioned geometry条件的阈值[2];
- 全局光束法平差后,检查退化的相机,这些相机是指异常视场角、畸变系数异常大,需要滤去,这就要求在光束法平差中焦距和畸变参数一起参与优化,同时为了解决像主点校正这个病态问题[4],对于未校正的相机,固定像主点为像片中心点。
再三角化:
- 与VisualSFM类似,在全局光束法平差之前,进行再三角化操作;
- 新增全局光束法平差之后的后三角化操作,针对光束法平差前由于估计的错误位姿导致的三角化失败的点以提高重建场景的完整性,这些点的观测值误差小于外点(outliers)滤除环节中对应阈值。实际软件测试时并没有发现后三角化。
本部分总结:不同于Bundler和VisualSFM非迭代进行光束法平差和外点(outliers)滤除的模式,COLMAP迭代进行光束法平差、再三角化、外点(outliers)滤除直到滤除的观测值数量及后三角化点的数量为减小,甚至为0。
本部分参考文献:
1 S. Agarwal, K. Mierle, and Others. Ceres Solver. http: //ceres-solver.org.
2 N. Snavely, S. Seitz, and R. Szeliski. Photo tourism: exploring photo collections in 3d. ACM TOG, 2006.
3 C. Wu. Towards linear-time incremental structure from motion. 3DV, 2013.
4 L. de Agapito, E. Hayman, and I. Reid. Self-calibration of a rotating camera with varying intrinsic parameters. BMVC, 1998.
Redundant View Mining
这部分的目的是利用ISFM的内在特性,对复杂稠密场景内的影像(相机)按照重叠度分组,每一组内使用一个相机模型表示,减小光束法平差计算量。
对于Internet photo collections,具有两个特点,一是高度不一致的视觉模式;二是这些无序数据集往往一块块聚集在一起,每一块中会有很多影像,分布不均匀。从数据集的这些特点出发,出现了不少提高光束法平差效率的策略[1][2][3]。[1]根据可视图案进行相机系统简化的预处理。COLMAP采用了[2]中相似的思想,将相机和空间三维点分成多个子地图,通过separator variables相连,转化为一个图割问题(graph cut);再此基础上,交叉进行固定相机和空间三维点,优化separator variables和固定separator variables,优化相机和空间三维点的操作。[3]中借鉴low-rank的概念,将空间接近的点对齐为一个单一点,在相机上增加high-rank约束。
在[2]的基础上,COLMAP做了如下的改进:
- 不再采用图割(graph cut)方法,利用ISFM内在特性:1)光束法平差只会影响新增的部分影像和相机[4];2)Internet
photo collections中通常相机分布不均匀,视角冗余。根据参数是否新增影像的影响分为两组,未受影响的部分按照重叠度等分为多组,具体后面进一步阐述,表示如下:
由于每一组内影像高度重叠,每一组G_r可只采用单一相机。受新增影像影响的单独分组进行后续的参数优化。所以如何确定哪些影像受影响哪些不受呢?文章中认为受影响的影像主要包括两类:1)在最新的一次重建场景拓展中新添加的;2)某些对应的观测值中重投影误差大于r个像素的占比大于阈值ϵ_r的影像。这部分优化策略用在四中global bundle adjustment之前,retriangulation之后。
接下来说一下前面提到的如何对未受影响的影像进行分组的流程。主要依据两个原则:1)同组内影像之间足够重叠、冗余[2];2)影像之间共视的空间三维点数量是影像间交互重叠程度的很好度量[1]。具体流程:
1.1.对于一个有三维空间点的场景,每一张图像描述为一个二进制向量,,对于其中第n项为1的话,说明在影像i中可见,否则不可见。因此,两张影像a、b关系可按如下公式确定:
1.2.按照的模先对影像降排序,得集合。初始化一个组,从中挑出排在第一位的向量对应的影像a,存入,此时再查找另一张影像b,如果且中的数量小于S,则存入影像b,否则存入一个新组。为了加快b的查找,参考空间最近邻的思想,找到a影像±𝛽度视角范围内的图像作为候选。 - 最终,每组的优化目标函数变为:
对比一般ISFM优化的目标函数:
可以看出每一个组内所有图像共享一个extrinsic group parameters Gr ∈ SE(3),优化过程中,每一个相机Pc矩阵固定,减少优化量,这是COLMAP的一大贡献。组的数量越多,表现越优,待优化的相机数量越少,使得原始直接优化方式的三次方时间复杂度降至非直接方式的线性复杂度。
本部分参考文献:
1 A. Kushal and S. Agarwal. Visibility based preconditioning for bundle adjustment. CVPR, 2012.
2 K. Ni, D. Steedly, and F. Dellaert. Out-of-core bundle adjustment for large-scale 3d reconstruction. ICCV, 2007.
3 L. Carlone, P. Fernandez Alcantarilla, H.-P. Chiu, Z. Kira, and F. Dellaert. Mining structure fragments for smart bundle adjustment. BMVC, 2014.
4 C. Wu. Towards linear-time incremental structure from motion. 3DV, 2013.
最后
以上为个人的理解,可能会存在一些理解偏差,欢迎指正~~