1 引    言

空间插值常用于将离散点数据转换为连续曲面，以便更为直观地考察数据要素的空间分布模式。在等值线绘图中，离散点空间插值是绘制等值线的关键基础，插值方法选取不仅影响要素值的空间分布形式，而且还影响着等值线绘制的视觉效果^[1,2,3,4,5]。现今多被采用的空间插值方法有反距离权重法、克里金法、样条法、三角剖分法、多项式法、趋势面法等。反距离权重算法简便、普适性强，但仅考虑插值点与样本点的空间距离,没有包含方向相关^[6,7]；克里金以空间统计学作为理论基础，使用变异函数度量样点间的空间相关性，可对插值误差进行逐点理论估计，但其变异函数需人为经验选定，且当变异函数为多个组合时，其计算量剧增^[8,9]；样条法更为适合高密度的样点内插，而对于稀疏有限的样点该方法多不适合^[10,11]；多项式和趋势面的空间插值更多依赖于插值要素固有的空间分布趋势，若样点布局遵从插值要素的空间变化，则可得到理想的插值效果，反之则难以保证插值的合理性^[12,13]；基于三角剖分的插值方法可获得高精度的插值曲面，在地表特征表达和地形等高线的绘制中经常使用^[14,15,16]，但该方法也存在类似于克里金的插值缺陷，即算法复杂计算量大，且其插值曲面的等值线棱角过于明显。简而言之，各种空间插值方式均有优缺，并不存在普适最优的插值方法。

在众多空间插值方法中，反距离权重因其插值原理易于理解、算法简单便于实现，且插值后能够保留原样点真值，常被作为离散点空间分析的传统方法之一。但反距离权重未考虑插值点与样本点之间的方向关系，当样点多集中于插值点某个方向时，该方向的样点总权重就会大大增加，而其他方向上的样点因权重过小而常被忽视，空间插值缺少方向均衡性，其插值合理性亦受到质疑^[17,18]。另一方面，反距离权重的等值线绘图通常有一弊病，即图中常会有孤立圆区存在，且易在样本极值点周围产生同心圆区(“牛眼”现象)^{[19,20,21,22]}，这严重影响了等值线图的审美效果。在大多情况下，为追求绘图美观不得不摒弃反距离权重而选用其他插值方式，不惜以数据失真来换取整幅等值线图的协调美观。本文旨在对反距离权重插值算法进行优化改进，减少或消除反距离权重等值线绘图中不合理的孤立圆和同心圆现象的存在。通过增加一个可反映插值样点方位的调和权重系数K，克服反距离权重插值中仅考虑插值点与样本点之间距离而未考虑它们之间方位的插值缺陷，由此提出了调和反距离权重(adjusted inverse distance weighted，AIDW)空间插值方式。 2 调和反距离权重插值 2.1 反距离权重算法

反距离权重(inverse distance weighted，IDW)算法于1968年由Shepard提出，1985年Watson 等将其应用于空间插值的等值线绘制^[23,24]，继而IDW算法被广泛应用于各行业领域的空间分析与制图^[25]。IDW算法公式通常被表示为

点击浏览公式

从IDW算法(式(1))可看出，插值点的估算值仅与插值样本距插值点的远近相关，并未考虑样本点的方位性(即样本点被表示为各向同性)。IDW插值的基本假设是样点在插值区呈均匀分布，当样点在各方向较均匀分布时，该插值算法十分可靠。然而在众多情况下，样点在各向分布并非均匀，甚至会出现样点集中于某一方向的现象，违背了基本假设，其插值合理性就难被保证。针对IDW这一插值局限，本文提出了调和反距离权重(AIDW)插值算法。

2.2 调和反距离权重算法 2.2.1 基本假设

