A Simple Adaptive DEM Void Filling Method Based on Raster Distance Transform

翻译

摘要——数字高程模型（DEM）在摄影测量和遥感领域是不可或缺的。在从多视角卫星图像生成DEM时，常常会在云覆盖或水体区域出现无值，导致空洞区域。现有的插值或融合方法在填充DEM时往往难以确保高程过渡的平滑性，而基于深度学习的方法则面临较大的内存和计算需求。为了解决这些问题，本研究提出了一种基于栅格距离变换的简单自适应DEM空洞填充方法，利用公开可用的全球DEM作为辅助数据，以重建完整的地形表面。首先，通过对空洞区域的掩模进行距离变换，创建初始权重网格。然后，通过在初始权重网格中设置适当的阈值，自适应确定DEM融合的缓冲区和权重图。最后，根据权重图将目标DEM与辅助DEM无缝融合。针对不同地形的比较实验表明，所提方法在无缝DEM空洞填充方面的有效性，根均方差（RMSE）从4.27米减少到1.47米，显著减轻了视觉上的阶梯效应，模拟实验也证明了该方法具有更高的地形保真度。

关键词——数字高程模型（DEM），栅格距离变换，空洞填充，权重网格。

解释

这段摘要介绍了数字高程模型（DEM）在遥感和摄影测量中的重要性，以及在生成过程中可能遇到的空洞问题，比如由于云或水体导致的数据缺失。现有的填充方法往往不能很好地处理这些空洞，尤其是在保证高程变化平滑方面。

为了解决这些问题，研究者们提出了一种新的方法：

基于栅格距离变换：利用这种技术，可以评估空洞区域与周围有效数据的距离，并创建一个初始的“权重网格”。
自适应确定融合区域：通过设置阈值，动态确定填充区域的缓冲区和权重，从而更精准地进行数据融合。
无缝融合：根据权重图将缺失数据的DEM与辅助的全球DEM进行融合。

通过不同地形的比较实验，该方法显著降低了高程过渡处的误差，从4.27米降至1.47米，改善了视觉效果，并提高了地形的保真度。

总体来说，这项研究为填补DEM中的空洞提供了一种有效且高效的解决方案，有助于提升地形数据的质量和应用效果。

翻译

I. 引言

数字高程模型（DEM）作为摄影测量和遥感领域的重要数字产品被广泛应用[1]。最近，大规模DEM生成利用了摄影测量技术，该技术通过密集匹配从多视角光学卫星图像中提取地形信息[2]。然而，光学卫星图像受到气象条件和地表特征多样性的显著影响。因此，图像中不可避免地包含被云层遮挡或缺乏纹理的区域，如水体和阴影，这些区域的匹配特别不稳定，导致出现空洞。

解释

这一段介绍了数字高程模型（DEM）的背景和重要性，特别是在摄影测量和遥感领域。随着技术的发展，利用多视角卫星图像生成DEM变得越来越普遍。然而，由于天气因素和地表特征的多样性，卫星图像往往会有被云遮挡或缺乏清晰纹理的区域。这些区域的高程数据匹配不稳定，容易产生缺失（空洞），这在进行高程模型构建时是一个主要挑战。

翻译

在生产过程中，填补这些空洞是至关重要的。数字高程模型（DEM）的空洞填充方法大致可分为三类：基于插值的方法、基于融合的方法和基于深度学习的方法。基于深度学习的方法能够从有效数据中表示全局信息并推导出填补空洞的趋势[3]，[4]，[5]，[6]。然而，这些方法在模型训练和预测时需要大量计算资源，可能会引入不必要的开销。基于插值的方法，如克里金插值、逆距离加权（IDW）和样条插值，通过利用相邻已知数据点的空间连续性来估计高程[7]，[8]，但这些方法仅适用于小空洞，插值高程与实际值之间的差异显著[9]，[10]。基于融合的方法通常利用公开可用的全球DEM，如航天飞机雷达地形测量（SRTM）DEM和ALOS世界3D-30米（AW3D30），作为辅助数据来填补较大的空洞[11]。消除由于数据采集时间和分辨率的变化而导致的高程突变，这是这些方法的一个关键问题，这种变化在接缝处表现为“阶梯”效应。

解释

这段文字讨论了在生成数字高程模型（DEM）时，如何填补因数据缺失而产生的空洞。空洞填充方法主要分为三类：

