基于三维点云的工程结构平面分割及表面检测

余润泽，单伽锃，李培振，杨红梅，王玉龙

(1.同济大学 结构防灾减灾工程系，上海，200092；2.广联达科技股份有限公司，北京，100193；3.广州地铁建设管理有限公司，广东 广州，510330)

工程结构的三维可视化、几何特征提取和可见物理表面检测对评价其建造质量与服役状态具有重要工程意义。近年来，随着激光扫描技术和图像三维重建方法的发展，物理对象的三维点云模型被广泛应用于结构检测等领域[1-3]，带来了数字化和智能化的发展趋势[4-5]。研究表明，基于图像的三维重建方法能够建立具备足够精度的点云模型[6-7]，例如KHALOO 等[8]设计的方法能够生成精度为0.1 mm 的密集点云，使得基于空间点云模型的工程结构几何特征识别与物理状态检测成为可能。

在工程结构检测中，基于几何体或特征面的结构准静态变形(例如倾斜、不均匀沉降等)和基于可见物理表面的结构损伤状态(例如墙体开裂、墙面剥落)都是重要的工程参数[9-10]。目前，传统结构状态测量与检测方法包括水准测量、全站仪测量等。传统测量方法本质上是针对离散点局部位置的测量，不同测点的选择都可能影响最后的测量结果，且无法充分映射工程结构的连续性物理体本质特征。而激光扫描测量与图像三维重建能建立高保真的虚拟三维数字模型，具备获得更全面结构状态信息的条件，且容易获得比离散点局部观测更准确和鲁棒的测量结果[11]。从工程结构三维物理实体出发，三维点云模型为工程结构状态的可视化与状态分析提供了全面的数据，处理与挖掘点云模型已成为结构数字化检测的重要内容之一。其中，基于工程分析需求的点云模型分析及其关键特征提取通常需要先高效分割几何体和特征面[12-13]。

目前，工程结构的点云分割方法已取得一系列的研究进展。例如，经典的随机采样一致性(random sampling consistency，RANSAC)算法已用于钢拱桥的虚拟预拼装[14]，改进的RANSAC 算法应用于复杂建筑物屋顶的平面分割[15]，区域成长法已应用于房屋群的屋顶分割[16]。面对不同的工程检测场景与个性化需求，点云分割算法得到了进一步改进与发展，例如针对室内场景的点云平面分割改进算法，针对砌体建筑外墙的构件分割分类算法[17]，针对建筑物立面特征提取的分割与聚类方法[18]，基于特征聚类方式的道路平面分割方法[19]。而在结构检测领域，针对检测工程结构准静态变形和可见物理损伤，主要采用手动方式分割平面，例如基于空间位置手动分割结构平面用于表面裂缝检测[20]、根据感兴趣区域手动提取基坑平面用于基坑变形检测[21-22]。考虑工程结构的体积巨大，以多层建筑为例，其点云模型往往具有百万量级的点云数据量，同时要人为处理点云噪点、周边环境等干扰因素，一定程度上限制了三维点云模型数据分析与信息挖掘的精度、效率与稳定性。因此，需要从工程结构检测的适用性角度出发，发展一类结构三维点云模型的快速分割方法与多元特征分析方法。本文提出了一种基于三维点云模型的工程结构特征面分割及表面检测方法，为工程结构构件层次物理状态可视化与数字化评估建立基础。

1 点云分割方法

1.1 总体框架

本文提出的点云分割与表面检测方法主要包括3个阶段，分别是整体点云的输入、关键点云的处理以及结构的状态估计，如图1所示。

图1 基于点云的结构平面分割及状态估计方法原理Fig.1 Principle of structural plane segmentation and state estimation method based on point cloud

1) 通过图像三维重建或者激光扫描直接生成初始点云模型。其中，基于图像的三维重建方法可利用开源的重建代码或使用国内外商用软件，例如大疆制图、Pix4Dmapper和ContextCapture等，均可实现基于图像的密集三维点云重建。对于生成的点云模型，可以通过软件CloudCompare、Open3D库或者Point Cloud Library编辑库实现点云的滤波，降噪，裁剪等预处理步骤，获取满足后续分割分析要求的点云模型。

