一种局部稀疏地面点云与已有DEM的融合方法

雷丽珍，林 超

(广东省国土资源技术中心，广东 广州 510075)

利用机载LiDAR点云数据快速完成高精度数字高程模型(DEM)的大规模生产，已成为近几年自然资源管理、城镇规划与建设等诸多领域的一个研究和应用热点[1-5]。机载LiDAR点云数据质量好坏对后期DEM生产具有重要影响。由于设备特点、测区地形地貌、植被覆盖、成本和天气等诸多因素影响[6]，实际获取的机载LiDAR点云数据存在地面点局部过于稀疏的问题，导致采用三角网内插算法构建DEM时，出现较严重的“三角面片化”问题，在一定程度上降低了成果精度[7]，影响DEM的制图与使用，因此需要对这些LiDAR点数据进行数据质量改善。而进行局部区域的重飞或补飞，要消耗较大的成本、时间，且后期获取的数据与前期获取的数据之间可能存在一定的误差，数据一致性与整合也会存在一定的问题[8]。针对以上难题，本文提出一种基于高斯核函数加权迭代插值算法的局部稀疏地面点云与已有DEM融合方法，以实现LiDAR点云生成DEM时局部地面点稀疏区域的地形特征保持及高程精度的改善。

1 研究方法

1.1 基本思路

LiDAR地面点云与已有DEM之间存在一定的高程差异，且高程差异为非刚性的偏差，不同位置的差异值不一致。本文提出一种内插算法，在LiDAR点云稀疏区域，以滤波后少量正确的地面点作为控制点，利用高斯核函数加权迭代插值方法改正低精度的DEM数据，对填补在LiDAR地面点较稀疏区域的DEM点云的高程值进行局部改善，使其高程数值与其邻近的LiDAR地面点的高程值保持一致，同时尽量保持原始DEM的地形形态特征[9]。

1.2 融合插值原理

地表形态通常具有连续性，在一定的空间邻域内具有一定的相似性。地理学第一定律所描述的“越邻近的事物越相关”，为空间数据插值等方法应用提供了基础的理论依据。空间插值方法是将不同粒度、不同尺度的空间数据进行有效融合的一种常用方法，它使用一组采样点的观测值的加权组合来确定待观测点的值。本文采用加权插值方法实现稀疏地面点云与已有DEM的数据融合，将稀疏地面点云作为高精度的控制点，通过插值原理估计控制点一定邻域范围内DEM高程的局部改正值，同时在对DEM高程进行纠正时尽量保持DEM的整体地形形态特征。其原理如图1所示。

如图1所示，估计p点处的高程改正值时，首先搜索其一定距离范围内的LiDAR地面点；然后，依据这些地面点云的真实高程与其估计高程值间的高程差异(即dH)，通过高斯核函数加权来估计p点位置的高程改正。计算公式如下

(1)

式中，p表示待估计高程改正值的DEM点；q表示LiDAR地面点云G中的一个点；d(q,p)表示q与p之间的欧氏距离。上述公式值中，权值wq的确定是影响高程改正的重要因素。本文采用待估点与控制点间距离的高斯核函数作为控制点的加权权重，即

(2)

式中，h为高斯核函数的带宽。

由于稀疏地面点云数据的空间分布密度不均，采用全局高斯核函数带宽参数很难获得比较理想的插值结果，为此，本文提出基于Delaunay三角网的高斯核函数带宽h的自适应确定方法[10]。首先采用稀疏地面点云构建Delaunay三角网，然后统计与每个地面点(如pi)相连的三角网边的平均长度(记为mi)和Delaunay三角网中所有边长的标准方差(记为s)，则对于点pi的高斯核加权核函数的带宽h计算为

h(pi)=mi+3s

(3)

进而，按照以下步骤迭代对已有DEM的高程进行局部改正：

(1) 根据已有的DEM数据，构建一定距离间隔(如2 m×2 m)的规则格网，通过双线性插值算法获得每个格网角点的高程值，将这些格网的角点保存为带高程的DEM点并输出。

(2) 针对每个DEM点，通过最邻近查询，搜索DEM点附近的LiDAR地面点，计算这些DEM点的平均高程作为LiDAR点的估计高程值，记为H1，LiDAR点原始高程值记为H0，计算两者高程差异

dH=H0-H1

(4)

(3) 统计所有LiDAR点的高程差异dH的平均值、最大差异值及大部分点的差异值。

