从点云到网格（二）VRIP介绍

VRIP（Volumetric Range Image Processing），顾名思义，是从深度图重建网格的一种方法。VRIP是Brian Curless和Marc Levoy在1996年提出来的方法，距今已经有20年的历史了，依然属于最好的方法之一。

VRIP的核心问题是

已知世界坐标系下，某物体表面(f)在不同视角下的深度图(hat{f_1},...,hat{f_K})，求(f)。这里隐含深度图在世界坐标系下的位姿是已知的。

许多三维测量方法，比如激光、TOF、结构光、双目视觉等，都可以得到深度图。因此这是一个非常有意义的问题。如何将这些深度图融合成为一个平滑的单一网格呢？这就是VRIP要解决的问题。
下面给出一个结构光成像的例子。左上图是将不同视角拍摄到的深度图（已转化为三维网格）匹配在一起后的情形。右上图是左图中某一部分的切面，可以看到很多层网格重叠在一起，噪声、匹配误差、采样率等都反映在这个局部中。左下图是VRIP的结果，右下图是VRIP融合的结果和深度图放在一起（棕色是深度图，蓝色是融合后的网格，只能看到一点点）。

VRIP的基本假设

VRIP最重要的假设是测量误差沿着光线传播方向（投影方向）并服从高斯分布。假设(hat{f_k}(u,v))是在第(k)个视角下，光线从传感器((u,v))位置沿着传播方向到达(f)的距离的测量值，(f_{k}(u,v))为其真实值。那么条件概率满足

(P(f_{k}(u,v)|hat{f_k}(u,v))=c_k(u,v)exp[-frac{1}{2}(frac{hat{f_k}(u,v)-f_k(u,v)}{sigma_k(u,v)})^2])

因此，VRIP算法和视角非常相关。

VRIP的模型

VRIP尝试从概率的角度来描述核心问题。对任意曲面(f)，(P(f|hat{f_1},...,hat{f_K}))为(f)的条件概率。那么求解核心问题转化为一个最大似然问题

(max_{f}P(f|hat{f_1},...,hat{f_K}))

经过一系列独立性假设，

(P(f|hat{f_1},...,hat{f_K})= prod_{k=1}^{K}prod_{i=1}^{M}prod_{j=1}^{N}P(f_{k}(i,j)|hat{f_k}(i,j)))

取对数，并转化为求和

(E(f)=-SigmaSigmaSigma log[P(f_{k}(i,j)|hat{f_k}(i,j))])

将离散情形转换为连续情形

(E(f)=-Sigmaintint_{A_k}log[P(f_{k}(u,v)|hat{f_k}(u,v))]dudv)

定义(d(u,v,f_k)=hat{f_k}(u,v)-f_k(u,v))，代入(P)的表达式，则有

(E(f)=frac{1}{2}Sigmaintint_{A_k}w_k(u,v)d_k(u,v,f_k)^2 du dv)

上式还是建立在各个视角下的局部坐标系（传感器坐标系）。将上式转换到世界坐标系下，有

(E(f)=frac{1}{2}Sigmaintint_{A}w_k(x,y,z)d_k(x,y,z)^2[mathbf{v_k}(x,y,z)cdot (frac{partial z}{partial x}) ,frac{partial z}{partial y},-1)] dx dy)

其中(z=f(x,y)) ，(mathbf{v_k})是第(k)个视角的投影方向，(d_k(x,y,z)=hat{f}(x,y,z,mathbf{v_k})-f(x,y)) ，为三维点((x,y,z))沿着第(k)个视角的投影方向与深度图对应点的测量值(hat{f})的距离，即Signed Distance Function。积分中点乘那项是从(uv)坐标系转换到(xyz)坐标系的(Jacobian)矩阵。
下图是SDF的示意图。注意靠近传感器的方向距离为正，远离传感器的方向距离为负。

下图是同一个三维点在两个视角下的SDF示意图。注意SDF是沿着视角方向（光线传播反方向）的。图中(d_1<0,d_2>0)。