2) 为了提高计算效率，在保证点云精度的条件下，使用Python 语言设计相关算法依次进行密集点云降采样、点云分割、同类点云聚类、点云上采样等步骤，提取结构点云模型各个目标特征面。

3) 根据获得的结构墙面点云，使用最小二乘法和回归分析等数学方法推导结构墙面的数理方程，使用Python 语言设计相关算法，测量结构各个墙面的整体倾斜与表面形状。

1.2 算法设计

点云可以通过激光扫描仪直接生成，也可以通过基于图像的三维重建生成。通常，原始点云模型数据庞大，内部噪点较多，难以直接用于工程结构的状态估计，点云预处理是为后续算法处理生成高质量的点云。点云预处理方法主要为滤波、降采样、裁剪、仿射变换等。预处理后的密集点云通常有数百万个点，直接处理数百万个点会耗费大量时间。因此，本文采用体素降采样算法来进行密集点云模型的降采样，其基本原理如图2所示。体素降采样算法是以体素网格中所有点的质心来近似代表该体素中所有的点。其中，算法的关键参数是体素尺寸V，体素尺寸与降采样后的点云规模呈反比。研究表明，参数合理的降采样算法能有效维持原有点云模型的主要空间特征。点云处理之后再把降采样的体素点还原成原始点云，实现不改变点云原始密度情况下的快速点云处理。

图2 三维点云数据的体素降采样与上采样处理方法原理Fig.2 Principle of Voxel Grid Method and Voxel Upsampling Method for 3D point cloud data

体素降采样之后，本文采用RANSAC 算法分割点云模型，其原理如图3所示。RANSAC通过随机采样点来拟合模型，根据拟合阈值确定几何模型内点，分割几何模型内点与外点。RANSAC 算法是一种随机参数估计算法，适用于在大量噪声点的情况下拟合提取直线、平面等具有明确数学模型的对象。RANSAC有2个关键参数，拟合阈值T和迭代次数K。

图3 随机采样一致(RANSAC)算法原理Fig.3 Principle of random sampling consensus(RANSAC) algorithm

使用RANSAC 算法处理点云，会使同一平面内不同墙面及噪点分割在一起。因此，需要对RANSAC 处理后的点云模型再次处理。本文采用基于密度的聚类算法(density-based spatial clustering of applications with noise，DBSCAN)分割目标特征面以及消除潜在噪点。DBSCAN 算法可以有效区分三维点云的高密度区域与低密度区域，并单独分割出各个高密度区域。DBSCAN 算法需要2个控制参数：点的邻域半径和形成密集区域所需点的最少数量，分别用E和M表示。

本文设计的点云处理算法融合了降采样、RANSAC 和DBSCAN 算法，能够有效分割结构各个墙面，其原理如图4所示。首先，使用体素降采样方法处理点云，每个体素中所有点的质心作为近似采样点，最近2个点间的距离约等于体素尺寸V。之后，使用RANSAC 算法分割出结构平面点。结构平面内点排列密集且有规律，结构平面噪点排列稀疏、随机。因此，使用DBSCAN 算法很容易区分不同墙面内点和噪点。在结构平面内点(DBSCAN 算法中核心点)的邻域内至少有4 个点，而结构平面噪点(DBSCAN算法中噪点)的邻域内最多有3 个点。设置DBSCAN 算法参数E=1.2V，M=4，便可以准确地提取与聚类不同墙面的结构平面内点。

图4 融合了降采样、RANSAC和DBSCAN算法的点云处理算法原理Fig.4 Principle of point cloud data processing algorithm integrating Voxel Grid, RANSAC, and DBSCAN method

1.3 状态估计

对于已经分割提取的各结构特征平面点云集，首先可以分析所在特征平面的空间姿态，例如平面整体倾斜。若结构平面平整如图5(a)所示，可以使用最小二乘法拟合同类结构平面离散点，推导拟合平面方程，并基于水平参考平面识别所在特征平面的倾角θ。若结构平面存在凹凸情况，如图5(b)所示，无法直观分析平面空间姿态时，假定该平面的拟合方程为：

图5 结构平面倾斜与表面形状测量原理Fig.5 Principle of structural plane inclination and surface shape measurement

