MVS Net:Depth Inference for UNstructured Multi-view Stereo
目录
论文地址: https://arxiv.org/abs/1804.02505
项目地址: https://github.com/YoYo000/MVSNet
一、Background(背景)
-
多视图立体匹配(Multi-view Stereo, MVS)是计算机领域中一个核心问题,其目的是通过不同视点拍摄的图像,恢复出真实的三维场景。
-
传统的三维重建使用手工设计的相似性度量指标和正则化方法计算场景的稠密对应关系。这些方法在非朗伯体表面、无弱纹理区域的场景可以达到很好的效果。但是在弱纹理区域和朗伯体表面,手工设计的相似性指标变得不可信,导致重建结果不完整。
二、Framework(网络架构)
-
针对传统的手工特征设计的特征鲁棒性不强的问题,作者将深度学习引入到三维重构中。
- 在2D图像上进行特征提取得到特征图。
- 通过可微分的单应变换,基于参考视图的相机视锥体构建3D代价体。
- 使用3D卷积对代价体进行正则化,回归得到初始深度图。
- 通过参考视图图像优化得到最后的估计的深度图。
1、Feature extraction(特征提取)
- 每个batch选定N张图片,其中一张为参考视图(Reference image),其余为源视图(Source images)。
- 对选定的 I 1 I_{1} I1(Reference image)和多张 [ I i ] 2 N [I_{i}]^{N}_{2} [Ii]2N(Source images)用CNN进行特征提取,经过两次strides为2的缩放,得到channel为32,长和宽分别为h/4和w/4的feature map;
- 为了减小参数量和提高计算效率,Feature Extraction Network权值共享。
- 经过CNN特征提取的特征图,具有较高的语义信息和上下文信息。
2、Differentiable Homography(可微分单应性变化)
- 因为 I 1 I_{1} I1和 [ I i ] 2 N [I_{i}]^{N}_{2} [Ii]2N的视角不一致,模型的目的是估计参考视图 I 1 I_{1} I1的深度图,需要将 [ I i ] 2 N [I_{i}]^{N}_{2} [Ii]2N经过可微分的单应变换warp到 I 1 I_{1} I1对应的相机坐标系内。
- 根据先验的深度范围信息,以
I
1
I_{1}
I1的主光轴为扫描方向,按照固定的最小深度间隔
Θ
s
c
a
l
e
\Theta_{scale}
Θscale,经过可微分的单应性变换,将
[
I
i
]
2
N
[I_{i}]^{N}_{2}
[Ii]2N提取出的特征体积从
Θ
m
i
n
\Theta_{min}
Θmin 映射到
Θ
m
a
x
\Theta_{max}
Θmax,得到
[
I
i
]
2
N
[I_{i}]^{N}_{2}
[Ii]2N的视锥体(源视图的特征体经过相机内外参和目标深度d,映射到目标视图中得到不同深度的特征体):
H i ( d ) = K i R i ( I − ( R 1 T ⋅ t 1 − R i T ⋅ t i ) d n T ⋅ R 1 ) R 1 T K 1 − 1 (1) \mathbf{H}_{i}(d)=K_{i} R_{i}\left(I-\frac{\left(R_{1}^{T} \cdot t_{1}-R_{\mathrm{i}}^{T} \cdot t_{i}\right)}{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1}\tag{1} Hi(d)=KiRi(I−d(R1T⋅t1−RiT⋅ti)nT⋅R1)R1TK1−1(1)
式中:
H i H_{i} Hi 代表可微分单应性变换 d代表深度
K i K_{i} Ki 代表第i个视图的内参矩阵
R i R_{i} Ri 代表第i个视图的旋转矩阵
t i t_{i} ti 代表第i个视图的平移矩阵
n 1 T n_{1}^{T} n1T 代表参考视图的主光轴方向
注:原文中的公式有误,此处公式与原文不同,具体推导在附录
3、Volume Cost(代价体)
![](https://i-blog.csdnimg.cn/blog_migrate/31a63d74260f8f9a5eae96dbe9ec6748.png)
- 为了使图像大小一致,将视锥体通过线性插值到[D, C, h/4, w/4]大小的特征体, N 张视图可形成 N 个特征体,这一步完成了从2D到3D的转换,后面基于3D的特征体进行运算。
- 在N个特征体的基础上,基于方差计算多个特征体的匹配代价,这样基于方差来计算匹配代价的方式优于基于均值计算匹配代价,因为方差更能体现不同视图之间的差异性信息。
C = M ( V 1 , ⋯ , V N ) = ∑ i = 1 N ( V i − V i ‾ ) 2 N (2) \mathbf{C}=\mathcal{M}\left(\mathbf{V}_{1}, \cdots, \mathbf{V}_{N}\right)=\frac{\sum_{i=1}^{N}\left(\mathbf{V}_{i}-\overline{\mathbf{V}_{i}}\right)^{2}}{N}\tag{2} C=M(V1,⋯,VN)=N∑i=1N(Vi−Vi)2(2)
4、Volume Cost Regularization(代价体正则化)
![](https://i-blog.csdnimg.cn/blog_migrate/8cbee2f64de9e7191fa71d998448efb0.png)
- 直接基于方差计算得到的**初始特征体(initial cost volume)**带有较多的噪声,这些噪声主要由视图间的遮挡等因素造成。
- 基于3D Convolutional Network 对 initial cost Volume 进行正则化,采用编码解码的原理,利用一个3D Unet网络的结构得到的特征体通过将C通道维度降为1,得到一个概率体(probability volume),其维度为[D, h/4, w/4],此时输出在D维度上做softMax,在D维度表示每个像素沿深度方向的概率,称之为概率体。
- 可以采用argmax的方式来求得该像素点的深度值(赢者通吃原则),但argmax是不可微的操作,这样的操作梯度不可反向传播,作者在这里采用了求期望的方式(soft argmin),估计每一个像素点的深度值。
D = ∑ d = d min d max d × P ( d ) (3) \mathbf{D}=\sum_{d=d_{\min }}^{d_{\max }} d \times \mathbf{P}(d)\tag{3} D=d=dmin∑dmaxd×P(d)(3)
5、Depth Map Refinement (深度图增强)
![](https://i-blog.csdnimg.cn/blog_migrate/04bb15c97a602ae65483c198e0e1f01e.png)
- 由于3D CNN在正则化过程中,较大的感受野造成重建深度图边界过平滑,作者将 I 1 I_{1} I1作为一个辅助项和Initial depth map拼接,送入一个带有残差的CNN模块中进行特征融合,得到一个包含更多边界信息的深度图。
6、Loss(损失函数)
![](https://i-blog.csdnimg.cn/blog_migrate/8e5551c87f07629be15667767fbd5c03.png)
- 使用真实深度图与所估计深度的
L
1
L_{1}
L1损失作为训练损失。计算损失时候不考虑背景的深度,背景深度不进行反向传播,代码实现中用mask将背景剔除。
Loss = ∑ p ∈ p valid ∥ d ( p ) − d ^ i ( p ) ∥ 1 ⏟ Loss 0 + λ ⋅ ∥ d ( p ) − d ^ r ( p ) ∥ 1 ⏟ Loss 1 (4) \text { Loss }=\sum_{p \in \mathbf{p}_{\text {valid }}} \underbrace{\left\|d(p)-\hat{d}_{i}(p)\right\|_{1}}_{\text {Loss } 0}+\lambda \cdot \underbrace{\left\|d(p)-\hat{d}_{r}(p)\right\|_{1}}_{\text {Loss } 1}\tag{4} Loss =p∈pvalid ∑Loss 0 d(p)−d^i(p) 1+λ⋅Loss 1 d(p)−d^r(p) 1(4)
三、post Processing
![](https://i-blog.csdnimg.cn/blog_migrate/89bc7693da1c6366c6f54a2d859e2ba6.png)
- 通过MVS Net得到参考图像每个像素的深度估计值,进行下一步的稠密点云计算之前,使用光度一致性与几何一致性作为条件对背景区域与遮挡区域的点进行进一步的剔除。
- 光照一致性条件:沿深度方向的概率分布可以反映深度估计的质量,3D CNN的正则化也有助于将概率分布调整为单峰分布,但由于错误匹配的存在,导致深度概率分布比较离散。因为深度假设是在相机视锥体内离散采样得到的,将四个最近深度的概率求和得到每个像素点最后的估计质量。文章中将估计质量概率值低于0.8的视为外点。
- 几何一致性条件: ∣ P p r e p r o j − P 1 ∣ < 1 |P_{preproj}-P_{1}|<1 ∣Ppreproj−P1∣<1, ∣ d p r e p r o j − d 1 ∣ / d 1 < 0.01 |d_{preproj}-d_{1}|/d_{1} < 0.01 ∣dpreproj−d1∣/d1<0.01,就称像素p_1 处的深度估计值d_1是两视图连续。
四、Result
![](https://i-blog.csdnimg.cn/blog_migrate/cfb708f4c2b78be7f9ca718bdce3ec27.png)
![](https://i-blog.csdnimg.cn/blog_migrate/f2ace3ec95b472f38af7b3343dcd5f88.png)
附录:可微分单应性变换的过程推导:
由对极几何的原理可知(世界坐标系到相机坐标系过程中采用先平移,再旋转),假设要将
i
i
i视图的特征warp到参考视图的
d
d
d 深度:
对于视图1:
d
⋅
p
1
=
K
1
⋅
p
c
1
=
K
1
R
1
(
P
w
−
C
1
)
(5)
d\cdot p_{1}=K_{1} \cdot p_{c_{1}}=K_{1} R_{1}\left(P_{w}-C_{1}\right)\tag{5}
d⋅p1=K1⋅pc1=K1R1(Pw−C1)(5)
对于视图
i
i
i:
Z
c
i
⋅
p
i
=
K
i
⋅
p
c
i
=
K
i
R
i
(
P
w
−
C
i
)
(6)
Z_{c_{i}} \cdot p_{i}=K_{i} \cdot p_{c_{i}}=K_{i} R_{i}\left(P_{w}-C_{i}\right)\tag{6}
Zci⋅pi=Ki⋅pci=KiRi(Pw−Ci)(6)
由投影的几何信息可知:
n
T
⋅
p
c
1
−
d
=
0
(7)
n^{T} \cdot p_{c_{1}}-d=0\tag{7}
nT⋅pc1−d=0(7)
其中:d为参考视图相机坐标系下的深度,
p
i
p_{i}
pi为
i
i
i 视图上的像素点坐标,
K
i
K_{i}
Ki 为
i
i
i 视图的内参矩阵,
p
c
i
p_{c_{i}}
pci 为
i
i
i 视图相机坐标系的坐标,
P
w
P_{w}
Pw 点为对应世界坐标系的坐标,
R
i
R_{i}
Ri 为
i
i
i 视图的旋转矩阵,
C
i
C_{i}
Ci 为
i
i
i 视图的平移矩阵,
n
T
n^{T}
nT 为参考视图
z
z
z 方向的方向向量;
联立式5和式6得:
Z
c
i
p
i
=
K
i
⋅
p
c
i
=
K
i
R
i
(
P
w
−
C
i
)
=
K
i
R
i
(
R
1
T
p
c
1
+
C
1
−
C
i
)
(8)
Z_{c_{i}} p_{i}=K_{i} \cdot p_{c_{i}}=K_{i} R_{i}\left(P_{w}-C_{i}\right)=K_{i} R_{i}\left(R_{1}^{T} p_{c_{1}}+C_{1}-C_{i}\right)\tag{8}
Zcipi=Ki⋅pci=KiRi(Pw−Ci)=KiRi(R1Tpc1+C1−Ci)(8)
将式7代入式8中,利用R矩阵为正交矩阵,且得:
Z
c
i
p
i
=
K
i
⋅
p
c
i
=
K
i
R
i
(
P
w
−
C
i
)
=
K
i
R
i
(
R
1
T
p
c
1
+
C
1
−
C
i
)
=
K
i
R
i
(
R
1
T
p
c
1
+
(
C
1
−
C
i
)
n
T
⋅
p
c
1
d
)
=
K
i
R
i
(
I
+
(
C
1
−
C
i
)
n
T
d
⋅
R
1
)
R
1
T
p
c
1
=
K
i
R
i
(
I
+
(
C
1
−
C
i
)
n
T
d
⋅
R
1
)
R
1
T
K
1
−
1
p
1
d
(9)
\begin{array}{l} Z_{c_{i}} p_{i}=K_{i} \cdot p_{c_{i}}=K_{i} R_{i}\left(P_{w}-C_{i}\right)=K_{i} R_{i}\left(R_{1}^{T} p_{c 1}+C_{1}-C_{i}\right) \\\\ =K_{i} R_{i}\left(R_{1}^{T} p_{c_{1}}+\frac{\left(C_{1}-C_{i}\right) n^{T} \cdot p_{c_{1}}}{d}\right) \\\\ =K_{i} R_{i}\left(I+\frac{\left(C_{1}-C_{i}\right) n^{T}}{d} \cdot R_{1}\right) R_{1}^{T} p_{c_{1}} \\\\ =K_{i} R_{i}\left(I+\frac{\left(C_{1}-C_{i}\right) n^{T}}{d} \cdot R_{1}\right) R_{1}^{T} K_{1}^{-1} p_{1} d \end{array}\tag{9}
Zcipi=Ki⋅pci=KiRi(Pw−Ci)=KiRi(R1Tpc1+C1−Ci)=KiRi(R1Tpc1+d(C1−Ci)nT⋅pc1)=KiRi(I+d(C1−Ci)nT⋅R1)R1Tpc1=KiRi(I+d(C1−Ci)nT⋅R1)R1TK1−1p1d(9)
即:
p
i
=
K
i
R
i
(
I
+
(
C
1
−
C
i
)
d
n
T
⋅
R
1
)
R
1
T
K
1
−
1
p
1
⋅
d
Z
c
i
(10)
p_{i}=K_{i} R_{i}\left(I+\frac{\left(C_{1}-C_{i}\right) }{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1} p_{1} \cdot \frac{d}{Z_{c_{i}}}\tag{10}
pi=KiRi(I+d(C1−Ci)nT⋅R1)R1TK1−1p1⋅Zcid(10)
由于
p
i
p_{i}
pi 与
p
1
p_{1}
p1 为齐次坐标,因此可以将常数
d
Z
c
i
\frac{d}{Z_{c_{i}}}
Zcid消除得:
p
i
=
K
i
R
i
(
I
+
(
C
1
−
C
i
)
d
n
T
⋅
R
1
)
R
1
T
K
1
−
1
p
1
(11)
p_{i}=K_{i} R_{i}\left(I+\frac{\left(C_{1}-C_{i}\right)}{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1} p_{1}\tag{11}
pi=KiRi(I+d(C1−Ci)nT⋅R1)R1TK1−1p1(11)
即:
p
i
=
H
(
d
)
⋅
p
1
(12)
p_{i}=H(d)\cdot p_{1}\tag{12}
pi=H(d)⋅p1(12)
其中
H
(
d
)
=
K
i
R
i
(
I
+
(
C
1
−
C
i
)
d
n
T
⋅
R
1
)
R
1
T
K
1
−
1
H(d)=K_{i} R_{i}\left(I+\frac{\left(C_{1}-C_{i}\right) }{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1}
H(d)=KiRi(I+d(C1−Ci)nT⋅R1)R1TK1−1
由
{
P
W
=
R
1
T
⋅
p
c
1
+
C
1
P
W
=
R
i
T
⋅
p
c
i
+
C
i
\left\{\begin{array}{l} P_{W}=R_{1}^{T} \cdot p_{c_{1}}+C_{1} \\ P_{W}=R_{\mathrm{i}}^{T} \cdot p_{ci}+C_{i} \end{array}\right.
{PW=R1T⋅pc1+C1PW=RiT⋅pci+Ci得:
C
1
−
C
i
=
R
i
T
⋅
p
c
i
−
R
1
T
⋅
p
c
1
(13)
C_{1}-C_{i}=R_{i}^{T} \cdot p_{c_{i}}-R_{1}^{T} \cdot p_{c_{1}}\tag{13}
C1−Ci=RiT⋅pci−R1T⋅pc1(13)
由
{
p
c
1
=
R
1
⋅
P
W
+
t
1
p
c
i
=
R
i
⋅
P
W
+
t
i
\left\{\begin{array}{l} p_{c_{1}}=R_{1} \cdot P_{W}+t_{1} \\ p_{c_{i}}=R_{i} \cdot P_{W}+t_{i} \end{array}\right.
{pc1=R1⋅PW+t1pci=Ri⋅PW+ti得:
R
i
T
⋅
p
c
i
−
R
1
T
⋅
p
c
1
=
R
i
T
⋅
t
i
−
R
1
T
⋅
t
1
(14)
R_{i}^{T} \cdot p_{c_{i}}-R_{1}^{T} \cdot p_{c_{1}} = R_{i}^{T} \cdot t_{i}-R_{1}^{T} \cdot t_{1}\tag{14}
RiT⋅pci−R1T⋅pc1=RiT⋅ti−R1T⋅t1(14)
将式13代入式14中得:
C
1
−
C
i
=
−
(
R
1
T
⋅
t
1
−
R
i
T
⋅
t
i
)
(15)
C_{1}-C_{i}=-\left(R_{1}^{T} \cdot t_{1}-R_{\mathrm{i}}^{T} \cdot t_{i}\right)\tag{15}
C1−Ci=−(R1T⋅t1−RiT⋅ti)(15)
将式15代入式11中得:
p
i
=
K
i
R
i
(
I
−
(
R
1
T
⋅
t
1
−
R
i
T
⋅
t
i
)
d
n
T
⋅
R
1
)
R
1
T
K
1
−
1
p
1
(16)
p_{i}=K_{i} R_{i}\left(I-\frac{\left(R_{1}^{T} \cdot t_{1}-R_{\mathrm{i}}^{T} \cdot t_{i}\right)}{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1} p_{1}\tag{16}
pi=KiRi(I−d(R1T⋅t1−RiT⋅ti)nT⋅R1)R1TK1−1p1(16)
即:
H
(
d
)
=
K
i
R
i
(
I
−
(
R
1
T
⋅
t
1
−
R
i
T
⋅
t
i
)
d
n
T
⋅
R
1
)
R
1
T
K
1
−
1
(17)
H(d) = K_{i} R_{i}\left(I-\frac{\left(R_{1}^{T} \cdot t_{1}-R_{\mathrm{i}}^{T} \cdot t_{i}\right)}{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1}\tag{17}
H(d)=KiRi(I−d(R1T⋅t1−RiT⋅ti)nT⋅R1)R1TK1−1(17)
由此得到
H
(
d
)
H(d)
H(d)的表达式,在代码中的实现:
由式16得:
Z
i
p
i
=
K
i
⋅
p
c
i
=
K
i
R
i
(
I
−
(
R
1
T
⋅
t
1
−
R
i
T
⋅
t
i
)
d
n
T
⋅
R
1
)
R
1
T
K
1
−
1
p
1
d
=
K
i
R
i
R
1
T
K
1
−
1
p
1
⋅
d
−
K
i
R
i
R
1
T
⋅
t
1
−
R
i
T
⋅
t
i
d
n
T
⋅
R
1
R
1
T
K
1
−
1
p
1
d
=
K
i
R
i
R
1
T
K
1
−
1
p
1
⋅
d
+
K
i
R
i
(
R
i
T
⋅
t
i
−
R
1
T
⋅
t
1
)
n
T
⋅
K
1
−
1
p
1
=
K
i
R
i
R
1
T
K
1
−
1
p
1
⋅
d
+
K
i
R
i
(
R
i
T
⋅
t
i
−
R
1
T
⋅
t
1
)
n
T
⋅
p
c
1
d
=
K
i
R
i
R
1
T
K
1
−
1
p
1
⋅
d
+
K
i
R
i
(
R
i
T
⋅
t
i
−
R
1
T
⋅
t
1
)
(18)
\begin{array}{l} Z_{i} p_{i}=K_{i} \cdot p_{c i}=K_{i} R_{i}\left(I-\frac{\left(R_{1}^{T} \cdot t_{1}-R_{\mathrm{i}}^{T} \cdot t_{i}\right)}{d} n^{T}\cdot R_{1}\right) R_{1}^{T} K_{1}^{-1} p_{1}d \\\\ =K_{i} R_{i} R_{1}^{T} K_{1}^{-1} p_{1}\cdot d-K_{i} R_{i}\frac{R_{1}^{T} \cdot t_{1}-R_{\mathrm{i}}^{T} \cdot t_{i}}{d} n^{T}\cdot R_{1}R_{1}^{T} K_{1}^{-1} p_{1}d\\\\ =K_{i} R_{i} R_{1}^{T} K_{1}^{-1} p_{1}\cdot d+K_{i} R_{i}(R_{i}^{T} \cdot t_{i}-R_{\mathrm{1}}^{T} \cdot t_{1})n^{T}\cdot K_{1}^{-1} p_{1}\\\\ =K_{i} R_{i} R_{1}^{T} K_{1}^{-1} p_{1}\cdot d+K_{i} R_{i}(R_{i}^{T} \cdot t_{i}-R_{\mathrm{1}}^{T} \cdot t_{1})n^{T}\cdot \frac{p_{c_{1}}}{d} \\\\ =K_{i} R_{i} R_{1}^{T} K_{1}^{-1} p_{1}\cdot d+K_{i} R_{i}(R_{i}^{T} \cdot t_{i}-R_{\mathrm{1}}^{T} \cdot t_{1}) \end{array}\tag{18}
Zipi=Ki⋅pci=KiRi(I−d(R1T⋅t1−RiT⋅ti)nT⋅R1)R1TK1−1p1d=KiRiR1TK1−1p1⋅d−KiRidR1T⋅t1−RiT⋅tinT⋅R1R1TK1−1p1d=KiRiR1TK1−1p1⋅d+KiRi(RiT⋅ti−R1T⋅t1)nT⋅K1−1p1=KiRiR1TK1−1p1⋅d+KiRi(RiT⋅ti−R1T⋅t1)nT⋅dpc1=KiRiR1TK1−1p1⋅d+KiRi(RiT⋅ti−R1T⋅t1)(18)