本 文 内 容 来 源 于《测绘学报》2026年第5期（审图号 GS京(2026)1499 号）
高精度激光点云配准驱动的毫米级地表形变检测方法
韦朋成1,2, 蒋贵宇1,2, 沈航毅1,2, 黄海峰1,2, 张溶玲3
1.防灾减灾湖北省重点实验室（三峡大学），湖北 宜昌 443002
2.三峡大学土木与建筑学院，湖北 宜昌 443002
3.香港理工大学土地测量及地理资讯学系，香港 999077
摘要
针对地质灾害监测中激光点云配准特征误匹配率高、稳健性差，以及海量离散点云中微小形变难以与植被干扰、环境噪声有效区分的难题，本文提出一种由粗到精的高精度激光点云配准驱动的毫米级地表形变检测方法。首先，利用图论构建可靠性度量剔除特征匹配粗差，并设计增强的GNC-Welsch稳健估计器，解决高粗差率下稳健估计目标函数的非凸优化难题；进而基于微观结构混合特征因子驱动实现高精度配准。在此基础上，提出精确提取与多维验证的形变检测策略：利用特征分析滤除植被干扰并提取形变候选区域，构建融合几何形态、统计分布及物理机理特征的多维验证框架，以剔除虚假形变。模拟试验表明，本文方法在不同形变量级下的均方根误差（RMSE）稳定在0.52~0.61 mm，在5 mm微小形变量下F1值达86.11%，8 mm形变量下F1值达95.39%，验证了该方法在毫米级形变检测中的有效性；真实边坡场景试验证实，本文方法能有效拒绝伪形变聚类，精准识别出15 mm的地表沉降与22.6 mm的地表隆起，验证了其在复杂野外环境下的稳健性与实用性。
关键词
点云配准; 稳健估计; 地表形变检测
基金项目
湖北省自然科学基金(2025ABF104); 防灾减灾湖北省重点实验室（三峡大学）开放基金(2025KJZ06)
作者简介
韦朋成（1993—），男，博士，讲师，研究方向为时空数据处理与三维重建。 E-mail： weipengcheng@ctgu.edu.cn
通信作者:张溶玲 E-mail：rongling.zhang@connect.polyu.hk
本文引用格式
韦朋成, 蒋贵宇, 沈航毅, 黄海峰, 张溶玲. 高精度激光点云配准驱动的毫米级地表形变检测方法[J]. 测绘学报, 2026, 55(5): 866-880 doi:10.11947/j.AGCS.2026.20250549
WEI Pengcheng, JIANG Guiyu, SHEN Hangyi, HUANG Haifeng, ZHANG Rongling. Millimeter-level surface deformation detection method based on high-precision laser point cloud registration[J]. Acta Geodaetica et Cartographica Sinica, 2026, 55(5): 866-880 doi:10.11947/j.AGCS.2026.20250549
阅读全文
http://xb.chinasmp.com/article/2026/1001-1595/1001-1595-2026-05-0866.shtml
地质灾害（如滑坡、崩塌等）具有隐蔽性强、突发性高和破坏性大的特点，对人民生命财产安全构成严重威胁[1]。实现对地质灾害隐患点的早期识别与全过程监测，是防灾减灾的关键环节。传统的接触式监测手段（如GNSS、全站仪、裂缝计等）虽然精度较高，但存在单点监测、布设困难、难以反映整体形变特征等局限性[1-2]。近年来，地面三维激光扫描（terrestrial laser scanning，TLS）技术凭借其非接触、高精度、高密度及全天候作业的优势，已成为地表形变监测的重要技术手段[3-6]。
文献[7]系统总结了TLS技术在工程变形监测中的应用进展，指出点云配准精度是制约形变检测能力的关键因素。文献[8]提出了基于加权最小二乘正则化的滑坡形变提取方法，验证了TLS技术在毫米级形变检测中的可行性。然而，在实际的大范围复杂野外场景监测中，基于TLS的微小形变检测仍面临诸多挑战。首先，多期点云的高精度配准是形变分析的前提[9-10]；靶标辅助配准精度高，但球形靶标的数量与布设位置对配准精度影响显著[11]，而野外地质灾害现场地形险峻、可达性差，靶标规范化布设往往难以保障。目前主流的点云自动配准方法包括粗配准和精配准两个阶段[12]。在粗配准阶段，基于特征的方法通过提取和匹配局部几何特征实现初始对齐。快速点特征直方图（fast point feature histogram，FPFH）因其计算效率高、特征区分度强而被广泛采用[13-14]。然而，由于地质场景的复杂性，特征匹配过程中往往存在大量误匹配（粗差），严重影响配准精度[15]。随机采样一致性（random sample consensus，RANSAC）算法虽然是处理粗差的经典方法[16]，但在高粗差率场景下计算效率较低且难以保证全局最优。针对高粗差率下的稳健估计问题，文献[17]提出了基于渐进非凸性（graduated non-convexity，GNC）的稳健估计框架，通过结合Black-Rangarajan对偶性将非线性损失函数转换为加权残差优化问题，在70%~80%粗差率下仍能获得准确的配准结果，为解决复杂场景下的点云配准问题提供了重要参考。在精配准阶段，迭代最近点（iterative closest point，ICP）算法及其变体是最常用的方法[18]。然而，传统ICP算法对初始位姿敏感，且在处理大规模点云时计算效率较低。基于法向量的点到面ICP方法因其收敛速度快、精度高而受到关注[19]，但在复杂地形条件下平面假设难以满足。尽管近年来深度学习配准方法发展迅速，但在跨时相、地貌变化大的场景下，其泛化能力仍有待验证[20]。此外，为了处理大规模的点云，通常需要对点云进行降采样处理，该步骤会造成精度损失，因此如何结合采样体素中的点云微观结构特征以进一步提升精配准精度，仍需深入研究。
在点云形变检测方面，多尺度模型到模型点云比较（multiscale model-to-model cloud comparison，M3C2）和点云到点云比较（cloud-to-cloud，C2C）等距离计算算法是目前应用最广泛的方法[21]。文献[22]提出M3C2-EP方法，通过误差传播实现了更准确的最小可检测变化量估计。文献[23]发展了基于对应关系驱动的平面M3C2方法，显著降低了形变量化的不确定性。文献[24]提出Patch-based M3C2方法，通过局部平面片段投影实现了更低不确定性和更高分辨率的形变分析。然而，现有M3C2类方法主要提供逐点形变量估计及其不确定性，缺乏对形变区域真实性的系统判别机制[25]：无法有效区分植被生长变化、系统噪声与真实土体形变，在植被覆盖的野外场景中容易产生大量伪检测。此外，M3C2依赖局部法向量方向进行点到点投影，当法向量估计不稳定时（如地形突变区域），其不确定性显著增大。
针对上述问题，本文提出一种高精度激光点云配准驱动的毫米级地表形变检测方法。本文的主要贡献如下：①构建了基于图节点与边双层可靠性的特征优选策略与增强的GNC-Welsch稳健估计器，有效抑制了复杂场景下的高粗差率问题；设计了同时优化平面参数因子与变换参数因子的微观结构驱动精配准方法，在保留体素内原始点云的前提下将配准精度提升至毫米级。②提出了一套融合几何形态、统计分布与物理机理的多维度形变验证框架，通过多类约束并行判别实现形变区域完整识别与真伪判别，弥补了现有点云形变检测方法普遍缺乏形变真实性系统判别的局限。
1 方法 1.1 高精度激光点云配准
将配准问题抽象为数学优化问题：给定待配准的两组点集
和
，假设通过特征点匹配方法得到N组同名点 且其精确对应，则可通过以下模型估计刚体变换参数，使同名点对满足pi=Rqi+t，i∈{1，2，…，N}
式中， 表示二阶范数。式（1）是一个经典的最小二乘目标函数[26]，旨在求解最优刚体变换参数 ）。
1.1.1 基于可靠性度量的特征匹配点粗差剔除
地质灾害场景（如滑坡、崩塌堆积体、岩质边坡等）局部几何特征丰富，通过特征点匹配能够提取大量候选同名点。然而，受地质形变、物体位移、数据缺失及描述算子精确性等因素的影响，候选同名点集中往往混有大量误匹配（粗差/外点）；在平坦地表、均质岩面等局部几何区分度较低的区域尤为显著。因此，需剔除几何特征不一致的点以降低粗差率。
。由于点云中点的相对位置在刚体变换下保持不变，在图 中对于由边 连接的节点 和 ，如果在图 中 对应的正确节点是 和 ，那么两个对应边的欧氏距离在无噪声情况下应当相等。假设内点的噪声容忍值是ϵ，得到如下约束
邻接矩阵A的构建方式如下：若对于节点对 和 满足几何一致性约束，则计算Aw（i，j）=wn（Δdij），其中，wn是指数函数；如果是外点，则值设置为0。通过计算每个节点的权重，得到节点的可靠性 。然后根据可靠性从高到低对节点进行排序，选择前Kopt个节点对组成优选同名点集 。
（2）基于匹配图边可靠性的优选方法。核心思想是把每条候选边作为“旋转轴”，用其他节点投票，通过这个旋转对齐的节点越多，说明这条边越可靠。本文通过设计约束函数 来评估匹配图中同名边的可靠性，并引入图论中的边-节点亲和向量 进行评价，约束函数 定义如下
式中， 表示图中任意两条边对齐后的同名边对 ； 表示以 为旋转轴、旋转角度θ后将点云对齐的旋转矩阵。除节点i、j外，其他同名点对 均可以计算一个 值，若满足 ，则将矩阵B中对应位置设置为1，否则记为0，则边的可靠性度量表示为
根据式（4）计算所有边的可靠度，选择最大可靠度的边所对应的同名点（B矩阵中值为1的元素）集合作为更新后的优选集 。
1.1.2 基于增强GNC-Welsch函数的粗配准参数稳健估计
上述过程获取的优选集 中仍可能存在误匹配。相比于直接使用奇异值分解（singular value decomposition，SVD）计算式（1）中的待求参数，稳健估计器提供了更为精确的结果，能够从包含大量误匹配的同名点集中估计出准确的刚性变换参数。与其他稳健估计器的损失函数如Huber、Geman-McClure和Cauchy相比，Welsch函数更适用于从包含大量极端粗差值的同名点对中估计配准参数，如图1所示。Welsch函数能够快速对大粗差值进行平滑处理，同时在小粗差范围内保持敏感度，这使得它在对噪声和粗差较多的同名点集的稳健性优化中表现尤为优异。
图1
图1Welsch函数与其他损失函数性能比较
Fig. 1Performances comparison of Welsch function with other loss functions
然而，将该稳健函数应用到点云配准的目标函数式（1）并进行线性求解是具有挑战的，因为这涉及严格的数学推导和证明。本文结合Welsch函数和GNC框架设计了一种增强的GNC-Welsch估计器。该估计器从凸优化逐渐过渡到非凸优化，避免了局部极小值，降低了对初始值的敏感性。本文基于式（1）设计的目标函数如下
(5)
式中，ρ（r）=1－exp（－r2/2σ2）为Welsch函数， 代表残差值。参数σ控制函数的形状（图2），其初始值由 中对应点的平均距离确定。为避免陷入局部极小，本文借鉴GNC框架的退火思想引入退火因子μ，按σ←σ/μ逐次衰减σ，使损失函数从近似凸过渡到强非凸。然而，由于式（5）的非凸性，直接求解仍然困难。本文参考Black-Rangarajan对偶性[27]将其转化为加权残差形式，引入缩放变量ω=r2/（2σ2），则原始的损失函数为 。由于ϕ（ω）为关于ω的凹函数，根据Fenchel对偶性，存在权重变量z∈（0，1]与惩罚函数Ψ（z），使得
(6)
图2不同σ值下的Welsch函数值
Fig. 2Welsch function values under different σvalues
最优值在z*=ϕ′（ω）=exp（－ω）处得到。将ω=r2/（2σ2）代入式（6），并由ρ（r）=ϕ（ω）的等价关系，即得ρ（r）的对偶形式
(7)
将刚体变换{R，t}、点对pi和qi转换为齐次坐标，表示为T、 和 ，残差记作 。将式（7）代入式（5），原非凸优化问题转化为关于T与权重{zi}的联合优化
(8)
由于目标函数关于{zi}可分离，且对T与{zi}分别条件最优时形式简单，采用交替最小化求解。在第k次迭代中，首先固定上一轮的权重 ，此时与T无关的项为常数，问题退化为加权最小二乘
(9)
将式（9）视为李代数 上的优化问题，采用Gauss-Newton法配合左扰动模型迭代求解，得到当前轮次的变换估计 。然后固定 ，目标函数关于各zi独立可分，求导并令其为零得到闭式解
(10)
该解自动满足约束z∈（0，1]。残差越大的点对应的权重越小，从而自适应抑制粗差。当相邻两次迭代满足 时算法收敛，输出粗配准变换矩阵参数 。
1.1.3 微观结构混合特征因子驱动的精配准
面向应用的点云配准方法常用体素采样降低数据量，但会带来精度损失。与之不同，本文在全局尺度上以采样点实现快速粗配准，同时保留体素内点并构建微观结构，从更精细层次实现精配准。核心思想是：同名特征点的邻域点应具有一致的几何特征，这一假设由物体表面的几何连续性与邻域局部一致性所决定。
粗配准的作用是稳健地获取两个点云之间初始的转换矩阵 ，利用该矩阵可以将原始同名点集 中的源点集转换为 ，根据条件 获取准确的同名点集 。本文根据同名点所在的体素构建体素内微观结构。如图3所示，图3（a）表示同名点 经 变换后所在的体素，由于粗配准的精度有限，导致体素内的点云未能完全重叠在一起（图3（b）），本文基于体素内微观结构的精配准实现两点云的完全重叠（图3（c））。
图3
图3点云粗配准到精配准结果
Fig. 3Point cloud results form coarse registration to fine registation
图4
图4基于平面特征的体素微观结构构建
Fig. 4Construction of voxel microstructure based on planar features
传统的点到面的优化目标函数通常最小化源点云中各点到目标点云平面的法向量距离，与之不同，本文提出同时优化平面参数因子和点云变换参数因子来实现更高效的配准。将 作为参考基准（对应变换为单位矩阵），目标是优化相对于 的 的变换矩阵。特征点在参考坐标系中表示为
(11)
优化的目标函数是每个平面特征点 到其对应平面的平均距离
(12)
式中，nK表示平面的单位法向量；
(13)
则目标函数可重新表达为式（14），二次型 是一个标准的特征值问题。式（13）中内层最小化的解是协方差矩阵COV中最小特征值λmin对应的特征向量umin
(14)
目标函数式（13）通过代入平面参数的解析解，简化为仅依赖于T的函数，从而降低了优化变量的维度
(15)
尽管该函数可以使用LM算法迭代求解，为提高求解效率，本文在LM算法中引入Anderson加速策略[28-29]，式（14）代表的优化迭代可视为寻找固定点的过程：令Tk+1=G（Tk），定义残差函数H（T）=G（T）－T，则当T收敛至不动点T*时，H（T*）=0。Anderson加速器通过利用历史m次迭代的信息构建更优的搜索方向。为保证稳定性，仅当加速后的目标函数值低于当前迭代值时才接受该更新，否则回退至标准LM迭代步；同时设置最大加速历史窗口数m≤5，防止数值不稳定。最终的优化结果 作为精配准的变换矩阵。
1.2 基于精确匹配的形变区域提取与多特征融合的形变区域验证
点云形变检测的核心挑战在于如何从海量离散点中准确提取真实形变区域，并有效区分真实形变与测量噪声、系统误差及环境干扰。本文提出了一种由精确提取到多维验证的两阶段方法：第1阶段通过植被过滤、精确匹配和聚类分析实现形变区域的准确提取；第2阶段构建融合几何形态、统计分布和物理机理的多维度验证框架，对提取的形变区域进行真实性判别。
1.2.1 基于PCA特征分析的植被点过滤
植被的生长、摆动或季节性变化会在点云对比中产生显著的垂直位移，是虚假形变的主要干扰源。本文采用两阶段滤除策略：首先使用PTD[30]地面滤波算法提取初始地面点集，滤除乔木和高矮灌木。由于PTD基于地面坡度模型，对低矮匍匐植被、贴地草丛等地表附着型植被的分离能力有限，易残留于地面点集中。然后，进一步采用主成分分析（PCA）在局部邻域尺度上提取几何特征，识别残留植被点：对于待判定的点pi，在半径r内搜索k个邻域点，构建邻域协方差矩阵描述其空间分布
(16)
对邻域协方差矩阵进行特征值分解得到λ1≥λ2≥λ3及对应特征向量。特征值的相对大小蕴含局部几何信息：λ1≫λ2≈λ3表示线状分布（树干），λ1≈λ2≫λ3表示面状分布（地面），λ1≈λ2≈λ3表示散乱体状分布（树冠）。基于特征值定义4个归一化几何描述符
(17)
式中，Lλ为线性度；Pλ为平面性；Cλ为曲率；Vλ为垂直性（最小特征向量Z分量绝对值，该值越大表示局部结构在垂直方向延展越小）；ε为数值稳定项。PTD滤波后残留的近地植被呈现两类典型几何特征：匍匐茎秆、水平根系等近地线性结构在垂直方向延展极小，表现为高线性度与高垂直性；灌木冠层及不规则植被表面呈散乱分布，表现为高曲率与低平面性。据此建立双准则判别规则
式中，lveg（pi）为点pi的植被指示函数，取值1表示该点为植被点，0表示非植被点。
1.2.2 同名点精确匹配与形变区域聚类策略
在滤除植被干扰后，需在两期点云间建立准确的点对点对应关系来分析形变。由于点的离散性，直接使用最近邻匹配容易产生误匹配。对于参考点云（形变后）中的点 ，首先通过KD-Tree在目标点云（形变前）中进行K近邻搜索，获取候选匹配点集 ；然后，在候选集中应用XY平面距离约束，筛选出满足 的点，形成精筛候选集；最后选择Z向距离最小的点作为最终匹配
匹配成功后计算垂直位移 作为形变量，并根据阈值判定类型：di≤τs为沉降（标记－1），di≥τu为隆起（标记+1），τs
经过匹配和标记后，形变点呈离散分布，既包含成片的真实形变区域，也包含测量噪声产生的孤立点。本文采用基于欧氏距离的区域生长算法提取空间连续的形变区域。聚类Cm定义为满足以下条件的最大连通点集：对于任意点pi，存在同聚类中的点pj使得 且形变标记相同。为过滤噪声，仅保留点数nm>nmin的聚类。需指出的是，聚类目的是空间区域划定而非坐标平滑，每个点的Z向位移量由点对点匹配独立计算，聚类过程不修改任何点的坐标或位移值。对每个有效聚类Cm，计算质心cm（几何中心）、平均形变量 （整体形变水平）及形变范围[dmin，m，dmax，m]（内部形变差异程度），为后续的验证分析提供基础数据。
1.2.3 多维度特征融合验证方法
聚类分析虽然提取了空间连续的形变区域，但仍需进一步验证这些聚类是否为真实形变。局部粗差和数据质量问题等都可能形成满足聚类条件但不符合规律的伪形变区域。本文构建了包含几何形态、统计分布和物理机理的多维特征验证框架，通过多角度交叉验证确保识别结果的可靠性。
（1）几何形态特征提取与约束。利用Graham扫描算法[31]计算聚类在XY平面的凸包Hm，凸包描述了形变区域的外轮廓。根据凸包顶点坐标，采用Shoelace公式计算面积： ，通过累加边长计算周长Pm。面积是最基本的形态参数，过小的面积通常对应噪声点聚集，过大的面积可能源于系统偏差或配准误差。为了进一步描述形状规则性，引入无量纲形状描述符紧致度 ，度量形状的复杂程度，圆形时达到最小值4π≈12.57，值越大表示边界越不规则；圆度 是紧致度的倒数归一化形式，圆形时为1，不规则形状时趋于0。
（2）统计分布特征分析与检验：从形变量的概率分布检验合理性。微观层面的形变是众多随机因素综合作用的结果，根据中心极限定理会呈现正态性。本文计算聚类内形变量的高阶统计矩：均值μm、方差 、偏度γ1和峰度γ2。偏度反映分布的对称性，正偏度表示存在较多的正向极值，负偏度则相反；峰度反映分布的尖锐程度，正峰度表示比正态分布更集中，负峰度表示更平坦。若聚类的偏度或峰度显著偏离正态分布，则表示可能存在异常数据。标准差比例Rσ，m=σm/（dmax，m－dmin，m）量化了数据的集中程度，该值过小表示几个极端值主导，过大表示分布过于离散，真实形变通常Rσ，m∈[0.2，0.5]。除了点值统计，空间梯度场分析能够揭示形变的空间连续性。对于聚类中的每个点，计算其与邻域点的形变梯度
该值反映了局部形变的陡峭程度。统计梯度场的均值 、标准差σG和最大值Gmax，m，这些参数描述了形变场的粗糙度和突变特征。真实形变由于受到土体连续性约束，其梯度场相对平缓且均匀，因此设定约束Gmax，m<5.0（避免过陡突变）。基于梯度定义平滑度因子Φsmooth=exp（－Gmax，m/τG），该值越接近1，表示形变场越平滑。此外，采用3σ准则检测离群值，计算离群值比例ηoutlier，m，该比例过高表明存在大量异常点，聚类可能不可靠。综合这些统计指标，可以从概率分布和空间结构两个角度判断形变的合理性。
（3）物理机理特征建模与检验：该部分将地质力学中的经验规律转化为可计算的数学模型。在均质连续介质假设条件下，地面沉降或隆起的形成遵循一定的物理规律：荷载作用下形变从加载中心向外传播，幅值逐渐衰减；形变场受土体连续性约束，呈现空间连续分布；在各向同性或近各向同性的地质条件下，形变呈径向对称或近对称分布[32]。基于此提取两类物理特征：①中心峰值特征检验形变的衰减规律。将聚类按照到质心距离划分为多个同心环带 }，计算各环带内形变量绝对值的平均值 ，真实形变应满足 ，即从中心向外单调递减。考虑到测量噪声，可以放宽为 。②空间连续性综合评估形变场的整体平滑性和一致性。真实的土体形变受到物质连续性约束，不会出现空间上的剧烈跳变。定义空间连续性指标 ，其中，第1项基于平均梯度，反映整体平滑度；第2项基于梯度标准差，反映梯度的均匀性。权重系数α、β根据经验设定。该指标越接近1，形变场越符合连续介质假设。
（4）综合置信度评估：综合上述特征，构建加权评分模型量化聚类可信度。几何特征评分fg基于形状描述符与参考值的接近程度，将紧致度和圆度归一化后平均；统计特征评分fs基于平滑度因子、离群值比例和梯度约束的满足程度；物理特征评分fp基于中心峰特征和空间连续性的综合表现。最终的综合置信度定义为3类评分的加权和
(21)
式中，wg、ws、wp为权重系数（wg+ws+wp=1）；fg、fs、fp为归一化评分函数。仅当置信度Confm≥Confmin且硬约束条件均满足时，聚类被判定为有效形变。
上述多类约束并行工作，而非依赖单一阈值进行判别。在集中形变场景下，统计约束与几何约束起主导作用；在缓形变场景下，单点位移量可能接近噪声水平，此时空间连续性约束通过识别相邻点一致性小位移的连续分布将其与随机噪声点区分。地表覆盖干扰已由预处理阶段的地面点提取与PCA植被过滤消除。并行设计保证了单一约束判别力下降时框架的整体稳定性。
2 试验与分析 2.1 模拟试验与分析 2.1.1 模拟试验设计与数据获取
在实际工程应用中，形变检测方法的性能评估面临真值获取困难的问题。现实场景中的地表形变缺乏高精度的独立测量基准，且形变过程具有不可重复性，使得方法精度的客观评估变得极其困难。为了准确评估形变检测的准确度，本文设计了一种形变模拟方法。首先，选取真实地形点云作为基准数据；然后，基于椭圆高斯衰减模型在基准点云上施加已知参数的人工形变，该模型通过椭圆方程确定形变影响范围，并采用非线性高斯衰减函数式（22）计算各点形变量
(22)
式中，dmax为形变中心的最大形变量；dmin为椭圆边缘的最小可检测形变量（试验设为3 mm）； 为归一化椭圆距离；α为衰减系数，控制形变分布的陡峭程度；εi为符合 正态分布的随机噪声项，模拟实际测量中的随机误差。考虑到两期点云采集存在测站位置偏移、仪器姿态差异等几何变化，通过施加已知的三维刚体变换P′=RP+t模拟几何差异，旋转矩阵R模拟仪器姿态变化，平移向量t模拟测站位置偏移；最后，运用本文方法对模拟数据进行形变检测，并与真值对比评估性能。需指出，该连续高斯衰减模型主要适用于土体压密沉降、缓慢滑移等渐变型形变，这也是TLS早期预警监测的主要目标形变类型。
模拟试验数据选用WHU-TLS Mountain公开数据集。该数据集由武汉大学测绘遥感信息工程全国重点实验室提供，采用Leica ScanStation C5三维激光扫描仪进行6次独立扫描采集，模拟试验采用第5站扫描（图5）。该数据集的场景为山地地形，无显著垂直建筑结构，主要包含未结构化岩石和草地等自然地物。
图5模拟试验数据
注：红框区域内为配准后显示的模拟形变区域。
Fig. 5Simulated experimental data
2.1.2 模拟试验分析
为全面、客观地评价本文方法的性能，本文从分类性能和形变检测精度方面进行评估。分类性能评估采用混淆矩阵构建4个基本评估指标：精确度（Precision）、召回率（Recall）、F1值及交并比（IoU），计算公式为
(23)
式中，TP、TN、FP、FN分别表示真正例、真负例、假正例和假负例的数量；IoU综合衡量检测区域与真实形变区域的空间重叠程度；A和B分别为方法检测区域和真实形变区域在XY平面上的投影多边形。这些指标能够同时惩罚漏检（FN）和误检（FP）。对于成功识别的形变点，需要进一步评估其位移检测精度。本文采用平均绝对误差（MAE）、均方根误差（RMSE）量化位移检测误差。
平均绝对误差反映位移检测的平均偏差水平，具体公式为
(24)
均方根误差强调大误差的影响，对异常值更加敏感，具体公式为
(25)
为评估本文方法在不同形变量级下的检测性能，设计多组试验。在保持其他参数不变的情况下，设置最大形变量dmax为[0.005，0.006，0.007，0.008，0.01，0.02，0.05，0.1]m。对比分析不同形变量级下各项评估指标的变化规律，以深入了解方法的适用范围和性能边界。
表1为本文方法在不同形变量级下的性能表现。从分类性能来看，所有测试场景下精确度均达到100%，说明理想情况下方法具有强抗误检能力，能有效避免将稳定区域误判为形变区域的问题。召回率呈现形变量级依赖特性，在0.005 m微小形变下为75.60%，随形变增大快速提升，0.008 m时达91.19%，更大形变量级下稳定在95%以上。这种现象的原因主要在于当形变幅度接近或略高于噪声水平时，部分真实的微小形变点被误判为噪声点而过滤。即便如此，0.005 m形变量下的F1值仍然达到86.11%，0.008 m形变量下F1值达到95.39%，验证了本文方法在毫米级范围内有效识别形变的能力。IoU指标与召回率在本试验中数值相同，这是由于精确度为100%（即FP=0），此时IoU=TP/（TP+FN）=Recall。在形变量≥6 mm时，IoU超过80%，说明检测区域与真实形变区域的空间重叠程度良好。
表1不同形变量参数下的方法性能评估结果
Tab. 1 Method performance evaluation results under different deformation magnitude parameters
dmax
/m
分类性能
形变检测精度
精确度/（%）
召回率/（%）
F1
值/（%）
IoU/（%）
MAE/mm
RMSE/mm
100
0.33
0.52
100
0.33
0.52
100
0.34
0.52
100
0.34
0.52
0.01
100
0.33
0.52
0.02
100
0.35
0.54
0.05
100
0.33
0.53
0.1
100
0.35
0.61
形变定量精度方面。MAE稳定在0.33~0.35 mm之间，RMSE维持在0.52~0.61 mm的范围内，这种误差高度一致性表明本文方法误差主要由系统性因素主导（源于点云配准的系统性偏差），但其量级稳定在亚毫米级别，能够充分支撑毫米级的形变检测需求。图6展示了在施加8 mm模拟形变下每个点模拟形变与本文方法检测形变的对比结果。
图6
图68 mm模拟形变下的结果对比
Fig. 6Result comparison under 8 mm simulated deformation
2.2 真实形变试验与分析
为验证方法在真实场景中的性能，在某大学校园边坡开展了受控形变试验。相比模拟试验，真实场景试验面临植被干扰、地表粗糙度差异等挑战。在真实地质环境中人工施加已知形变，可以在接近实际应用条件下全面评估本文方法的稳健性和实用性[33-34]。
2.2.1 可控形变施加与数据采集
试验区域为校园内典型风化岩质边坡，地表覆盖层主要由风化岩屑、碎石和土壤混合组成，植被以低矮灌木和草本植物为主，地形起伏不平，具有典型的自然边坡地质特征。为模拟真实地质灾害中的沉降和隆起两类典型形变模式，设计了两种受控形变方案：①采用压实模板法模拟地表沉降过程。如图7（a）所示选取尺寸为40 cm×30 cm×3.0 cm的刚性木板作为压实工具，在选定的监测区域内，通过均匀施加垂直压力使地表土层产生压密沉降（由于土质较硬，实际产生形变应在1.0~2.0 cm之间）。②通过放置物体法模拟地表隆起现象。如图7（c）所示，在监测区域内放置一块纸板（尺寸约15 cm×12 cm×2.5 cm），隆起区域边界相对清晰，便于验证本文方法对局部隆起形变的检测能力。
图7
图7试验现场与形变布设
Fig. 7Experimental sites and deformation setup
试验采用天宝X7三维激光扫描仪进行数据采集，扫描距离为8~12 m，单站扫描设置为标准模式，分辨率为10 Mpts，单次扫描耗时约为2 min。点云平均点间距为5~8 mm，能够满足毫米级形变检测精度需求。采用两期扫描方案：第1期为形变前的参考点云，第2期为施加形变后的变形点云。两期扫描均在不同测站位置完成。对原始点云进行预处理，包括噪声滤波和地面点提取等步骤，得到第1期参考点云点数为4 406 631个，第2期变形点云点数为4 529 193个。通过高精度点云配准方法实现两期数据的精确对齐（图8），随后进行形变区域提取与多维度验证。
图8
图8真实扫描点云配准
Fig. 8Real scanned point clouds registration
2.2.2 试验结果与分析
根据试验场景特点和预期形变量级，本文方法参数设置如下：点匹配阶段采用K=20的邻近搜索，XY平面距离阈值为10 mm；沉降检测阈值为－5 mm，隆起检测阈值为+5 mm；聚类分析采用30 mm欧氏距离容差，最小聚类点数为500。验证框架关键参数：最小面积为0.05 m2，最大紧致度为25.0，梯度最大值为5.0，最小置信度为60%。
空间聚类分析在沉降点集中提取3个聚类，在隆起点集中提取两个聚类，共5个初始聚类。对这5个聚类进行几何形态、统计分布和物理机理三维特征的定量分析，为简洁起见，本文仅列出各维度的代表性指标：几何形态特征（面积、紧致度）、统计分布特征（位移标准差比例、梯度最大值）、物理机理特征（空间连续性）及综合置信度，验证结果见表2。
表2真实边坡试验形变聚类验证结果
Tab. 2 Quantitative validation results of deformation clusters in real slope experiment
聚类ID
点数
平均位移/mm
面积/m 2
紧致度
位移标准差比例
梯度最大值
空间连续性
综合置信度/（%）
验证结果
S-0
895
－15.0
18.3
3.6
63.6
通过
S-1
707
－15.1
19.7
3.9
63.2
通过
S-2
646
－168.0
0
拒绝
U-0
810
17.3
21.5
6.8
54.3
拒绝
U-1
566
22.6
17.6
4.2
66.4
通过
验证通过的聚类分析：①沉降聚类S-0和S-1分别包含895和707个点，平均沉降量约为15 mm，对应试验中两处压实模板的施加位置。几何形态特征方面，两者面积分别为0.112和0.089 m2（二者均大于0.05 m2阈值），紧致度分别为18.3和19.7（二者均小于25.0阈值），形状规则；统计分布特征方面，两者位移标准差比例分别为0.182和0.143（远小于0.5阈值），梯度最大值分别为3.6和3.9（二者均小于5.0阈值），表明位移分布均匀且过渡平滑；物理机理特征方面，两者的空间连续性分别为0.101和0.124（二者均大于0.01阈值），符合土体连续介质要求；综合置信度分别为63.6%和63.2%，均超过60%阈值，被判定为有效沉降区域。②隆起聚类U-1包含566个点，平均隆起量为22.6 mm，位于纸板放置位置。几何形态特征面积为0.073 m2，紧致度为17.6；统计分布特征梯度最大值为4.2，物理机理特征空间连续性为0.150（所有聚类中最高），综合置信度为66.4%，成功识别为真实隆起形变。
被拒绝的聚类分析：①沉降聚类S-2检测到平均沉降量高达168.0 mm，远超试验施加的形变量级（约15 mm）。几何形态特征显示其面积仅为0.026 m2，低于最小面积阈值（0.05 m2），在第一轮验证即被拒绝。经现场核查，该聚类位于测站附近，异常位移由未滤除的行人扫描点引起。②隆起聚类U-0虽然平均隆起量17.3 mm在合理范围，但统计分布特征显示其梯度最大值达到6.8，显著超过阈值（5.0），表明存在位移的剧烈空间变化，违反土体连续形变规律。物理机理特征显示其空间连续性仅为0.068，综合置信度为54.3%，未达标，被拒绝。现场分析发现，该聚类位于植被区，检测到的“隆起”实际是植被变化引起的伪形变。
图9为形变检测结果对比。有效沉降区域准确覆盖了两处压实位置；有效隆起区域集中在纸板位置，空间范围与其平面尺寸基本吻合。
图9
图9形变检测结果对比
Fig. 9Comparison of deformation detection results
所有试验在Intel Core Ultra 9 285K处理器（24核）、96 GB内存、Windows 11、MSVC 2022（Release x64，OpenMP并行）平台上完成。表3为全流程各阶段计算耗时，总耗时为31.68 s，能够满足地质灾害监测中非实时处理的需求。本文提出的粗配准、精配准和多维验证3个核心阶段合计仅耗时0.38 s，占比不足2%；预处理（14.52 s）和形变计算（16.77 s）为主要耗时环节，后续可通过GPU并行化进一步压缩处理时间。
表3全流程各阶段计算耗时
Tab. 3 Computational times for each stage of the complete processing pipeline
统计指标
耗时/s
占比/（%）
预处理
45.8
粗配准
0.19
0.6
精配准
0.05
0.2
形变计算
52.9
多维验证
0.14
0.4
全流程总计
100
3 讨论 3.1 消融试验与分析
为验证各模块的有效性，本文设计了5组消融试验，定量结果见表4。表中，M1采用传统配准基线（FPFH+RANSAC+ICP）结合标准C2C欧氏距离度量，代表最基础的点云形变检测流程；M2在相同配准下将度量替换为标准M3C2法向投影距离；M3在M2基础上引入本文多维验证框架；M4将配准替换为本文方法，保留标准M3C2度量但不加验证；M5为本文完整方法。配准重叠率定义为全分辨率点云中点对面距离小于1 cm的配准点占比，重叠区均值P2Pl为这些重叠点的平均点到面距离，两者共同反映配准质量；初始形变点数为阈值判别后的形变候选数量；有效聚类数和有效形变点数为经多维验证后保留的最终检测结果。图10为5组消融试验的形变检测结果（蓝色为沉降，红色为隆起，深灰色为稳定区域）。
表4消融试验定量结果
Tab. 4 Quantitative results of ablation experiments
方法
配准重叠率/（%）
重叠区均值P2Pl/mm
初始形变点数
有效聚类数
有效形变点数
M1
5.06
3 115 409
M2
5.06
3 076 770
M3
5.06
3 076 770
35
40 218
M4
1.73
383 734
M5
1.73
129 011
3
2168
图10
图105组消融试验的形变检测结果
Fig. 10Deformation detection results of five ablation experiment groups
M1与M2采用相同的传统基线配准，配准重叠率仅47.21%，重叠区均值P2Pl达5.06 mm，超过半数配准点的点对面距离大于1 cm。在此配准精度下，M1产生3 115 409个初始形变点，M2替换为M3C2法向投影后初始形变点数为3 076 770个，两者差异甚微。图10中M1与M2的场景均被蓝色和红色的伪形变大面积覆盖，几乎不可见稳定区域。由于本试验将±5 mm作为形变判别阈值（扫描仪标称精度的2~3倍），而传统基线配准的重叠区均值P2Pl（5.06 mm）已达到该阈值水平，配准残差本身即被大量误判为形变响应，度量方式的改进无法从根本上降低虚假检测。
M4仅将配准替换为本文方法，配准重叠率提升至81.87%，重叠区均值P2Pl降至1.73 mm，初始形变点数从3 076 770个骤降至383 734个。然而图10中M4仍存在大量零散的伪形变响应，分析表明其383 734个初始形变点中45.7%集中在5~8 mm的近阈值区间，且约41.6%分布于PTD地面滤波后残留的植被区域，标准M3C2的最近邻匹配无法区分植被变化与地形形变。M5针对上述问题，通过PCA植被过滤剔除残余植被点、法向一致性约束排除不可靠匹配，将初始形变点数压缩至129 011个；经多维验证后最终输出3个空间集中的有效聚类，图10（e）中仅在3处局部区域呈现集中的形变标记。
M3与M5同样包含验证框架，但两者的输入质量截然不同。M3建立在传统基线配准上（配准重叠率47.21%），验证后仍保留35个有效聚类（40 218个点），是M5的11.7倍。图10中M3虽较M1、M2有所改善，但仍残留大量分散的伪形变区域。这表明当配准残差与检测阈值处于同一量级时，伪形变聚类在空间分布和统计特征上与真实形变相似，仅靠验证框架无法有效区分。
综上，配准精度是整个检测链路的核心瓶颈。本文配准方案将配准重叠率从47.21%提升至81.87%，重叠区均值P2Pl从5.06 mm降至1.73 mm，为下游检测提供了可靠基础；在此基础上进行形变检测，本文方法通过植被过滤与法向约束进一步消除非地形噪声，验证框架精确筛选出空间连续且统计可靠的形变区域，最终从3 115 409个初始形变点中定位3个有效形变聚类（2168个点）。3个模块协同工作，缺一不可。
3.2 稳定区基准噪声验证
真实场景下的形变检测精度受配准残差与传感器测量噪声的共同影响，有必要对其进行实测量化。基于2.2节真实边坡试验的两期点云数据，在重叠区域内自动选取远离已知形变聚类质心（距离≥2 m）的2 m×2 m稳定地面区域，对该区域内所有有效匹配点对统计Z向差值与点-平面距离分布，以综合反映配准残差与传感器测量噪声的叠加效果，结果见表5。
表52 m×2 m稳定参考区域多期扫描Z向差值与点-平面距离统计
Tab. 5 Z -direction difference and point-to-plane distance statistics of multi-epoch scans in the 2 m×2 m stable reference region
有效样本点数
Z
向差值均值/mm
Z
向差值标准差/mm
95%置信区间/mm
点-平面距离均值/mm
点-平面距离标准差/mm
429 801
－0.099
±0.723
±1.417
±0.752
由表5可知，Z向差值均值为－0.099 mm，表明两期点云间的整体配准结果无明显系统性偏移趋势；Z向差值标准差σ=±0.723 mm与点-平面距离均值0.707 mm数值相近，两种独立精度指标相互印证，说明配准残差处于亚毫米量级。95%置信区间为±1.417 mm，基于3σ准则，本文方法的统计检测下限约为2.2 mm，具备毫米级形变的识别能力。试验中检测到的15 mm沉降对应显著性达20.7σ，22.6 mm隆起对应显著性达31.3σ，二者均远超噪声水平，形变检测结论具有统计可靠性。
3.3 局限性与展望
本文方法在现有试验条件下取得了良好效果，但仍存在若干局限性有待改进。在参数设置方面，验证框架中的置信度阈值、梯度约束等关键参数目前依赖先验经验设定，不同地质场景（软土、硬岩、松散堆积体）下最优参数存在差异，后续将研究面向场景自适应的参数估计机制以降低人工调参依赖。在精度验证方面，后续研究将引入全站仪或高精度GNSS对形变量进行独立测量，以更严格地量化本文方法在真实场景下的误差分布特征。在形变验证的完整性方面，本文模拟试验基于高斯衰减模型生成了空间连续递减的渐变型形变，验证了本文方法对连续梯度形变的检测能力；真实边坡试验在植被干扰、地表粗糙度差异等复杂条件下验证了本文方法的稳健性。但受控形变区域内位移相对均匀、边界过渡较陡，与自然蠕滑形变的连续梯度分布和渐变边界特征仍有差距，后续将结合实际滑坡监测数据进一步验证。在监测能力方面，本文仅开展双期点云对比分析，后续将引入Kalman滤波或LSTM模型融合多期时序信息以追踪累积形变演化；并探索TLS与UAV[35]、InSAR的多源联合验证机制，构建多尺度形变监测体系。
4 结论
针对地质灾害监测中激光点云配准精度受限和形变检测抗干扰能力弱等问题，本文提出了一种高精度激光点云配准驱动的毫米级地表形变检测方法。首先基于匹配图节点和边可靠性度量的特征点粗差剔除方法，显著降低特征匹配粗差率；然后构建了增强GNC-Welsch函数的配准参数稳健估计方法，有效解决了局部极小值问题，在高粗差率条件下仍能获得可靠的粗配准参数。设计了微观结构混合特征因子驱动的精配准方法，保留体素内点云并构建微观结构，将配准精度提升至毫米级。建立了融合几何形态、统计分布和物理机理的多维度形变验证框架，有效区分真实形变与测量噪声。模拟试验表明，本文方法在所有测试场景下精确度均达100%，5 mm形变量下F1值达86.11%，8 mm形变量下达95.39%，平均绝对误差稳定在0.33~0.35 mm，均方根误差稳定在0.52~0.61 mm，充分验证了毫米级形变检测能力。真实边坡监测试验表明，本文方法能够准确识别15 mm级沉降和22.6 mm级隆起形变，多维度验证框架有效剔除了植被变化和异常点等引起的伪形变。本文方法为复杂场景下地质灾害的毫米级形变监测提供了技术支撑。未来研究将把本文方法推广至更大范围的实际工程场景，并开展系统性验证。
来源：智绘科服
热门跟贴