CN112799088A - 一种基于三维激光雷达的全自动隧道湿喷状态检测方法 - Google Patents

Abstract

本发明公开了一种基于三维激光雷达的全自动隧道湿喷状态检测方法。首先利用三维激光雷达采集在建隧道的表面点云信息，然后对点云数据做三阶段预处理：利用三次投影完成隧道点云位姿矫正；结合随机抽样一致性算法和欧式聚类算法完成隧道表面点云的提取，滤除障碍物及地面点云数据；结合自适应阈值法和分段三次埃尔米特插值算法对点云数据进行补偿。最后建立隧道岩石表面湿喷状态检测模型，结合高斯‑牛顿迭代法求解模型的解。本发明可以实现对椭圆形隧道的隧道表面湿喷状态的自动检测，并计算隧道表面每一点的待湿喷深度，整个过程无需人工干预，可应用于隧道湿喷混凝土的状态检测及待湿喷深度求解。

Description

一种基于三维激光雷达的全自动隧道湿喷状态检测方法

技术领域

本发明涉及一种隧道混凝土施工领域，特别涉及一种基于三维激光雷达的全自动隧道湿喷状态检测方法。

背景技术

隧道施工的初期支护主要有干喷和湿喷混凝土技术，在湿喷混凝土技术中，往往采用工人凭经验判断是否湿喷到位，因此实现隧道湿喷自动化的前提条件是对隧道湿喷状态进行有效检测并实时感知隧道表面待湿喷深度。然而，隧道环境恶劣复杂，光照条件差，能见度低，采用传统的图像处理往往达不到预期的效果。因此需要重新选择传感器，不仅如此，由于隧道的施工进度和安全隐患问题，可供研究人员出入隧道的时间极短，且隧道内不仅布置有施工设备而且还有施工人员走动，不可避免地给数据采集造成一定的影响，往往会直接导致采集的数据出现倾斜、丢失等现象。

为解决以上问题，常用的传感器设备采用LiDar、深度相机、全站仪等。目前，对已建隧道应用激光雷达进行缺陷检测成为国内外研究的热点，但对在建隧道的研究和指导还尚少。为实现隧道湿喷的自动化，及时感知隧道表面的湿喷状态和待湿喷深度，必须对在建隧道建模分析。

发明内容

由于在建隧道是一个高度粗糙、噪声随机分布且带有不确定性的非线性系统，直接对在建隧道进行整体建模分析显然会耗费大量的资源，且难以保证隧道表面湿喷状态和待湿喷深度的实时性。本发明针对隧道表面湿喷状态的检测和待湿喷深度的计算等问题，在充分理解隧道湿喷混凝土的工作原理和隧道表面湿喷前后的映射关系的基础上，提出一种能够对在建隧道表面湿喷状态快速检测并快速计算隧道表面待湿喷深度的模型，将其运用到在建隧道中，最终达到对在建隧道湿喷状态实时检测的目的。

为了实现上述技术目的，本发明的技术方案是，

一种基于三维激光雷达的全自动隧道湿喷状态检测方法，包括以下步骤：

S1：采集隧道点云数据，通过三次投影算法对倾斜的隧道点云进行矫正；

S2：利用随机抽样一致性算法滤除隧道点云中的地面数据，再利用欧式聚类算法实现对隧道表面点云的提取；

S3：采用自适应阈值算法和分段三次埃米尔特插值模型相结合的方法对隧道表面信息丢失的部分做针对性数据补偿；

S4：由已湿喷隧道表面的点云信息，结合待湿喷隧道表面与完成湿喷后的隧道表面之间的映射关系，建立隧道表面湿喷状态检测模型，并利用高斯-牛顿迭代法求解模型。

所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，所述步骤S1包括以下步骤：

S101：首先为激光雷达建立笛卡尔空间坐标系，激光雷达在扫描隧道表面的过程中保持高度不变，采集隧道面点云在激光雷达的球坐标下的数据，然后将该球坐标下的点云数据转换到激光雷达笛卡尔空间坐标系，转换公式为

其中r_i是隧道第i点距离激光雷达的实测数据，(x_i,y_i,z_i)是隧道第i点在笛卡尔空间坐标，ω为该点与坐标原点连线的垂直角度，α表示该连线的水平旋转角度。

