An automatic extraction and segmentation algorithm of underground cable point cloud

随着城市化进程的不断深化，城市架空输电线路通道资源日益紧张，地下电缆取代了传统架空线路成为了城市电网网架的主要支柱，其稳定运行是保证电网安全可靠供电的重要一环。然而地下电缆位于地下电缆隧道中，隧道内部线路复杂，各线路管线交错，在实际运行维护中，运维人员往往难以凭借传统2维计算机辅助设计(computer aided design，CAD)图纸对各线路电缆进行有效识别与管理，需要电缆隧道内电缆布设的3维分布信息辅助检修维护。同时，地下电缆隧道内电缆布设的3维信息是实现地下电缆智能管理的必要基础数据，如何快速、精确获取地下电缆布设3维信息，实现地下电缆3维重建成为目前地下电缆运维领域的研究热点。

3维激光扫描技术能够快捷、精确地获取地下电缆隧道中的3维点云，成为目前地下电缆隧道3维重建的主流技术，国内已有学者开展了基于激光点云在地下电缆3维重建中的相关研究^[1-9]。在已有地下电缆隧道点云3维重建的相关文献中，2019年，ZHAO基于电缆点云数据，采用改进标架的Sweep算法实现对电缆隧道内电缆形状走向进行模拟，进而得到地下电缆网络模型^[10]；同年，ZHAO基于点云边界点提取算法实现了地下电缆工井墙壁点云提取，并结合工井墙壁空洞提取电缆孔洞，采用一种混合型的拓扑模型结构对地下工井电缆网络进行构建^[11]。上述文献中都将提取分割后电缆点集数据作为输入数据，并未实现电缆点集从地下电缆隧道点云的自动提取与分割，电缆点仍需要人工提取。

为解决目前地下电缆隧道点云中电缆支架与电缆点分离困难、地下电缆点云需要人工提取的问题，本文作者提出一种基于区域圆柱面拟合与抗差自适应Kalman滤波的地下电缆点云自动提取分割算法。在本文算法中，首先基于电缆局部呈圆柱的形状特征，将区域点云进行圆柱面拟合，将圆柱面拟合的单位权中误差以及拟合半径作为评价标准来确定区域点云是否为电缆部分，提取出区域电缆的中轴线与半径；其次，区域电缆中轴线为起点与走向，采用抗差自适应Kalman滤波算法实现对电缆区域轴线进行延长估计，通过抗差估计以及自适应因子调整权阵的方式抑制其它类别点、噪点以及因遮挡导致缺点对电缆中轴线滤波估计的影响，得到整条电缆中心轴线实现单条电缆点集分割。实验证明，该方法有效提升了地下电缆提取的准确性和可靠性，在地线电缆点云存在遮挡、噪点情况下仍能进行有效提取分割。

地下电缆隧道内部一般由电缆隧道墙面、电缆以及电缆支架组成，其中电缆在空间分布上是一个弯曲的圆柱体，从较小局部范围上可视为标准圆柱体，其点云符合标准圆柱面分布特征，如图 1所示。

Figure 1.  Distribution characteristics of underground cables in point cloud

根据电缆局部点云符合标准圆柱面分布的形状特征，本文作者提出一种通过拟合区域点云所在圆柱面确定电缆中轴线段的区域电缆中轴点提取算法，具体理论如下所述^[12]。

设区域点云的3维点集为{P₁，P₂，…，P_m}，点总数为m。由于圆柱体仅需要中轴线上任意一点(x_l，y_l，z_l)、中轴线走向参数(a，b，c)以及圆柱体半径R这7个参数即可唯一确定，因此本文中圆柱面曲面方程为：

由于中轴线任意一点(x_l，y_l，z_l)与走向参数(a，b，c)并不独立，存在无穷组解，故增加约束条件：(1)a²+b²+c²=1；(2)z_l为区域点云z坐标均值，求解上述情况的唯一解，因此圆柱面拟合的误差方程为：

式中，v_I为区域点云中第I个点(x_I，y_I，z_I)到拟合圆柱面的距离残差，z_l，0为z坐标初始值，固定为区域点云z坐标均值，对(2)式求偏导进行球线性化，可得到系数矩阵B与L：

式中，(x_l，0，y_l，0，z_l，0)为中轴线上任意一点初始值，设置为区域点云重心；(a₀，b₀，c₀)为中轴线走向初始值，设置为地下电缆隧道走向；R₀为电缆半径初始值，设置为0.1 m。采用最小二乘法进行参数估计，参数矩阵X=(B^TPB)^－1B^TPL，残差矩阵v=BX－L，单位权中误差为$\sigma_0=\sqrt{\boldsymbol{v}^{\mathrm{T}} \boldsymbol{P} \boldsymbol{v} /(m-5)}$，其中参数矩阵X=[x_l，y_l，a，b，R]^T，权阵P为单位矩阵。

