基于实测点云数据的隧道分割及形变检测方法

2021-05-18刘树亚徐剑敏鲁德宁

隧道建设(中英文) 2021年4期

刘树亚，徐剑敏，鲁德宁

(1.深圳市地铁集团有限公司，广东深圳 518040；2.南京航空航天大学，江苏南京 210016)

0 引言

隧道结构变形检测是一项长期任务，城市经济、建设的发展对地铁隧道形变检测的高效性、准确性提出了更高的要求。地铁隧道因交通设施的穿越等近接工程对周围土层产生的应力场改变，容易发生受力失衡情况，产生裂缝、形变等，从而危害其结构的安全。同时，隧道供工作人员进行检测的天窗时间极短。传统的结构检测，如采用测距仪、水准仪、全站仪等测量仪器测定待定点的变形值，不仅需要消耗大量人力、物力，而且采集的数据量有限，分析结果准确度不高[1-4]。随着隧道结构复杂程度逐渐增加以及对状态评估的可靠性要求不断提高，无论是在精度还是在效率方面都给传统的检测方法带来了巨大的挑战。

近年来，三维激光技术快速发展，并越来越广泛地应用于隧道检测领域。除隧道整体结构建模外[5]，国内外学者对三维激光扫描仪检测隧道形变也进行了大量研究[6-10]。与传统方法相比，三维激光扫描技术具有精度高、效率高、采集信息丰富及非接触性等特征，使得基于三维扫描技术的隧道形变检测方法具有极大的优势。

受地下不均匀压力的影响，隧道容易发生不规则形变。而目前的基于三维点云数据的隧道形变检测往往将隧道断面假定为圆或椭圆[11]，这种假定在隧道形变分析过程中引入了近似误差，从而降低了隧道检测的精度和准确性。因此，亟需一种准确的隧道形变检测方法对隧道进行结构安全评估。

本文提出一种新的基于实测点云数据的隧道形变检测方法。首先，基于最小二乘优化的3次B样条拟合算法提取隧道轴线，提高轴线提取准确率；其次，改进RANSAC算法进行隧道断面圆拟合，提高拟合效率并保证拟合效果；然后，基于隧道的分层结构，结合隧道点云数据几何特征及反射率信息对隧道进行分割，提取各隧道管片块；最后，通过对分割后隧道各管片块进行单独拟合，准确获取隧道断面全角度的形变数据，从而全面评估隧道结构健康安全状态。

1 隧道形变检测方法

采用三维激光扫描仪获取地铁隧道的原始扫描数据，并基于地铁隧道设计轴线，提取隧道实测轴线。基于实测轴线，对隧道进行分割，提取各个管片块，最后设计形变分析方法实现隧道形变检测。算法流程如图1所示。算法主要由3个核心部分组成。

图1 算法流程图

1)提取隧道的实测轴线。以隧道原始扫描点云数据和设计轴线为输入，通过优化设计轴线获取实测轴线。具体地，以设计轴线为初始轴线，以等间距提取一系列隧道截面，并应用圆拟合算法拟合断面圆，得到断面圆心。获取到断面圆心点集后，采用基于最小二乘优化的3次B样条拟合算法提取新的轴线，并应用迭代优化算法得到最终的隧道实测轴线。

2)断面分割及管片块提取。基于提取的实测轴线，将点云数据进行展开，生成反射率影像。进一步地，应用边界线检测并利用几何特征生成几何权重图，进行二维边界线提取，然后将二维边界线的像素坐标映射至三维点云数据空间，实现隧道点云数据分割。

3)隧道断面形变检测。通过对分割后的隧道管片块进行单独拟合，准确获取隧道断面全角度的形变数据，从而全面评估隧道结构健康安全状态。

1.1 隧道轴线提取

隧道轴线反映了隧道的整体走向与姿态，是隧道断面提取和形变检测的基础。由于隧道长期受到外力的挤压，整体结构发生形变，导致隧道实测轴线相对设计轴线发生偏移，从而对后续的隧道形变检测产生影响。从已有的隧道三维扫描点云数据中，准确提取全局实测轴线，对隧道的状态分析与评估有着重要意义[12]。

隧道轴线是一条空间连续曲线，简单地将其表示为一系列直线段的轴线拟合算法无法实现对轴线的准确提取，同时这种方法在大规模隧道处理情况下，工作量十分巨大。为此，基于设计轴线优化算法，对隧道实测轴线进行提取。

1.1.1 隧道断面去噪及中心点提取

