一种基于Fréchet距离的断裂等高线内插算法

，， ，，

(中国测绘科学研究院，北京 100830)

等高线是描述地球表面地理信息最直观的方法之一，在地图制图中有着重要的应用。等高线内插是指在已知高程的两条等高线之间，通过内插方法自动生成与已知等高线相一致的新等高线[1]，其在制图综合、三维重建等领域具有广泛的应用，因此如何准确地内插等高线具有重要的意义。

目前等高线数据的获取方式分为两种：遥感影像中获取[2-3]和实地测量高程获取。通过遥感数据提取的等高线一般都是连续且封闭的。但由于遥感数据分辨率所限，通过这种方式提取的等高线的等高距较大，不能满足需求。除此之外，还有一部分通过实地测量的符合精度的等高线，但在实际的高程测量中，由于遇到建筑物不能实地测量等情况，绘制出的等高线是断裂的，因此需要通过遥感数据提取的等高线人作内插补充。

国内外学者对等高线内插进行了广泛研究，为解决等高线内插提供了思路和算法。目前内插等高线有如下几种方法：规则网格内插、约束三角网内插、辅助线内插[4-7]等。对于马鞍型地形，文献[8]提出使用“搭桥法”的插值方法解决，即将多条同一高程的等高线用线连接起来，组成一条等高线，再用上述两条等高线内插等高线的算法插入所需要的等高线。但上述算法不适用于目前有部分断裂等高线存在的情况，由于有一部分断裂等高线的存在，在插值时需要考虑已有的等高线。对于马鞍型地形的处理，传统方法不能满足当前插值的需要。如图1 所示，若用“搭桥法”连接，则会出现与已有等高线不协调的情况。

综上，目前的插值算法不能很好地解决与已有等高线的协调关系，因此本文提出一套基于相似度的自动内插算法，通过等高线的相似度对参考等高线进行选择，再进行插值，以满足与已有等高线的协调关系。

1 相似性度量

由于空间距离无法考虑较长等高线的形状，因此引入相似性度量。线状要素可以看作是一组空间点集，线状要素的几何相似性可以看作空间点集的相似性[9]。空间点集相似运用了两种重要的空间点集距离：Hausdorff距离[10]和Fréchet距离[11]。但由于Hausdorff距离没有考虑点的顺序，局限性较大[12-13]，因此有可能形状差异很大的两个几何形状，其Hausdorff距离却很小，因此本文选用Fréchet距离作为相似度测量的方法。在计算时，真实的Fréchet距离即为线的连续Fréchet距离。但由于线要素是离散化的，无法计算线的连续Fréchet距离，因此，可以使用离散Fréchet距离对等高线相似度进行度量。离散Fréchet距离定义如下[14]：

给定一个有n个至高点的多边形链P=〈p1,p2,…,pn〉，一个沿着P的k步，分割P的至高点成为k个不相交的非空子集{Pi}i=1,…,k(1≤k≤n)，使得Pi=〈pni-1+1,…,pni〉和0=n0

给定两个多边形链A=〈a1,a2,…,an〉，B=〈b1,b2,…,bn〉，一个沿着A和B的组合步为：一个沿着A的k步{Ai}i=1，…，k和一个沿着B的k步{Bi}i=1，…，k，使得1≤i≤k，或|Ai|=1，或|Bi|=1(表示Ai，Bi恰好包含一个至高点)。

一个沿着链A和B的组合步W={(Ai，Bi)}的花费为

则链A和B间的离散Fréchet距离为

2 等高线插值算法

等高线插值首先需要对首曲线和计曲线数据作节点加密，对于有多个插值参考线的情况，需要通过相似性度量选择正确的参考等高线，最后根据参考等高线进行插值，从而对断裂的等高线进行补充。插值算法的具体流程如下：

(1) 对首曲线和计曲线数据作节点加密。

(2) 判断是否有多个参考等高线，若有，则根据Fréchet距离选择相似性高的参考等高线。

(3) 用“最近点”的方式，通过上步选取的参考等高线对断裂等高线进行插值。

2.1 节点加密

离散Fréchet距离衡量等高线的相似性时，若部分等高线的节点之间距离过大，线上的点较稀疏，在计算时会与真实Fréchet距离相差较大，不能准确地衡量相似度。图2所示有两条曲线〈a,b,c〉和〈d,e〉，它们之间连续的Fréchet距离就是从d到线段ab的距离，如图2(a)中的虚线所示。

但由于线ef中节点较稀疏，计算得到的离散Fréchet距离为d(b,e)，远大于真实的Fréchet距离。如图2(b)中，对图2(a)中线ef中节点进行加密，则离散Fréchet距离为d(b，h)，更接近真实的Fréchet距离，可以较为准确地判断相似度。当曲线上的样本点足够多，离散Fréchet距离会无限接近连续Fréchet距离，但也会增加运算量。因此需要对等高线的节点进行适当加密以保证判断的准确性。

2.2 基于相似度选择参考等高线

对于鞍部等特殊地形中等高线参考线不确定的情况，需要基于相似度对参考等高线进行选择，以保证插值的准确性。根据上文中所提到的线状要素相似度的计算，基于离散Fréchet距离选择与断裂等高线有较高相似度的计曲线作为参考线。

