面状要素渐变变换利用渐变思想，使面状同名要素在由大比例尺向小比例尺变化的过程中，实现要素形状的平滑、连续渐变且最大程度保留同名要素主体特征，弱化细节特征。

面状要素渐变变换首先需要对同名的面状要素建立两个数据集，即找出同名要素的特征点，并通过数据集形式存储；其次，对两个数据集建立之间的对应关系，即确定匹配算法对特征点进行匹配；最后，以一定方式的平滑内插，在匹配的特征点间获得内插数据，形成介于始末状态间的多尺度面状要素。

矢量面状数据的渐变变换研究中特征匹配通常选择面状要素的视觉特征点^[3]或者边界弧度较大的顶点^[4]作为特征点，并通过特定算法将同名要素的特征点匹配成功。

这种取特征点的方式，图形层能较好体现出面状图形的特征。对于地理要素未能顾及空间数据的地理特征，造成体现地理特性的特征点未选取，导致渐变时主体特征的弱化；匹配的算法考虑到同名要素在大小比例尺下的两重表达存在几何细节差异，会形成两个不等长序列，故采用Latecki等^[6]提出的最优子序列双射法。对两个不等长序列构造有向无环图，求出最短路径，实现不等长序列的点匹配。该算法能较好解决特征点匹配问题，但由于面状要素边界分割序列首尾闭合，会存在起点选择的失误，造成匹配误差，需要进行算法改进。

结合地图综合原理与图形匹配算法可知，空间数据的地理特征为特征匹配, 受地图综合的影响，同名要素在大、小比例尺下的两重表达存在几何细节差异，但其主体形状结构特征相似。多边形骨架线作为二维形状特征的一维表达，必然存在类似的变化趋势。鉴于这一特征，顾及地理特征的特征匹配算法改进有两个关键点。一方面是提取同名面状要素的骨架线，对于同名要素提取骨架线，Delaunay三角网对面状要素地理特征有较好的反应，连接三角网内三角形的几何特征点，便可获得面状要素的骨架线。骨架线如何最大程度反映面状要素的特征，关键在三角网的构造方式和网内三角形的几何点的选择；另一方面，利用骨架线端点完成面状要素的匹配。在文献[6]提出的最优子序列双射法上进行优化，减少匹配误差。

基于以上讨论，本文提出一种基于骨架线匹配的面状同名地理要素渐变算法。首先，利用Delaunay三角网及相应规则提取不同比例尺下同名要素的骨架线特征；然后，对两组骨架线特征形成的特征向量进行定量描述，并利用改进最优子序方法进行全局最优匹配，得到相匹配的多边形边界特征点集合；最后，根据特征点对多边形进行剖分，在对应节点之间进行线性插值获得连续渐变的数据多尺度表达效果。

利用文献[11]算法提取面状要素的三角网，定义邻近三角形和邻近边，将三角网中的三角形分为Ⅰ类、Ⅱ类、Ⅲ类和孤立三角形, 并根据三角形的类别，提取骨架线(见图 1)。

图  1  三角形分类及骨架线提取

Figure 1.  Classification and Extraction

为便于后续的匹配，需确定骨架线的根结点并对骨架线进行结构化组织生成骨架树^[12]。仅有一条骨架线端点的骨架点称为骨架端点；三条骨架线交汇的骨架点称为骨架分支点；其他骨架点为连接点。根结点应为骨架分支点且尽可能处于面状要素的中心，离各端点的距离尽可能接近，处于骨架线最平衡位置以方便匹配的要求。据已有的文献，根结点有不同的定义方法^[13-14]，但存在出现多个根结点或选择根结点时计算复杂的问题。本文提出一种新的方法, 记录每个骨架端点到该分支点经过骨架长度的总和，其总和数值越小，表示该骨架分支点作为骨架分支点的可能性越大。选择总和最小值相应的骨架分支点，定义其为根结点。

根结点确定后，遍历骨架线上的所有骨架端点，对骨架端点进行描述可得到表示面状要素形状及其骨架特征的骨架端点特征信息。定义从根结点到骨架端点间的最短路径映射到原来骨架上所形成的一段骨架为骨架路径。为便于匹配将拓扑关系抽象为特征向量，设每个骨架端点的特征向量为F，对每个骨架端点的骨架路径进行M等分，得到M-1个等分骨架点，寻找每个等分骨架点与面状要素边界的最短距离r，作为F向量的前M-1维，每一维表示该点所通过区域的最小缓冲区，前M-1维整体表示骨架路径影响的面状要素区域；骨架路径的长度h作为第M维，表示骨架路径的影响长度。两组维度形成特征向量F =[r₁ r₂…r_M-1 h]表示形状的骨架长度及其通过形状的区域，如图 2所示。

图  2  骨架路径

Figure 2.  Path of Skeleton