(4) 利用高斯核函数加权插值方法，根据LiDAR点云平面位置、dH数值和DEM点的平面位置计算出DEM点位置处的高程改正项correctH，将原始DEM高程数值加上改正后的高差作为改正后DEM格网点的高程值[11]。

(5) 迭代执行以上步骤，直到DEM点的高程值收敛(或改变很微小)时停止。

1.3 技术流程

本文方法在LAS点云格式下进行运算，技术流程如图2所示。选取需处理的稀疏地面点云区域的数据，与转换为LAS格式的DEM点进行融合处理生成插值点，与LiDAR点云地面点叠加内插生成DEM数据。

1.4 算法实现

基于以上原理开发设计了点云数据融合工具软件，界面设计如图3所示。

点云数据融合工具主要包括以下功能：① 数据输入：选择需要进行处理的LiDAR点云数据文件和DEM点云数据文件；② 参数设置：设置点搜索的邻域半径，LiDAR点按类别读取及高程改正的平滑参数等；③ 结果输出与保存；④ 运行状态提示，反馈中间计算的运行情况。

2 试验结果

试验采样了某测区内3.9 km×1.5 km范围内的点云数据作为原始待融合的LiDAR点云数据。该区域LiDAR数据中一些区域地面点较密，也有一些区域地面点非常稀少，这些区域并非平地，过于稀少的LiDAR地面点会影响这些区域的DEM构建。为此，采样同一区域已有的较低精度DEM点云数据进行LiDAR地面点的填补与改善。

测试结果的局部区域显示结果如图4所示。

从测试结果可以发现，该算法可以有效地对原始稀疏LiDAR地面点进行填补和改善，在对稀疏LiDAR点进行插值的同时，也尽可能地保持了原始地形走势和总体的地貌形态特征。根据插值后的点云生产的DEM如图5所示。

通过比较融合前后地面点生产DEM可以发现，基于原始地面点生产的DEM在地面点稀少的局部区域存在较严重的“面片化”现象，地形过渡不自然、不流畅；而经过插值后的点云在这些区域的地形更自然，贴近自然的地貌形态特征。

3 精度及可靠性验证

为了综合验证本文方法插值的精度及可靠性，选取了某区域作为验证区(如图6(b)所示)。该区域范围为700 m×700 m，数据获取时设备穿透性高，分类后地面点分布密集，点云数量为4 223 559，统计密度为5.2点/平方米，地面点数量为215 440，该测区实际检验的点云高程中误差为0.2 m，通过人工剔除120 133个(占55.8%)地面点模拟实际中可能出现的稀疏点云状态(如图6(c)、图6(d)所示)，并与该区域已有低精度DEM(如图6(a)所示)数据进行本文融合算法处理，并生成2 m格网DEM成果(如图6(e)、图6(f)所示)。

由图6可以看出，已有低精度DEM与基于点云生成的高精度DEM在山体地形走势上基本一致，但高精度DEM在细节上更为丰富和准确，模拟稀疏点云剔除部分关键地面点后局部山脊及山谷明显失真，通过高斯核函数加权内插后，失真情况在视觉上有明显改善。从剔除的地面点中均匀选取1881个模型关键点作为本次验证的检查点，对以上数据进行误差统计。检查点如图6(g)所示，主要分布在模拟稀疏点云区域，且数量较多，对试验中各个DEM数据进行误差计算分析。统计后中误差情况见表1，对各区间残差分布的百分比对比情况分析如图7所示。

表1 中误差统计 m

利用高斯核函数加权内插融合后的DEM的中误差为1.796 m，在精度上相对于融合前的稀疏地面点情况有一定改善，且误差呈正态分布，与稀疏地面点DEM误差分布曲线相比，内插后粗差有明显减小，DEM在弱精度区域的可靠性显著提升，鉴于对点云细微地物完全删除后难以通过内插恢复，误差在±2 m区间内的变化差异不大。

4 结 语

本文针对机载LiDAR点云数据中地面点稀少问题，提出了基于已有低精度DEM数据与LiDAR稀疏地面点进行融合的技术方法，并通过试验对融合效果、精度及可靠性进行了分析验证。结果表明，使用融合后的点云数据生产的DEM产品，可以减少稀疏地面点云生成的DEM中存在的“三角面片化”问题，并提高了DEM的精度[12-13]，尤其是在存在较大误差的区域，对实际生产中出现的地面点稀疏区域进行局部融合，能够改善DEM的整体质量。