顾及点云曲率的快速点云表面模型重建算法

A Fast Surface Reconstruction Algorithm Considering the Curvature of Point Cloud

摘要: 针对现行点云表面模型重建算法的效率问题，提出一种顾及点云曲率的快速点云表面模型重建算法。首先将点云Delaunay三角化，然后利用点云曲率删减可视信息，用剩余可视信息构建图割问题后，求解图割问题得到点云表面模型。实验结果表明，本文算法能得到完整度高，细节丰富的表面模型，重建速度快。

Abstract: We describe a fast surface reconstruction algorithm considering the curvature of point cloud from a set of merged range scans. Our key contribution is improving the efficiency of the algorithm by deleting part of visual information. First, Delaunay edges are added to the point cloud to construct Delaunay structure. Then, part of visual information is deleted base of curvature of point cloud, and a graph-cuts problem is established based on the remaining visual information. Finally, a surface model is obtained by solving the graph-cuts problem. We tested our method on several publicly available sets of range scans. The experimental results show that the method can efficiently reconstruct high-quality surface model with rich details and high integrity.

1. 引言

逼真的三维模型数据正逐渐成为城市管理和服务的重要数据源 [1]。 三维重建的方法主要包括手工建模 [2] 、基于三维激光扫描数据建模 [3] 和基于影像建模 [4]。 近年来，基于影像的三维重建由于其高效便捷、成本低廉等优势成为研究热点。而利用影像多视匹配获取的密集点云进行点云表面模型重建是三维重建的核心问题之一 [5]。 点云表面模型重建主要面临的问题有弱纹理区域的重建、误匹配的影响、噪声的影响、数据冗余问题、精度问题、效率问题等。

国内外已有大量学者进行了点云表面模型重建的研究。Kazhdan等提出了基于隐式曲面方程的泊松点云表面模型重建算法，采用隐函数拟合三维点云，得到水密表面模型，然而算法易受噪声影响，细节容易丢失 [6]。 Rocchini提出了Marching Cubes算法，在每一个体素中抽取相应的等值面，但是该方法对照相机位置有严格要求，且体素大小直接影响三维重建的效果 [7]。 Calakli等提出了Smooth Signed Distance表面模型重建，用变分方法得到表面模型，算法也易受噪声干扰 [8]。 Labatut提出了利用图割问题构建表面模型的算法，鲁棒性强，不易受噪声干扰，但算法时间代价高 [9]。 Michal Jancosek提出了一种用可视信息进行增强的表面模型重建算法，对弱纹理区域重建的准确性高 [10]。 这些方法各有优点与特色，针对算法鲁棒性、弱纹理区域重建、消除噪声、数据冗余问题、精度问题、效率问题等各个方面做出了研究与改进，但每种算法都不能面面俱到，也不能完全解决问题，尤其是效率问题。因此，提出一种快速点云表面模型重建算法迫在眉睫。

本文提出一种顾及点云曲率的快速点云表面模型重建算法，利用图割算法进行点云表面模型重建。首先对点云进行三维Delaunay三角剖分，将空间分割为多个四面体 [11]。 然后利用点云曲率删减可视信息，用剩余可视信息对Delaunay结构中的四面体加权，构建图割问题。最后求解图割问题，将四面体分为S和T两个子集，提取分割边界得到点云表面。实验结果证明，本文算法效率高，重建模型细节丰富，完整度高。

2. 快速点云表面模型重建算法

本文算法分为Delaunay三角化、基于点云曲率删减可视信息、利用可视信息加权构建图割问题、图割求解四步。

2.1. Delaunay三角化

给三维点云添加互不相交的边，将点云所在空间分割为多个四面体，这个过程叫做三角剖分。如果点云的一个三角剖分只包含Delaunay边，该三角剖分叫做Delaunay三角剖分。Delaunay三角剖分具有空球属性，即任意四面体的外接球中仅包含点集中的四个点。任意点集的Delaunay三角剖分结构是唯一的。Delaunay结构分割出的四面体即待分割的目标群体。四面体是图割算法网格图中的顶点，四面体之间的邻近关系形成网格图的边。

2.2. 基于点云曲率删减可视信息

