一种基于密集点云数据构建道路DEM的方法

本发明属于高精度数字高程模型构建的技术领域，涉及顾及道路的形态特征和密集点云数据特征来实现道路dem构建，以保证所构建得道路dem具有较高得形态精度和高程精度。

背景技术：

数字高程模型(digitalelevationmodel，简称dem)是通过有限的地形高程数据实现对地面地形的数字化模拟，它是用一组有序数值阵列形式表示地面高程的一种实体地面模型。dem是国家基础地理数据库的核心内容，在交通、军事、水利、农业、环境、资源管理、规划与旅游等领域都有广泛的应用。

当前主要的dem构建方法适用于连续变化的自然地形，然而在城市等人类活动痕迹明显的区域，各种人类工程造就了许多非连续的地形，例如道路、堤坝、坡面等，其中道路作为城市区域地形的骨架部分，其高精度地形建模受到相关领域的重视。与自然起伏地形不同，道路具有显著的形态特征，例如从形态的角度看一般道路边界较为清晰且形状规则，从高程变化方面看道路一般是纵向在一定坡度范围内起伏延伸，而在横向高程变化基本可以忽略。因此，在道路dem构建时，必需要顾及道路这类几何信息与语义信息共生的特征，保证构建的道路dem在形态上能符合实际道路形态特征。

在道路dem构建领域，以往较少考虑道路的形态特征，主要以高程点、特征线等地形数据为基础，采用克里格、样条等插值方法进行空间内插得到道路dem。后来有学者意识到道路形态与其他区域的明显差异，试图采用道路边线作为约束信息来区分道路与其邻接区域，但是由于在道路内部dem构建时并没有考虑道路的特殊形态，构建结果改善并不明显。近期有学者提出一种基于形态与语义特征的道路dem建模方法，该方法以大比例尺地形图数据为数据源，以保证道路建模结果的形态特征为主要目标，对高程点数据进行筛选，然后依赖空间插值技术实现道路dem构建。该方法构建的道路dem形态上符合实际道路形态，但是由于采用的数据过于稀疏，难以保证所构建道路dem的高程精度，此外，大比例尺地形图获取代价高昂，也决定了该方法难以推广应用。

技术实现要素：

为了解决上述现有技术的问题，可以将道路表面的形态特征与可用建模数据的密集特征相结合，一方面要强调以观测数据为主，以空间插值技术为辅来计算道路dem高程值；另一方面，也需要对道路dem做适当的平滑处理以保证道路正确的形态。因此，本发明提供一种基于密集点云数据构建道路dem的方法，在保持高程精度的同时能够对道路形态进行很好地表达。

本发明采用的技术方案如下：

一种基于密集点云数据构建道路dem的方法，包括以下步骤：首先对道路进行子单元划分，并基于密集点云数据计算各个道路子单元的平均高程值；然后，根据设置好的道路最大起伏阈值对各个道路子单元平均高程值进行修改；接着，生成建模区域dem格网，与道路子单元叠加；最后，遍历dem格网，基于格网中心点与道路子单元的空间位置关系以及道路子单元的高程值计算dem格网的高程值。

进一步地，所述方法具体包括如下步骤：

s1，划分道路子单元，即按照一定步长，将面状道路数据划分成一系列单元格，其中若存在道路相交的情况，则将相交路口作为一个子单元；子单元划分完成后，生成道路子单元中心点；

s2，基于密集点云数据，计算每一个道路子单元中的高程平均值，将该值赋给道路子单元中心点；

s3，采用道格拉斯算法，对道路的各个子单元的高程平均值进行计算调整，使得相邻道路子单元之间的高程变化率不高于预先设置的阈值；

s4，生成建模区域格网dem格网的中心格网点，并与道路子单元数据进行空间叠加；

s5，根据dem格网的中心格网点与其两侧道路子单元中心点的位置关系，按照两道路子单元之间高程线性变化的原理，计算dem中心格网点的高程值；

s6，依次完成所有dem格网中心点的高程值计算，并按照栅格数据规范输出计算结果。

进一步地，所述步骤s3的具体过程为：

s31，将当前所处理的道路中，位于道路头a和道路尾b处道路子单元的高程值标记为不可编辑，即该处的高程值是确定的；

s32，以当前所处理的道路为主干道，判断该道路是否与其他道路相交，如果有相交道路，则将相交道路口形成的道路子单元高程值标记为不可编辑，即道路口处的高程值是确定的；