式中：a，b，c和d是平面方程的系数；(x，y，z)代表平面上各个点的坐标。

可设定垂直于水平面的结构数值参考平面方程为：

计算结构平面各离散点到竖直参考平面的距离，即可识别该平面的表面形状及其空间分布。针对分割点云平面的空间姿态和平面形状分析，再结合工程结构检测需求，可拓展识别该平面的相对变形。

2 工程实例

2.1 工程概况

本文以广州地铁某在建车站项目为实例，对其进行平面分割与表面测量。此项目基坑长为177.9 m，宽为23.5 m，现场试验时开挖深度约为9 m，基坑的整体与局部状态如图6(a)和6(b)所示。本研究使用大疆精灵4 RTK 无人机对基坑拍摄了1 753张实景照片。同时，使用全站仪测量19个地面控制点(ground control point，GCP)。每张图像的分辨率为5 472 像素×3 648 像素。图6(c)展示了大疆制图软件重建的基坑三维点云模型，其中使用GCP控制点云模型质量，GCP在x，y和z方向的均方根误差(RMSE)分别为6.18，10.38 和35.68 mm，整体均方根误差为37.67 mm，各个GCP 的误差如表1所示。生成的点云有44 508 258个点，文件容量约为1.41 GB，包括基坑壁、基坑的梁柱和周边的环境。为了便于分割基坑平面，即分割提取基坑壁，需要对基坑点云进行预处理。采用降噪滤波和仿射变换等方式预处理点云，裁剪掉周边环境与基坑内的部分梁和柱，以便之后的点云处理。预处理后的点云有15 589 178个点，文件容量约为0.49 GB，点云模型如图6(d)所示。

表1 点云模型三维重建中各个GCP误差Table 1 GCP errors in 3D reconstruction of point cloud model mm

本文首先手动分割了结构点云模型，把各个墙面的分类结果作为分割的真实点集，据此来分析算法自动分割的准确性。交并比(intersection over union，IoU)被用来作为分析各个平面的分割准确性的指标。第i面墙的分割交并比Ri与点云模型所有墙面的整体分割交并比R的计算公式分别为：

式中：Ai为第i面墙的真实点的集合；Bi为算法自动分割的第i面墙的点集合；n为点云模型的墙面数量。

2.2 算法处理

设置算法参数拟合阈值T=0.22 m，体素尺寸V=0.2 m，迭代次数K=50 000，E=1.2V，M=4，使用本文提出的融合降采样，RANSAC 和DBSCAN 的算法处理基坑点云，得到点云模型所有墙面的整体分割交并比R为0.852，所分割的所有墙面点总数为10 572 564 个，效果图如图7(a)所示。由于基坑壁凹凸不平，分割时基坑壁上方横梁和下方地面的部分点会和基坑壁分割到一起。此时，可根据结构平面点在各高度处的点密度对基坑壁精细分割。精细分割原理如图8所示，基坑壁面上的点云较均匀，而横梁和地面与基坑壁相接处会有较大的点云密度，发生点云密度突变。因此，根据点云沿高度方向的密度分割处理该基坑，得到最终的效果图如图7(b)所示，精细分割处理后，点云模型所有墙面的整体分割交并比R为0.889，比初始分割交并比提高了4.36%，所有墙面点总数为9 434 826 个，比初始分割的点数减少10.76%，精细分割剔除了横梁与地面上的噪点，展示了完整的墙面。

图7 点云分割算法初始效果与精细分割效果对比Fig.7 Comparison of initial effect and fine segmentation effect of point cloud segmentation algorithm

图8 基坑壁平面精细分割原理Fig.8 Principle of fine plane segmentation for foundation walls

各个基坑壁初始分割交并比与精细分割的分割交并比如表2所示。基坑壁的11个墙面中有6个墙面上存在横梁与地面上的噪点。在精细处理后，各个墙面的分割交并比都超过0.7，其中基坑壁4获得了最高的分割交并比为0.935。在精细分割处理后，基坑壁4 减少315 402 个噪点，占基坑壁初始分割点数的15.97%，其分割交并比相比初始分割交并比提升11.10%。部分基坑壁在初始分割之后没有大量噪点，所以精细处理不会对这类基坑壁起作用。在精细处理之后，其他有噪点的基坑壁分割交并比至少提升4%，点数至少减少9%。因此，精细分割能有效提高基坑壁分割的准确性。

