基于椭圆拟合的隧道点云数据去噪方法

李先帅，武斌

基于椭圆拟合的隧道点云数据去噪方法

李先帅，武斌

(天津城建大学 计算机与信息工程学院，天津 300384)

针对地铁隧道这一狭长密闭环境，激光雷达相较于传统隧道形变检测方法具有显著优势。由于受环境以及隧道中安装的大量设备的影响，数据中会包含大量噪声点，无法得到隧道完整内表壁数据，从而影响到数据在形变检测及三维建模等方面的应用。结合最小二乘法和拉格朗日乘数法对隧道点云模型进行椭圆拟合，并将其应用于圆形隧道点云滤波，利用各区域数据对隧道形变进行检测。对滤波结果的精度评定及其有效性分析表明，此算法针对隧道点云数据的滤波具有良好效果及稳定性。

椭圆拟合；最小二乘法；拉格朗日乘数法；点云滤波；形变检测

随着城市交通的发展，交通的安全问题成为人们日益关注的焦点问题。地铁作为城市交通的重要组成部分，对其日常的安全检测不容忽视。地铁隧道因受到地质运动以及人为活动的影响，会发生一定程度的形变，当形变超过一定限度便会造成巨大的安全事故，相较于传统的全站仪对隧道形变的检测方式，激光雷达[1−2]因具有精度高、速度快、操作简单等优势被广泛应用于精密工程测量领域。然而由于隧道内部环境复杂，有大量支架，电缆，供配电装备、信号设备、照明设备等，以及隧道管片的安装偏差，导致激光雷达检测得到的点云数据中包含大量非隧道壁点(噪声点)，对后期形变分析、隧道建模等造成干扰。一些常用的滤波算法无法满足精度要求，因此如何在保留隧道内壁点云数据的同时最大限度地剔除噪声点成为了急需解决的关键问题。DU等[3]提出一种利用高密度点云数据进行隧道变形监测方法，通过K-Nearest算法计算每个横截面的缓冲区，并通过投影方法快速构造初始横截面点集的方式将截面去噪。KANG等[4]通过RANSAC算法和曲线拟合确定地铁隧道的中心轴, 然后使用RANSAC优化提取的垂直截面，以滤除噪声。XU等[5]基于隧道结构的对称性和圆形形状，在消除大部分噪声同时避免许多低效的迭代，提高了检测效率。程效军等[6]介绍了基于中轴线拟合的方式，通过双向投影提取隧道点云中轴线并设定距离阈值进行滤波。Delaloye等[7]通过分析LiDAR扫描相关的精度问题，完成了收敛测量技术的灵敏度测试，证明了椭圆拟合的方法在隧道形变检测中的有效性。Walton等[8]通过优化椭圆拟合算法，提高了椭圆拟合算法对于形变检测的能力，并将其应用于圆形隧道及竖井的形变检测中。ZHU等[9−11]都采用了椭圆拟合的方式进行点云数据滤波，分别基于最小二乘法的不同优化进行椭圆拟合，通过计算标准差，剔除掉误差较大的点，取得了较好的滤波效果。张立朔等[12]提出一种以残差１范最小为平差准则的椭圆拟合方法方式，通过分析圆心角，将隧道曲线分割为5段圆弧进行分段去噪。本文通过移动LiDAR采集点云数据，针对椭圆隧道模型点云分布特点，提出一种结合最小二乘法和拉格朗日乘数法进行椭圆拟合的隧道点云去噪方法，本算法通过提高曲线拟合精度从而在保证点云数据完整的基础上提高去噪精度，经过滤波结果分析，证明本文算法具有较好的有效性和鲁棒性，得到的完整的隧道内壁点云滤波结果可为形变分析、断面拟合等提供高质量的点云数据。

1 隧道点云去噪方法