在对线的相似性的判断中，用全部断裂等高线作相似度判断是不准确的，因此需要对参考等高线作截断。如图3(a)所示，等高线1有A、B两个端点，则在等高线2上的截断线为ab。由于弧段AB与ab不相似，因此无法找到插值的参考等高线。

因此，选用断裂处附近的一系列点而不是整条线作为相似度判定的依据比较合适。基于离散Fréchet距离的计算方法，采用“部分”Fréchet距离计算方法，截取部分线作为相似度的判断依据，而截取的长度则是判断准确与否的关键。截取的线太长容易将不匹配的部分代入，求得的Fréchet距离变大。若截取的线太短则会失去Fréchet判断的准确性。如图3(b)所示，等高线1的部分与内部环相似，若取等高线1中AB作截取，则截取的等高线过短，容易造成截取点数量不足而判断错误，若截取部分为AC，则截取过长，形状不同造成相似度判断错误。因此需要选择合适数量的点对等高线进行截取。

本文提出的解决方法是先选取相对少的点作判断，避免选取的点过多，截取过长而导致相似度判断错误。同时，为避免线的长度过短造成的Fréchet距离判断不准确，可以在选取后对线的长度进行限制，若线的长度过短，则可以增加少量点数达到长度限制后再进行相似性的判断。

图4所示为完整等高线A、B、C，断裂等高线D。D的端点为d1和d10，取部分线d1d5在B上截取，截取后的线为b1b4，在C上截取，截取后的线为c1c2。分别计算d1d5与b1b4、c1c2之间的离散Fréchet距离。经过计算得到d1d5与b1b4的离散Fréchet距离为d(d1,b1)，d1d5与c1c2的离散Fréchet距离为d(d1,c2)。由于d(d1,b1)

2.3 等高线的内插

根据上节中得到参考等高线，进而对断裂等高线进行插值。由于断裂等高线的存在，插值时需要考虑参考等高线与内部已存在断裂等高线的方向。延伸方向不同的等高线，插值的方向不同。判断两条等高线的延伸方向是否一致的方法如下：断裂等高线B的两个断裂点为P1、P2，其中P1、P2分别对应断裂线的首、尾节点，对应在完整等高线C上距离最近的两个点Q1和Q2，计算得到等高线C上的序号Pi1，Pi2。如图5所示，可得Pi1为5，Pi2为11。

同样用此方法获得P1、P2在等高线B上临近的点P3、P4在等高线2的序号Pi3为4，Pi4为12。若Pi1=Pi3或Pi2=Pi4，则在等高线上对相应的点继续取下一个点，直到Pi1≠Pi3或Pi2≠Pi4。则若Pi1Pi4且Pi1>Pi2时，等高线B和等高线C延伸方向相同，则插值顺序为Pi2、Pi2+1、…、Pi1-1，Pi1。反之，若Pi1>Pi3、Pi2Pi3、Pi2

计算连接顺序后，即可通过点的坐标和这两条等高线的高程进行插值，得到插值后的线。具体算法如下：

由上文中获知插值顺序和序号，对序号点进行插值，最后连接已有等高线的两个端点。设这两条等高线分别为等高线C和等高线A，设上步中求得的点顺序为P1、P1+1、…、Pm，则从P1开始对等高线A上节点计算最近点R1，直到最后一个点为止。求得等高线3上的点序为R1、R2、…、Rm，并据此在等高线2、3间进行插值。连接的公式为

式中，hc为当前断裂点的高程值；h1为等高线1的高程值；h2为等高线2的高程值；xn为点n的横坐标；yn为点n的纵坐标。

3 试验与分析

将本文的修复算法搭载到中国测绘科学研究院自主研发的缩编软件平台——WJ-Ⅲ无级地图工作站上，开发环境为Microsoft Visual Studio 2010和C++编程语言，试验测试环境为单台PC机，系统版本为64位Windows 7操作系统，CPU为Intel Core i7-6700，主频3.40 GHz，内存(RAM)4 GB。

试验数据为江苏某地区的断裂等高线数据，分为计曲线和首曲线两个图层，计曲线完整，首曲线部分断裂。计曲线的基本等高距为25 m，首曲线的基本等高距为5 m。试验结果见表1。

表1 试验结果汇总

图6为试验结果，其中粗实线代表完整计曲线，细实线代表断裂首曲线，虚线代表插值后的补充等高线。

从试验结果可以看出，应用相似度的等高线插值算法，其插值的正确比率在90%以上，精确度较高。尤其是在鞍部的插值，相较于空间距离度量，使用相似性度量的方法提升了正确率。基于自动插值的结果，再辅以少量的人工作业，即可将等高线插值完整。

4 结 语

本文对目前数据的等高线插值流程进行了讨论。在准确性上，结合相似性判断，减少了鞍部等复杂地形的插值错误，并通过试验证明了方法的有效性和合理性。但在对等高线进行自动插值的过程中，还需要对部分情况进行人工处理，在今后的研究中可针对这个问题进行改进。