S102：根据激光雷达坐标系下的点云位置坐标，利用RANSAC算法取地面点云数据点集合；

S103：由步骤S102所得的地面点云数据点集合，分别投影到激光雷达的笛卡尔坐标系的三个平面：YOZ面、XOZ面、XOY面，再分别对投影区域进行边缘点检测与曲线拟合，得到三条中线，并对中线分别用一次函数做直线拟合，得到拟合方程f(y,z)、f(x,z)、f(x,y)，最后计算得到直线斜率k₁、k₂、k₃；

S104:由步骤S103所得的直线斜率k₁、k₂、k₃，利用罗德里格旋转公式对点云数据完成旋转矫正。

所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，所述步骤S102中的RANSAC算法包括以下步骤：

1)任意选取原始点云中3个数据点，判断其是否共线，共线则重新选取3个数据点，否则计算3个数据点所形成的平面P的平面方程；

2)计算原始点云中每个点到P平面的距离d；

3)设定阈值t，当d>t时被认为是无效数据，标记为无效点，反之则是有效数据，标记为局内点，标记完所有数据点后统计局内点的个数；

4)设置迭代次数k次来重复1-3步，在迭代形成的k个平面中选择局内点个数最多的平面，即为点云地面的拟合平面。

所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，所述步骤S104中的利用罗德里格旋转公式对点云数据完成旋转矫正包括以下步骤：

1)根据直线斜率得到倾斜锐角θ_i＝arctan(k_i),i＝1,2,3；

2)根据罗德里格旋转公式，即向量v绕单位向量ω旋转θ后的向量v′为：

v′＝(1-cosθ_i)(ω·v)ω+vcosθ_i+(ω×v)sinθ_i

则将隧道点云的每一点均视为一个向量v，将x轴的单位向量视为ω，并取θ₁代入后得到向量v′，然后再将向量v′、y轴的单位向量和θ₂代入并计算，再将结果、z轴的单位向量和θ₃代入，即得到隧道点云依次绕x轴、y轴、z轴的单位向量旋转θ₁、θ₂、θ₃矫正后的隧道点云。

所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，所述步骤S2中提取隧道表面点云的方法为：

S201：经过步骤S1对点云完成矫正后，再次利用RANSAC算法重新拟合地面点云数据并从隧道点云中剔除拟合的地面点云数据，得到除去地面点云数据的点集{P}；

S202：对S201中得到的点集{P}采用改进的欧式聚类算法得到隧道表面点云点集{P_s}，改进的欧式聚类算法步骤如下：

1)将点集{P}投影在激光雷达的笛卡尔空间坐标系的XOZ面，并用椭圆方程拟合得到拟合的焦点坐标和长轴长；

2)对点集{P}构建KD树，在点集{P}中任取一点P_i，并根据KD树搜索点P_i的k个近邻点{p_i,i＝1,...k}；

3)求解近邻点{p_i,i＝1,...k}到椭圆焦点的距离和D，计算该距离与椭圆长轴长的差值并与阈值t'比较，将差值小于阈值的点标记为隧道面有效点，否则标记为噪声点；

4)重复步骤1)-3)，直到所有的点均完成标记，剔除噪声点得到隧道面点云数据。

所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，所述的步骤S3中针对性数据补偿的步骤包括：

S301：根据线激光雷达上各激光收发器之间的安装间距以及自身的扫描角度范围，得到第i个激光收发器的角度范围为θ_ilow≤θ_i≤θ_ihigh；

S302：根据S301中角度扫描阈值θ_ilow、θ_ihigh，并结合激光收发器的安装位置，计算第i个激光收发器形成的点云数据的纵坐标y_i＝rtan(θ_i)，得到纵坐标的分布范围为y_ilow≤y_i≤y_ihigh，遍历隧道面点云数据点，根据每个点的纵坐标来匹配到具有相应纵坐标范围的激光收发器，分别得到每个激光收发器的扫描数据{S_i}；

S303：将S302中得到的单线扫描数据利用分段三次埃米尔特插值算法对丢失的部分进行数据补偿，得到补偿点数据{S_i′,i＝1...16}，数据补偿的过程如下：

