电离层效应是卫星导航系统最严重的误差源之一，高精度电离层延迟信息对提高单频用户时延修正精度和双频/多频用户精密定位收敛速度具有重要作用^[1-4]。获取整个空间高精度电离层延迟信息一般需要通过插值方法计算。电离层总电子含量（total electron content，TEC）插值是指基于区域内有限个电离层穿刺点（ionospheric pierce point，IPP）的TEC值，采用一定的插值算法获取区域内任意未知点的电离层TEC，它是电离层研究的重要内容^{[1, 5]}。利用全球卫星导航系统（global navigation satellite system, GNSS）双频电离层TEC观测量构建格网电离层、缺失格网点数据补插、大尺度格网数据的区域加密等都离不开电离层插值^[5-7]。反距离加权（inverse distance weighting，IDW）插值方法是一种应用广泛、原理简单、运算量小的空间插值方法，已广泛应用气象、水文、地质、地磁、重力等领域^[8-10]。但是，普通IDW插值方法在电离层TEC插值中，由于精度较低导致应用受限。本文通过分析IDW插值方法的插值参数，构造一种考虑经纬度方向异性的IDW插值方法，提高这种简单实用方法的电离层TEC插值精度。

插值参数是构成插值方法的基本元素，插值精度取决于插值参数对要素空间变异性和相关性的反映^[8]。IDW插值方法的插值参数包括插值范围和插值核函数。其中，插值范围是以待插值点为中心、一定半径的邻域范围；插值核函数是反映两点间空间关系及其对待插值点贡献的函数^[11]。本文考虑到电离层TEC分布的经纬度方向异性，引入经纬度方向异性调节因子r，设计了包含等权在内的6种电离层距离计算方案，选用一组全球电离层格网（global ionospheric map, GIM）产品，分析不同方案的插值精度，获取最优方案的IDW插值。之后，针对南北纬度20°电离层活动剧烈区域，采用连续12 a“两分两至”日前后GIM数据分析最优方案的插值精度提升情况，验证电离层TEC的经纬度方向异性特点及本文改进IDW方案的正确性；采用2014年太阳活动高年“两分两至”日数据，分析IDW插值精度与地方时的相关性。

1 电离层经纬度方向异性的IDW插值函数 1.1 基本原理

IDW插值方法也称为距离乘方倒数法，基本原理是利用区域内有限个电离层IPP Z(l_i，b_i)(i = 1，2…n)的电离层观测值T_i，求解区域内任意未知点Z(l₀，b₀)的电离层T₀，公式为：