s33，将道路子单元中心点由三维空间投影到二维空间，其中道路头a的x坐标为0，y坐标为其高程值；其余单元格中心点的x坐标为距道路头a的距离，距离值为沿道路延伸方向计算得到，y坐标为其高程值；

s34，连接道路头a和道路尾b两个子单元的几何中心点对应的投影点，分别标记为a1和b1，根据其三维坐标确定直线a1b1的表达式；

s35，遍历投影点a1、b1之间的未标记投影点，确定一个点c1，使得点c1到直线a1b1的距离最远；

s36，计算直线a1c1和c1b1的高程变化率，与设置好的阈值比较，若高程变化率大于设置阈值则修改点c1的高程值，使得计算的高程变化率等于设置阈值，然后将原道路划分为新的两段道路，即直线a1c1和b1c1，分别以点a1与c1、b1与c1作为新道路两端点，重复步骤s34-s35；若高程变化率小于设置阈值，则点c1处高程值不变，并将点c1标记为已处理，再执行步骤s35；

s37，待所有投影点都标记完成，将投影点高程值赋给原始道路子单元中心点，完成道路子单元中心点高程值的更新。

进一步地，所述步骤s5的具体过程为：

s51，过每个道路子单元的中心点，生成道路中心线的垂线；

s52，依次遍历与当前处理道路叠加的dem格网中心点，确定其两侧的中心线的垂线，分别标记为l1和l2；

s53，计算待计算高程值的dem格网中心点距垂线l1和l2的距离，并分别标记为d1和d2；

s54，基于d1和d2以及l1、l2所对应的子单元中心点的高程值v1和v2，按照线性变化原理，计算当前dem格网中心点的高程值，计算公式为：

s55，重复执行步骤s52-s54，直至与当前道路相交的所有dem格网中心点高程值计算完成。

本发明方法基于密集点云数据，包含道路子单元高程值确定与dem格网中心点高程值计算两个关键技术环节，能够有效保证道路dem建模的高程精度与形态特征。本发明的技术特点及有益效果：

(1)本发明提出的道路dem构建方法充分顾及了道路的形态特征，摈弃了传统的dem插值方法，依据道路面在不同方向上高程变异特征设计了dem格网高程值的计算方法，能够非常好的保证dem构建结果的形态特征。

(2)本发明提出的道路dem构建方法以密集点云数据为高程信息源，能够最大限度的保证dem结果地高程精度，符合当前数据采集与建模地发展趋势。

附图说明

图1为本发明方法的流程图；

图2为道路子单元划分示意图；

图3为道路子单元中的密集高程点示意图；

图4为道路子单元中心点示意图；

图5为道路子单元中心点投影到二维坐标系下示意图；

图6为调整道路子单元中心点高程值示意图；

图7为道路子单元与拟构建dem格网单元叠加示意图；

图8为道路子单元内生成的道路中心线垂线的示意图；

图9为dem格网中心点与道路子单元空间位置关系示意图；

图10为dem格网中心点高程值计算示意图；

图11为实施例中的道路图；

图12为实施例中由道路生成的道路子单元图；

图13为实施例中本发明方法构建的道路dem；

图14为实施例中反距离加权方法构建的道路dem；

图15为实施例中本发明方法结果对应的坡度图；

图16为实施例中反距离加权方法结果对应的坡度图。

具体实施方式

以下结合附图对本发明的原理和特征进行描述，所举实例只用于解释本发明，并非用于限定本发明的范围。

本发明方法首先将道路划分为道路子单元，基于密集点云数据，计算每个子单元平均高程值，并确定道路首、道路尾和道路口处高程值；再根据设置好的道路最大起伏阈值对各个道路子单元平均高程值进行修改；接着，生成建模区域dem格网，与道路子单元叠加；最后，遍历dem格网，基于格网中心点与道路子单元的空间位置关系以及道路子单元的高程值计算dem格网的高程值。

如图1所示，本实施例基于密集点云数据构建道路dem方法的具体步骤如下：

步骤1：道路子单元划分。将道路面数据按照一定步长划分为一系列单元格，其中若存在道路相交的情况，则将相交道路作为一个子单元(图2)。

步骤2：基于密集点云数据，根据单元格内点的高程(图3)计算每一个道路子单元的高程平均值，然后将计算的平均值赋给每个道路子单元的几何中心点(图4)。

步骤3：确定道路头和道路尾的高程值。将当前处理的道路中，位于道路头和道路尾处道路子单元的高程值标记为不可编辑，即该处的高程值是确定的。