1)对{S_i′,i＝1...16}的x坐标值做一阶差分运算，若|x_k+1-x_k|≥0.005则判定区间[x_k，x_k+1)存在数据丢失，取这两点的函数值及导数值分别为y_k、y_k+1、m_k、m_k+1；

2)则三次埃米尔特插值多项式H₃(x)满足

3)采用基函数方法，构造满足步骤2)的H₃(x)方程为：

H₃(x)＝α_k(x)y_k+α_k+1(x)y_k+1+β_k(x)m_k+β_k+1(x)m_k+1，其中α_k(x),α_k+1(x),β_k(x),β_k+1(x)是关于节点x_k及x_k+1的三次埃尔米特插值基函数，它们分别满足以下条件

4)结合步骤2和步骤3，计算H₃(x)的α_k(x),α_k+1(x),β_k(x),β_k+1(x)得到H₃(x)的表达式，并对数据缺失区间[x_k,x_k+1]以预设值为间隔进行插值，得到一系列数据补偿点{S_i′,i＝1...16}；

S304：取补偿点数据{S_i′,i＝1...16}与原点云数据{S_i,i＝1...16}的并集，得到最终数据补偿后的隧道表面点云。

所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，所述步骤S4中建立隧道湿喷状态检测模型包括以下步骤：

S401：对完成湿喷的隧道面点云几何中心点和待湿喷的隧道面点云中心点通过x方向和y方向的平移使得两个中心点处于同一竖直线，并在同一笛卡尔空间坐标系中展现；

S402：在步骤S401的坐标系下，任取隧道横截面，分别得到已完成湿喷的隧道面点云和待湿喷的隧道面点云与横截面的交线L_sprayed、L_spraying；

S403：利用高斯-牛顿迭代法求解L_sprayed的曲线方程和几何中心O_c的位置，求解步骤包括：

1)设L_sprayed的曲线方程为

其中u为待估计的回归系数，n是L_sprayed的点云点数，x_i为L_sprayed的第i个点的实际横坐标值，误差项ε_i～N(0,σ²)，表示ε为一个以0为期望值，以σ²为方差的正态分布随机变量，N表示正态分布；

2)设

为待估计的回归系数u＝(u₀,u₁,…,u_p-1)^T的初始解，其中p表示回归系数的个数，T为转置符号；将曲线方程中f(x_i,u)在g⁽⁰⁾处做泰勒展开，并略去二阶及二阶以上的偏导数项，得到

3)计算曲线在该点的函数估计值

残差项为y_i-f(x_i,g⁽⁰⁾)，残差平方和

4)令

并用矩阵形式表示误差有Y⁽⁰⁾≈D⁽⁰⁾b⁽⁰⁾+ε，修正估计因子为b⁽⁰⁾，则第二次迭代值g⁽¹⁾为g⁽¹⁾＝g⁽⁰⁾+b⁽⁰⁾；

5)根据经验设置允许误差率K，当

时，得到曲线拟合参数的解，结束迭代，否则重复步骤2-4；

S404：取L_spraying上任意一点p，根据步骤S403中L_sprayed的曲线方程，求解直线O_cp与曲线L_sprayed的交点p′，则p点待湿喷深度为|pp′|，即为该位置隧道面的待湿喷深度d。

一种电子设备，包括：

一个或多个处理器；

存储装置，用于存储一个或多个程序，

当所述一个或多个程序被所述一个或多个处理器执行，使得所述一个或多个处理器实现如前述的方法。

一种计算机可读介质，其上存储有计算机程序，所述程序被处理器执行时实现如前述的方法。

本发明的技术效果在于，创新性地关注在建隧道湿喷混凝土的过程，充分考虑隧道施工的实时性、隧道表面的粗糙性、非线性等特征，通过三次投影算法实现对隧道点云位姿的自动校正，并结合RANSAC算法和欧式聚类算法提取隧道表面点云信息，同时，为解决数据采集过程中出现的信息丢失等问题，采用自适应阈值算法和分段三次埃米尔特插值模型对隧道数据进行补偿，进而得到隧道湿喷检测模型所需要的模型数据。最后融合已完成湿喷的隧道表面点云数据和待湿喷的隧道表面点云数据，建立隧道表面湿喷状态检测模型，并结合高斯-牛顿迭代法求解隧道表面待湿喷深度，为保证隧道湿喷混凝土的施工质量提供量化指标，同时为实现隧道湿喷混凝土自动化提供数据支撑。

