电离层作为日地空间环境中的重要一环，既是电波传播的媒介，又是耦合磁层与中性大气的纽带.其受到太阳辐射、地磁扰动、大气活动等多种因素作用，呈现出复杂的气候学和天气学变化特征，从而会对卫星导航、雷达成像、短波通信等涉及电波传播和电磁环境的技术手段产生严重影响.而电离层总电子含量(Total Electron Content, TEC)是描述电离层的重要参量之一，其与电离层引起的电波传播附加时延成正比；故而分析电离层形态结构、研究电离层扰动特性、消除电离层延迟误差，亟需获取及时、准确、可靠的电离层TEC时空分布和变化信息.为了解决电离层对无线电通讯、卫星导航影响的关键问题，许多研究针对电离层TEC进行了全面监测、实时反演、并最终建立了区域或者全球性的模型(e.g., Mao et el., 2008; Ghoddousi-Fard et al., 2011; Liu et al., 2011; A et al., 2012; Wan et el., 2012; Chen et al., 2015; Mukhtarov et al., 2013).然而，在面向电离层天气业务应用需求的TEC现报领域，传统电离层经验模型和理论模型仍然有一定的局限性.电离层经验模型是以大量电离层观测数据为基础而进行统计分析所得，其描述的是电离层气候学意义上的统计平均特征，在电离层平静期的TEC现报效果较好，而在暴时则有较大误差.电离层理论模型主要是通过数值模拟等手段，解算电离层等离子体物理化学方程所得，其能从物理层面反映电离层动力学特征和变化规律；但理论模型运算时间较长，再加上受边界条件随机误差和其自身稳定性与完备性方面的限制，故而通常多用于研究电离层演化机理和结构特征，而较难直接应用于电离层TEC业务化现报当中.因此，如何准确把握电离层TEC的分布特征和变化规律，建立起既包含内在物理，又反映真实观测，同时满足空间天气业务应用需求的电离层TEC现报模型及系统，提供高精度准实时的电离层误差修正信息，成为了电离层研究领域的重点内容之一.

近十年来全球导航卫星系统(Global Navigation Satellite System, GNSS)的空前发展和地基观测网的迅速增加，为实现区域甚至全球的高精度电离层TEC监测和现报提供了有效的手段.当前国际GNSS服务组织(International GNSS Service，IGS)每天能收集到全球数百个GNSS接收站的数据.在此基础上，共有5家分析中心利用IGS台站网格提供全球电离层TEC地图产品(Global Ionospheric Maps, GIMs)，它们分别是：欧洲定轨中心(CODE)，美国喷气动力实验室(JPL)，欧洲空间局(ESA)，西班牙加泰罗尼亚理工大学(UPC)，武汉大学与中国科学院(WHU-CAS)(Hernández-Pajares et al., 2009).其中以JPL和CODE的产品较为通用，JPL采用三角网格插值和卡尔曼滤波相结合的方法来进行全球TEC网格点现报(Mannucci et al., 1998).CODE和ESA建模的原理基本一致，即采用球谐函数模型，在日固地磁参考系下用多阶次的球谐函数进行展开，TEC值在时间域内采用分段线性函数进行表达(Schaer, 1999; Feltens, 2007; Feltens and Schaer, 1998).WHU-CAS的GIM产品采取广义三角级数逐基准站建立局部电离层TEC模型，再采用球谐函数外推构建全球TEC格网(Li et al., 2015).上述TEC产品都是全球覆盖的，其时空网格分辨率基本一致(经度：5°，纬度：2.5°，间隔：30 min(CAS)、1 h(CODE)、2 h(JPL\ESA\UPC))，对于重点区域的电离层高精度现报需求而言仍有待进一步提升.