单位权中误差σ₀表示区域点云中各点到拟合圆柱面距离的期望值，可用来表征区域点云符合圆柱面分布的程度。为了准确识别电缆点云，将电缆半径范围设置为[0.02 m，0.25 m]，单位权中误差要求为小于0.01 m，因此界定区域点云是否为电缆的判定公式如下：

当区域点云圆柱面拟合结果满足(5)式时，将区域点云投影至拟合直线上，获取投影后直线上两侧端点p_A和p_B，将p_A与p_B两点连线视为电缆初始中轴线段。

地下电缆可视为由多个标准圆柱体连接而成，相邻局部区域圆柱体的中心轴线连续，走向有细微差别。根据以上特征，作者在电缆中轴线延长估计中引入Kalman滤波理论，通过圆柱面拟合确定起始电缆区域与电缆半径，将起始圆柱体中轴线走向视为电缆走向，按步长对初始圆柱面中轴线进行延迟，作为其相邻圆柱体中轴线参数预测值。以预测电缆中轴线端为中心，以2倍电缆半径为半径构建缓冲区提取点云，利用区域点云分布信息对预测中轴线进行修正得到本段电缆中轴线参数滤波值，具体原理示意如图 2所示。

Figure 2.  Schematic diagram of extension of central axis section of cable

Kalman滤波是一种根据状态向量对下一时段、下一测点进行预计估计，在得到下一时段实际观测值后对预计值进行修正的滤波方法。本文中将观测噪声与动态噪声均视为零均值白噪声，Kalman滤波的观测方程与状态方程如下^[13-15]：

式中，B_k、L_k为第k段区域圆柱面拟合的B、L矩阵(详见(3)式、(4)式)；Φ_k，k－1为状态变化向量，默认第k-1段圆柱面与第k段圆柱面的中心轴线方向与半径一致，第k段圆柱面参数(x_l，k，y_l，k，z_l，k)固定为图 2中绿色箭头的起点位置；v_k为观测值残差矩阵；u_k为预测参数残差矩阵；X_k和X_k－1分别为k和k-1段区域圆柱面拟合参数矩阵。

Kalman滤波算法本质上是一种结合预测信息与观测信息的最小二乘估计，求解v_k^TP_vkv_k+u_k^TP_uku_k最小情况下的最优解，P_vk为第k段区域点云观测权阵，默认为单位矩阵，P_uk为第k段预测参数权阵。

由于电缆隧道中电缆与电缆支架、井壁存在接触，在提取局域点云时，不可避免将其余类别点以及部分噪点纳入到滤波计算中，此时观测值权阵P_vk仍为单位矩阵，将点云中各点视为等权，将导致滤波结果受到其余类别点以及部分噪点的影响，偏离实际值。其次，由于点云采集时观测视角受限，电缆往往会被电缆支架遮挡导致区域电缆点云缺失，导致中轴线延伸提取区域点云时点云为空，此时算法将直接将预测值作为滤波值继续向下一区域进行延伸。由于每一个相邻圆柱体的中轴线方向有细微差距，连续多个区域点云为空后，滤波预测值将逐渐偏离电缆中轴线，此时仍将预测值权阵P_uk保持不变不符合实际，同时单条电缆无限延伸也并不合理。

针对以上问题，本文中采用抗差自适应Kalman进行参数求解。基于抗差估计理论，使用IGGⅢ权因子公式(见下式)进行选权迭代，计算观测值等价权阵P_vk，通过降低其余类别点以及部分噪点权重达到抑制其对滤波估计的影响。

式中，常数K₀=1.5，K₁=2.5；标准化残差 $\bar{e}_{v_{k, j}}=v_{k, j} / \sigma_{v_{k, j}}$，v_k，j为残差矩阵v_k的第j个残差值，σ_vk，j为v_k，j的中误差；ω_j为权调整因子。抗差估计是在迭代计算过程中对观测值进行验后方差评价，按(7)式计算权因子，令P_vk[j]=ω_jP_vk[j]调整权重后并重新平差，直至|X_k(i)－X_k(i－1)| < ε结束(i为迭代计算次数)，P_vk为迭代结束时调整后的权阵，ε为常数阈值，取值为0.001。

其次构造自适应因子a_k，在区域点云为空后自适应降低预测值权重，以平衡观测值与预测值之间偏差，自适应因子a_k取值如下式所示^[16-18]：

