Accurate Multiple View 3D Reconstruction using Patch-Based Stereo for Large-Scale Scenes——论文阅读笔记

论文阅读笔记：Accurate Multiple View 3D Reconstruction using Patch-Based Stereo for Large-Scale Scenes

• • 一、立体对选择
  • 二、深度图的计算
  • 三、深度图细化
  • 四、深度图融合

一、立体对选择

对于图像集中的每个图像，我们需要为其选择一个参考图像以进行立体计算。立体图像对的选择不仅关系到立体匹配的准确性，而且关系到最终的MVS结果。立体对选择对于车载摄像头采集的数据或受控于环境的摄像头采集的数据（Middlebury数据）来说是一项相对容易的任务，但针对无序图像就比较困难了。好的候选参考图像应具有与目标图像相似的观察方向，并且具有合适的基线，既不能太短而降低重建精度，也不能太长而导致具有较少的场景公共覆盖。

假设我们有n个图像，对其中的第i个图像，我们计算相机i和j的主视图方向之间的角度θ_{ij}(j=1,\ldots,n)。如果使用SFM算法校准相机姿态，则生成一组稀疏3D点，使用相机i和j的可见点和相机中心之间角度的平均值可以计算出更好的θ_{ij}。除了θ_{ij}，我们还计算了另一个参数d_{ij}(j=1,\ldots,n)，对于每个图像i，d_{ij}表示摄像机i和j的光学中心之间的距离。然后，对于满足5˚<θ_{ij}<60˚的图像，我们计算其d_{ij}的中值\bar{d}，并移除d_{ij}>2\bar{d}或d_{ij}<0.05\bar{d}的图像。在这些计算之后，如果剩余图像的数目小于k_1，则将这些图像视为图像i的相邻图像，表示为N(i)。否则，剩余图像按照θ_{ij}·d_{ij}按升序排序，前k_1张图像组成相邻图像N(i)（在本文中，我们设置k_1=10）。最后，选择N(i)中θ_{ij}·d_{ij}最小的图像作为第i幅图像的参考图像，形成立体对。

二、深度图的计算

对于每个符合条件的立体对，我们按照1中的想法计算深度贴图。其核心思想是，对输入图像中的每个像素，我们试图找到一个与参考图像具有最小聚合匹配代价的平面，如图{\rm{1}}所示。支撑平面f本质上是场景曲面的局部切平面，在相关摄像机的坐标系中，它由一个3D点和它的法线n_i表示。

对于图像集中的第i个输入图像I_i，给定其参考图像为I_j，以及相关的相机参数{\lbrace}K_i，R_i，C_i{\rbrace}和{\lbrace}K_j，R_j，C_j{\rbrace}，我们首先将I_i中的每个像素p分配给一个随机三维平面，其中K是相机内参，R是旋转矩阵，C是相机中心。假设p的齐次坐标为：p=\left[ \begin{array}{c} u\\ v \\ 1 \end{array} \right]
3D点X_i必须与像素p在同一条直线上，我们在深度范围\lambda\in[\lambda_{min},\lambda_{max}]中为\lambda选择一个随机，则X_i在C_i的坐标中计算为：X_i={\lambda}K^{-1}_ip
\color{blue}{这个公式就是将像素坐标转化为世界坐标，转换公式为：}Z_c\left[ \begin{array}{c} u\\ v\\ 1 \end{array} \right]=\left[ \begin{array}{c} \frac{1}{d_x} & 0 & u_0\\ 0 & \frac{1}{d_y} & v_0\\ 0& 0& 0 \end{array} \right] \left[ \begin{array}{c} f& 0 & 0&0\\ 0 & f& 0&0\\ 0& 0&1&0 \end{array} \right] \left[ \begin{array}{c} R& T\\ \vec{0} &1\\ \end{array} \right] \left[ \begin{array}{c} X_w\\ Y_w\\ Z_w\\ 1 \end{array} \right]= \left[ \begin{array}{c} f_x& 0&u_0&0\\ 0 &f_y&v_0&0\\ 0&0&1&0 \end{array} \right] \left[ \begin{array}{c} R& T\\ \vec{0} &1\\ \end{array} \right] \left[ \begin{array}{c} X_w\\ Y_w\\ Z_w\\ 1 \end{array} \right]
\color{blue}{其中}\left[ \begin{array}{c} f_x& 0&u_0&0\\ 0 &f_y&v_0&0\\ 0&0&1&0 \end{array} \right]\color{blue}{为相机的内参数矩阵，}\left[ \begin{array}{c} R& T\\ \vec{0} &1\\ \end{array} \right]\color{blue}{为外参数矩阵，}Z_c=\lambda\color{blue}{由于世界坐标系原点设置在相机坐标系原点，故外参数矩阵为零}