因此，许多国家和地区都构建了精度更高、实时性更好的全球或区域性TEC地图.其中一些产品采取了网格法，即利用较高密度的GNSS站网布局，假设固定网格内TEC短时不变的原则进行解算.例如，美国麻省理工大学的GPS自动处理系统(MAPGPS)利用不限于IGS的上千个GPS台站构建了全球TEC地图，其分辨率可达1°×1°×5 min，并通过Madrigal系统进行在线发布(Rideout and Coster, 2006).日本情报通信研究机构(NICT)旗下的空间天气预报中心同样基于日本高密度的GNSS台站，构建了分辨率最高可达0.5°×0.5°×10 min的地区性TEC地图产品(Miyake and Jin, 2010).网格法的优点在于其结果严格基于实测数据，并没有采用插值或者函数扩展，所以精度相对较高；而缺点在于其TEC地图覆盖并不完整，在没有台站的沙漠或海洋区域会有留空，所以其更适用于范围较小而台站密度较高的区域.另外一些TEC产品则是在一定数量的地基GNSS台站数据基础上通过函数展开或者插值得到的.例如，澳大利亚的空间天气服务中心和加拿大自然资源部的地学中心都采用球冠谐函数展开的方法构建了区域TEC地图产品(Bouya et al., 2010; Ghoddousi-Fard et al., 2011).我国的中国科学院地质与地球物理研究所利用经验基函数模式算法，利用实测数据估算对应系数，构建了覆盖中国区域的TEC现报系统(万卫星等, 2007).德国航天局(DLR)采用多项式模型对TEC进行泰勒展开，建立了欧洲TEC现报地图(Jakowski et al., 2011).比利时皇家天文台利用样条插值方法同样给出了欧洲区域TEC现报地图(Bergeot et al., 2014).总体而言，利用函数展开或者插值得到的TEC地图具有完整的覆盖性，但其准确性受所用函数和插值方法的制约.

随着近年来电离层数据同化理论的蓬勃发展，利用数据同化方法，建立既包含内在物理过程，又反映真实观测，同时满足人们对电离层空间天气现报要求的模型，成为了电离层研究和应用领域的重点内容之一.数据同化在气象学和海洋学领域已经得到了广泛使用，并积累了诸如最优估计理论、卡尔曼滤波、变分同化等一系列算法.其基本原理是在基于物理机制的背景模型上加入时空不规则分布的观测资料，将观测值通过同化过程融合到背景模型中，使模型与数据相互匹配以达到更精确的分析结果.在电离层研究领域，美国国防部制定的综合大学优先研究计划(Multidisciplinary University Research Initiative, MURI)中，支持了两个团队研究全球电离层数据同化.其中之一是美国犹他州立大学研究的全球电离层观测同化模型(The Utah State University Global Assimilation of Ionospheric Measurements, USU GAIM).USU GAIM采用了其自行开发的电离层预测模型(Ionosphere forecast Model, IFM)和Gauss-Markov卡尔曼滤波技术，对TEC数据、电子密度数据、以及UV数据进行近实时的同化(Scherliess et al., 2004, 2006; Schunk et al., 2004, 2005).MURI支持的另外一方是喷气动力实验室和南加州大学共同开发的全球同化电离层模型(Jet Propulsion Laboratory/University of Southern California Global Assimilative Ionospheric Model, JPL/USC GAIM).JPL/USC GAIM基于一个全球时变电离层和等离子体层动力学模型，采用四维变分法同化测高仪、GPS卫星、电场、沉降粒子在内的各项观测资料(Pi et al., 2003; Wang et al., 2004).除此之外，Bust等(2004, 2007)利用三维变分法，发展了一个电离层数据三维同化模型(Ionospheric Data Assimilation Three-Dimensional, IDA3D).IDA3D可以利用经验模型或理论模型，同化测高仪、GPS接收机、卫星观测等多种手段的观测资料.Fuller-Rowell等(2006)利用国际参考电离层(International Reference Ionosphere, IRI)作为背景模型，使用卡尔曼滤波算法和经验正交函数的结合来进行数据同化，开发了覆盖美国区域的USTEC模型，并进一步将其发展为覆盖北美范围的NATEC(North American TEC).近来有越来越多的研究人员利用测高仪、非相干散射雷达、COSMIC/FORMOSAT-3掩星观测等多源电离层数据进行同化研究(e.g., Angling, 2008; Yue et al., 2007, 2011, 2012a, 2012b; Galkin et al., 2012; Lee et al., 2012; Gardner et al., 2014; Lin et al., 2015).