附图说明

图1为本发明的流程示意图。

图2为本发明实施例中采集到的隧道点云图。

图3为本发明实施例中完成数据补偿的隧道表面点云图。

图4为本发明实施例中隧道湿喷状态检测算法的模型示意图。

图5为本发明应用在在建隧道的隧道表面各点待湿喷深度分布示意图。

具体实施方式

下面结合具体实施例，并参照附图，对本发明作进一步详细说明。

本实施例的步骤流程如图1所示。这里以图2所示的一帧隧道点云图为例，本实例所述的基于三维激光雷达的全自动隧道湿喷状态检测方法包括以下步骤：

S1：如附图2所示，原始的隧道点云数据由于设备安装或路况颠簸导致数据倾斜，本实施例根据隧道点云数据完成自动位姿矫正，采用的是对地面点云数据进行三次投影，分别得到点云相对于x轴、y轴和z轴的相对偏移角度，从而对倾斜的隧道点云实现自动矫正，具体步骤如下：

S101：首先为激光雷达建立笛卡尔空间坐标系并安装，激光雷达在扫描隧道表面的过程中保持高度不变，充分采集隧道面点云在激光雷达内置三维球坐标下的数据，由于激光雷达封装的数据包为水平旋转角度和距离参量，即以球坐标的形式对数据进行封装，故接下来将该球坐标下的点云数据转换到激光雷达笛卡尔空间坐标系，转换公式为

其中r_i是隧道第i点距离激光雷达的实测数据，(x_i,y_i,z_i)是隧道第i点在笛卡尔空间坐标，ω为该点与坐标原点连线的垂直角度，α表示该连线的水平旋转角度。

S102：根据激光雷达坐标系下的点云位置坐标，利用RANSAC算法取地面点云数据点集合，RANSAC算法步骤如下：

1)若原始点云共n个数据点，任意选取其中3个点，判断其是否共线，共线则重新选取3个数据点，若不共线则计算其对应的平面P的平面方程。

2)计算原始点云中每个点到P平面的距离d。

3)设定阈值t，当d>t时被认为是无效数据，标记为无效点，反之则是有效数据，标记为局内点，标记完所有数据点后统计局内点的个数。其中阈值是根据路面起伏情况估计阈值范围，然后进行实际实验确定的，一般来说取值范围为0.4±0.05m，效果最佳，在具体实施实施时可根据具体情况调整。

4)重复1-3步，根据经验设置迭代次数k次，在k个平面中选择局内点个数最多的平面，即为点云地面的拟合平面。

S103：由步骤S102所得的地面点云数据点集合，分别投影激光雷达的笛卡尔坐标系的三个平面：YOZ面、XOZ面、XOY面，分别对投影区域进行边缘点检测与曲线拟合，得到三条中线，并对中线分别用一次函数做直线拟合，得到拟合方程f(y,z)、f(x,z)、f(x,y)，并根据计算直线斜率k₁、k₂、k₃；

S104:由步骤S103所得的直线斜率k₁、k₂、k₃，利用罗德里格旋转公式对点云数据完成旋转矫正，矫正过程如下：

1)根据直线斜率得到倾斜锐角θ_i＝arctan(k_i),i＝1,2,3；

2)根据罗德里格旋转公式，即向量v绕单位向量ω旋转θ后的向量v′通过下式计算：v′＝(1-cosθ)(ω·v)ω+vcosθ+(ω×v)sinθ，将隧道点云依次绕x轴、y轴、z轴旋转θ₁、θ₂、θ₃后得到矫正后的隧道点云。即将隧道点云的每一点均视为一个向量v，将x轴的单位向量视为ω，并取θ₁代入后得到向量v′，然后再将向量v′、y轴的单位向量和θ₂代入并计算，再将结果、z轴的单位向量和θ₃代入，即得到隧道点云依次绕x轴、y轴、z轴的单位向量旋转θ₁、θ₂、θ₃矫正后的隧道点云。

S2：然后根据隧道表面点云的提取算法，首先利用随机抽样一致性算法对隧道点云滤除地面数据，再利用改进的欧式聚类算法实现对隧道表面点云的提取，具体步骤如下：