基于影像的三维重建算法中，可视信息指每个三维点由哪些摄像机产生 [12] ，三维点与相应摄像机的摄影中心连线构成一段光线。用于三维重建的影像必须有一定重叠率，重叠率不满足最低要求会产生航摄漏洞，影响重建效果。对重叠率的要求导致每个三维点由多个摄像机产生，每个三维点都有多条光线。据统计，基于无人机影像的三维重建，光线数量是三维点数量的20~40倍。光线数量庞大导致点云表面模型重建算法时间代价高，也造成数据冗余。

因此，本文通过点云曲率对光线进行删减，用剩余光线进行可视信息加权部分。光线删减方法详见第三章。

2.3. 利用可视信息加权

利用可视信息构建图割问题，需要依次遍历所有光线，根据每条光线与Delaunay三角剖分结构的相交结果对四面体加权。单条光线的加权过程如图1所示。加权公式为公式(1)。

$w=\alpha \times \left(1.0-e^{\frac{-d^2}{2\times \sigma }}\right)$ (1)

公式(1)中α为光线权值，σ为平滑系数。

Figure 1. Weighting process of single ray

图1. 单条光线加权流程

2.4. 图割求解

可视信息构建出一个有向图G(V,E)， $V=\left\{v1,v2,v3,\cdots ,vn\right\}$ 为顶点集合，E为有向边集合。额外增加两个特殊顶点，源点s和终点t。运用最大流最小割算法 [13] 将顶点集合分割为两个子集S和T。S集合中的顶点与源点s相连，T集合中的顶点与终点t相连。 $S\cap T=\varnothing$， $S\cup T=V$。割集的能量方程为

$C\left(S,T\right)={\displaystyle \underset{Vp\in S\backslash \left\{s\right\}\wedge Vq\in T\backslash \left\{t\right\}}{\sum }{w}_{pg}}+{\displaystyle \underset{Vp\in S\backslash \left\{s\right\}}{\sum }{w}_{pt}}+{\displaystyle \underset{Vp\in T\backslash \left\{t\right\}}{\sum }{w}_{ps}}$ (2)

公式(2)中第一项为顶点之间的连接边，衡量顶点的不连续性。第二项与第三项衡量顶点连接特殊顶点须耗费的代价值。使得能量方程最小的割集为图割结果。分割边界对应的三角面是模型表面。

3. 基于点云曲率的光线删减

3.1. 单点光线删减

图割构网算法最终将点云所在空间划分为“内”和“外”两类 [14]。 若摄影中心C拍摄到点P，则给线段CP所在空间提供“外”空间支持，使其更倾向于被划分为“外”。如图2(a)所示，球形物体、墙壁、地面所在空间是“内”空间，即存在实体的空间。其余空间(被光线穿过的空间)是不存在实体的“外”空间。点云所在空间被分割为多个四面体，四面体最终被分类为S和T两类，对应两种空间。点云表面是内空间和外空间的分界，即两类四面体的分割边界。

Figure 2. Full space and free space

图2. “内”、“外”空间示意图

图2(a)、图2(b)中均有一个球形物体、一面墙壁和地面。红、绿、蓝三种颜色分别标注出A、B、C三个传感器对应的光线。图2(a)均匀删减光线后得到图2(b)，删减光线后的“内”、“外”空间未改变，则空间分割边界同样未改变。另外，论文 [9] 证明图割构网算法具有鲁棒性，在光线欠采样、分布不均、没有准确传感器位置的情况下，依然可以正确重构出三维网。因此，均匀删减部分光线后，依然可以重建正确的点云表面。

每个三维点有多条光线。删减光线时，遍历所有三维点，每个点随机删除部分光线，以保证光线删减均匀。

3.2. 基于点云曲率的光线删减

过曲面上某个点具有无穷个正交曲率，其中存在一条曲线使得该曲线的曲率为极大，这个曲率为极大值Kmax，垂直于极大曲率面的曲率为极小值Kmin。这两个曲率属性为主曲率，代表法曲率的极值。Kmax和Kmin的平均值称为曲面的平均曲率，刻画了曲面的弯曲程度。