$ \left\{ {\begin{array}{*{20}{l}} {{T_0} = \frac{{\sum\limits_{i = 1}^n {\frac{1}{{d_i^\mu }}} {T_i}}}{{\sum\limits_{i = 1}^n {\frac{1}{{d_i^\mu }}} }}}\\ {{d_i} = \sqrt {{{\left( {{l_i} - {l_0}} \right)}^2} + {{\left( {{b_i} - {b_0}} \right)}^2}} } \end{array}} \right. $

(1)

式中，l_i和b_i为IPP经纬度；d_i为采样点和待插值点间距离；μ为衰减指数。

IDW插值方法忽略了电离层TEC在经纬度方向的变化梯度。文献[12]表明，中纬地区磁平静期间，纬度方向梯度大于经度方向，磁暴期间梯度差异更加明显。因此，考虑电离层TEC经纬度方向梯度，本文对采样点与待插值点间距离函数d进行修正，增加r，公式为：

$ {{d}_{i}}=\sqrt{{{({{l}_{i}}-{{l}_{0}})}^{2}}+r{{({{b}_{i}}-{{b}_{0}})}^{2}}} $

(2)

IDW插值的基础是采样点和待插值点间的空间相关性，如果插值区域各向均匀，则相关性只与距离有关。由于电离层TEC经纬度方向变化梯度不同，纬度方向电离层TEC变化大，空间相关性低；而经度方向电离层TEC变化小，空间相关性高。对于与待插值点相同距离的样本点，显然当调节因子r = 1时，即为普通IDW插值；当r > 1时，纬度方向距离增大，也即该方向贡献减小；当r < 1时，纬度方向距离减小，也即该方向贡献增大。以与待插值点距离为10°的所有样本点为例，分析取值1.5和0.5时距离插值因子1/d的变化，如图 1所示。

图 1中半圆蓝线为与待插值点（0, 0）距离为10°的样本点分布位置，红色星号为r取值1.5时距离插值因子，对比r取值1.0（黑色加号），纬度方向样本点对待插值点的贡献随距离增大而减小（下凹）；粉色倒立三角形为r取值0.5时距离插值因子，对比r取值1.0，纬度方向样本点对待插值点的贡献随距离增大而增大（上凸）。当r ≠ 1时，已知点对待插值点贡献由等权直线分布变为非等权分布。

1.2 电离层IDW插值参数优化选取

由式（1）可知，插值参数d是IDW插值方法的基本元素，决定着插值精度。通过优化IDW插值参数有望改善其插值效果。

插值范围决定着采样点数量、分布，进而影响插值精度。文献[13]认为GNSS单站电离层IPP的覆盖范围约为30°×30°，选择10°×15°和20°×20°插值范围进行插值分析，如图 2所示。

由插值原理可知，插值核函数主要包括距离函数d和距离衰减指数μ两部分。

1）距离函数d。主要考虑经纬度方向异性调节因子r，设l、b为经纬度均值，T为电离层TEC均值。设计以下5种计算方案^{[5, 14-15]}。

（1）取r=1。

（2）考虑电离层TEC在经纬度方向变化梯度，计算公式为：

$ r=\sum\limits_{{}}^{i\ne j}{\frac{\left( {{l}_{i}}-{{l}_{j}} \right)}{\left( {{T}_{i}}-{{T}_{j}} \right)}/}\sum\limits_{{}}^{i\ne j}{\frac{\left( {{b}_{i}}-{{b}_{j}} \right)}{\left( {{T}_{i}}-{{T}_{j}} \right)}} $

(3)

（3）考虑电离层TEC与经纬度相关系数，计算公式为：

$ r=\frac{\sum\limits_{i=1}^{n}{\left( {{l}_{i}}-\bar{l} \right)\left( {{T}_{i}}-\bar{T} \right)}}{\sqrt{\sum\limits_{i=1}^{n}{{{({{l}_{i}}-\bar{l})}^{2}}}}\sqrt{\sum\limits_{i=1}^{n}{{{({{T}_{i}}-\bar{T})}^{2}}}}}/\frac{\sum\limits_{i=1}^{n}{\left( {{b}_{i}}-\bar{b} \right)\left( {{T}_{i}}-\bar{T} \right)}}{\sqrt{\sum\limits_{i=1}^{n}{{{({{b}_{i}}-\bar{b})}^{2}}}}\sqrt{\sum\limits_{i=1}^{n}{{{({{T}_{i}}-\bar{T})}^{2}}}}} $

(4)

（4）考虑电离层TEC与待插值点地理纬度φ的关系，计算公式为：

$ r = {\rm{sin}}\varphi $

(5)

（5）参考Kriging插值中电离层距离概念，取经验值r = 0.2。

2）距离衰减指数μ。随着μ增大，样本点对待插值点影响敏感程度增加，表现为距离近样本点影响显著增加，反之则显著降低，一般不宜采用过大或过小值。文献[16]认为取1或2比较合适，文献[17]认为取2能取得更好实验效果，文献[18]则认为取0~5之间插值效果最优，综上，本文将μ取值范围设为2~5。

2 TEC插值实测分析 2.1 方法验证 2.1.1 实验数据

选择太阳活动高年2014年，年积日（day of year, DOY）001天通用协调时(universal time coordinated，UTC)2时，国际GNSS服务组织（International GNSS Service，IGS）发布的全球电离层格网（global ionospheric map，GIM）产品（产品精度2.5~3.5总电离子含量单位（total electron content unit, TECU），时间分辨率2 h，空间分辨率2.5°×5°）作为样本，选择北半球纬度0°~ 60°、间隔20°，经度0°~ 360°、间隔15°所有格网点为待插值点，分析插值精度^[19-20]。样本点选取范围是以待插值点为中心10°×15°和20°×20°的矩形区域，距离衰减指数μ取值2~5，距离函数d_i为§1.2的5种方案。计算待插值点插值结果与已知垂直总电子含量（vertical total electron content，VTEC）之差d_vteci，统计同一纬度上所有格网点d_vteci的均方根误差（root mean square，RMS）值，验证不同方案插值精度

2.1.2 5种方案比较分析

1）插值精度分析。不同插值范围、不同插值函数的插值精度如图 3~6所示：（1）高纬度区域插值精度高于中低纬度区域，纬度20°区域插值误差较大；（2）10°×15°插值范围的插值精度高于20°×20°插值范围；（3）考虑经纬度方向变化梯度，插值结果总体好于其他方案；（4）随着距离衰减指数μ增大，插值精度呈提高趋势，不同距离计算函数插值精度趋于相近。

2）精度提升分析。针对10°×15°插值范围，分析各种方案插值精度提升(见表 1~4），总体上，各种方案相比方案（1）插值精度均有提升；采用方案（2）时，插值精度提升最高，平均提高20.57%；在纬度0°，方案（4）提升最高，达到15.90%，中高纬度区域方案最优，提升约为28%。

进一步分析：（1）IDW插值方法，样本点距离待插值点距离越近，相关性越大，调节衰减指数μ可以有效控制距离较远样本点对待插值点的影响；（2）电离层有明显的经纬度方向异性，考虑经纬度方向变化梯度的插值精度较高；（3）纬度20°区域属于电离层活动剧烈区域，插值误差较大。

2.1.3 最优方案精度分析

§2.1.2分析表明，方案(2)并不是最优的，如在纬度0°时，方案(4)提升更高，即是r取值0时，插值效果最优。分析认为：（1）低纬度区域电离层TEC变化剧烈，有限样本不足以反映电离层TEC经纬度方向异性；（2）样本点与待插值点距离相等时，经度方向相关性高于纬度方向。为此，本节设计一种综合方案（方案（6）），该方案计算距离函数d_i时，r取值为方案(2)~(5)的最小值。从表 5~6可以看出：相对于方案(1)，平均提升23.62%，相对于方案(2)，平均提升20.57%；μ取值2~3时，精度提升较大。

综上分析：（1）方案(6)计算r的插值精度优于其他5种方案；（2）当其他插值参数相同时，插值范围10°×15°结果好于20°×20°，即样本点充足时，可以适当缩小样本点选取范围或者通过增大距离衰减指数μ 来降低距离较远样本点的影响；（3）电离层TEC经度方向相关性高于纬度方向，样本点充足时可以适当选用经度差更小的点或通过调节因子r来降低纬度差较大的样本点的影响。

2.2 长期插值精度分析 2.2.1 与太阳活动周期有关的插值分析

选用2004—2015年共12 a“两分两至”日前后UTC 6时GIM产品（对应经度120°E的地方时为14时），插值范围10°×15°，μ 取值3，采用方案(6)的IDW插值方法，分析电离层活动剧烈区域（南北纬度20°）的插值精度，并与方案(1)插值结果对比。从图 7~10和表 7~8可以看出:（1）电离层活动剧烈区域，方案(6)相对于方案(1)精度提升约25%，总体插值精度约为2.0~3.0 TECU；（2）春分日和秋分日前后电离层活动相对剧烈，电离层插值精度相对较差，采用方案(6)精度提升明显；（3）太阳活动剧烈年份，IDW插值精度相对较差，但方案(6)相对于方案(1)在不同年份精度提升稳定，认为是GIM样本精度所致。

2.2.2 与地方时有关的插值分析

选用2014年太阳活动高年“两分两至”日GIM产品，采用方案(6)IDW插值方法，分析电离层活动剧烈区域（南北纬20°），经度范围90°E~120°E（间隔10°，中心经度105°，对应地方时UTC 7时），统计同一时刻所有待插值点插值结果与已知值之差d_vteci的RMS值，分析随地方时变化情况。插值范围为10°×15°，μ取值为3，比较方案(6)和方案(1)的插值精度。从图 11~14可以看出：（1）方案(1)和方案(6)插值精度均随时间变化明显，在地方时14 h（UTC 7时）后3~5 h（冬至日时间滞后），插值精度明显变差；（2）方案(6)相对于方案(1)在午后电离层活动剧烈时，插值精度提升明显，RMS最大不超过4.0 TECU。

3 结语

本文从IDW插值方法插值范围、插值核函数等入手，考虑电离层TEC分布经纬度方向异性，引入经纬度方向异性调节因子r ，设计多种电离层距离计算方案，选用一组GIM数据分析不同方案插值精度，获取最优方案IDW插值，之后，采用2004—2015年部分GIM数据，分析最优方案IDW插值精度，得到一些有益结论，具体如下：

1）插值范围影响电离层TEC插值精度。基于GIM数据插值，当其他插值参数相同时，插值范围10°×15°结果好于20°×20°，插值范围对待插值点的影响可以通过调整样本选取范围或距离衰减指数μ来调节。

2）电离层TEC经纬度方向异性影响插值精度。经度方向相关性高于纬度方向，当样本点充足时，可以适当选用经度差更小的样本点或通过调节因子r来降低纬度差较大的样本点的影响，提高插值精度。

3）综合方案（方案(6)）能较好地反映经纬度方向异性特点，插值精度优于其他方案。采用该方案的IDW插值分析长期插值精度，电离层活动剧烈区域（南北纬20°）连续12 a“两分两至”日前后的GIM数据插值，最优方案比普通IDW插值精度提升约25%；2014年太阳活动高年“两分两至”日GIM数据插值，在地方时14时后3~5 h电离层活动剧烈时，最优方案插值精度提升明显，插值误差RMS最大不超过4.0 TECU。