S201：经过步骤S1对点云完成矫正后，再次利用步骤S102中RANSAC算法重新拟合地面点云数据并从隧道点云中剔除拟合的地面点云数据，得到除去地面点云数据的点集{P}。

S202：对S201中得到的点集{P}采用改进的欧式聚类算法得到隧道表面点云点集{P_s}，改进的欧式聚类算法步骤如下：

1)对点集{P}投影在激光雷达的笛卡尔空间坐标系的XOZ面，并用椭圆方程拟合得到拟合的焦点坐标和长轴长。

2)对点集{P}构建KD树，在点集{P}中任取一点P_i，并根据KD树搜索点P_i的k个近邻点{p_i,i＝1,...k}。

3)求解近邻点{p_i,i＝1,...k}到椭圆焦点的距离和D，计算该距离与椭圆长轴长的差值，根据经验设置阈值t'，将差值小于阈值的点标记为隧道面有效点，否则标记为噪声点。

4)重复步骤1-3，直到所有的点均完成标记，剔除噪声点得到隧道面点云数据，如附图3所示。

S3：由S2中提取的隧道表面点云，分析隧道表面点云是否完整，并对存在数据丢失的部分进行补偿，补偿数据如附图3中深色部分所示。

S301：线激光雷达的点云数据具有均匀分布的特点，参考线激光雷达的内置参数，第i个激光收发器的角度范围为θ_ilow≤θ_i≤θ_ihigh。本实施例中以16线激光雷达为例，即安装有16个激光收发器的激光雷达，激光雷达在垂直方向的角度范围是-15°～+15°，每个激光发射器间隔2°均匀分布。

S302：根据S301中角度扫描阈值θ_ilow、θ_ihigh，并结合激光雷达的安装位置，计算第i个激光收发器形成的点云数据的纵坐标y_i＝rtan(θ_i)，其中r是隧道某点距离雷达的实测数据，得到纵坐标的分布范围为y_ilow≤y_i≤y_ihigh，遍历隧道面点云数据点，根据该点的纵坐标进行分割，即该点的纵坐标值与哪个激光收发器的纵坐标分布范围一致，就是属于哪个收发器，分别得到每个激光收发器的扫描数据{S_i,i＝1...16}。

S303：将S302中得到的单线扫描数据利用分段三次埃米尔特插值算法对丢失的部分进行数据补偿，得到补偿点数据{S_i′,i＝1...16}，数据补偿的过程如下：

1)对{S_i′,i＝1...16}的x坐标值做一阶差分运算，若|x_k+1-x_k|≥0.005则判定区间[x_k，x_k+1)存在数据丢失。这是由于点云数据非常致密，若是点云不存在数据丢失，则差分值的平均值范围应为0.0012m～0.0027m，所以若差分值大于0.005，则可以判定该区间存在数据丢失。取这两点的函数值及导数值分别为y_k、y_k+1、m_k、m_k+1；

2)则三次埃米尔特插值多项式H₃(x)满足

3)采用基函数方法，构造满足步骤2的H₃(x)方程为H₃(x)＝α_k(x)y_k+α_k+1(x)y_k+1+β_k(x)m_k+β_k+1(x)m_k+1，其中α_k(x),α_k+1(x),β_k(x),β_k+1(x)是关于节点x_k及x_k+1的三次埃尔米特插值基函数，它们分别满足以下条件

4)结合步骤2和步骤3，计算H₃(x)的α_k(x),α_k+1(x),β_k(x),β_k+1(x)得到H₃(x)的表达式，并对数据缺失区间[x_k,x_k+1]以预设值为间隔进行插值，得到一系列数据补偿点{S_i′,i＝1...16}。

S304：取补偿点数据{S_i′,i＝1…16}与原点云数据{S_i,i＝1...16}的并集，得到最终数据补偿后的隧道表面点云，如附图3深色部分所示。

S4：由已湿喷隧道表面的点云信息，分析待湿喷隧道表面与完成湿喷后的隧道表面之间的映射关系，建立隧道表面湿喷状态检测模型，如附图4所示，并利用高斯-牛顿迭代法求解模型的解，具体步骤如下：