AIDW算法与IDW基本类似，关键在于增加了可反映插值点与样本点方位关系的调和权重系数K，其基本假设是：距插值点近的样本点，对其后方的样本点有遮蔽效应，当两样本点与插值点的连线夹角α<360°/n(n为插值搜索邻域内的样点个数)时，遮蔽效应存在，当α≥360°/n时，遮蔽效应消失。在AIDW插值过程中，受遮蔽影响的样本点，其插值权重将被削弱，削弱的程度取决于该样点K值的大小(参见2.2.2节)。

按上述假设：图 1(a)所示的5个样点在方向上均匀地分布在插值点(中心点)周围，任意两样点与插值点的连线夹角均大于或等于72°(即α≥360°/5)，即认为该5个样点间相互不存在遮蔽效应；在图 1(c)中，任意两样点与插值点的连线夹角均小于72°，即认为距插值点的近样点，对其后的样点均具有遮蔽效应；在大多情况下，样点在插值点周围的分布应类似图 1(b)，既不像图 1(a)均匀分布，也不像图 1(c)集中分布。图 1(b)中Z₁、Z₃对任一样点均无遮蔽，Z₂对Z₄、Z₅有遮蔽，Z₄对Z₅也有遮蔽。

图 1 插值搜索邻域中样点散布示例(5个样点) Fig. 1 Examples of sample points scattering in searching neighborhood(five points)

图选项

将IDW传统的算法思想与本文的基本假设结合，提出了AIDW算法，该算法可被表达为

点击浏览公式

点击浏览公式

图 2 sinpθ随样点位置的变化示例 Fig. 2 Examples of sinpθ value changing with sample point moving

图选项

图 2示例解释了样点Z_i所受遮蔽效应与其位置的关系以及sin^pθ对这种关系的表达能力。在Z_i以固定a角逐渐移近插值点Z_o的过程中(图 2(a))，θ角和sin^pθ值在不断增大，Z_j对Z_i的遮蔽效应逐步降低；当Z_i移动到距插值点的距离等于Z_j距插值点距离时(即两样点已无先后，Z_j、Z_i与插值点构成了等腰三角形)，θ角为90°，sin^pθ=1，此时表明Z_j已不会对Z_i产生遮蔽影响。另一种情况，当Z_i保持与插值点的距离不变，但逐渐向Z_j背后靠近时(图 2(b))，显然样点Z_j对Z_i的遮蔽效应在逐步增大，这一过程中θ角和sin^pθ值也在不断变小；当Z_i移动到Z_j与插值点连线的正后方，此时θ角变为0°，sin^pθ=0，表明Z_j已完全遮蔽了Z_i对插值点的影响。

若将样点Z_i的变化位置看作是其他样点的分布，则2解释了插值样点在两种特例分布下θ和sin^pθ的极端变化情况。即两样点与插值点的距离趋近于等同，θ→90°，sin^pθ→1；某样点位置趋向于另一样点(与插值点连线)正后方，θ→0°，sin^pθ→0。通常情况下，有遮蔽效应的两样点，它们的θ角多介于0°～90°，sin^pθ值多介于0～1。

2.3 AIDW的特点分析

由于AIDW综合考虑了样点距离和样点方位对插值点的共同影响，因而其插值过程(3)和插值结果(4)较IDW更趋于合理。从图 3(a)中可看出，插值样点在各向上的分布，较难满足样点均匀布局这一IDW插值基本假设，此种情况下IDW的插值可靠性难以保证。AIDW插值考虑了样点间的方位影响，通过计算各样点间的遮蔽效应，削减受遮蔽样点的插值权重。按AIDW算法，在图 3(b)中因Z₁对Z₆、Z₃对Z₇和Z₈、Z₄对Z₇有遮蔽影响，这些受遮蔽样点的插值权重被削减，Z₁₀、Z₁₁、Z₁₂分别被Z₄、Z₃、Z₇完全遮蔽，它们的插值权重降至为0。依照式(2)和式(3)，最终插值点估算值的计算式为

点击浏览公式

图 3 IDW和AIDW两种插值的样点处理方式 Fig. 3 Treating mode of sample points of the IDW and AIDW interpolations

图选项