然后，我们在摄影机C_i的球坐标中随机指定平面的法线，如下所示：n_i=\left[ \begin{array}{c} \cos{\theta}\sin{\phi}\\ \sin{\theta}\sin{\phi}\\ \cos{\phi} \end{array} \right]
其中，{\theta}是[0°,360°]中的任意角度，{\phi}是[0°,60°]内的任意角度。这些范围设置源于一个简单的假设，即当\rm{patch}法线与摄影机C_i坐标系的z轴之间的角度低于某个阈值（在本文中，我们将该阈值设置为60˚）时，面片在图像I_i中可见。

\color{blue}{如上图所示，}n_z=R\cos\phi,n_x=R\sin\phi\cos\theta,n_y=R\sin\phi\sin\theta，\color{blue}{其中}R=1。
上述随机初始化过程很可能为图像中的每个场景平面找到至少一个良好的猜测，尤其是对于场景平面包含更多像素的高分辨率图像，这意味着高分辨图像比低分辨率图像有更多的猜测。我们应该注意，一旦计算了图像I_i的depthmap，我们就可以在计算I_i的参考图像I_j的{depthmap}时改进纯粹的随机初始化过程，在计算I_j深度图时，I_i深度图中每个像素的深度和\rm{patch}法线可以作为I_j的初始估计值，I_j中I_i和I_j之间没有映射的像素仍然使用随机初始化。通过这种方式，我们可以为I_j中的每个像素分配一个更好的平面，因为平面对于立体对I_i和I_j是一致的。根据2，给定两个摄像机的投影矩阵P=[I_{3×3}|0_3]]和P^{'}=[R|t]，以及由{\pi}^{T}X=0和{\pi}=(V^{T},1)^{T}定义的平面，则该平面推导出来的单应性矩阵H为：
H=R-tV^{T}
这里，I_{3×3}是3×3单位矩阵，0_3是零矢量，表示选择的世界坐标与相机P重合。本文中，图像对的相机参数为{\lbrace}K_i，R_i，C_i{\rbrace}和{\lbrace}K_j，R_j，C_j{\rbrace}，平面f_p={\lbrace}X_i,n_i{\rbrace}在相机C_i坐标系中定义。因此，投影矩阵和平面参数可以转换为标准形式（将世界坐标系原点置于C_i），如下所示:
P_i=K_i[I_{3×3}|0_3],P_j=K_j[R_jR^{-1}_{i}|R_j(C_i-C_j)], V^{T}=-\frac{n^{T}_i}{n^{T}_{i}X_i}
根据公式H=R-tV^{T}，相机\tilde{P_i}=[I_{3×3}|0_3]和\tilde{P_j}=[R_jR^{-1}_i|R_j(C_i-C_j)]的单应性矩阵为：\tilde{H_{ij}}=R_jR^{-1}_i+\frac{R_j(C_i-C_j)n^{T}_i}{n^{T}_{i}X_i}
将变换K_i和K_j应用于图像，我们得到摄像机P_i=K_i\tilde{P_i}，P_j=K_j\tilde{P_j}，得到的诱导单应性矩阵为：H_{ij}=K_j(R_jR^{-1}_j+\frac{R_j(C_i-C_j)n^{T}_i}{n^{T}_{i}X_i})K^{-1}_i我们设置了一个以像素p为中心的方形窗口B，其中B=w×w（在本文中，我们设置了w=7个像素）。对于B中的每个像素q，我们使用单应映射H_{ij}(q)在参考图像I_j中找到其对应的像素。然后，像素p的聚合匹配成本m(p，fp)用1减去q和H_{ij}(q)之间的归一化互相关（NCC）分数，
m(p,f_p)=1-\frac{\sum_{q\in{B}}(q-\bar{q})(H_{ij}(q)-\bar{H_{ij}(q)})}{\sqrt{\sum_{q\in{B}}(q-\bar{q})^2\sum_{q\in{B}}(H_{ij}(q)-\bar{H_{ij}(q)})^2}}
请注意，一些更复杂、更稳健的聚合技术，如345，可用于生成比NCC更可靠的结果。然而，高分辨率图像可以提供比低分辨率图像更可靠的匹配，并且简单的NCC足以可靠地测量光度一致性。此外，在深度图细化步骤中可以去除NCC生成的大部分不可靠像素，使得NCC的最终重建结果与其他复杂聚集方法的重建结果几乎相同。因此，在本文中，我们使用简单NCC作为聚合匹配成本，这与6相同。
初始化之后，图像I_i中的每个像素与一个3D平面相关联。然后在n_2次迭代中处理I_i中的每一个像素，以细化平面。在奇数迭代中，我们从左上角的像素开始，按行顺序遍历，直到到达右下角的像素。在偶数迭代中，我们反过来，从右下角遍历到到左上角像素，同样是按行顺序。本文将平面细化的次数设为k_2=3。
在每次迭代中，每个像素都有两个操作，空间传播和随机分配。空间传播用于将相邻像素的平面与当前像素的平面进行比较和传播。在奇数迭代中，相邻像素是左边、上面和左上的相邻像素，在偶数迭代中是右、下和右下的相邻像素。令p_N表示当前像素p的相邻像素，f_{p_N}为p_N的平面，我们使用等式m(p,f_p)=1-\frac{\sum_{q\in{B}}(q-\bar{q})(H_{ij}(q)-\bar{H_{ij}(q)})}{\sqrt{\sum_{q\in{B}}(q-\bar{q})\sum_{q\in{B}}(H_{ij}(q)-\bar{H_{ij}(q)})^2}}中的匹配代价来检查条件m(p，f_{p_N})，如果满足这个条件，我们认为f_{p_N}是比p的当前平面f_p更好的选择，并令f_p=f_{p_N}，即空间传播。这种空间传播过程依赖于一个事实，即相邻像素很可能与当前像素具有相似的3D平面，尤其是对于高分辨率图像。理论上，，即使在初始化的时候只有一个好的猜测，在第一次和第二次空间传播之后也足以将该平面传播到该区域的其他像素上。
对于每个像素p，在空间传播之后，我们使用随机分配进一步细化平面f_p。随机分配的目的是通过测试多个随机平面参数来进一步降低等式m(p,f_p)=1-\frac{\sum_{q\in{B}}(q-\bar{q})(H_{ij}(q)-\bar{H_{ij}(q)})}{\sqrt{\sum_{q\in{B}}(q-\bar{q})\sum_{q\in{B}}(H_{ij}(q)-\bar{H_{ij}(q)})^2}}中的匹配成本。给定一个范围{\lbrace\Delta\lambda,\Delta\theta,\Delta\phi\rbrace}，我们

1. 在范围\lambda^{'}\in(\phi-\Delta\phi,\phi+\Delta\phi)，\theta^{'}\in(\theta-\Delta\theta,\theta+\Delta\theta),\phi^{'}\in(\phi-\Delta\phi,\phi+\Delta\phi)中随机选择一个平面参数\lbrace\lambda^{'},\theta^{'},\phi^{'}\rbrace。
2. 使用公式X_i={\lambda}K^{-1}_ip和公式n_i=\left[ \begin{array}{c} \cos{\theta}\sin{\phi}\\ \sin{\theta}\sin{\phi}\\ \cos{\phi} \end{array} \right]重新计算新的平面f^{'}_p=\lbrace{X^{'}_i},n^{'}_i\rbrace
3. 如果m(p,f^{'}_p),则令f_p=f^{'}_p,\lambda=\lambda^{'},\theta=\theta^{'},\phi=\phi^{'}。
4. 将范围\lbrace\Delta\lambda,\Delta\theta,\Delta\phi\rbrace减半。
5. 重复第一步

上述过程重复k_3次。在本文中，我们将初始范围和重复时间设置为：\Delta\lambda=\frac{\lambda_max-\lambda_min}{4},\Deltaθ=90˚,\Deltaφ=15˚,k_3=6。此随机分配过程逐渐缩小搜索范围，以便捕获深度和法线细节。

空间传播和随机分配思想已成功应用于patchmatch立体匹配方法1和混合递归匹配（HRM）方法7。本文将1的思想扩展到使用高分辨率图像的大规模场景的多视图立体。所提出方法的新颖之处在于以适当的方式修改pathmatch立体算法1，使其对于大规模MVS问题更加强大和高效。所提出的方法与1中的方法之间的主要区别在于，论文1是将平面定义在图像坐标系，而本文是将坐标系定义在相机坐标系，因为本文没有对立体对进行校正。此外，通过充分利用多幅高分辨率图像，我们对1进行了三次简化，以减少计算开销。首先，我们的聚合匹配成本m(p，f_p)是一个简单的归一化互相关，而不是1中更复杂的自适应支持权重版本，因为NCC对于高分辨率图像足够可靠，并且剩余的误差可以在下面的深度图细化步骤中消除。其次，与1中的空间和视图传播相比，我们只使用空间传播，因为我们只计算I_i上的dethmap，而不是[32]中的I_i和I_j上的深度贴图都要计算。第三，[32]中的方法包含另一个后处理步骤，该步骤通过左/右一致性检查来进行遮挡处理，并填充无效像素。我们的方法没有这个过程，因为下面的深度贴图细化步骤可以达到类似的效果。在空间传播和随机分配过程之后，我们移除深度图中聚合匹配代价高于某个阈值\tau_1的不可靠点（在本文中，我们设置\tau_1=0.3）。

三、深度图细化

由于深度误差，原始的depthmap可能在公共区域上彼此不完全一致，因此执行细化过程以加强相邻视图的一致性。对于图像I_i中的每个点p，我们使用其深度\lambda和相机参数将其反向投影到3D，如下所示：
X=\lambda{R^{T}_i}K^{-1}_ip+C_i
其中p是等式1中定义的齐次坐标，X是世界坐标中的3D点。然后，我们将X投影到的I_i的相邻图像N(i)。假设N_k是N(i)中的第k个相邻图像，我们将d(X，N_k)表示X相对于摄像机N_k的深度，\lambda(X，N_k)表示将X投影到N_k的深度图上，投影处的深度值。如果\lambda(X,N_k)与d(X，N_k)足够接近，即：\frac{|d(X,N_k)-\lambda(X,N_k)|}{\lambda(X,N_k)}<\tau_2，其中τ_2是一个阈值（在本文中，我们设置τ_2=0.01），我们就说X在I_i和N_k中是一致的。如果X与N(i)中至少k_4张图像是一致的（在本文中，我们设置k_4=2），它就被视为可靠的场景点，并且它在I_i的深度图中对应的像素p被保留，否则它被移除。经过上述细化过程后，大多数错误都可以被删除，从而在每个视图中生成一个相对干净的深度图。

四、深度图融合

细化后，可以合并所有深度图以表示场景。然而，直接合并深度贴图可能包含大量冗余，因为不同的深度图可能覆盖场景的相同区域，尤其是对于相邻图像。为了消除这些冗余，通过相邻深度图测试进一步减少深度贴图。如图4所示，对于相机C_i深度图中的每个像素，我们使用公式7将其反向投影到3D 中记为X，并将X重新投影到C_i的相邻摄像头。如果X相对于相邻摄像机的深度小于相邻摄像机深度图中X的投影位置的深度值，如图4中的摄像机N_1和N_2，我们说该相邻摄像机中X的投影被遮挡，并从该相邻摄像机的深度图中移除。如果这两个深度值足够接近，如图4中的相机N_4，我们说该相邻相机中的X投影与X是相同的点，这是一个冗余，将其从该相邻相机的深度贴图中移除。

最后，将所有深度图反向投影到3D中，并合到单个点云中。最终的点云通常非常密集，尤其是在使用高分辨率图像时。如果我们想使其稀疏，我们可以简单地在深度图中的稀疏位置返回投影点。例如，仅使用在深度图中图像位置为(2n，2n)的像素点，这将使点云的大小大约减小到使用所有点的大小的四分之一。这为我们提供了一种根据内存和存储限制控制点云大小的方法。

---

1. http://www.bmva.org/bmvc/2011/proceedings/paper14/paper14.pdf ↩︎ ↩︎ ↩︎ ↩︎ ↩︎ ↩︎ ↩︎ ↩︎ ↩︎

2. http://www.r-5.org/files/books/computers/algo-list/image-processing/vision/Richard_Hartley_Andrew_Zisserman-Multiple_View_Geometry_in_Computer_Vision-EN.pdf ↩︎

3. http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.99.9581&rep=rep1&type=pdf ↩︎

4. http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.160.5145&rep=rep1&type=pdf ↩︎

5. https://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.219.3390&rep=rep1&type=pdf ↩︎

6. http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.228.340&rep=rep1&type=pdf ↩︎

7. https://web.archive.org/web/20170810233955id_/http://www.item.ntnu.no/fag/ttm4142/host/PensumArtikler/StereoMatchingAtzpadin04.pdf ↩︎