S401：在步骤S101中，激光雷达保持相对高度不变，因此对完成湿喷的隧道面点云几何中心点和待湿喷的隧道面点云中心点通过x方向和y方向的平移使得两个中心点处于同一竖直线，并在同一笛卡尔空间坐标系中展现；

S402：在步骤S401的坐标系下，任取隧道横截面，分别得到已完成湿喷的隧道面点云和待湿喷的隧道面点云与横截面的交线L_sprayed、L_spraying；

S403：利用高斯-牛顿迭代法求解L_sprayed的曲线方程和几何中心O_c的位置，具体求解步骤如下：

1)设L_sprayed的曲线方程为

其中u为待估计的回归系数，n是L_sprayed的点云点数，x_i为L_sprayed的第i个点的实际横坐标值，误差项ε_i～N(0,σ²)，表示ε为一个以0为期望值，以σ²为方差的正态分布随机变量，N表示正态分布。

2)设

为待估计的回归系数u＝(u₀,u₁,…,u_p-1)^T的初始解，其中p表示回归系数的个数，T为转置符号。将曲线方程中f(x_i,u)在g⁽⁰⁾处做泰勒展开，并略去二阶及二阶以上的偏导数项，得到

3)计算曲线在该点的函数估计值

残差项为y_i-f(x_i,g⁽⁰⁾)，残差平方和

4)令

并用矩阵形式表示误差有Y⁽⁰⁾≈D⁽⁰⁾b⁽⁰⁾+ε，修正估计因子为b⁽⁰⁾，则第二次迭代值g⁽¹⁾为g⁽¹⁾＝g⁽⁰⁾+b⁽⁰⁾；

5)根据经验设置允许误差率K，当

时，得到曲线拟合参数的解，结束迭代，否则重复步骤2)-4)。

S404：取L_spraying上任意一点p，根据步骤S403中L_sprayed的曲线方程，求解直线O_cp与曲线L_sprayed的交点p′，则p点待湿喷深度为|pp′|，如附图4所示。对隧道面点云求解各位置隧道面的待湿喷深度，如附图5所示。

根据本发明的实施例，本发明还提供了一种电子设备和一种计算机可读介质。

其中电子设备，包括：

一个或多个处理器；

存储装置，用于存储一个或多个程序，

当所述一个或多个程序被所述一个或多个处理器执行，使得所述一个或多个处理器实现前述的方法。

具体使用中，用户能够通过作为终端设备的电子设备并基于网络来与同样作为电子设备的服务器进行交互，实现接收或发送消息等功能。终端设备一般是设有显示装置、基于人机界面来使用的各种电子设备，包括但不限于智能手机、平板电脑、笔记本电脑和台式电脑等。其中终端设备上根据需要可安装各种具体的应用软件，包括但不限于网页浏览器软件、即时通信软件、社交平台软件、购物软件等。

服务器是用于提供各种服务的网络服务端，如对收到的从终端设备传输过来的隧道点云数据提供相应计算服务的后台服务器。以实现对接收到的隧道点云数据进行计算并得到椭圆形隧道的隧道表面湿喷状态的自动检测数据，并将最终的检测数据结果返回至终端设备。

本实施例所提供的状态检测方法一般由服务器执行，在实际运用中，在满足必要条件下，终端设备亦可直接执行状态检测。

类似的，本发明的计算机可读介质，其上存储有计算机程序，所述程序被处理器执行时实现本发明实施例的一种基于三维激光雷达的全自动隧道湿喷状态检测方法。

Claims (9)

1.一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，包括以下步骤：

S1：采集隧道点云数据，通过三次投影算法对倾斜的隧道点云进行矫正；

S2：利用随机抽样一致性算法滤除隧道点云中的地面数据，再利用欧式聚类算法实现对隧道表面点云的提取；

S3：采用自适应阈值算法和分段三次埃米尔特插值模型相结合的方法对隧道表面信息丢失的部分做针对性数据补偿；

S4：由已湿喷隧道表面的点云信息，结合待湿喷隧道表面与完成湿喷后的隧道表面之间的映射关系，建立隧道表面湿喷状态检测模型，并利用高斯-牛顿迭代法求解模型。

2.根据权利要求1所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，所述步骤S1包括以下步骤：