表2 点云分割算法初始效果与精细分割效果交并比Table 2 Intersection union ratio of initial effect and fine segmentation effect of point cloud segmentation algorithm

分割出该基坑的11 个基坑壁之后，根据本文提出的方法，对分割出的每个基坑壁按照式(1)和(2)设置垂直参考平面，可计算出基坑壁表面各点到各个参考平面的距离，如图9(a)所示，其中正值代表基坑壁向基坑内部突出，负值代表基坑壁向基坑内部凹陷。从图9可以看出，各基坑壁上的点距各参考平面的距离大多在-0.2～0.2 m之间，距参考平面最远的点大多在各基坑壁的边缘转角处。为了更好地观察每一个基坑壁的具体情况，本文提取了11 个基坑壁，如图9(b)所示。单独提取了每个基坑壁的点云模型并计算了其表面形状。选择了基坑壁4 中的区域A，基坑壁7 中的区域B，基坑壁11 中的区域C 分别如图9(c)，9(d)和9(e)所示。基坑壁点云的凹陷与突出在基坑壁表面形状的计算中都得到了体现，证明了本文方法对基坑壁表面形状计算的准确性。

图9 基坑壁表面测量效果整体及局部展示Fig.9 Overall and local display of surface measurement for foundation walls

3 参数分析

本文进一步对比算法关键参数对算法自动分割效果的影响。点云模型整体的分割交并比R主要受算法关键参数T，V和K的影响。算法参数T，V和K对分割效果的影响如图10 所示，其中，控制算法参数V和K分别为0.2 m和50 000次。

图10 不同算法参数对自动分割效果的影响Fig.10 Influence of different algorithm parameters on performance of automatic segmentation

1) 当算法参数T较小时，如0.07 m，会导致欠分割，基坑壁上有较多缺失，不能完整地分割出基坑壁，如图10(a)所示，分割效果交并比R为0.701。

2) 当算法参数T较大时，如0.45 m，基坑壁不能被优先提取，以致基坑壁上下有较大区域的缺失，如图10(c)所示，分割效果交并比R为0.797。

3) 只有当算法参数T适中时，如0.22 m，能够有最优的分割交并比R为0.889，同时有最佳的分割效果，如图10(b)所示。

图11 所示为不同算法参数对自动分割交并比的影响。算法参数T代表了平面分割时的拟合阈值，由图11(a)可见，当阈值在一定的合适范围内(如0.13～0.32 m)，算法分割的交并比R均超过0.840，即根据分割效果或者经验选择合适的算法参数T能够获得较好的分割效果。算法参数V反映了点云处理时模型的点密度，较大的密度意味着更多的信息和更多的计算时间。由图11(b)可见：当点云密度足够大(如V小于0.35 m)时，算法能够获得较好的分割交并比，R均大于0.880，更大密度的点云即更小的V不能明显提高分割的精确性。由图11(c)可见：当迭代次数够大(如K大于10 000次)时，算法能够获得较大的分割交并比，R均大于0.875，更多的迭代次数对自动分割的交并比R提升有限。

图11 不同算法参数对自动分割交并比的影响Fig.11 Influence of different algorithm parameters on intersection union ratio of automatic segmentation

4 结论

1) 以工程结构作为研究对象，引入机器学习方法，提出一种基于三维点云的结构平面快速分割及表面检测方法，提升工程结构点云研究中平面分割的分析效率，实现工程结构点云模型的快速处理。

2) 针对某在建地铁车站深基坑实例，实现了无人机图像驱动的基坑整体点云高精度重建、点云模型基坑壁分割提取、基坑壁表面形状识别和平面分割精确度分析。

3) 设计的算法效果主要受其参数影响。其中，算法参数T的影响最大，当T为0.22 m时有最佳的分割效果及最大的交并比，为0.889；更小的算法参数V和更大的算法参数K有利于算法的效果提升，但会消耗更多的计算时间。