基于插值的方法：利用已知数据点之间的空间关系来估算空洞中的高程，例如克里金插值和逆距离加权。这些方法适合小空洞，但在处理较大空洞时效果不佳，插值结果与真实值相差较大。
基于融合的方法：利用现有的全球DEM（如SRTM和AW3D30）作为辅助数据来填补大空洞。这种方法可能会在数据接缝处出现“阶梯”效应，即高程突然变化，这是由于不同数据源的采集时间和分辨率差异引起的。
基于深度学习的方法：这些方法可以有效提取全局信息，填补空洞，但需要大量计算资源，可能导致额外开销。

总体来说，选择合适的填补方法对于提高DEM的完整性和准确性至关重要。

翻译

为了应对传统空洞填充过程中接缝处的阶梯效应，最简单的策略是填充和羽化，但这可能导致地形信息的丢失。增量表面填充（DSF）方法[12]，[13]及其改进版本，如基于TIN的增量表面方法（TIN-DSF）[14]和层次增量表面插值方法[15]，显著减少了不同DEM之间的系统性垂直差异。然而，由于空洞的复杂几何形状和分辨率差异，这些方法无法处理接缝。基于多尺度分解（MSD）[16]和变分模型[17]的DEM融合方法在频域中操作，未考虑空洞的空间性和几何形状，导致无缝整合效果不佳。目前，权重融合方法在DEM空洞填充中被广泛使用[18]。然而，权重融合区域的选择和权重的确定与空洞的几何形状相关，使得这些方法的通用性和适用性较低。

解释

这段文字讨论了填补DEM空洞时常见的阶梯效应及其解决方案。虽然简单的填充和羽化方法可以缓解阶梯效应，但可能会导致地形信息的损失。

一些更先进的方法，如增量表面填充（DSF）及其改进版本（如TIN-DSF），能够有效减少不同DEM之间的高程差异，但由于空洞的复杂形状和不同分辨率，它们在处理接缝问题上仍然存在困难。

此外，基于多尺度分解和变分模型的融合方法主要在频域中工作，但没有考虑空洞的空间特性和几何形状，这使得无缝整合效果不理想。

目前，权重融合方法被广泛使用，但这些方法在选择权重区域和确定权重时受到空洞几何形状的影响，导致其适用性有限。因此，寻找更通用和有效的填补策略仍然是一个挑战。

翻译

因此，本研究提出了一种基于栅格距离变换的简单自适应DEM空洞填充方法。栅格距离变换是地图代数中的一种基本理论，广泛应用于图像分割[19]和水文分析[20]。本研究提议利用这一理论融合多源栅格DEM。所提方法根据距离变换自适应地确定栅格数据融合的缓冲区和权重。利用由无云图像或同一区域的开源DEM构建的DEM作为辅助数据，所提方法能够无缝填补DEM空洞，显著减轻接缝处的阶梯效应，同时保持原始地形信息。该方法通过将向量问题重新表述为栅格问题，消除了显式跟踪拼接线以确定融合权重的必要性。其数学基础是地图代数中简单统一的操作，确保了低复杂性和高计算效率，同时显著提高了通用性。本研究的结构如下：第二节介绍所提方法的详细实现；第三节展示与几种代表性的DEM空洞填充方法的比较实验结果；第四节总结本研究。

解释

本段文字介绍了一种新的DEM空洞填充方法，利用栅格距离变换来处理多源数据。该方法可以自适应地确定填补空洞所需的区域和权重，使用无云图像或开源DEM作为参考数据，从而减少填补过程中可能出现的阶梯效应，保持地形的真实样貌。

重要的是，这种方法避免了复杂的拼接线跟踪，简化了计算过程，提高了效率和适用性。文中还概述了研究的结构，包括实现细节、实验结果和总结部分。

翻译

II. 方法论

基于栅格距离变换的DEM空洞填充方法由三个关键组件组成：1）权重网格的创建和更新；2）缓冲区设置和融合权重计算；3）加权融合。

通俗解释

本节介绍了新提出的DEM空洞填充方法的三个主要部分：

权重网格的创建和更新：首先建立一个用于填补空洞的权重网格，并在填充过程中不断更新这些权重。
缓冲区设置和融合权重计算：根据距离变换确定缓冲区，并计算填补空洞所需的融合权重。
加权融合：根据计算得到的权重，将辅助数据与原始DEM进行融合，实现空洞的填补。

这些步骤共同确保了填补过程的有效性和精确性。