式中，Z为插值点(中心点)估算值;Z₁—Z₉为样点观测值；d₁—d₉为样点与插值点的欧氏距离；p是幂指数；θ角如图 3(b)所示。由此看来，从某种意义上可以说AIDW算法可使插值样点的分布在各向上更趋于均匀。

为了便于理解AIDW和IDW两种插值在等值线绘制中的各自表现，本文用一个简单特例对其加以表述。在图 4中紫色样点与插值点(黑色点)呈等距离分布，距离为1；浅橙色样点与插值点也呈等距离分布，距离为2，且浅橙色样点位于插值点与紫色样点连线的正后方，各样点的观测值如图 4所示。按IDW算法(p=1)，插值点的估算值为50；按AIDW(p=1)算法，插值点的估算值为20(浅橙色样点完全被遮蔽，仅需考虑紫色样点插值权重)。从样点数据值的分布可看出(图 4，等值线为30和40)，很大程度上AIDW要比IDW的插值结果更趋于合理。IDW算法使低值点中出现了高值，其等值线绘制也就增大了“牛眼”的出现几率；而AIDW插值结果与临近数据值接近，很好地避免了“牛眼”现象(图 4(b))。有关AIDW和IDW的插值误差及其等值线绘制的实例对比分析，参见第4节

图 4 IDW和AIDW插值方式的等值线绘制特例 Fig. 4 Typical examples of contour plotting for the IDW and AIDW interpolation methods

图选项