通过激光雷达获取隧道三维点云数据，首先提取隧道中轴线，然后沿中轴线方向等间距划分提取隧道截面，将点云坐标转化到载体坐标系中，最后利用本文算法对各区域的数据进行滤波，去除非隧道内壁噪声点，从而得到完整隧道内壁数据。

1.1 隧道中轴线提取及点云断面提取

由于提取到的点云数据量较大，进行直接拟合影响计算效率，所以滤波前必须进行区域划分，其中中轴线的走向准确的反应了隧道的姿态，所以首先依据中轴线进行区域划分。综合考虑隧道点云特点，本文采用双向投影法提取隧道中轴线，对数据分别投影至(水平)面和(竖直)面，其中轴方向为隧道走向，轴方向为竖直向上，，和构成右手坐标系，如图1所示。

沿中轴线方向进行点云断面提取，计算中轴线上任意一点的法平面，其中，Pn和Pn+1为中轴线上2点，其间隔为d，过两点的法平面ln，ln+1可实现对区域点云的分割，同理可实现对整体隧道椭圆柱面的点云划分，如图2所示。

1.2 载体坐标系坐标获取

由于隧道走向并非直线，导致点云分割后的断面随隧道走向确定方向，无法直接应用于连续断面点云滤波，需经过载体坐标系坐标转化，本文以断面切片中心点为坐标原点重新构建坐标系，平台前进方向为坐标系轴，轴竖直向上，，和构成右手坐标系，如图3所示，经坐标转化后隧道点云坐标可转化为相对于当前坐标系的坐标，用于之后点云数据滤波。

图3 载体坐标系定义

1.3 椭圆隧道柱面模型点云滤波

经过椭圆隧道模型中轴线提取、点云断面提取及坐标重建后的点云切片数据为椭圆柱面分布且与轴平行的离散点集，坐标和构成二维椭圆平面上的点，其中含有大量噪声点，需使用点云滤波算法对噪声点进行滤除，从而得到隧道内壁完整数据。

目前在椭圆拟合中应用最广泛的最小二乘 法[13]，通过计算数据点到椭圆曲线的距离，并最小化误差的平方和来寻找数据的最佳函数匹配，以下是以代数距离作为距离度量的最小二乘法椭圆拟合方法。

空间中，椭圆的一般方程为：

若要椭圆函数满足当前数据点需要满足目标函数：

其中，为当前界面数据点个数，取函数极值，令：

由此可得线性方程组：=，通过求解此方程组可得椭圆方程参数，设定拟合误差阈值，并淘汰数据点中距离大于阈值的数据点，经过多次迭代，当方程参数不再变化，则认为完成对点云数据的滤波。

在最小二乘法拟合椭圆滤波的实验过程中，数据表明滤波结果会有误差，误差主要来源于2方面，一方面来自椭圆曲线拟合过程中产生的误差，此误差相对于点云较复杂以及非完整截面较为明显；另一方面来自噪声剔除过程中标准差阈值选择时的误差。本文算法针对第一种情况进行改进，由于椭圆曲线属于圆锥曲线，而圆锥曲线根据判别式的不同，还包含了椭圆、双曲线、抛物线以及各种退化情形。从而导致迭代过程中由于数据点复杂程度不同，单纯的求取距离平方和的最小值，并不能保证每一次都能拟合出椭圆函数曲线，当迭代过程中出现一次非椭圆曲线拟合时，会导致计算的拟合误差均值并不准确，以至于出现噪声点无法滤除，而非噪声点被滤除的情况，影响算法的稳定性及滤波精度。针对此问题，本文通过添加约束条件：

确保迭代过程中的每一次拟合都为椭圆曲线，由于添加了约束条件，所以在计算目标函数(3)时，通过引入拉格朗日乘数法[14]对其进行求解，获得等式组，然后求解等式组得到最优的椭圆曲线。算法原理如下：

1) 由椭圆一般方程(2)，令：

椭圆方程可表示为：

令：

拟合椭圆的最优化问题可转化为优化目标：