翻译

A. 权重网格的创建和更新

为了自适应地确定后续融合的权重，基础步骤是根据目标DEM中的空洞创建一个新的权重网格 ( p )。初始化时，目标DEM中对应空洞的网格单元被赋予零值，表示缺少数据，而其他单元则接收一个大的占位值 ( D_0 )，表示它们与空洞区域的初始“距离”。对于目标DEM中具有高程值 ( z_{i,j} ) 的网格单元 ( (i, j) )，权重网格中该单元的初始值 ( p_{i,j} ) 被给定为：

$p_{i,j} = \begin{cases} 0.0, & z_{i,j} = \text{无值} \\ D_0, & \text{否则} \end{cases}$

其中 ( D_0 ) 是DEM中任意两个单元之间的最大距离，网格DEM的对角线长度可以选作 ( D_0 ) 的值。

通俗解释

这一部分描述了如何创建和更新权重网格，以便为后续的空洞填充提供必要的权重。在初始化过程中，空洞区域的网格单元被标记为零，表示没有数据，而其他单元则被赋予一个较大的值（( D_0 )），表示它们与空洞的距离。

具体来说，权重网格中每个单元的初始值根据目标DEM的高程值来设定。如果该单元没有高程值（即是空洞），权重为零；否则，权重则为 ( D_0 )，这个值可以代表整个DEM网格中可能的最大距离。这一过程为后续的融合步骤打下基础。

翻译

在初始化后，接下来进行一个顺序扫描和更新过程，从网格的左下角开始，向上和向右移动。在此步骤中，更新每个单元的权重，考虑到左侧和下方的八个相邻单元及当前单元。根据相邻单元与中心单元之间的欧几里得距离，重新计算每个单元的权重 ( p’_{i,j} )，公式如下：

$p'_{i,j} = \min \left( p_{i,j}, p_{i,j-1} + 1, p_{i-1,j-2} + \sqrt{5}, p_{i-1,j-1} + \sqrt{2}, p_{i-2,j-1} + \sqrt{5}, p_{i-1,j} + 1, p_{i-2,j+1} + \sqrt{5}, p_{i-1,j+1} + \sqrt{2}, p_{i-1,j+2} + \sqrt{5} \right) .$

通俗解释

这一部分描述了如何更新权重网格中的每个单元的权重值。具体步骤如下：

顺序扫描：从网格的左下角开始，逐步向上和向右移动，依次处理每个单元。
考虑邻接单元：在更新权重时，会考虑当前单元及其相邻的八个单元的权重。这个过程可以确保新权重是基于周围单元的影响。
计算新权重：使用给定的公式 ( p’_{i,j} )，计算当前单元的新权重。这个计算过程使用了相邻单元的权重值，并根据它们与当前单元的距离进行了调整：
- 选择当前单元的原始权重 ( p_{i,j} ) 和相邻单元权重的组合，取最小值来更新 ( p’_{i,j} )。
- 其中，权重的调整使用了不同的距离系数，例如，对于相邻的上下单元加上1，对于对角线单元加上 ( \sqrt{2} ) 等，这些值是根据欧几里得距离的计算而来。

表达式

$p'_{i,j} = \min(p_{i,j}, p_{i,j-1} + 1, p_{i-1,j-2} + \sqrt{5}, p_{i-1,j-1} + \sqrt{2}, p_{i-2,j-1} + \sqrt{5}, p_{i-1,j} + 1, p_{i-2,j+1} + \sqrt{5}, p_{i-1,j+1} + \sqrt{2}, p_{i-1,j+2} + \sqrt{5})$

用于计算网格中每个单元的新权重 ( p’_{i,j} )。这个公式的含义如下：

解释

最小值选择：该公式的核心是取一组权重中的最小值。通过比较当前单元及其相邻单元的权重，找到最佳的权重值。
当前单元的权重：
- ( p_{i,j} )：当前单元的原始权重。
相邻单元的权重：
- ( p_{i,j-1} + 1 )：左边相邻单元的权重加1，表示与当前单元的距离为1。
- ( p_{i-1,j-2} + \sqrt{5} )：上方偏左的相邻单元权重加上 ( \sqrt{5} )，表示该单元与当前单元的对角距离。
- ( p_{i-1,j-1} + \sqrt{2} )：上方偏右的相邻单元权重加上 ( \sqrt{2} )，表示对角距离。
- ( p_{i-2,j-1} + \sqrt{5} )：上方的相邻单元权重加上 ( \sqrt{5} )。
- ( p_{i-1,j} + 1 )：上方相邻单元权重加1。
- ( p_{i-2,j+1} + \sqrt{5} )：上方偏右的相邻单元权重加上 ( \sqrt{5} )。
- ( p_{i-1,j+1} + \sqrt{2} )：下方偏右的相邻单元权重加上 ( \sqrt{2} )。
- ( p_{i-1,j+2} + \sqrt{5} )：下方偏右的相邻单元权重加上 ( \sqrt{5} )。