AIDW的插值程序可分为插值前准备和插值计算两个过程。插值前准备主要是用于搜索合适的插值样点，并为下一步的插值计算提供d_i和f_ij值；插值计算过程主要是求算反映样点遮蔽程度的sin^pθ_ij值，并结合d_i、Z_i值求算插值点的Z值(图 5)。插值搜索邻域的大小以格点数(k×k)表示，m是搜索邻域内的样点数，n是插值所需的样点数，d、f分别为样点与插值点的欧氏距离和两样点间的欧氏距离，i、j、u、v 均为插值样点的序号(i=1，2,…,n；有1≤j≤i-1、uij、θ_ij、Z_i与式(4)中的定义相同。

图 5 AIDW的插值程序设计流程(n>3) Fig. 5 Flow chart of the AIDW interpolation program (n>3)

图选项

sin^pθ_ij的求算涉及到两次运用三角余弦定理，一次是判别条件cosα_ij的求算，另一次是cosθ_ij的求算。为减少程序运算过程，若插值搜索邻域固定样点数的设定大于3个(即n>3)，则两样点是否存在遮蔽的判别条件，程序直接使用cosα_ij与cos(360°/n)对比，而不使用α_ij 角度与360°/n对比，以便去除对α_ij 角度的求解。另一方面，当出现某样点i被另一样点j完全遮蔽(即cosα_ij=1)时，该样点i与其他样点的遮蔽效应，程序采取了判别赋值方式。因已有sin^pθ_ij=0，故i样点与其他样点的遮蔽关系可不再考虑，方便起见可直接赋值sin^pθ_ui=1和sin^pθ_iv=1(usin^pθ求算。 4 AIDW空间插值实例分析

本文对浙江省69个气象台站的某年降水量观测数据，分别利用IDW和AIDW插值算法进行空间插值，进而绘制降水量等值线图。插值实例中，IDW和AIDW算法的p取值为2，插值搜索邻域的固定样点数为5，插值范围为浙江全省，插值空间分辨率即格点大小设置为1km×1km。IDW插值及等值线绘制由ArcGIS软件的ArcMap 9.3完成，AIDW插值和等值线绘制由文中的插值程序及文献[27—28]中的等值线绘制方式完成。

浙江省年降水量IDW和AIDW两种空间插值的误差分析见表 1。AIDW插值的平均误差(ME)、平均绝对误差(MAE)、平均相对误差(MRE)以及均方根误差(RMSE)，均小于IDW插值方式，尤其是ME和RMSE这两种误差更是明显低于IDW方式，这表明了AIDW相比于IDW，其获得的空间插值结果更趋于合理。

从图 6可看出(等值线基值为1000mm雨量，间隔为50mm雨量)：IDW和AIDW两者的等值线图总体分布趋势基本一致(图 6(a))，但在局部区域却存在显著差异。区域放大图显示：在浙西开化、浙西北安吉和长兴、浙东南温岭、乐清、温州、瑞安及平阳等区，IDW绘图中出现的同心圆或孤立圆现象(图 6(b))，在AIDW绘图中(图 6(c))得以消除。

图 6 IDW和AIDW插值方式的等值线绘制实例对比 Fig. 6 Comparison of contour plotting cases of the IDW and AIDW interpolation methods

图选项

在浙西开化区，随插值点距开化站(高雨量站)的距离变远，IDW插值的开化站权重会逐渐减小，而常山、江山、衢州和淳安等低雨量站的权重相应增加，故越远离开化站，其降水插值结果越小，从而在开化站周围形成了同心圆等值线(图 6(b))。在AIDW插值过程中，由于考虑了站点方位遮蔽效应，当对开化站左侧区域插值时，开化站对周围其他站产生了不同程度的遮蔽，且其他站之间也存在遮蔽，由此常山、江山、衢州、淳安等低雨量站的权重被大幅削减，故开化站左侧区的插值未被低雨量站越拉越低，从而避免了同 心圆等值线出现(图 6(c))。与开化情况相类似的还有浙西北的安吉和浙东南的温岭。在浙西北的长兴区，IDW的等值线图中有一个孤立圆(图 6(b))，这是因为湖州和德清两个高雨量站的存在，使得长兴站周围的插值，随距长兴站越远，其插值结果越高，由此等值线绘制也不得不形成一个孤立圆区；而在AIDW插值过程中，当对长兴站左侧区插值时，湖州和德清两个高雨量站因受长兴站的遮蔽以及湖州、安吉两站对德清站的遮蔽，湖州和德清两站的插值权重(尤其是德清站)明显被减小，因此长兴站左侧区的插值结果未被明显拉高，等值线绘制时此区也就与其他区域合并在了一起，并未出现孤立圆现象(图 6(c))。与长兴情况相类似的还有浙西南乐清、温州、瑞安及平阳等区。

值得注意的是：AIDW插值算法并未消除IDW绘图中浙西龙游的孤立圆以及浙西北德清和杭州的同心圆现象(如图 6(b)和图 6(c))。仔细查看可发现，这些站周围的各站点几乎呈均匀布局(上下左右各1个)，并不像开化和长兴那样集中于一侧，这也意味着AIDW插值算法即可保留IDW插值思想，又可减小样点分布不均IDW的插值缺陷。简而言之，AIDW相对于IDW，其插值结果更趋于合理，等值线绘制也更接近于人工作业，如开化、安吉、长兴、温岭、乐清、温州、瑞安及平阳等区的等值线走向判识。

5 结 论

AIDW算法依据样点与插值点的连线夹角a判别样点间是否有遮蔽效应，采用两样点连线与其中线(过插值点)夹角的sin^pθ值表征样点受遮蔽程度，以方位权重系数调和距离权重，插值过程考虑样点距离和样点方位共同影响，其空间插值结果更趋于合理。基于AIDW空间插值曲面的等值线绘制，可有效减少IDW等值线绘制中孤立圆和同心圆的出现几率，使等值线绘图更具科学美观性。

因AIDW算法涉及cosα和cosθ值求算，故其空间插值运算速度较IDW插值略慢。在浙江省69气象站点的年降水量空间插值实例中，AIDW和IDW插值运算时间均不足1s，但当插值的空间分辨率由1km×1km提高到500m×500m，插值搜索邻域固定样点数设为15个时，浙江全省的AIDW插值运算时间为11s，而IDW插值运算时间为6s。随插值格点总数以及插值搜索邻域固定样点数的增加，AIDW插值运算速度会略慢于IDW插值。