式中，$\Delta \widetilde{\boldsymbol{X}}_k=\left\|\widetilde{\boldsymbol{X}}_k-\overline{\boldsymbol{X}}_k\right\| / \sqrt{\operatorname{tr}\left(\sum \overline{\boldsymbol{X}}_k\right)}, \widetilde{\boldsymbol{X}}_k$为第k段圆柱体参数抗差值，X_k=Φ_k，k－1X_k为第k段圆柱体参数预报值，c₀、c₁为常数，分别选取1.5和8。在Kalman滤波过程中，通过对参数抗差值与预报值之间进行比较，当抗差值与预报值之间发生较大偏差时，认为电缆中轴线因多个局部区域点云缺失已偏离实际状态，通过(8)式构造自适应因子降低预报值在滤波估计中的权重。

综上所述，本文中采用抗差自适应Kalman滤波，求解v_k^TP_vkv_k+a_ku_kTP_uku_k最小情况下的最优解，进而得到状态参数的抗差自适应滤波解为：

设地下电缆隧道点云为{p₁，p₂，…，p_n}，基于区域圆柱面拟合与抗差自适应Kalman滤波的地下电缆点云自动提取分割算法具体步骤如下：(1)令i=1；(2)获取点p_i在0.5 m临域内的点集Φ_i，点集Φ_i内各点按第1.1节中所述进行圆柱面拟合，若拟合结果满足(5)式限制条件时，将点集Φ_i所在区域视为单条电缆起始区域，计算起始中轴线段两侧端点p_A、p_B，进入步骤(3)；否则令i=i+1，重复步骤(2)，直至i>n；(3)令k=1，断线计数D=0；(4)以单次中轴线延迟为例，以向量p_Ap_B为中轴线方向参数，p_B为中心线上一点，构造第k段圆柱体参数预报值X_k; 以预测中轴线为中心、0.2 m为步长、2倍电缆半径为半径构建缓冲区提取点云ψ_k; (5)若D>10，认为本条电缆已到末端，本条电缆分割结束，令i=i+1，返回步骤(2)，直至i>n；反之进入步骤(6)；(6) 若点云ψ_k中点数少于6个，以p_B为起点，按向量p_Ap_B方向延伸0.2 m得到点p_C，令p_A=p_B，p_B=p_C，D=D+1，k=k+1，返回步骤(4)；反之，进入步骤(7);(7)按(6)式构造观测方程与状态方程，按(7)式计算抗差等价权阵P_vk，用(8)式计算自适应降权因子a_k；若a_k=0，则认为预测值与抗差值差距过大，本条电缆分割结束，令i=i+1，返回步骤(2)，直至i>n；0反之进入步骤(8)；(8)按(9)式计算第k段圆柱体参数滤波值$\hat{\boldsymbol{X}}_k$，以p_B为起点，按$\hat{\boldsymbol{X}}_k$中的方向参数延迟.2 m得到点p_C，令p_A=p_B，p_B=p_C；将等价权阵P_vk中权重为1对应的点加入到单条电缆点云中，同时将电缆点中点序号大于i的点从地下电缆隧道点云移除，令k=k+1，D=0，返回步骤(4)。

详细算法流程如图 3所示。

Figure 3.  Flow chart of underground cable segmentation and extraction algorithm

为验证本文中算法实际应用效果，采用北科天绘UA-1000地面扫描仪对3段地下电缆隧道进行激光扫描，获取地下电缆隧道点云进行电缆点云分割提取实例分析，其中扫描仪测距精度为5 mm~8 mm@100 m(每100 m差距时误差在5 mm~8 mm)，测角精度为0.001°；算法运行硬件环境为英特尔i7-6820HQ处理器、32 GB内存；软件环境为Windows10 64位操作系统、VS2013与PCL1.80。为便于直观展示，将地下电缆隧道点云遮挡部分进行分割，实例点云数据如图 4所示。

Figure 4.  Example of cable tunnel original point cloud

图 5中彩色部分为电缆提取分割算法提取出的电缆点云，分割出的单条点云采用不同颜色着色。如图 5所示，本文中算法实现了从地下电缆隧道中自动提取电缆点云，并能够将布设的多个电缆分割成单条电缆，帮助电缆隧道运维人员快速获取电缆布设位置信息，对于电缆隧道3维重建与智能管控具有重要意义。

Figure 5.  Cable point cloud extraction effect of underground cable tunnel

为了进一步阐述本文中算法的优越性，借助实例数据从初始区域电缆中轴线确定以及电缆中轴线延伸两方面进行进一步阐述。