全部骨架端点特征向量中第M维最大值所对应的骨架端点为起点，其中大比例尺的起点记作l₁，小比例尺记作s₁，其他的骨架端点按顺时针排序，分别表示为l_j、s_i，其中i={2, 3…m}，j={2, 3…n}, m、n为大比例尺和小比例尺骨架线端点数。分别对两个同名骨架线的每个骨架端点进行特征向量计算，得到整个骨架特征向量的集合，设大比例尺下的特征向量集合为F(l_j)=[r′_j1 r′_j2…r′_j(M-1) h′_j]，小比例尺下为F(s_i)=[r_i1 r_i2…r_i(M-1)h_i]。

由于几何表达详细程度的差异，特征点数目必然不同，从而导致特征向量集合基数不一致，为获得两特征向量集合的最优匹配，需要排除异常点，形成一对一的匹配。为达到弹性匹配的目的，引入最优子序列双射^[6]算法。该算法通过两个骨架端点所表达范围之间的相似关系进行骨架端点匹配。定义两个骨架端点的非相似度为d，d值越小则相似度越大，反之，相似度越小。算式为：

式中，d(s_i, l_j)表示小比例尺中第i个端点和大比例尺中第j个端点的相似度。

计算集合中任意骨架端点与另一个集合中所有骨架端点的非相似度，代入矩阵D₁为：

分析矩阵D₁中的数值，选择每行最小元素，它对应的行列为两个待匹配集合中最相似的骨架端点序号，这是最为简单的骨架端点匹配方法，但弹性匹配跳过端点没有受到制约，有一定不可控性。为实现最优的骨架端点的匹配，引入动态规划思想。利用动态规划思想对式(2)的矩阵做出进一步改进，设定控制因素以限定跳跃的范围，选择一组匹配方案。

为防止起点选择失误所造成的匹配误差，将目标序列延长一倍，使其循环寻找最大相似度。增加矩阵列数为2n后，得到非相似性矩阵D为：

对式(2)中的距离矩阵进行有向无环图构造，设定跳跃权重J，并根据J值计算有向无环图各边的权重：

式中, α为可自定义参数; (a, b)、(p, q)分别表示非相似矩阵D中第a行第b列的元素和第p行第q列的元素。通过比较矩阵位于第1行的第1~(n-1)列，终点位于第m行的各种情况的路径，记录每一个路径中经过的元素值和边的权重总值，选择其中总和最小值所对应的路径。通过该路径所经过的元素相对的行列值可得到两个待匹配集合中最相似的骨架端点序号，这些匹配的骨架端点形成特征点序列用于剖分同名要素边界。该方式一定程度上控制因跳过骨架端点匹配所造成的误差，且小比例尺下的骨架端点必能获得大比例尺下与其匹配的端点，匹配效果更为优化。

为了实现渐变的效果，最直接的方式就是线性插值算法^[4]：在一个对象上找到一个点使其到初始点的相对距离与另一个对象的每一个顶点到初始点的距离相等，将找到的点插入到两个对象中。定义一个对象的调节参数为g，则另一个调节参数为1-g。不同的g值变化形成不同尺度的渐变效果。根据§2.2中的匹配算法，形成的特征点序列可将始末多边形剖分形成个数相同且对应的线段。通过分段线性插值的方法可提高精确度，得到更好的渐变变换结果。首先对线段进行归一化，记录各节点在总线段的相对位置比率；其次，同名线段根据对方各节点的比率增加相应节点，使得同名线段节点数相同；最后，在对应节点之间进行线性插值获得中间序列的线段，形成介于始末状态间的多尺度面状要素，实现同名要素的渐变变换。

实验环境为DELL OptiPlex 3020MT台式机(Intel(R) Core(TM) i5-4590 CPU、8 GB内存)、Microsoft Win7系统，算法开发语言为C#、集成开发环境为VS2012，全部算法基于Arc Engine10.1开发包完成。实验数据选取三组两个不同比例尺下具有光滑边界特征的图形：湖泊数据H₁、H₂，坑塘水面数据Q₁、Q₂以及耕地数据W₁、W₂，如图 3所示。数据H₁、Q₁、W₁表示小比例尺1：S下概略的图形信息，数据H₂、Q₂、W₂表示大比例尺1：L下细节特征详细的图形。本文以第一组湖泊数据图形H₁和H₂，进行详细说明。

图  3  小、大比例尺下的数据表达

Figure 3.  Representation of Polygons at Small and Large Scales

对同名湖泊要素进行约束Delaunay三角网的构建，并按照三角形分类连接得到图形的骨架线。在实际应用中，大比例尺下的要素的骨架端点密集，且矢量数据点抖动会导致产生过多的Ⅲ类三角形，产生较小分枝。分析发现产生较小分枝的三角形为Ⅰ类三角形且三个顶点为边界上的顺序点。考虑匹配效率，对大比例尺下骨架线分枝进行删减。分析发现分枝由三个顶点在边界上的小面积Ⅰ类三角形导致，故当其面积小于预设阈值ε(ε为大比例尺下多边形面积的千分之一)时，将其忽略并删除其连接的骨架线。处理得到化简后骨架线，如图 4所示，其中左图为大比例尺的骨架线示意图，右图为小比例尺示意图。

图  4  化简后湖泊数据的骨架线