于是，问题转换为

求(，使得z=f(x,y)，使得)(E(f))最小。

求解这个问题涉及到很多偏导数、方向导数的知识，本人也并没有完全搞明白。但问题的解却是出乎意料的简单

若(z=f(x,y))满足(Sigma w_k(x,y,z)d_k(x,y,z)=0) ，则(z)是最优解。

VRIP算法流程

定义(D(mathbf{x})=frac{Sigma w_k(mathbf{x})d_k(mathbf{x})}{Sigma w_k(mathbf{x})})，那么(D(mathbf{x})=0)就是我们要重构的三维表面。实际应用中，(D(mathbf{x}))可看作一个三维体数据（volume），(D(mathbf{x})=0)通过提取(D)的零等值面即可得到。因此，算法先增量构建(D)，然后通过等值面提取方法得到三维网格。例如，Marching Cube就是一种高效提取等值面的算法，而且非常适合在GPU上实现。

需要指出的是，作者在论文中采用了TSDF（Truncated signed distance function），即在一条光线上，只考虑测量值(hat{z})附近一定范围(hat{z}pm delta z)内的体素。一方面是因为一条光线可能会穿过物体不止一次。另一方面这样也可以减小搜索范围，加速算法。

VRIP算法的框架如下：

/* 初始化 */ 
对每个三维体素，设其权重为0。
/* 深度图融合 */ 
对每个深度图 {
  /* 准备 */ 
  网格化深度图；
  计算每个点的权重；
  /* 更新体素 */
  对该视角下FOV中的体素 {
    沿投影方向找到深度网格中的对应点；
    计算其沿投影方向的SDF；
    插值得到其权重；
    更新这个体素的权重和SDF。
  }
}
/* 表面提取 */ 
提取零等值面。

着重强调一点，深度图的权重需要尽量准确，特别是（1）对噪声大的点能够赋予较小的权重，（2）对于法向量和视角方向角度比较大的三维点，可以降低其权重（曲率大，采样率不够，重建误差大）。

下图是两个视角下的TSDF的融合的示意图。

下图是多个视角下TSDF的融合（两颗真实的牙齿）

VRIP的优缺点

VRIP的优点主要有

• 它是一定意义下的最大似然解。这保证了解的精度。
• 它是一个增量方法，每一次得到新的深度图后，可简单快速地加入到TSDF中。
• 适合并行处理，可用GPU加速

VRIP的主要缺点有

• VRIP生成的网格会附加一定的平滑效果，在存在噪声和匹配误差的情况下有时不能重建出细微的结构。
• 如果深度图存在匹配误差，VRIP并不能消除这些误差。这些误差会反映在融合后的网格中（分层、噪声等）。

VRIP的加速

VRIP算法的时间复杂度比较高。作者在CPU端做了很多优化和加速工作。假设深度图平行于(xy)平面，且投影方向为正交投影，那么沿(z)方向的所有体素在深度图上的投影相同，因此其TSDF可以在(z)方向上简单计算得到，而且权重相等，不必要重复运算。作者的主要思路是，将不同视角下的深度图和坐标系，通过仿射变换和重采样，映射到相对标准的位置（深度图平行于(xy)平面），从而减少投影和权重的计算量。
如下图所示，（a）某个视角下的深度图及正交投影方向，投影方向和Voxel slices有一个夹角。（b）通过一个仿射变换，将投影方向变换为与Voxel slices相垂直。一般而言，仿射变换后深度图和Voxel slices仍然有夹角，因此还要将它在平行于slices的平面上重采样。（c）将计算得到的TSDF变换回原坐标系下的距离。

VRIP在SLAM中的应用

VRIP算法在RGBD SLAM中有着广泛的应用。Kinect Fusion的作者并没有采用CPU加速的方法，而是将VRIP算法移植到GPU上，利用Ray casing算法做到了网格的实时显示，效果非常好。后来的Kintinuous、Elastic Fusion、Dynamic Fusion也都采用类似的架构去生成网格。