步骤4：确定道路口处高程值。以当前所处理的道路为主干道，判断该道路是否与其他道路相交，如果有相交道路，则将相交道路口形成的道路子单元高程值标记为不可编辑，即道路口处的高程值是确定的。

步骤5：将道路单元中心点从三维空间投影到二维空间中，其中道路头(a)投影到二维空间中为a1，其x坐标为0，y坐标为a的高程值，其他道路单元中心点投影到二维平面中，其x坐标为中心点到道路头中心点a的距离，y坐标为其高程值(图5)。

步骤6：连接投影后的首(a1)、尾(b1)点，即原道路对应的道路头和道路尾，由a1、b1两点确定直线a1b1的表达式。

步骤7：遍历a1、b1之间的未标记的道路子单元中心点对应的投影点，并计算每个投影点到直线a1b1的距离，并选择出距离最远的子单元中心点对应的投影点c1(图6)，同时计算出最大距离为habmax；最大距离计算公式如下：

式中：a为直线a1b1的y坐标变化值即a＝yb-ya，b为直线a1b1的x坐标变化值即b＝xb-xa，c的表达式为：c＝(yb-ya)*xb-(xb-xa)*yb。

步骤8：计算直线a1c1和c1b1的高程变化率kac和kbc，与设置好的阈值比较，若高程变化率大于该阈值则修改c1高程值，使得计算的高程变化率等于设置好的阈值，然后将原道路划分为新的两段道路，即直线a1c1和b1c1，分别以点a1与c1、b1与c1作为新道路两端点，重复步骤6-8；若高程变化率小于该阈值，则c1处高程值不变，并将点c1标记为已处理，重新执行步骤7。

步骤9：待所有投影点都标记完成，将投影点的高程值赋给原始道路对应子单元的中心点，完成道路子单元中心点高程值的更新。

步骤10：在研究区域生成格网dem中心格网点，格网尺寸小于道路子单元，并与道路子单元数据进行叠加(图7)。

步骤11：基于道路中心线，过每个道路子单元的中心点，生成道路中心线的垂线(图8)。

步骤12：依次遍历与当前处理道路叠加的dem格网中心点，若dem格网完全位于一个道路子单元内，则将子单元几何中心点高程值赋给格网中心点(图9)。

步骤13：若dem格网横跨两个道路子单元，确定所在子单元中心点(v1和v2)两侧的道路中心线的垂线，分别标记为l1和l2。过待计算高程值的dem格网中心点生成l1、l2的垂线，计算格网中心点距l1和l2的距离，标记为d1，d2(图10)。

步骤14：基于d1，d2以及v1、v2的高程值，按照线性变化原理，计算当前dem格网中心点的高程值，计算公式为：

步骤16：重复执行步骤12到步骤14，直至与当前道路相交的所有dem格网中心点高程值计算完成，并按照栅格数据规范输出计算结果。

实施例

为了验证本发明方法构建dem的效果，本实例选取了位于安徽省滁州市的滁州学院内部部分道路(图11)为例进行dem构建测试，所采用的密集高程点数据通过无人机搭载高程信息采集仪器获取。应用本发明方法的dem构建实施步骤如下：

步骤1：按照5m的步长，将面状道路数据划分成一系列单元格，其中若存在道路相交的情况，则将相交路口作为一个子单元。(图12)

步骤2：基于密集点云数据，计算每一个道路子单元中的高程平均值。

步骤3：对道路各个子单元的高程平均值进行计算调整时，设置高程变化率阈值为0.05，对超于该阈值的高程点高程值进行修改，使得计算的高程变化率等于设置好的阈值。

步骤4：生成格网dem中心格网点，并与道路子单元数据进行叠加。

步骤5：根据dem中心格网点位置与道路子单元高程值，计算各个dem中心格网点的高程值。

步骤6：构建像元大小为1米的dem(图13)，作为对比，采用传统空间插值方法(以反距离加权插值方法为例)进行道路dem插值，结果见图14。能够看出，两种方法构建的道路dem在高程变化方面基本一致，本发明方法高程变化略小于传统方法。

基于本发明结果与反距离加权方法结果生成坡度图(图15，图16)，可以看出本发明方法得到的坡度基本上在道路延伸方向上变化，符合道路形态特征，而传统方法坡度变化明显更复杂，不符合道路形态特征。

综合以上分析，可以看出本发明方法构建的道路dem很好地保持了道路的形态特征，适用于城市等区域高精度dem构建。