根据平均曲率进行点云排序，利用排序结果将点云分类为特征点、中间点和平面点。特征点曲率最大，中间点次之，平面点曲率最小。图3将特征点显示为红色，中间点显示为蓝色、平面点显示为白色。图3(a)，图3(b)，图3(c)，图3(d)中的点云分别为Stanford 3D扫描库中的blade，Lund University提供的Alcatraz courtyard、House、Uoft。图3表明，平均曲率较大的点勾勒出模型大致轮廓，可以较多展现模型细节。平均曲率较小的点一般分布在平面上。平均曲率可以代表点的重要程度，而且平均曲率计算简单，时间代价小。因此本算法根据点的平均曲率决定单点光线的删减比例：特征点光线全部保留，中间点和平面点删减部分光线，平面点的光线删减比例高于中间点。

Figure 3. Mean curvature of point cloud

图3. 点云平均曲率显示

4. 实验结果分析

4.1. 数据准备

以VS2015为平台，采用C++语言实现本算法，用CGAL开源库 [15] 实现Delaunay结构构建，PCL开源库 [16] 实现点云曲率计算。计算机硬件环境为Intel Core i3 CPU 2.4 GHz。用Colmap开源代码对空三解算后的影像进行多视密集匹配，生成实验需要的三维密集点云。

4.2. 实验结果

用Lund University, Carl Olsson提供的House、Alcatraz courtyard、Uoft公开数据集进行实验。House数据集包含12张影像，AC数据集包含133张影像，Uoft数据集包含77张影像。分别用本文方法、OpenMVS表面重建方法和Poisson表面重建方法 [6] 对密集匹配得到的点云进行点云表面模型重建。本文算法得到的重建结果如图4所示。

Figure 4. Surface model of the point cloud

图4. 点云表面模型

图4中的模型分别为Alcatraz courtyard、house、Uoft。对比图3中相应点云，可以看出，本文算法重建出的点云表面模型完整还原了真实地物的几何结构信息，重建结果符合真实状况，较好的呈现出阶梯、门窗等表面细节。

在模型局部细节处对三种重建算法进行对比，对比结果如图5和图6所示。

Figure 5. Comparison of vegetation reconstruction results

图5. 植被重建结果对比

Figure 6. Comparison of building reconstruction results

图6. 建筑物重建结果对比

图5由上至下依次为本文方法、OpenMVS方法、Poisson方法和原图，结果表明，本文算法重建出的点云表面模型较为完整的还原了植物的表面轮廓，OpenMVS方法重建出的植被在红圈处有一个空洞，Poisson方法未能正确重建出植物表面，重建出的表面有多处缺失。图6由左至右依次为本文方法、OpenMVS方法、Poisson方法和原图。本文算法和OpenMVS方法都正确重建出建筑物的窗户轮廓，Poisson重建的窗户轮廓存在扭曲变形。

Table 1. System resulting data of standard experiment

表1. 标准试验系统结果数据

统计三种方法的运行时长，得到表1。表1结果表明，本文方法建模效率高，在三种方法中用时最短。

实验表明，Poisson表面重建算法重建出的模型具有水密性的封闭特征、良好的几何表面特性和细节特性，算法速度较快。但是Poisson表面重建算法抗噪性较差，而且处理非密闭点云时会产生多余面片。OpenMVS方法是一种基于Delaunay三角化的构网方法，鲁棒性强，但曲面重建耗时较长，不适合处理海量点云。本文提出的方法鲁棒性强，建模效率高，且能获得较为完整、符合真实状况的模型，在高效率的同时兼顾模型精度。

5. 结论与展望

本文提出了一种顾及点云曲率的快速点云表面模型重建算法，首先将点云Delaunay三角化，然后利用点云曲率将点云分类，根据类别删减可视信息，利用剩余可视信息构建图割问题。最后，求解图割问题得到点云表面。实验结果表明，本文方法重建效率高，重建结果符合真实状况，重建出的地物细节丰富，完整度高。

本文在已有研究基础上，仍需在以下方面进一步深入研究：1) 如何修正模型表面的拓扑连接错误。2) 如何对复杂、质量缺陷的点云进行表面模型重建。

致谢

本文工作在郭丙轩教授的指导和帮助下完成。感谢Lund University提供的数据集。

基金项目

国家自然科学基金重大研究计划重点支持项目(91638203)；

国家重点研发计划(2016YFB0502200，2018YFB0505401)；

国土资源部城市土地资源监测与仿真重点实验室开放基金资助课题(KF-2018-03-052)。

NOTES

^*通讯作者。