总结

这个公式通过考虑当前单元及其八个相邻单元的权重，能够有效地更新当前单元的权重 ( p’_{i,j} )。最小值的选择确保了新权重反映了与周围单元的关系，有助于在空洞填充时保持地形的连续性和一致性。

在更新之前和之后的权重网格如图2所示，其中黑色区域表示目标数字高程模型（DEM）中的空缺，而白色区域表示有效值。通过距离变换，较小的权重会逐渐从空缺区域向外扩展，形成图2(b)中的灰色区域。在这个区域中，权重值从0开始，逐渐增加到最大值 ( D_0 )，并且是从空缺的边界开始计算的。

解释

黑色区域（空缺）：这些是数据缺失的部分，通常是由于云层或水体等因素造成的。
白色区域（有效值）：这些区域包含有效的高程数据，表示可以用来填充空缺的参考数据。
灰色区域：这是通过距离变换生成的区域，权重值从空缺边界向外逐渐增大。这种设计使得离空缺区域越近的单元权重越小，离空缺越远的单元权重越大。

重要性

这种权重分布方式有助于在填充过程中保持更高的地形一致性。由于权重值的渐变，算法可以更自然地将有效数据与空缺区域结合，减少因填充而导致的高程不连续性，进而降低视觉上的“台阶效应”。这种方法使得填充结果更加平滑，同时保留原有地形的特征。

B. 缓冲区设置与融合权重图计算

缓冲区的划定在高程融合中至关重要，以确保填充空缺与完整区域之间的平滑过渡。设定一个阈值 ( T ) 来确定更新后的权重网格中的缓冲区，其中条件 ( 0 < p’'_{i,j} \leq T ) 成立，如图3(a)中的蓝色区域所示。

基于此，每个目标DEM单元的融合权重 ( p^f_{i,j} ) 可以通过以下公式计算：

$p^f_{i,j} = \frac{p''_{i,j}}{T}$

解释

缓冲区的重要性：缓冲区用于创建填充区域与有效区域之间的过渡，使得高程变化更平滑，减少突兀的高程跳跃。
阈值 ( T )：这个阈值决定了缓冲区的范围。较大的 ( T ) 值意味着更平滑的过渡，类似于正射影像拼接中的羽化距离。
设置推荐：根据文献[21]的建议，羽化距离应大于20个像素，因此在本研究中，选择将 ( T ) 设置为32，以平衡过渡的平滑性和地形保真度。

总结

通过设定合适的缓冲区和融合权重计算，研究能够有效地融合空缺区域与有效区域，增强DEM的整体一致性和准确性。这样的设计在实际应用中能够显著改善地形重建的质量。

融合高程公式解释

公式：

$Z^{fusion}_{i,j} = p^f_{i,j} \times Z^{tar}_{i,j} + (1 - p^f_{i,j}) \times Z^{aux}_{i,j}$

通俗解释

什么是 (Z^{fusion}_{i,j})：
- 这是我们要得到的最终高程值，用于填补DEM中的空缺。
组成部分：
- (p^f_{i,j})：这是一个权重，代表目标高程（即原始DEM高程）在最终计算中的重要性。这个值在0到1之间，值越高，说明我们越依赖目标DEM的高程。
- (Z^{tar}_{i,j})：这是目标DEM中的原始高程值，也就是我们想要填补的区域的实际高程。
- (Z^{aux}_{i,j})：这是辅助DEM中的高程值，用于提供参考。
计算过程：
- 加权平均：我们通过权重 (p^f_{i,j}) 和 (1 - p^f_{i,j}) 来平衡这两个高程值。
  - 如果 (p^f_{i,j}) 很高（接近1），最终的高程值会更倾向于目标DEM的高程。
  - 如果 (p^f_{i,j}) 很低（接近0），最终的高程值会更多地依赖于辅助DEM的高程。