电离层有很强的地域特性，中国扇区所覆盖的地磁纬度较低，且包含了电离层赤道异常峰的北驼峰；这导致了该区域有着较大的梯度变化和相对复杂的电动力学过程，加之该区域的电离层观测资料相对缺乏，从而使得以IRI为代表的电离层经验模型在应用于这一区域时有较明显的偏差.因此利用该区域的电离层实测资料，结合背景模型，基于数据同化算法，建立起更适用中国扇区的同化模型是一项很有必要的工作.例如，一些研究人员利用部分地基GNSS数据和IRI模型，采取卡尔曼滤波算法(乐新安等，2010；A et al., 2015b)或三维变分方法(Yu et al., 2014；A et al., 2016)尝试构建了中国区域电离层同化模型.本文将在上述研究的基础上，依据中国及周边区域的地基GNSS数据和国际参考电离层IRI模型，结合三维变分与Gauss-Markov卡尔曼滤波两者在算法上的优势，进行数据同化，开发中国区域电离层数据同化现报系统，生成覆盖中国及周边区域的电离层TEC和GPS单频接收机定位误差的格点地图，并在中国科学院空间环境预报中心网上(http://sepc.ac.cn/TEC_chn.php)进行准实时发布，每15 min进行更新，为电离层空间环境用户提供现报服务.下文将针对该电离层数据同化系统进行介绍，文章第1节将介绍GNSS台站布局、观测数据处理、以及背景模型；第2节介绍三维变分数据同化方法；第3节给出同化现报结果和分析；最后是总结与展望.

同化系统共使用了50个左右的GNSS台站数据，主要来自于中国科学院空间环境监测网(SEMnet)和国际GNSS服务组织(IGS).中国科学院空间环境监测网集成了八家单位空间环境研究领域的优势力量，能够提供太阳、宇宙线、地磁、电离层和中高层大气等多种空间环境要素的监测信息，有12个台站配置有电离层GNSS接收机.其中隶属于中国科学院地质与地球物理所的有四个站：漠河站、北京站、武汉站、三亚站；隶属于中国科学院空间中心的有八个站：福州站、厦门站、广州站、南宁站、富克站、昆明站、张掖站、北京(空间中心)站.此外IGS在中国及周边也有38个台站数据可使用.图 1给出了上述GNSS台站的分布示意图.

图 1

Fig. 1

图 1 同化现报系统所用的中国及周边区域GNSS台站分布示意图 Fig. 1 Distribution of GNSS stations over China and adjacent areas for assimilation now-casting system

电离层TEC的解算可以通过GPS双频L波段信号(f_L1=1575.42 MHz, f_L2=1227.60 MHz)在电离层斜向传播路径上的伪距和载波相位的联合演算得到(Mannucci et al., 1998).采用GPS测算的TEC，其精度高、连续性好，故而成为电离层形态学和电离层建模研究的极佳数据来源.在TEC解算过程当中，一个关键任务是剔除卫星与接收机内部的硬件偏差，从而获得准确的TEC值.本系统中采用了一种自适应网格(Self-Calibration of Pseudo-Range Errors, SCORE)方法(Biship et al., 1996; Zhang et al., 2009).其核心思想是将电离层视为没有厚度的薄层球壳，位于其质心所在高度，卫星-接收机连线和电离层球壳的交点称为电离层穿刺点(IPP).将电离层按照地方时和纬度划分成足够小的均匀网格，并假设每个网格点内部的垂直TEC相等.如果两条连线的电离层穿刺点落在同一网格内，则可以将其表述为含有TEC和硬件偏差的等式.SCORE方法的基本理论可以用公式进行如下表达：VTEC_αi(θ_αi, τ_αi)表示GPS卫星α和接收机i连线路径上的电离层穿刺点垂直TEC，θ表示地理纬度，τ表示地方时；同理VTEC_βj(θ_βj, τ_βj)表示GPS卫星β和接收机j连线路径上的电离层垂直TEC.如果两条连线的电离层穿刺点落在同一网格点内，则有

(1)

上式可改写为

(2)

其中B_α和B_β为接收机的硬件偏差；B_i和B_j为卫星的硬件偏差；z是电离层穿刺点处的太阳天顶角.最终50个台站的数据构成以硬件偏差为未知量的超定方程组，可以用奇异值分解法来进行解算.关于SCORE方法的详细描述和结果验证，读者可参考A等(2015a)和Zhang等(2010, 2014).

同化系统所采用的背景模型是国际参考电离层模型IRI(Bilitza, 2001; Bilitza and Reinisch, 2008; Bilitza et al., 2014).IRI的数据来源基于地面、卫星、火箭探测，能够提供电子、离子等一系列电离层参量的密度和温度等表现其全球分布形态和时空变化的信息.IRI的电子密度是通过对包括F₂区、F₁区、E区的峰值电子密度的锚定点的扩展来实现的，其顶部电子密度由NeQuick模型给出(Radicella and Leitinger, 2001; Coisson et al., 2006).本系统中所采用的是最新版本的IRI-2016.

数据同化的基本原理是将观测数据融合到基于物理规律的模型中，产生既逼近真实状态又包括内在物理的高可靠度的结果；同化过程可以被简单的看作是观测资料和背景模型基于误差协方差矩阵的最优结合.同化方法根据其理论基础主要可以分为两类：一类是基于统计学理论的估计方法，如最优插值、卡尔曼滤波、扩展卡尔曼滤波、集合卡尔曼滤波等方法；另一类是基于变分理论的方法，如三维变分、四维变分法等.目前三维变分法是国际上主要数值预报中心进行资料同化所采取的一般方法，在本文中也采用这种方法进行数据同化.三维变分法的基本概念是通过迭代一个代价函数(或称泛函)，产生一个分析时刻的电离层真实状态的最优估计.泛函一般表示为观测值和背景值的加权和，相应的权重与对应的观测误差协方差和背景误差协方差有关，其通常的表达式为

(3)

其中x指的是系统的状态向量，即同化网格点的电子密度；x_b是背景场模型向量，由IRI提供；P为背景误差协方差矩阵，P=E[e_be_b^T], e_b=x_b－x.y是观测向量，即STEC数据；H是观测算子，使得模式向量向观测向量转化，即GNSS信号传播的卫星-接收机连线在电离层网格当中的长度；R是观测误差协方差矩阵，R=E[e_oe_o^T], e_o=y－Hx.变分问题就是使得泛函取极小值时的系统分析状态的解，这个解代表了电离层真实状态的一个最大似然(最小方差)估计；而为了求泛函的极小值，一般需要利用其梯度信息：

(4)

当上式等于0时，即可求解出系统状态向量的最优估计值.

在同化系统的连续递推运算当中，主要包括时间更新和测量更新两步.时间更新是指对下一时刻的状态向量和误差协方差矩阵进行先验估计，而测量更新是指利用观测数据对状态向量和误差协方差矩阵进行同化修正.由于背景场是IRI经验模型，不具备前向预报的功能，我们在这里采用Gauss-Markov卡尔曼滤波来进行时间更新，它假设观测对背景场的修正随时间的变化指数衰减，用公式可以表示为

(5)

(6)

其中下标b表示背景值，下标f表示预报值，t_n为当前时刻，t_n+1为下一时刻.Δt为时间步长，本系统中为15 min，τ表示电离层时间尺度，此处取为30 min.

在得到了先验估计和观测数据后，就可以进行测量更新，用公式表示为

(7)

(8)

此处下标a表示同化分析场.通过上述两步的递推反复，就可以使得同化系统持续运行更新.

在三维变分数据同化中，观测误差和背景误差协方差矩阵决定了观测资料和背景场的权重分配，对误差协方差的描述决定着同化效果的优劣.误差协方差(Covariance)可以被分为两部分：第一部分表征的是观测(背景)和电离层真实状态的偏差(Variance)，第二部分表征的是空间中两点之间的相关性(Correlation).在前人的研究结果当中(e.g., Bust et al., 2004; Wang et al., 2004; Mandrake et al., 2005; Yue et al., 2007, 2011)，观测值通常被认为是不存在相关性的，且观测误差是无偏的，并与观测值的平方成正比；背景误差也与背景值的平方成正比，其相关性呈现高斯分布；则观测误差协方差矩阵和背景误差协方差矩阵可以分别表示为

(9)

(10)

其中C_R和C_P是可调节的比例系数，可以通过对观测和背景进行误差统计而予以设定.d_ij和L_ij分别是第i个和第j个格点间的地球大圆弧长以及电离层水平方向相关距离；z是格点高度，L_z是电离层垂直方向相关距离.关于电离层相关距离的数学表达式，水平方向可以表示为格点i与格点j之间方位角α的函数：

(11)

其中L_θ和L_φ分别对应于子午向和纬向的电离层相关距离，可以表示为随地方时变化的影响因子乘以一个随地磁纬度变化的函数；电离层垂直方向的相关距离可以用电离层标高来表示，上述参量相应的取值如图 2所示.关于同化流程、误差协方差矩阵和电离层相关距离的详细描述，读者可参考A等(2016).

图 2

Fig. 2

图 2 (a) 电离层相关距离影响因子随地方时的变化；(b)垂直相关距离随高度变化；(c，d)子午向和纬向电离层相关距离随地磁纬度变化示意图 Fig. 2 (a) Diurnal variation of impact factor; (b) Vertical ionosphere correlation length with respect to altitude; (c, d) Meridional and zonal ionosphere correlation length with respect to magnetic latitude

本系统选择中国及周边区域进行数据同化，经度范围为70°E—140°E，纬度范围15°N—55°N，高度范围100~1000 km，划定的网格点间隔为2°×2°×50 km，图 3给出了网格点划分和一个同化时段内观测射线分布的示意图.采用上述三维变分与Gauss-Markov卡尔曼滤波结合的算法，同化融合50个GNSS台站的斜向TEC数据和背景IRI模型，最后积分网格点电子密度可以生成覆盖中国及周边区域的电离层TEC格点化准实时现报地图.由于电离层延迟误差和TEC成正比，和接收机频率成反比，同化产品同时叠加给出了GPS单频接收机的延迟误差.相比于IGS分析中心的全球TEC地图，本产品能够更有效地反映电离层赤道异常东亚扇区的北驼峰区域和中纬区域的变化态势.为了开展电离层环境服务保障工作，TEC和接收机延迟误差的格点数据产品会发送给定制用户，用于工程任务和科学研究；地图产品则在中国科学院空间环境预报中心的网站上(http://sepc.ac.cn/TEC_chn.php)发布, 每15 min进行更新.图 4给出了网页上发布的TEC地图现报示意图.

图 3

Fig. 3

图 3 网格点划分和观测射线在一个同化时段内的分布示意图 Fig. 3 The grid division and the distribution of selected rays for certain assimilation period

图 4

Fig. 4

图 4 网页发布的电离层TEC和GPS单频接收机误差的同化现报地图 Fig. 4 On-line data assimilation now-casting map of ionospheric TEC and approximate L1 position error

为了验证同化效果，图 5给出了一天24小时数据同化前后的TEC结果与第三方实测的电离层TEC比较.此处实测值由MIT Madrigal Database提供.图中上方的两幅子图给出了同化前后TEC残差的直方图分布.可以看出数据同化前的残差分布明显偏离中心，与观测值有较大偏差，其平均误差为13.93 TECU，标准差为11.71 TECU；而数据同化后的TEC残差更接近于无偏的高斯分布，平均误差下降至0.21 TECU，标准差下降为3.20 TECU.图中下方的两幅子图给出了同化前后TEC的散点对比，可以看到相关系数由同化前的0.7761上升至同化后的0.931，同化结果与实测值的一致性较好.上述分析说明数据同化能够有效地将观测数据融合到背景场中，获得更为合理准确的现报结果.关于同化结果与多源电离层数据进一步的比较验证，读者可参阅A等(2015b, 2016)中的相关内容.

图 5

Fig. 5

图 5 数据同化前后的TEC残差直方图分析(上图)和散点对比图(下图)，以MIT Madrigal TEC为参考 Fig. 5 Histogram of the TEC residuals (top) and scatter points comparison (bottom) before and after data assimilation; Verification is implemented on the basis of MIT Madrigal TEC

本文介绍了我国首个基于数据同化算法的电离层TEC现报系统的构建机理.该系统利用中国科学院空间环境监测网和IGS的50个地基GNSS台站的TEC数据作为观测值，使用国际参考电离层(IRI)作为背景场，并采用三维变分与Gauss-Markov卡尔曼滤波相结合的算法进行背景场和观测值的数据同化，生成覆盖中国及周边区域(15°N—55°N，70°E—140°E)的电离层TEC和GPS单频接收机定位误差的格点化准实时现报地图，并将产品在中国科学院空间环境预报中心(http://sepc.ac.cn/TEC_chn.php)网上发布，每15 min进行更新.该系统在数据解算、同化算法等方面均有一定创新，已用于中国及周边区域的电离层环境实时监测，可为卫星导航、雷达成像、短波通信等科学研究和工程应用提供相对及时、准确、有效的电离层TEC和误差修正信息.

目前的观测数据整体上基本覆盖了我国及周边区域，然而台站分布仍然存在不均匀的情况，例如在华南沿海区域较为密集，而在中西部地区则分布较为稀疏.未来为了改善我国电离层TEC现报的精度，进一步提升地图网格点分辨率，需要适当增加中西部地区的台站数量.为了弥补这点，本同化系统在设计时也预留了未来增加台站的开放性原则，对于新的GNSS台站，只要其观测数据是典型的RINEX格式，系统都可以将其兼容进入同化TEC地图的计算，以提供更精准的电离层空间环境服务.同化系统计划在未来业务化、常规化引入中国大陆构造环境监测网络(陆态网络)的GNSS数据.陆态网络是以GNSS为主，辅以多种空间观测技术，主要用于监测中国大陆地壳运动、重力场形态及变化、大气圈对流层水汽含量变化及电离层离子浓度的变化的地球科学综合观测网(李强等，2012).陆态网络在全国共包括260余个GNSS基准站，利用这样稠密的监测网有助于给出更高精度的电离层空间天气同化现报产品，对于电离层高分辨率动态监测、分析暴时扰动形态提供了有效的手段.此处尝试以2015年3月17—18日大磁暴期间的中国及周边区域电离层扰动为例，给出了使用陆态网络GNSS数据的TEC同化结果图.图 6给出了2015年3月16—20日的行星际太阳风和磁场参数(ACE卫星所测)，以及环电流Dst指数.本次磁暴事件是由15号爆发的日冕物质抛射所致，可以看出在3月17日世界时4点左右，激波到达L1点，太阳风速度出现陡然上升，行星际磁场南向分量(IMF Bz)随之出现多次南北翻转，并在11点之后保持长时间南向，Dst指数在磁暴急始后出现大幅下降，并在23时左右达到最低值-223 nT.本次磁暴事件给中国区域上空电离层造成了较大的扰动，其中主要的扰动时段出现在3月18日磁暴恢复相期间.图 7给出了从3月17日下午到3月18日全天，中国及周边区域电离层的相对TEC(减去磁暴前五日静日平均)的变化态势，这是在同化系统中引入了陆态网络GNSS数据所得的结果.从图中可以看出，在17日下午，部分中高纬度区域出现了电离层TEC的正相扰动，而低纬地区则有TEC负相扰动的特征.随着Dst指数在17日23时达到最低值，磁暴进入恢复相，中高纬度区域也呈现出电离层负暴的倾向.在18日，整个中国区域出现了普遍的电离层负暴，中低纬度的负暴效应尤为明显，DTEC最低达到了-48 TECU.强烈的电离层负暴在18日晚些时候才逐渐消退，有些地方甚至持续到19日.对比图 4和图 7可以看出，在增加了陆态网络数据驱动之后，同化系统可以给出更为精细化的电离层结构变化，为监测分析我国上空电离层TEC的时空变化特性，开展电离层空间环境观测研究和服务保障提供更好的支持.

图 6

Fig. 6

图 6 2015年3月16—20日期间，行星际磁场南向分量Bz、太阳风速度、太阳风动压、环电流指数Dst Fig. 6 Interplanetary Magnetic Field (IMF) Bz component, solar wind speed, dynamic pressure, and Dst index during March 16—20, 2015

图 7

Fig. 7

图 7 2015年3月17日下午和18日全天，利用数据同化得到的中国及周边区域电离层DTEC变化 Fig. 7 Variation of ionospheric DTEC derived by data assimilation over China and adjacent areas during March 17—18, 2015

未来随着同化系统可用的驱动数据的增加，在电离层现报中会涉及到超大规模误差协方差矩阵的存储、求逆等运算问题，目前可设想的是采用稀疏矩阵和SVD分解的方法来予以解决.此外，还需要持续优化算法，提升同化系统的运行效率和结果的可靠性.未来为了进一步发展具有现报预报功能的电离层同化系统，不排除在系统中引入理论模型作为背景场，使用同化结果作为初始值，利用理论模型解决前向预报问题等.上述工作都在规划和开展当中，为提升电离层同化系统性能和满足相关空间环境应用做进一步的改进.

陆态网络的GNSS数据由中国地震局地壳运动监测工程研究中心提供；ACE卫星数据由美国航天局Goddard空间飞行中心提供；Dst指数由京都世界地磁数据中心提供；部分TEC验证数据由MIT Madrigal Database提供.在此一并表示感谢.