S101：首先为激光雷达建立笛卡尔空间坐标系，激光雷达在扫描隧道表面的过程中保持高度不变，采集隧道面点云在激光雷达的球坐标下的数据，然后将该球坐标下的点云数据转换到激光雷达笛卡尔空间坐标系，转换公式为

其中r_i是隧道第i点距离激光雷达的实测数据，(x_i,y_i,z_i)是隧道第i点在笛卡尔空间坐标，ω为该点与坐标原点连线的垂直角度，α表示该连线的水平旋转角度。

S102：根据激光雷达坐标系下的点云位置坐标，利用RANSAC算法取地面点云数据点集合；

S103：由步骤S102所得的地面点云数据点集合，分别投影到激光雷达的笛卡尔坐标系的三个平面：YOZ面、XOZ面、XOY面，再分别对投影区域进行边缘点检测与曲线拟合，得到三条中线，并对中线分别用一次函数做直线拟合，得到拟合方程f(y,z)、f(x,z)、f(x,y)，最后计算得到直线斜率k₁、k₂、k₃；

S104:由步骤S103所得的直线斜率k₁、k₂、k₃，利用罗德里格旋转公式对点云数据完成旋转矫正。

3.根据权利要求2所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，所述步骤S102中的RANSAC算法包括以下步骤：

1)任意选取原始点云中3个数据点，判断其是否共线，共线则重新选取3个数据点，否则计算3个数据点所形成的平面P的平面方程；

2)计算原始点云中每个点到P平面的距离d；

3)设定阈值t，当d>t时被认为是无效数据，标记为无效点，反之则是有效数据，标记为局内点，标记完所有数据点后统计局内点的个数；

4)设置迭代次数k次来重复1-3步，在迭代形成的k个平面中选择局内点个数最多的平面，即为点云地面的拟合平面。

4.根据权利要求2所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，所述步骤S104中的利用罗德里格旋转公式对点云数据完成旋转矫正包括以下步骤：

1)根据直线斜率得到倾斜锐角θ_i＝arctan(k_i),i＝1,2,3；

2)根据罗德里格旋转公式，即向量v绕单位向量ω旋转θ后的向量v′为：

v′＝(1-cosθ_i)(ω·v)ω+vcosθ_i+(ω×v)sinθ_i

则将隧道点云的每一点均视为一个向量v，将x轴的单位向量视为ω，并取θ₁代入后得到向量v′，然后再将向量v′、y轴的单位向量和θ₂代入并计算，再将结果、z轴的单位向量和θ₃代入，即得到隧道点云依次绕x轴、y轴、z轴的单位向量旋转θ₁、θ₂、θ₃矫正后的隧道点云。

5.根据权利要求1所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，所述步骤S2中提取隧道表面点云的方法为：

S201：经过步骤S1对点云完成矫正后，再次利用RANSAC算法重新拟合地面点云数据并从隧道点云中剔除拟合的地面点云数据，得到除去地面点云数据的点集{P}；

S202：对S201中得到的点集{P}采用改进的欧式聚类算法得到隧道表面点云点集{P_s}，改进的欧式聚类算法步骤如下：

1)将点集{P}投影在激光雷达的笛卡尔空间坐标系的XOZ面，并用椭圆方程拟合得到拟合的焦点坐标和长轴长；

2)对点集{P}构建KD树，在点集{P}中任取一点P_i，并根据KD树搜索点P_i的k个近邻点{p_i,i＝1,...k}；

3)求解近邻点{p_i,i＝1,...k}到椭圆焦点的距离和D，计算该距离与椭圆长轴长的差值并与阈值t'比较，将差值小于阈值的点标记为隧道面有效点，否则标记为噪声点；

4)重复步骤1)-3)，直到所有的点均完成标记，剔除噪声点得到隧道面点云数据。

6.根据权利要求1所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，所述的步骤S3中针对性数据补偿的步骤包括：

S301：根据线激光雷达上各激光收发器之间的安装间距以及自身的扫描角度范围，得到第i个激光收发器的角度范围为θ_ilow≤θ_i≤θ_ihigh；