其中：

2) 引入拉格朗日因子，构造拉格朗日函数：

对其求导数得：

令：

则：

通过求取特征向量可得6个，由于满足：

要使得：

则：

便可筛选出一个。

因为>0，所以计算此时的特征向量及缩放系数，可得到拟合曲线的方程解。根据拟合曲线方程可进行隧道点云滤波，具体步骤如下：

1) 由本文算法可得当前数据点云拟合曲线方程参数，，，，，。

2) 计算数据点在拟合中的误差以及平均误差：

其中为数据点编号，为当前截面中点云数据个数。

5) 确定最终椭圆曲线方程参数。

2 实验验证及分析

2.1 实验数据

为验证本文算法在椭圆柱面模型隧道中的滤波效果，采用三维激光扫描仪，选取一处地铁隧道进行隧道三维点云数据采集。选取20 m长隧道作为数据采集区域，此区域共采集859 291个数据点，点云数据中噪声点主要包含轨道、线缆、照明设施及管道附着的支架等。本文算法在MATLAB 2014a环境下进行点云滤波实验。提取的隧道点云数据如图4所示。

(a) 侧视图；(b) 正视图

2.2 中轴线提取及点云切片生成

采用本文1.1中的方法，通过计算边缘点均值来进行隧道点云中轴线提取，中轴线长20 m，切片间隔设置为0.2 m，并将切片进行坐标转换，图中展示其中连续提取的10个截面，如图5所示。

图5 隧道点云切片

2.3 滤波结果分析对比

为证明本文算法的可行性及滤波精确度，在提取的10个隧道截面中随机选取一个截面，首先剔除粗差点(残差大于0.1 m)，然后使用本文算法进行拟合中误差分别为，2，3时滤波并进行滤波效果对比，从滤波效果图中可以看出，算法在拟合中误差为时滤波效果最明显，但是会造成较多数据的缺失，造成过度滤波，在拟合中误差为2时既能保留完整隧道壁数据，同时对噪声点进行滤除，拟合中误差为3时对于一些距离较近的噪声点滤除效果不理想，仍有较多噪声点存在，由此可以证明，其噪声分布符合统计学特征，假设其标准差为高斯分布[15]，可通过本文算法随噪声进行滤除，如图6所示。

(a) σ滤除结果；(b) 2σ滤除结果；(c) 3σ滤除结果

为证明本文算法对于滤波精度的提升，选取同一截面分别使用RANSAC算法、最小二乘法、法向偏差法与本文算法进行滤波效果对比，如表1所示。

表1 不同算法滤波效果对比

法向偏差法4 811点2 085点 本文算法5 481点2 503点

表中统计结果表明，在隧道点云数据滤波中，针对断面(1)，RANSAC算法在大尺度噪声的滤除中有较好的效果，且较好地保留了非噪声点，但对于距离隧道内壁较近的小尺度噪声并不能很好的滤除，而最小二乘法相较于RANSAC算法在小尺度噪声滤除中有显著优势，本文算法相较于最小二乘法，同以2倍中误差(2)作为去噪准则，去噪精度提升约2%，在去噪效果图中可以观察出，在最小二乘法滤波过程中会出现个别噪声点无法滤除，而误将非噪声点滤除的情况。同样被应用于点云数据滤波的法向偏差法是利用点云法向量与隧道轴线方向的偏差提取出噪声点，在滤波效果对比图中可以看出，其对于噪声点的抑制效果与本文算法相近，但由于法向量对噪声的敏感特性导致算法在滤除噪声点的同时丢失了大量的隧道壁点，根据当前数据无法准确计算隧道形变量，而本文算法在抑制噪声点的同时对隧道壁点进行了较好的保留。