Figure 4.  Simplified Skeletons of Lake Polygons

根据§2.2和§2.3定义，确定两者根结点并计算骨架端点的特征向量。根据§2提到的算法得到匹配后的特征点序列，如图 5所示。由于匹配算法自身的设定，小比例尺地图要素骨架树的每一个端点在大比例尺上都有与之对应的特征点。

图  5  特征点匹配

Figure 5.  Matching of Feature Points

由图 5(a)可以看出，在合理情况下，点A会和点C进行匹配，但实际却是点A和点B匹配(图中绿线所示)，其原因在于B点特征向量第M维的长度更接近A点第M维长度。为避免出现此种结果，在匹配算法完成后，计算每一骨架端点与其相连的分支节点之间的角度作为端点角度，当两匹端点的角度差大于60°，且大比例尺端点的兄弟端点与小比例尺下匹配端点的端点夹角更小的情况下，以大比例尺的兄弟端点代替原有端点与小比例尺的端点进行匹配，进一步减少匹配的误差。通过检查代替后，点A成功与点C匹配。另外两组实验匹配结果如图 5(b)、5(c)所示。

匹配后的特征点序列将同名要素面状边界剖分为相对应的若干曲线段，对这些曲线段进行常规线性插值获得渐变变换效果。定义H₂的调节参数为g，则H₁的调节参数为1-g。线性插值形成中间比例尺尺度1：X，X由式(6)计算：

实验中设g分别等于0.1…0.9，间隔为0.1，小比例尺为1：50 000，大比例尺为1：10 000。实验完成，得到介于两个比例尺间的湖泊面状要素。同理得到其他两个实验数据的多尺度表达，如图 6所示。

图  6  本文算法的渐变效果图

Figure 6.  Morphing Effects of this Text Algorithm

同时，对同名湖泊数据H₁、H₂，坑塘水面数据Q₁、Q₂以及耕地数据W₁、W₂采用常规的线性内插，即§2.3中的线性内插方法，进行渐变变换，同样得到一组介于两个比例尺间的面状要素，实验效果如图 7所示。两组实验内插算法相同，差别在于是否进行边界分段内插，对比图 6和图 7可以比较两种算法在渐变整体效果上的优劣。

图  7  线性插值方法的渐变效果

Figure 7.  Morphing Effects of Linear Interpolation

选取g=0.3，0.4，0.5，0.6，0.7的中间尺度渐变图形，将局部放大进行对比，如图 8所示。图 8中矩形标注处，线性插值算法出现明显的锯齿，渐变不平稳；本文算法图形渐变的效果较好，中间一直保持光滑平稳渐变，变换结果与手工综合类似。渐变过程中，大比例尺下图形的细节被慢慢忽略，整体向概略图靠近，渐变取得良好的效果。圆圈标注处，线性插值算法出现异形刺尖，而本文算法未出现。

图  8  对比实验放大效果图

Figure 8.  Magnified Morphing Effects of the Contrast Experiment

为进一步验证本文算法的平稳渐变，以H₁、H₂为例，对两种算法的中间状态多边形进行内角计算。为方便验证，内角范围统一转化为0°~180°，定义小于10°的内角为异形刺尖，对比结果如表 1。线性插值出现五个刺尖，而本文算法只出现一个。对比分析发现，本文算法使渐变效果更为平滑，接近于手工渐变。实验证明本文算法在光滑边界的同名要素渐变方面具有可行性。

在效率方面，本文算法考虑空间特征的对应关系，比常规线性内插增加了匹配部分，具体匹配过程利用动态规划思想，以空间换取时间效率。

通过时间频度T(n)和时间复杂度O(n)进行评价。匹配大致分为三个部分：①骨架线生成，最坏情况下的时间复杂度为O(n²)；②骨架线特征向量生成，T(n)=N·D，其中N为边界结点，G为骨架端点，时间复杂度为O(n²)；③骨架匹配，利用动态规划时间复杂度为O(n²)。表 2为不同数据情况下的时间代价，其中N为多边形边界结点数，G为骨架端点数，U为骨架线生成时间，V为特征向量生成时间，P为矩阵生成和匹配时间，O为分段线性插值时间，本文算法总时间为T₁，直接线性插值时间为T₂。

根据表 2及时间复杂度，可知本文算法时间消耗主要与边界结点数量和骨架端点数量有关。通过删除微小骨架线分枝，可提高算法效率。对线段剖分分段插值，使运行时间较常规线性插值有所减少。

本文算法从空间数据的地理特征入手，通过三角网的构建，提取出更能反映地理特征的骨架线，并采用基于动态规划改进的最优子序双射算法匹配骨架端点，使端点的匹配正确率进一步提升，内插结果多边形呈现平滑、渐变特征。较以前关于渐变的研究，本文算法在特征匹配过程中考虑到空间数据的地理特征，根据这些特征进行骨架匹配和分段内插，内插结果在保留数据主体特征的基础上有效地避免面状要素的自相交问题，且渐变图形边界连续光滑，可用于空间数据的连续多尺度表达。