地下电缆隧道点云是由隧道井壁、电缆支架以及电缆3个类别3维点组合而成，为验证基于圆柱面拟合的初始区域电缆中轴线段提取的方法正确性，本文中从上述3个实例分类后点云中随机提取3个类别点各20个，提取各点对应的邻域点云利用圆柱面拟合计算半径与单位权中误差如表 1所示。其中符合初始区域电缆中轴线段确定条件(见(5)式)要求的邻域点云采用加粗字体标记。

从表 1中可以看出，隧道井壁区域点云由于呈平面或大弧度曲面分布，拟合计算出的圆柱体半径到达几十米以上；电缆支架区域点云由于不符合圆柱面分布特征，圆柱面拟合出的曲面无法完全贴合，单位权中误差在0.02 m以上，均无法通过初始区域电缆中轴线段确定条件检验。

表 1中电缆本体区域点云中第10号、第19号点云并未通过初始区域电缆中轴线段确定条件的检验，为了探究其原因，作者查看了其点云分布，如图 6所示。从图 6中可以看出，上述两个电缆点邻域点云均由于噪点以及其它类别点存在的影响，圆柱面拟合的单位权中误差将增大，通过初始区域电缆中轴线段确定条件中对单位权中误差的限制，避免将存在噪点以及其它类别点的邻域点云作为初始区域电缆点云，保持对初始电缆中轴线的稳健估计。

Figure 6.  Cable point neighborhood point cloud without determining conditions through the central axis of the cable in the initial area

在电缆中轴线延迟估计过程中，电缆点云不可避免存在噪点、与电缆支架相连以及因遮挡导致缺点现象，导致Kalman滤波估计发生偏差。本文作者针对上述情况，在Kalman滤波中引入抗差等加权与自适应因子，以保证在电缆中轴线延长过程中出现上述情况下仍能正确估计。为验证抗差自适应Kalman滤波在电缆中轴线延长估计中的应用效果，从3个实例数据中截取对应情况点云，并采用通用Kalman滤波进行对比分析，具体效果如图 7和图 8所示。

Figure 7.  Classification effect of cable central axis extension estimation based on normal Kalman filter in each case

Figure 8.  Classification effect of cable central axis extension estimation based on robust adaptive Kalman filter in each case

从图 7与图 8对比可以看出，通用Kalman滤波与抗差自适应Kalman滤波在电缆点云正常或缺点处无严重弯曲的情况下中轴线延迟均能够正确估计。然而通用Kalman滤波面对电缆点云含有噪点、与电缆支架相交时，滤波结果将收到非电缆点影响发生偏差，面对缺点处存在较大弯曲时，由于滤波预测值权重并未按照实际情况进行调整，滤波结果与实际电缆产生较大偏差，导致原本单条电缆分割为两条。与此相反，本文作者采用的基于抗差自适应Kalman滤波的电缆中轴线估计在面对存在噪点、与电缆支架相连以及因遮挡导致缺点现象仍能够保证正确延长估计，抗差估计能够有效抑制噪声点、其它类别点对延长估计的影响，自适应Kalman滤波能够在缺点处保证延长估计连续、正确，从而确保地下电缆点提取分割正确。

提出了一种基于区域圆柱面拟合与抗差自适应Kalman滤波的地下电缆点云自动提取分割算法，为解决地下电缆隧道中地下电缆点云的自动提取分割提供了一种新算法。

(1) 提出了基于圆柱面拟合的初始区域电缆中轴线段提取方法，基于电缆局部呈圆柱的形状特征，将区域点云进行圆柱面拟合，将圆柱面拟合的单位权中误差作为评价标准来确定区域点云是否为电缆部分，提取出区域电缆的中轴线与半径。该方法能够有效区分隧道井壁、电缆支架与地下电缆，并过滤受噪点与其它类别点影响的电缆区域，选取可靠性高的电缆区域作为初始区域。

(2) 提出了基于抗差自适应Kalman滤波的电缆中轴线延长估计算法，采用抗差自适应Kalman滤波算法实现对电缆区域轴线进行延长估计，通过抗差估计以及自适应因子调整权阵的方式抑制其它类别点、噪点以及因遮挡导致缺点对电缆中轴线滤波估计的影响，得到整条电缆中心轴线，实现单条电缆点集分割。

(3) 经实例点云测试证明，本文作者提出的地下电缆自动提取分割算法能够较好地从地下电缆隧道点云中提取出电缆点集，在地下电缆点云存在遮挡、噪点情况仍能进行有效提取分割，对地下电缆3维重建与智能管理具有重要意义。