针对断面(2)，由于此断面为非完整断面，传统最小二乘法无法根据当前数据点拟合出椭圆曲线，而拟合出抛物线、双曲线或其他函数模型，导致噪声点无法完全去除，其滤波精度与RANSAC算法相近，只能滤除掉一些大尺度噪声，而本文算法由于添加约束条件，依然可以根据非完整断面数据点拟合出椭圆曲线，并根据当前曲线完成滤波，具有较好的稳定性，去噪精度相较于传统最小二乘法提升约18%，效果提升较为明显。在法向偏差法的滤波结果中，同样因其法向量对噪声的敏感特性，导致了非噪声点的丢失，因此，在隧道形变分析、模型重建等领域，本文算法可以提供更高质量的数据支持。

(a) 断面长半轴提取；(b) 断面短半轴提取；(c) 断面偏心率

2.4 形变分析

本次采样隧道衬砌外径为6.2 m，内径5.5 m，通过本文介绍的提取点云切片的方法，在选取的20 m隧道点云数据中，以间隔0.2 m连续提取100个切片进行隧道形变分析，经本文算法进行滤波处理后的点云切片提取轴长及偏心率，如图7所示，结果表明长、短半轴均稳定在2.75 m，偏心率均值为0.004 4，共剔除噪声点549 946个，且本区域内隧道形变量均在最大形变限度内。

3 结论

1) 根据隧道点云噪声分布特点，提出一种结合最小二乘法和拉格朗日乘数法进行椭圆拟合的滤波算法，保证了椭圆拟合的稳定性及滤波精度。

2)相较于被广泛应用于曲线拟合的其他算法，本算法在非完整断面或者数据点信息较为复杂的情况下，均可得到较为有效拟合结果，具有较好的稳定性，且对于噪声具有较好的鲁棒性。

3) 随着三维扫描技术的发展及其在地铁隧道中的应用，此滤波算法既可以应用于隧道施工质量检查、竣工验收等，也可为隧道建模，形变检测等应用提供高效的数据支持。

[1] 张小红. 机载激光雷达测量技术理论与方法[M]. 武汉:武汉大学出版社, 2007. ZHANG Xiaohong. Theory and method of airborne lidar measurement technology[M]. Wuhan: Wuhan University Press, 2007.

[2] 秦海明, 王成, 习晓环, 等. 机载激光雷达测深技术与应用研究进展[J]. 遥感技术与应用, 2016, 31(4): 617−624. QIN Haiming, WANG Cheng, XI Xiaohuan, et al. Development of airborne laser bathymetric technology and applications[J]. Remote Sensing Technology and Application, 2016, 31(4): 617−624.

[3] DU L, ZHONG R, SUN H, et al. Automatic monitoring of tunnel deformation based on high density point clouds data[J]. ISPRS-International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, 2017, XLII−2/W7: 353−360.

[4] KANG Z, TUO L, Zlatanova S. Continuously deformation monitoring of subway tunnel based on terrestrial point clouds[J]. ISPRS - International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, 2012, XXXIX−B5: 199−203.

[5] XU Xiangyang, YANG Hao, Neumann I. Time-efficient filtering method for three-dimensional point clouds data of tunnel structures[J]. Advances in Mechanical Engineering, 2018, 10(5): 168781401877315.

[6] 程效军, 贾东峰, 刘燕萍, 等. 基于中轴线的隧道点云去噪算法[J]. 同济大学学报(自然科学版), 2015, 43(8): 1239−1245. CHENG Xiaojun, JIA Dongfeng, LIU Yanping, et al. Tunnel point cloud denoising algorithm based on centerline[J]. Journal of Tongji University (Natural Science), 2015, 43(8): 1239−1245.

[7] Delaloye D, Diederichs M S, Walton G, et al. Sensitivity testing of the newly developed elliptical fitting method for the measurement of convergence in tunnels and shafts[J]. Rock Mechanics and Rock Engineering, 2015, 48(2): 651−667.