以隧道原始扫描数据Tp和设计轴线A0为输入，通过对A0轴线进行迭代优化，准确提取隧道的实测轴线。具体地，首先基于设计轴线A0等间距地提取一系列隧道断面数据集S={s1,s2,…,si}，其中si为提取的第i个隧道断面。由于隧道点云数据中包含电缆、管线、螺栓等附件设施，对任意断面sk应用改进的RANSAC算法(random sample consensus，随机采样一致)[13]拟合断面圆，去除断面sk中的噪声和异常值。传统RANSAC圆拟合方法采用随机采样策略，在给定的数据集中随机选取3个点确定1个候选圆Cd，并基于设定的阈值带计算Cd的局内点(落入阈值带内的点)个数。通过迭代选取一系列候选圆后，选择局内点最多的候选圆Cm作为拟合结果。由于隧道点云数据规模较大，密度较高，应用传统RANSAC算法效率较低。为此，对RANSAC算法进行改进。与传统的随机采样(采样概率一致)不同，本文是基于隧道断面中每一点的密度值，重新计算每一点被采样的概率，即每个点的采样概率不同。由于密度越小的点越有可能是噪声点或离群点，因此本文算法通过降低该类点的采样概率，提高算法运行效率，同时保证拟合效果。本文定义待拟合的圆方程Fa,b,r(x,y)为

(x-a)2+(y-b)2=r2。

(1)

式中：r为改进的RANSAC拟合算法候选圆的半径；(a,b)为改进的RANSAC拟合算法候选圆的圆心。

针对任意断面数据sk，首先进行参数初始化，包括拟合次数N=0，拟合参数p=0.99，初始最大迭代次数I=∞，初始局内点集Sinlier=∅，初始候选圆得分S(C0)=0，阈值t=0.5。改进的RANSAC算法的步骤如下：

1)断面数据sk中每个点的密度计算。

2)断面数据sk中每个点的采样概率计算。sk中第i个点的采样概率

(2)

(3)

式(2)—(3)中：Dki为sk中第i个点的密度；Dmax、Dmin为sk中点的最大、最小密度值；|B|为以Dki为中心建立的邻域B中点的个数；r′为邻域B的半径。

3)候选圆选取。本文改进的RANSAC算法并非采取随机采样策略，而是基于2)中计算的采样概率进行采样，并计算对应的候选圆Cd参数(ad,yd,rd)。

4)局内点集计算。统计位于候选圆Cd阈值带内的点集Sinlier中每个点的密度，作为该候选圆的得分，即

(4)

式中：w为Sinlier内的点数；Dti为Sinlier内第i个点的密度值。

6)算法迭代停止。若拟合次数NI，则算法停止。

1.1.2 基于迭代优化的隧道轴线拟合

采用上述方法，最终得到一系列的断面圆的圆心坐标{c1,c2,…,cm}。基于该断面圆心点集，采用3次B样条的最小二乘拟合方法[14]进行轴线拟合。隧道轴线提取如图2所示。

A0绿色虚线表示设计轴线；A1红色虚线表示第1次迭代优化后的新轴线。

为了保证拟合的曲线准确表达隧道中轴线，应用迭代算法重复上述过程，即基于新轴线进行等间距提取断面、断面圆拟合及3次B样条曲线最小二乘拟合算法。迭代过程直至3次B样条曲线拟合趋于收敛，即断面圆的圆心坐标{c1,c2,…,cm}的误差小于一定阈值，即

(5)

式中ci和ci′分别表示2次迭代拟合的第i个断面的圆心坐标。

若E(c)<0.03，则认为拟合曲线趋于收敛。最终，经过迭代优化，可提取获得准确的隧道中轴线。

1.2 隧道分割

盾构隧道由多个单环管片拼装而成，而每环管片由若干个预制混凝土管片块铆接而成。地铁隧道由于长期受到不均匀地下压力的影响，容易在管片块接缝位置发生不规则形变。但管片块是整块灌注而成，发生形变的可能性较小[15-17]。在隧道形变检测中，由于地下不均匀压力引起的不规则形变，直接应用圆或椭圆拟合容易引入较大的拟合误差。只有精确定位各个管片块的位置，并分别进行拟合，才能得到精确的隧道断面的形变检测结果。为此，结合点云数据与影像信息对隧道进行分割。

1.2.1 隧道反射率影像生成

基于隧道实测轴线(1.1章节隧道轴线提取)，由于三维激光扫描技术获取的点云数据不仅包含三维坐标信息，同时也包含反射率信息，即激光雷达脉冲回波强度，拟将隧道点云数据展开至参数平面[5]，并进一步生成反射率影像。给定隧道点云数据P={p1,p2,…,pn}，对P中任意一点pi=(xi,yi,zi)检索提取的隧道轴线A1中的最近点qi，并建立Frenet局部坐标系，将pi坐标从全局坐标系转换至局部坐标系下，得到pi′=(xi′,yi′,zi′)。基于该局部坐标系，设计如下转换方程，将三维数据点pi转换至二维参数平面，其中(ui,vi)表示二维参数平面上与pi对应的点的坐标。

ui=xi。

(6)

(7)

式中：D为隧道设计直径，6 m；sign(xi′)为指示函数，xi′为正时，sign(xi′)为正，反之亦然；N为轴线A1在qi点的单位法向量。