S302：根据S301中角度扫描阈值θ_ilow、θ_ihigh，并结合激光收发器的安装位置，计算第i个激光收发器形成的点云数据的纵坐标y_i＝rtan(θ_i)，得到纵坐标的分布范围为y_ilow≤y_i≤y_ihigh，遍历隧道面点云数据点，根据每个点的纵坐标来匹配到具有相应纵坐标范围的激光收发器，分别得到每个激光收发器的扫描数据{S_i}；

S303：将S302中得到的单线扫描数据利用分段三次埃米尔特插值算法对丢失的部分进行数据补偿，得到补偿点数据{S_i′,i＝1...16}，数据补偿的过程如下：

1)对{S_i′,i＝1...16}的x坐标值做一阶差分运算，若|x_k+1-x_k|≥0.005则判定区间[x_k，x_k+1)存在数据丢失，取这两点的函数值及导数值分别为y_k、y_k+1、m_k、m_k+1；

2)则三次埃米尔特插值多项式H₃(x)满足

3)采用基函数方法，构造满足步骤2)的H₃(x)方程为：

H₃(x)＝α_k(x)y_k+α_k+1(x)y_k+1+β_k(x)m_k+β_k+1(x)m_k+1，其中α_k(x),α_k+1(x),β_k(x),β_k+1(x)是

关于节点x_k及x_k+1的三次埃尔米特插值基函数，它们分别满足以下条件

4)结合步骤2和步骤3，计算H₃(x)的α_k(x),α_k+1(x),β_k(x),β_k+1(x)得到H₃(x)的表达式，并对数据缺失区间[x_k,x_k+1]以预设值为间隔进行插值，得到一系列数据补偿点{S_i′,i＝1...16}；

S304：取补偿点数据{S_i′,i＝1...16}与原点云数据{S_i,i＝1...16}的并集，得到最终数据补偿后的隧道表面点云。

7.根据权利要求1所述的一种基于三维激光雷达的全自动隧道湿喷状态检测方法，其特征在于，所述步骤S4中建立隧道湿喷状态检测模型包括以下步骤：

S401：对完成湿喷的隧道面点云几何中心点和待湿喷的隧道面点云中心点通过x方向和y方向的平移使得两个中心点处于同一竖直线，并在同一笛卡尔空间坐标系中展现；

S402：在步骤S401的坐标系下，任取隧道横截面，分别得到已完成湿喷的隧道面点云和待湿喷的隧道面点云与横截面的交线L_sprayed、L_spraying；

S403：利用高斯-牛顿迭代法求解L_sprayed的曲线方程和几何中心O_c的位置，求解步骤包括：

1)设L_sprayed的曲线方程为

其中u为待估计的回归系数，n是L_sprayed的点云点数，x_i为L_sprayed的第i个点的实际横坐标值，误差项ε_i～N(0,σ²)，表示ε为一个以0为期望值，以σ²为方差的正态分布随机变量，N表示正态分布；

2)设

为待估计的回归系数u＝(u₀,u₁,…,u_p-1)^T的初始解，其中p表示回归系数的个数，T为转置符号；将曲线方程中f(x_i,u)在g⁽⁰⁾处做泰勒展开，并略去二阶及二阶以上的偏导数项，得到

3)计算曲线在该点的函数估计值

残差项为y_i-f(x_i,g⁽⁰⁾)，残差平方和

4)令

并用矩阵形式表示误差有Y⁽⁰⁾≈D⁽⁰⁾b⁽⁰⁾+ε，修正估计因子为b⁽⁰⁾，则第二次迭代值g⁽¹⁾为g⁽¹⁾＝g⁽⁰⁾+b⁽⁰⁾；

5)根据经验设置允许误差率K，当

时，得到曲线拟合参数的解，结束迭代，否则重复步骤2)-4)；

S404：取L_spraying上任意一点p，根据步骤S403中L_sprayed的曲线方程，求解直线O_cp与曲线L_sprayed的交点p′，则p点待湿喷深度为|pp′|，即为该位置隧道面的待湿喷深度d。

8.一种电子设备，其特征在于，包括：

一个或多个处理器；

存储装置，用于存储一个或多个程序，

当所述一个或多个程序被所述一个或多个处理器执行，使得所述一个或多个处理器实现如权利要求1-7中任一所述的方法。

9.一种计算机可读介质，其上存储有计算机程序，其特征在于，所述程序被处理器执行时实现如权利要求1-7中任一所述的方法。

Images