[8] Walton G, Delaloye D, Diederichs M S. Development of an elliptical fitting algorithm to improve change detection capabilities with applications for deformation monitoring in circular tunnels and shafts[J]. Tunnelling and Underground Space Technology, 2014, 43: 336−349.

[9] ZHU Ningning, JIA Yonghong, LUO Lun. Tunnel point cloud filtering method based on elliptic cylindrical model[J]. ISPRS-International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, 2016, XLI−B1: 735−740.

[10] 吕亚磊, 李永强, 董亚涵, 等. 基于椭圆约束的隧道点云滤波[J]. 测绘工程, 2019, 28(5): 44−50. LÜ Yalei, LI Yongqiang, DONG Yahan, et al. Tunnel point cloud filtering based on ellipserestriction[J]. Engineering of Surveying and Mapping, 2019, 28(5): 44−50.

[11] 黄帆, 李维涛, 侯阳飞, 等. 激光点云的隧道数据处理及形变分析[J]. 测绘科学, 2019, 44(5): 132−137. HUANG Fan, LI Weitao, HOU Yangfei, et al. Tunnel data processing and deformation analysis study based on laser point cloud[J]. Science of Surveying and Mapping, 2019, 44(5): 132−137.

[12] 张立朔, 程效军. 基于激光点云的隧道形变分析方法[J]. 中国激光, 2018, 45(4): 0404004. ZHANG Lishuo, CHENG Xiaojun. Tunnel deformation analysis based on lidar points[J]. Chinese Journal of Lasers, 2018, 45(4): 0404004.

[13] 陈若珠, 孙岳. 基于最小二乘法的椭圆拟合改进算法研究[J]. 工业仪表与自动化装置, 2017(2): 35−38, 46. CHEN Ruozhu, SUN Yue. The study of an improved randomized algorithm for detecting ellipses based on least square approach[J]. Industrial Instrumentation & Automation, 2017(2): 35−38, 46.

[14] Magaril-Il’yaev G G, Tikhomirov V M. Newton’s method, differential equations, and the Lagrangian principle for necessary extremum conditions[J]. Proceedings of the Steklov Institute of Mathematics, 2008, 262(1): 149−169.

[15] 黄晓杰, 陈宇磊, 邵跃堂, 等. 基于激光雷达的地铁隧道形变检测方法[J]. 城市轨道交通研究, 2019, 22(11): 47−51. HUANG Xiaojie, CHEN Yulei, SHAO Yuetang, et al. Research on metro tunnel deformation detection based on laser radar[J]. Urban Mass Transit, 2019, 22(11): 47−51.

Denoising method of tunnel point cloud data based on ellipse fitting

LI Xianshuai, WU Bin

(School of Computer and Information Engineering, Tianjin Chengjian University, Tianjin 300384, China)

In view of the narrow and closed environment of subway tunnels, lidar has significant advantages over traditional tunnel deformation detection methods. Due to the influence of the environment and the large number of equipment installed in the tunnel, the data will contain a lot of noise points, and the data of the complete inner wall of the tunnel cannot be obtained, thus affecting the application of the data in deformation detection and 3D modeling. This paper combined the least square method and Lagrange multiplier method to ellipse fit the tunnel point cloud model, and applied it to the circular tunnel point cloud filtering, and used the data of each region to detect the tunnel deformation. The accuracy evaluation and effectiveness analysis of the filtering results show that the algorithm has good filtering effect and stability for the tunnel point cloud data.

ellipse fitting; least square method; Lagrange multiplier method; point cloud filtering; deformation detection

U456.3+3

A

1672 − 7029(2021)03 − 0703 − 07

10.19713/j.cnki.43−1423/u.T20200403

2020−05−12

国家高技术研究发展计划(863计划)资助项目(2015BAF09B02-3)

武斌(1966−)，男，天津人，教授，从事嵌入式系统及FPGA逻辑前段研究和智能信息处理等方向的研究；E−mail：wubin@tjuci.edu.cn；

(编辑 涂鹏)