通过上述公式，将三维点云数据转换至u、v坐标系下的二维参数平面R，通过将R进行网格划分(设置分辨率，为8 387×l，其中l取决于扫描隧道长度)，并将每个三维点的反射率信息填充至对应的像素网格中，作为该像素的像素值。当同一个网格中有多个数据点时，取它们的平均反射率进行填充。通过该算法，可获得隧道反射率影像。

1.2.2 隧道反射率影像分割

获取到隧道反射率影像后，采用Canny图像边缘检测算法[18]提取隧道反射率影像中的所有边界特征。由于管片环的边界线均为垂直线段，进一步采用霍夫变换[19]，通过设定斜率区间(试验中设定为[80,100])，可提取反射率影像中的垂直边界。根据得到的垂直边界，即可将隧道反射率影像分割为1组单环管片子影像。隧道管片环分割如图3所示。

(a)扫描点云数据 (b)隧道影像 (c)影像分割

在单环管片子影像中，由于部分边界线较短且受噪声等因素影响，存在特征不明显、检测不准确的问题，这导致无法准确全面地检测出所有边界线特征。为了突出影像中的管片块边界线，利用隧道管缝处的几何特征代替反射率特征生成影像。管缝特征法线如图4所示。在隧道断面点云数据中，由于管片块管缝处凹陷，导致断面点云数据在相应位置产生同样的变化，从而引起该处点云数据的法线异于其他位置。

(a)断面点云数据 (b)管缝特征法线

进一步地，发现点云数据法线的变化也有可能是隧道附件如管线等结构引起，而管线处点云数据和管缝处点云数据的不同在于前者位于轮廓内侧，而后者位于轮廓外侧。管缝、管线特征区分如图5所示。

红框表示断面点云数据的管线特征；绿框表示管缝特征。

基于上述发现，为隧道点云数据中每一个点pe定义几何权重g(p)，从而突出管缝特征：

g(pe)=C·H(|pc|,r)·|sin[arccos(np,pc)]|。

(8)

式中：C为权重参数，取255；H(|pc|,r)为指示函数，用于区分轮廓内外侧，若|pc|>r，则该点位于轮廓外侧，H(|pc|,r)=1，反之则位于轮廓内侧，H(|pc|,r)=-1，从而区分管缝和管线等附件特征，其中r为改进的RANSAC拟合算法的断面圆的半径；np为pe点法向量；pc为断面圆心点c到点pe的方向向量(pe-c)。

通过上述定义，可准确区分管缝特征。采用同样的点云数据展开方法[5]，生成几何权重影像图。几何权重影像生成如图6所示。可以看出：几何权重影像中管片块边界线清晰可见，识别准确率将大大提升。

(a)单环管片点云数据 (b)管缝特征检测 (c)几何权重影像

基于生成的几何权重影像，采用Canny算法及霍夫变换，可准确提取出管片块边界特征。

1.2.3 隧道点云数据分割

基于隧道影像与点云数据的对应关系[5]，将隧道反射率影像中的边界像素映射至三维点云数据空间，从而完成对隧道点云数据的分割。隧道点云数据分割如图7所示。

(a)几何权重影像 (b)边界检测结果数据 (c)点云数据分割

1.3 隧道形变分析方法

基于隧道点云数据分割结果，对隧道断面进行形变检测。首先，根据隧道实测轴线上一点(xc,yc,zc)提取隧道断面sc，其中包含n个圆弧段{sc1,sc2,…,scn}，即隧道管片块。由于隧道长期受到不均匀地下压力的影响，其断面通常为不规则形状。为区别误差较大的圆或椭圆拟合，提出一种更加准确的断面拟合和分析方法。具体地，基于提取的n个圆弧段{sc1,sc2,…，scn}，分别对其进行圆拟合，获取拟合半径{rc1,rc2,…，rcn}及圆心{cc1,cc2,…，ccn}。基于断面中心点(xc,yc,zc)，将断面sc与理论设计断面进行套合(理论半径为5.4 m)，并计算不同角度下的形变量dα。隧道形变分析如图8所示。由于该形变分析方法是依据盾构隧道设计规律对隧道进行分段分析，相比于常用的断面分析方法，该方法可更加精确地计算隧道形变，反映隧道的真实状态。

(a)隧道断面点云数据 (b)断面拟合 (c)形变计算

2 试验结果及分析

选取某市地铁1号线某段竣工隧道进行试验。试验采用Leica ScanStation P20三维激光扫描仪获取点云数据。

通过三维激光扫描仪获取隧道点云数据后，应用本文算法进行隧道点云数据分割。为了更好地验证本文算法的效果，与经典边界线检测算法LSD(line segment detection)[20]进行对比(见图9)。对比结果表明，LSD算法可以检测出较为明显的垂直边界线，但是很容易漏检一些不明显的管片块边界线，同时LSD算法将隧道附件结构如管线等也误检为边界特征，这将降低隧道分割和后续的形变检测的准确度。相比之下，本文算法可以准确地检测到所有边界线特征，且由于引入了几何权重，突出边界特征，算法不会受到管线等附件结构的干扰。分割结果表明，本文的隧道分割方法能够准确提取隧道各个管片块。