Towards Linear-Time Incremental Structure from Motion¶
论文链接
原文：Michela Farenzena, Andrea Fusiello, Riccardo Gherardi. Towards Linear-Time Incremental Structure from Motion
Abstract¶
增量式运动恢复结构 (SFM) 的时间复杂度通常被认为是相对于相机数量的 \(O(n^4)\) 。由于捆绑调整 (BA) 最近通过预条件共轭梯度 (PCG) 得到了显著的改善，因此，增量 SfM 速度有多快是值得重新审视的。我们引入了一种新的 BA 策略，它在速度和精度之间提供了良好的平衡。通过算法分析和大量实验表明，在包括 BA 在内的许多主要步骤上，增量 SfM 只需要 \(O(n)\) 时间。我们的方法通过定期重新三角量测最初未能三角量测的特征匹配来保持较高的精度。我们在不同设置下的大型照片集和长视频序列上测试了我们的算法，结果表明我们的算法为大规模重建提供了最先进的性能。
1 Introduction¶
运动结构 (SFM) 已经成功地用于重建越来越大的不受控制的照片集 【8, 2, 5, 4, 7】。尽管可以将大的照片集合减少到用于重建的图标图像【8, 5】 或骨架图像 【15, 2】 的子集，但是增量 SfM 的公知的 \(O(n^4)\) 成本对于大的场景分量仍然很高。最近，用离散优化和连续优化相结合的方法证明了一种 \(O(n^3)\) 方法 【4】。然而，实现低复杂度和高可伸缩性仍然是 SfM 的关键。
本文从以下几个方面展示了我们为进一步提高 SfM 的效率所做的努力：
- 我们提出了一种抢占式特征匹配方法 (preemptive feature matching)，在恢复足够好匹配的同时，图像匹配对减少了 95%。
- 我们分析了共轭梯度捆绑调整方法 (conjugate gradient bundle adjustment) 的时间复杂度。通过理论分析和实验验证，证明捆绑调整在实际应用中所需时间为 \((O(n))\)。
- 我们展示了增量 SfM 的许多子步骤，包括 BA 和点滤波。当我们使用新的 BA 策略后，只需要 \((O(n))\) 的时间。
- 没有牺牲任何时间复杂度，我们引入了重新三角量测 (re-triangulation) 的步骤，用来处理没有显著闭环的累计漂移的问题。
本文剩余部分安排如下：第二章介绍本文的研究背景和相关工作。第三章首先介绍了抢占式特征匹配。第四章分析了捆绑调整的时间复杂度，第五节提出了新的 SfM 算法，第六节和第七节给出了实验和结论。
2 Related Work¶
在典型的增量式 SfM 系统 (如【14】) 中，首先在两幅图像之间成功匹配特征的基础上估计双视角重建，然后通过从良好的双视图重建结构初始化，重复添加匹配图像，三角化特征匹配以及捆绑调整结构和运动来重建 3D 模型。对于 n 幅图像，这种增量式 SfM 算法的时间复杂度通常为 \(O(n^4)\)，这种高的计算代价阻碍了这种简单的增量式 SfM 算法在大型照片集上的应用，
大规模数据往往包含大量冗余信息，这使得许多 3D 重建的计算可以近似以利于速度，例如：(1) 对于 n 幅图像，这样的增量 SfM 算法的时间复杂度通常是 \(O(n^4)\)，这阻碍了这种简单的增量式 SfM 算法在大规模照片收集上的应用，这使得 3D 重建的许多计算可以进行近似以利于提升速度，例如：
Image Matching
图像匹配而不是所有图像彼此匹配，Agarwal 等人提出。【2】首先使用词汇树识别 【11】 为每幅图像确定少量候选，然后使用近似最近邻法对特征进行匹配。结果表明，图像匹配的二重近似仍然保持了足够的特征匹配。Frahm 等人【5】利用近似的 GPS 标签，仅将图像与近处的图像进行匹配。本文提出了一种抢占式特征匹配算法，进一步提高了匹配速度。
Bundle Adjustment
由于 BA 已经利用了非线性优化问题的线性近似，所以通常不需要求解精确的下降步骤。当前算法通过使用预条件共轭梯度求解线性系统达到显著的提速。同样，没有必要为每个新图像运行完整的 BA，或者总是等待 BA/PCG 收敛到非常小的残差。在这篇文章中，我们证明了大规模不受控制的照片集的捆绑调整需要线性时间。
Scene Graph
大型照片集通常包含足够多的图像来实现高质量的重建。对于高成本的 SfM，可以通过减少图像数量来提高效率。[8，5]使用图标图像作为场景图的主骨架，而 [15，2] 从密集的场景图中提取骨架图。在实践中，对大规模 SfM 的其他类型的改进应该与场景图简化结合使用。然而，为了突破增量 SfM 的极限，我们考虑了在不简化场景图的情况下构造单个连通部件。
与增量 SfM 相比，其他工作试图避免这种贪婪的方式。Gherard 等人 【6】 通过平衡分枝和合并，提出了一种层次化的 SfM 算法，由于需要的大型模型基数较少，降低了算法的时间复杂度。Sinha 等人 【13】 从消失点恢复相对相机旋转，并将 SfM 转换为高效的 3D 模型合并。最近，Crandall et al.【4】利用基于 GPS 的初始化方法，将 SfM 建模为全局 MRF 优化。这项工作是对增量式 SfM 的重新研究，仍然有它的一些局限性，如初始化影响完备性，但不依赖于额外的校准、GPS 标签或消失点检测。
Notations
我们用 n、p 和 q 分别表示重建过程中摄像机、点和观测值的数量。假设每个图像都有有限数量的特征，我们有 p= \((O(n))\) 和 q= \((O(n))\) 。由于本文考虑的是单个连通的场景图，n 也被称为滥用记号的输入图像的数目。
3 Preemtive Feature Matching¶
图像匹配是 SfM 最耗时的步骤之一。对于输入的图像，完全成对匹配需要 \(O(n^2)\) 个时间，但是对于大数据集计算匹配的子集 (例如，通过使用词汇树【2】) 是很好的，并且总体计算量也可以减少 \(O(n)\) 。此外，图像匹配可以很容易地并行化到多台机器 【2】 或多个 GPU【5】上，以获得更大的加速比。然而，特征匹配仍然是一个瓶颈，因为典型的图像包含数千个特征，这是我们试图处理的方面。
由于大型照片收集中视角的多样性，大多数图像对不匹配 (对于我们实验的大型数据集，75%−98%)。如果我们能够鲁棒有效地识别出好的配对，就可以节省大量的匹配时间。通过利用不变特征的尺度，我们提出了以下抢占式特征匹配的方法。
- 将每幅图像的特征按比例递减排序。这是一次 \(O(n)\)预处理。
- 生成需要匹配的配对列表，要么使用所有配对，要么选择一个子集 ([2, 5])。
- 对于每个图像对 (并行)，执行以下操作：
- 匹配两个图像的前 h 个特征。
- 如果来自子集的匹配数目小于 \(t_h\) ，则返回并跳过下一步。
- 进行常规匹配和几何估计。
其中 h 是子集大小的参数，\(t_h\) 是预期匹配数的阈值。子集和全集的特征匹配使用相同的最近邻算法和距离比检验 【10】，并要求匹配的特征是相互最近的。
设 \(k_1\) 和 \(k_2\) 是两幅图像的特征数， \(k=max(k_1, k_2)\) 。设 \(m_p(h)\) 和 \(m_i(h)\) 是当使用多达 h 个顶级尺度特征时推定的和内嵌匹配的数目。我们将特征匹配的成品率定义为：
\(Y_p(h)=\dfrac{m_p(h)}{h}\) 和 \(Y_i(h)=\dfrac{m_i(h)}{h}\) 。
我们感兴趣的是特征子集的产率如何与最终的产率 \(Y_i(k)\) 相对应。如图 1(a) 所示，即使对于很小的 h，如 100， \(Y_p(h)\) 和 \(Y_i(h)\)的分布也与 \(Y_i(k)\) 密切相关。图 1(b) 展示了，在顶级尺度特征中匹配的可能性与其他特征相当。这意味着顶级子集大约有 h/k 次机会保持匹配，远高于随机抽样 h 个特征的 \(h^2/k_1k_2\) 次机会。此外，大多数的匹配时间减少到 \(h^2/k_1k_2\) 的倍数，同时有更好的缓存。在本文中，我们综合考虑效率和鲁棒性，选择 h=100。
顶层特征匹配良好的原因有以下几个
-
少量的顶层特征可以覆盖相对较大的尺度范围，这是由于高斯层上特征数量的减少。
-
特征匹配是结构化的，使得一幅图像中的大尺度特征经常与另一幅图像中的大尺度特征相匹配。两个匹配特征之间的尺度变化是由摄像机运动和场景结构共同决定的，因此多个特征匹配的尺度变化并不是相互独立的。图 1(c) 显示了我们对特征尺度变化的统计。由于视点变化较小或场景深度结构良好，图像对的特征尺度变化通常较小。基于同样的原因，我们使用多达 8192 个特征进行常规特征匹配，这对大多数场景来说已经足够了，虽然大型照片集合包含冗余的视图和特征，但我们的抢占式特征匹配允许将最大的努力放在更有可能匹配的特征上。
4 How fast is bundle adjustment?¶
捆绑调整 (BA) 是结构参数和运动参数的联合非线性优化，Levenberg-MarQuardt(LM)是其中的一种方法。近年来，预条件共轭梯度 (PCG)[1，3，16] 使大规模束平差的性能有了明显的提高，因此值得重新研究 BA 和 SfM 的时间复杂度。
设 x 是参数向量，\(f(x)\) 是 3D 重建的重投影误差向量。我们要解决的最优化问题是非线性最小二乘问题： \(x^*=argmin_x\left| f(x) \right|^2\) 。设 \(J\) 是 \(f(x)\) 的雅可比，\(\Lambda\) 是用于正则化的非负对角向量，然后 LM 重复对以下线性系统求解：
\((J^TJ+\Lambda)\delta=-J^Tf\)
如果 \(\left| f(x+\delta) \right|<\left| f(x) \right|\)，升级参数 \(x\leftarrow x+\delta\)x。矩阵 \(H_\Lambda=J^TJ+\Lambda\) 被称为增广 Hessian 矩阵。高斯消去法通常用于首先获得摄像机参数的降阶线性系统，称为舒尔补。
Hessian 矩阵需要 \(O(q)=O(n)\) 空间和时间来计算，而 Schur 补需要 \((O(n^2))\) 空间和时间来计算。利用 Cholesky 分解求解线性方程组需要\(O(n^3)\) 或 \(O(p^3)\) 的立方时间，由于摄像机数量远小于点的数量，Schur 补法可以将分解时间减少一个很大的常数因子，该算法的一个令人印象深刻的例子是 Lourakis 和 Argyros 的 SBA【9】 在 Photo Tourism【14】 中使用。
对于共轭梯度法，在多次 CG 迭代中的主要运算是矩阵 - 向量乘法，其时间复杂度由所涉及矩阵的大小决定。通过使用仅仅 \(O(n)\) 空间 Hessian 矩阵，且避免舒尔补，CG 迭代可以达到 \(O(n)\) 时间复杂度。最近，多核捆绑调整 【16】 更进一步，使用 Hessian 矩阵和 Schur 补的隐式乘法，只需要构造 O(N)空间的 Jacobian 矩阵。在这项工作中，我们使用了多核捆绑调整的 GPU 版本。图 2(a)显示了我们的捆绑调整问题中 CG 迭代的时间，它显示了 CG 迭代的时间 \(T_{cg}\) 与 n 之间的线性关系。
在每个 LM 步骤中，PCG 需要 \(O(\sqrt{K})\) 次迭代才能精确地求解线性方程组 【12】，其中κ是线性方程组的条件数。使用好的预条件算子可以得到较小的条件数，例如 block - 雅可比预条件算子[1，16] 证明了它的快速收敛速度。在我们的实验中，PCG 平均使用 20 次迭代来求解线性系统。
令人惊讶的是，假设一个 BA 中有 \(O(1)\) 个 CG 迭代，则捆绑调整的时间复杂度已经达到 \(O(n)\)。虽然 CG/LM 迭代的实际次数取决于输入问题的难度，但我们从大量不同问题大小的 BA 收集的统计数据很好地支持了 CG/LM 迭代的 O(1) 假设。图 2(b) 给出了 BA 使用的 LM 迭代的分布，其中平均值为 37%，93% 的 BA 在少于 100 次 LM 迭代后收敛。图 2(c) 显示了 total CG 迭代的分布 (由 BA 的所有 LM 步骤使用)，它们也很小。在实践中，我们选择每个 BA 最多使用 100 次 LM 迭代，每个 LM 最多使用 100 个 CG 迭代，这保证了大多数问题的良好收敛性。
5 Incremental Structure from Motion¶
这一部分介绍了我们的 SFM 的设计，它实际上具有线性运行时间。捆绑调整可以在 \(O(n)\) 时间完成，打开了推动增量重建时间近似于 \(O(n)\) 的机会。如图 3 所示，我们的算法在每次迭代中添加一幅图像，然后运行完全 BA 或部分 BA。在 BA 和滤波之后，我们在继续下一次迭代和重新三角测量步骤之间做出另一个选择。
5.1 How fast is incremental SfM?¶
相机和 3D 点在重建过程中通常会很快稳定下来，因此没有必要在每次迭代中优化所有的相机和点参数。典型的策略是在添加固定数量的相机 (如α) 后执行全 BAs，所有全 BAs 的累积时间为
当使用 PCG 时。与基于 Cholesky 因式分解的 BA 相比，这已经是 \((O(n^4))\) 时间的显著减少。随着模型变得更大和更稳定，跳过更昂贵的 BA 是可以承受的，所以我们想要研究在不损失精度的情况下我们还能走多远。
5.2 Re-triangulation(RT)¶
由于摄像机相对位姿的累积误差，增量式 SfM 存在漂移问题，两个摄像机位姿之间的约束由它们的三角特征匹配来提供。由于初始估计的姿态，甚至在部分 BA 之后的姿态可能不够准确，一些正确的特征匹配可能无法对某些三角剖分阈值和滤波阈值进行三角剖分。正确特征匹配的累积损失是造成漂移的主要原因之一。
为了解决这个问题，我们建议在增量 SfM 过程中对失败的特征匹配 (带延迟) 进行重新三角剖分 (RT)。两个摄像机之间可能存在的不良相对姿态的一个很好的标志是它们的公共点和它们的特征匹配之间的比率很低，我们称之为欠重构。为了等到姿势变得更加稳定，我们在几何序列(例如，当模型大小增加 25% 时，r'=25%) 下重新三角化未重构的对。为了获得更多的点数，我们还提高了 RT 过程中重投影错误的阈值。在对特征匹配进行重新三角剖分后，我们运行 FullBA 和点滤波来改进重建。由于几何序列，每个 RT 步需要 \(O(n)\) 个时间，并且累加到相同的 \(O(n)\) 个时间。
提出的 RT 步骤类似于闭环法，但是只在检测到循环时才处理漂移，通过寻找重构不足的对，我们的方法能够在没有显式的循环检测的情况下减少漂移误差，前提是有足够的特征匹配。实际上，RT 步骤更通用，因为它适用于任何匹配图像之间的相对位置，而且它还使得环路检测变得更加容易。图 4 显示了我们在一个 4K 图像循环上使用 RT 的增量 SfM，它使用 RT 正确地重建了长循环，而没有显式关闭循环。
6 Experiment¶
我们将我们的算法应用于五个不同大小的数据集。罗马中部数据集和艺术四边形数据集是从 【4】 的作者那里获得的，分别包含 32768 幅图像和 6514 幅图像。Loop 数据集是一个街区的 4342 帧高分辨率视频序列。圣彼得大教堂和罗马竞技场的数据集分别有 1275 张和 1164 张图像。我们在配备 Intel Xenon 5680 3.33 GHz CPU(24 核)、12 GB RAM 和 NVIDIA GTX 480 GPU 的 PC 上运行所有架构。
6.1 Feature Matching¶
为了在尽可能大的模型上进行实验，我们首先尝试匹配足够的图像对。完整的逐对数据集在 St. Peter’sBasilica and Colosseum 数据集上进行计算。对于 Arts Quad 和 Loop 数据集，根据 GPS 将每个图像与附近的图像进行匹配。对于 Central Rome 数据集，应用了 h=100 和 \((t_h)\) =4 的抢先匹配。
然后，我们使用满足抢占式匹配 h=100 和不同的\((t_h)\) 的图像匹配子集运行增量 SfM。表 1 显示了特征匹配的统计数据和最大重建模型的重建摄像机数量。通过先发制人匹配，我们能够为罗马重建 15065 台相机的大型 SFM 模型，并为其他数据集重建完整的模型。抢占式匹配能够显著减少常规匹配的图像对数，并且 STILL 保留了很大一部分正确的特征匹配，例如，对于圣彼得匹配，用 6% 的图像对获得了 43% 的特征匹配，值得注意的是，前 100 个特征的匹配速度非常快。我们的系统使用 24 个线程，平均速度为每秒 73K 对。
当使用较大的阈值时，抢占匹配具有较小的丢失弱链路的机会。当 \((t_h)\)=2 时，对所有数据集重建完整的模型。然而，对于 \((t_h)\)=4，我们发现艺术广场中的一座建筑由于遮挡而丢失，竞技场模型从内部和外部分开。我们认为，未来应该探索一种更具适应性的抢先匹配方案。
6.2 Incremental SfM¶
我们使用相同的设置对五个数据集的所有计算特征匹配运行我们的算法。表 2 总结了 r=5% 和 r‘=25% 的实验数据和时间。图 4、图 5 和图 7 显示了我们的 SfM 模型的屏幕截图。我们的方法高效地重建了大的、精确的、完整的模型，并且具有很高的点密度。特别地，两个 1K 图像数据集在少于 10 分钟进行重构，且 15K 罗马的相机模型在 1.67 小时中计算。各种设置下的其他结果可以在表 3 和补充材料中找到。
表 2 中的时序显示 BAS(全部 + 部分) 占据了重建时间的大部分，因此 BAS 的时间复杂度接近于增量 SFM 的时间复杂度。我们记录重建的时间和每一步的累积时间。如图 6 所示，重建时间（和累积的全 BAs 和部分 BAs）随着模型的增大大致呈线性增长，这验证了我们对时间复杂性的分析。Loop 重建具有较高的斜率，这主要是由于其较高的特征匹配数。
6.3 重建质量和速度¶
图 4、图 5 和图 7 展示了可与其他方法相媲美的高质量重建，我们也已经重建比在两个数据集上 DISCO 更多的相机。特别值得一提的是，图 5 显示了我们的 15065 相机型号的罗马中心在航拍图像上的正确覆盖。环路重构表明，RT 步长在没有显式闭环的情况下，正确地处理了漂移问题。首先通过 RT 将重建推离局部极小值，然后通过后续的 BAS 进行改进。该方法的鲁棒性和准确性得益于 BA 和 RT 的混合策略。特别值得一提的是，图 5 显示了我们的 15065 相机型号的罗马中心在航拍图像上的正确覆盖。环路重构表明，RT 步长在没有显式闭环的情况下，正确地处理了漂移问题。首先通过 RT 将重建推离局部极小值，然后通过后续的 BAS 进行改进。该方法的鲁棒性和准确性得益于 BA 和 RT 的混合策略。
对于 ARTICS 四边形数据集，我们的重建 (r=5%，r'=25%) 包含了 348 幅图像中的 261 幅，这些图像的地面真实 GPS 位置由 【4】 提供。我们使用 RANSAC 估计了 261 个摄像机位置与它们的欧几里得坐标之间的 3D 相似变换。经过最好的变换，我们的 3D 模型的平均误差为 2.5 米，中位误差为 0.89 米，比文献 【4】 报告的 1.16 米的误差要小。
7 Conclusion and Future Work¶
参考文献¶
References
[1] S. Agarwal, N. Snavely, S. Seitz, and R. Szeliski. Bundle adjustment in the large. In ECCV, pages II: 29–42, 2010.
[2] S. Agarwal, N. Snavely, I. Simon, S. M. Seitz, and R. Szeliski. Building Rome in a day. In ICCV, 2009.
[3] M. Byrod and K. Astrom. Conjugate gradient bundle adjustment. In ECCV, pages II: 114–127, 2010.
[4] D. Crandall, A. Owens, N. Snavely, and D. P. Huttenlocher. Discrete-continuous optimization for large-scale structure from motion. In CVPR, 2011.
[5] J. Frahm, P. Fite Georgel, D. Gallup, T. Johnson, R. Raguram, C. Wu, Y. Jen, E. Dunn, B. Clipp, S. Lazebnik, and M. Pollefeys. Building Rome on a cloudless day. In ECCV, pages IV: 368–381, 2010.
[6] R. Gherardi, M. Farenzena, and A. Fusiello. Improving the efficiency of hierarchical structure-and-motion. In CVPR, pages 1594–1600, 2010.
[7] A. Kushal, B. Self, Y. Furukawa, C. Hernandez, D. Gallup, B. Curless, and S. Seitz. Photo tours. In 3DimPVT, 2012.
[8] X. Li, C. Wu, C. Zach, S. Lazebnik, and J. Frahm. Modeling and recognition of landmark image collections using iconic scene graphs. In ECCV, 2008.
[9] M. A. Lourakis and A. Argyros. SBA: A Software Package for Generic Sparse Bundle Adjustment. ACM Trans. Math. Software, 36(1):1–30, 2009.
[10] D. G. Lowe. Distinctive image features from scale-invariant keypoints. IJCV, 60:91–110, 2004.
[11] D. Nister and H. Stewenius. Scalable recognition with a vocabulary tree. In CVPR, pages 2161–2168, 2006.
[12] J. R. Shewchuk. An introduction to the conjugate gradient method without the agonizing pain, 1994.
[13] S. N. Sinha, D. Steedly, and R. Szeliski. A multi-stage linear approach to structure from motion. In ECCV RMLE workshop, 2010.
[14] N. Snavely, S. Seitz, and R. Szeliski. Photo tourism: exploring photo collections in 3D. In SIGGRAPH, pages 835–846, 2006.
[15] N. Snavely, S. M. Seitz, and R. Szeliski. Skeletal graphs for efficient structure from motion. In CVPR, 2008.
[16] C. Wu, S. Agarwal, B. Curless, and S. M. Seitz. Multicore bundle adjustment. In CVPR, 2011.