短波通信选频中的电离层同化短期预报方法

汤 军 姬生云 王 健 王先义

(1．海司信息化部，北京 100841；2．中国电波传播研究所，山东 青岛 266107)

引 言

电离层短期预报技术被广泛应用于高频频率管理，与可用频率预测、频率分配等方法成为短波链路通信选频不可或缺的重要组成部分，为短波通信选频提供重要支撑.当前，电离层垂直探测获得的foE、foF2、M(3000)F2等电离层特性参数是短波链路通信选频重要的基础数据[1-2]；其中，foF2是短波可用频率预测最主要的参数，对foF2预测的准确性直接影响了短波频率预报的准确性，进而影响了短波频率资源划分、资源管理的可靠性.目前采用的方法主要有自相关函数法[3-4]，多元线性回归法[5]，人工神经网络法[6-7]，等效太阳黑子数法[8]， 卡尔曼滤波法[9-10]、 相似日法[11-12]、暴时电离层预报方法[13]等等.上述方法各有优势，如自相关法在低于3 h的预报时间内能够取得高的预测精度，相似日在1～24 h情况下均取得较高的预测精度，且稳定性高.然而这些方法都需要较长时间的数据积累以及一些其他的参数支撑，例如自相关方法需要20～30 d的数据，神经网络法需要5～7 d的电离层观测数据以及地磁K指数和Ap指数等信息，而相似日法则需要预先利用自相关函数判断相似性，再进行预测，所需数据量与自相关法接近.本文将卡尔曼同化技术引入电离层参数短期预报中，基于电离层的日变化规律利用同化技术，利用1～2 d的电离层观测数据，实现电离层特征参数foF2的1～24 h的短期预报.

1 基于同化技术的短期预报方法

1.1 数据同化方法

数据同化方法的出现解决了两类问题：一是如何把观测和模式所带来的两种不同但又“互补”的信息融合起来[14]，二是如何利用同化思想进一步提高资料分析质量和数据预报的准确性[15].方法有早期的统计订正法[16](Successive Correction Method，SCM)、最优插值法[17](Optimal Interpolation，OI)，以及现在较为成熟的变分同化法[18](Variational Analysis)，卡尔曼滤波(Kalman Filter，KF)和拓展卡尔曼滤波(Extended Kalman Filter，EnKF)[19]等等.

1.2 Kalman滤波同化原理

鉴于KF方法已经成熟，且易于工程实现[20]的特点，本文选用了该方法进行数据同化分析.KF最早由Kalman于1960年首先提出，其基本思想是：首先进行模式状态的预报，在此基础上，根据观测数据对模式状态进行修正.Kalman滤波是一种有严格数学理论基础的预报—校正统计优化方法，这种优化的特点是能够不断地将观测资料同化到动态系统中，以数值模式为动力约束条件，在使观测和模式结果均方误差达到最小的条件下，得到基于以前所有观测和当前观测的系统变量的最优线性估计.而且，该方法具有递推特性，可以借助前一时刻的滤波结果推出当前时刻的状态估计值，从而大大减少计算量和存储量.

由Kalman滤波同化基本思路可知：随着模式状态预报的持续进行和新的观测数据的陆续输入，同化过程可不断向前推进.

Kalman滤波同化预报方程可表示为

xa=xb+K(x0-Hxb).

(1)

式中：xb是模式预报向量(即背景场)；xo是观测向量(即观测场)；H是观测算子，即模式预报向量与观测向量转换因子；K是增益矩阵，可表示为

K=PbHT(HPbHT+R)-1.

(2)

式中：R是观测误差的协方差矩阵；Pb是背景场误差协方差. 通过这两个参数就可以求出分析误差协方差矩阵Pa.

Pa=Pb-PbHT(HPbHT+R)-1HPb.

(3)

简单地讲，同化预报的状态xa是通过背景场xb和观测增长(观测场与背景场的差异)xo-Hxb之间的加权平均来估计的，权重是K，它起到用观测增长在相对附近的模式预报点上来调整背景场的作用.一般而言，离观测点越近的模式预报点要比远离观测点的模式预报点受到的影响要大.

1.3 foF2短期预报

利用上述思路，基于Kalman滤波同化理论的foF2预报过程如图1所示.

图1 foF2短期预报流程图

1) 读取观测数据：读取预报前8 h的数据作为观测场xo；预报前9～32 h(观测场前24 h)的数据为初始猜测场；读取中国参考电离层预报所得的24 h月中值数据作为背景场xb.

2) 选择特征尺度：特征尺度l是计算背景场相关函数的一个重要参数，其大小直接影响到了预报准确性.l的选择需要根据实际的情况，通过对前期预报结果进行比较，选择一个误差相对较小的数值，作为最终确定的特征尺度.

3) 计算误差协方差：分别计算背景场和观测场的误差协方差Pb和观测误差的协方差矩阵R.

4) 计算背景场的相关矩阵：利用欧氏距离公式计算背景场各个时刻数据之间的距离d，利用高斯相关函数计算相关矩阵.

(i,j= 1,…,4)，

(4)

(5)

式中，fi和fj分别为i和j时刻的foF2，单位为MHz.

5) 计算增益矩阵：其中关键是计算观测算子H，由于电离层的日变化规律特性，同时背景场数据与观测场数据都是foF2，变量参数也一致，因此，观测算子是“完美的”，即

H=E⟹H(i,j)=1.

(6)

令b=r·Pb，进而可得

(i=1,…,8;j=1,…,24);

(7)

(8)

至此，可以根据式(2)计算卡尔曼系数矩阵，即增益矩阵.

6) 同化分析：因为观测算子是完美的，所以有Hxb=xb，根据式(1)，获得同化后的预报结果.在计算增益矩阵时，部分参数都已经获得，但在计算分析误差协方差时还有一个未知量HPb，该变量的计算方法如下

(i=1,…,8;j=1,…,24),

(9)

就可以利用式(3)计算分析误差协方差矩阵.

7) 判断同化过程是否结束：设定同化过程结束条件，约束条件有两个：一是设定一个合理的分析误差协方差，每次同化过程结束时都判断这个条件是否满足，满足则终止同化，不满足则继续同化；二是设定一个最多的同化次数，当超出这个同化次数时，自动终止同化预报过程,将同化的结果作为预报结果进行输出，如果两个条件都不满足则将同化结果作为初始猜测场，跳入第一步，重新开始同化，直到两个条件满足一个为止.

2 预报方法评估

因电离层受太阳周期活动、季节等影响明显，故选取了太阳活动高年、中期、低年三个典型时期进行分析，如表1所列.

表1 太阳活动期选择

选取中国境内北京、长春、广州、海南、拉萨、兰州、满洲里、青岛、乌鲁木齐和重庆电离层观测站点数据.基于这三组实测数据利用Kalman同化技术对未来1～24 h的foF2进行预报.

特征尺度的差异会改变预报精度，图2显示了特征尺度从0.001到100变化，不同年份预报精度的相对误差.

图2 特征尺度对预报结果的影响

从图2中可以看出，特征尺度对预报结果的影响较大，以2000年为例，预报误差接近5%.因此，为了获得较好的预报精度，2000年的特征尺度l=10，2003年的特征尺度l=10，2008年的特征尺度l=0.001.

根据上面确定的特征尺度，分别对太阳活动高年、中期以及低年三个不同月份进行1～24 h的预报，并与实测结果进行比较，如图3、图4和图5所示：为了便于观察分析，所有的预报都是从北京时间0点开始.

图3 太阳活动高年预报结果与实测数据对比图

图4 太阳活动中期预报结果与实测数据对比图

图5 太阳活动低年预报结果与实测数据对比图

对上述不同太阳活动时期的1～24 h的预报结果进行均方根误差统计，如表2所列，均方根误差的统计公式为

(10)

式中，Δfi为每个时刻预测值与实测值之间的绝对误差，单位为MHz.

表2 均方根误差分析表(单位MHz)

对上述不同太阳活动时期的1～24 h的预报结果进行相对误差统计，如表3所列，相对误差的统计公式为

(11)

式中，fi为每个时刻的实测值.

表3 相对误差的统计分析

通过统计，太阳活动高年预报的均方根误差最大，均值为0.82 MHz，中期次之，均值为0.70 MHz，低年最小，均值为0.50 MHz.结合表2可以看出，在2000年，太阳活动高年的平均均方根误差最大，拉萨、海南和长春三站预报的均方根误差尤其高，最高达到了1.74 MHz，其余各站则较小，满洲里站最小为0.29 MHz；在太阳活动低年，广州和海南站预报的均方根误差较大，分别为0.99 MHz和1.27 MHz，其余各站误差基本上都小于0.4 MHz；太阳活动中期预报的均方根误差较为平稳，除了海南、重庆和拉萨站外，其他各站的误差都在0.5 MHz左右.

由表3可以看出，海南、拉萨和广州站的预报精度稍差，相对误差在15%左右，北京、青岛和乌鲁木齐的预报精度较高，相对误差在8%左右，其他站相对误差大体上10%左右，总体上的预报精度达到了10%.综合比较均方根误差和相对误差可以发现，中高纬度地区的预报准确性要高于中低纬度地区.

下面，对不同预报尺度下的预报误差分析，预报尺度间隔分别为1～24 h，25～48 h以及49～72 h，如表4所示.

表4 不同预报尺度的误差分析统计

从表4可以看出，均方根误差和相对误差的变化具有一致性，前24 h的预报误差最小，25～48 h的预报误差次之，49～72 h的预报误差最大，即随着预报时间的增长，预报准确性会下降.1～24 h预报的均方根误差在0.57～0.72 MHz之间，相对误差在9%～12%之间.根据相对误差和预报结果的稳定性，本文推荐的预报时间为1～24 h.与文献[11]介绍的相似日分析方法相比，本方法对1～24 h的预报相对误差为10%，与相似日方法的相对误差基本相同.

3 结 论

本文以探测数据为基础，利用Kalman滤波同化方法，提出了一种1～24 h的电离层参量foF2短期预报方法.并利用该方法对太阳活动高、低、中期年分不同季节不同地区预报结果进行分析.分析结果显示了该方法在不同太阳活动周期时的预报准确度相近，对中高纬度地区的预报准确性较高，如青岛、北京和满洲里等地，同时，预报准确度随着预报尺度的增加而下降，例如1～24 h的预报准确度优于25～48 h，该方法对数据积累的时间要求低，更利于实现短波频率预报工程化，为短波频率管理提供可靠的支撑[1-2]，同时可为相关领域的扩展应用提供支撑[21].

通过对数据误差进行分析，该方法还有需要改进的地方，如中低纬地区的预报准确性较低.因此，借助于逐步成熟的集合Kalman同化以及4维变分同化等技术，提升中低纬度电离层短期预报准确性是下一阶段研究工作的重点.

[1] ITU. P. 533-9 Method for the Prediction of the Performance of HF Circuits[S]. Geneva: ITU-R, 2008.

[2] ITU. P.1240-1 Methods of Basic MUF, Operational MUF and Ray-path Prediction[S]. Geneva: ITU-R, 2008.

[3] 刘瑞源, 刘顺林, 徐中华, 等. 自相关分析法在中国电离层短期预报中的应用[J]. 科学通报, 2002, 50(24): 2781-2785.

[4] 冯 静, 柳 文, 焦培南等. 电离层特征参量的自相关原理插值方法[J]. 空间科学学报, 2009, 29(2): 195-201.

FENG Jing, LIU Wen, JIAO Peinan, et al. Autocorrelation method for interpolation of ionospheric characteristic parameters[J]. Chinese Journal of Space Science, 2009, 29(2): 195-201. (in Chinese)

[5] MARIN D, G MIRO A V MIKHAILOV. A method forfoF2short-term prediction, 4th COST 251 Workshop proceedings [C]// 4th COST 251 Workshop proceedings, Madeira, 1999: 214-222.

[6] 陈 春, 吴振森, 赵振维, 等. 基于神经网络技术的foF2短期预报方法[J]. 电波科学学报, 2008, 23(4): 708-712.

CHEN Chun, WU Zhensen, ZHAO zhenwei, et al. A short-termfoF2forecasting method based on neural network techniques[J]. Chinese Journal of Radio Science, 2008, 23(4): 708-712. (in Chinese)

[7] 陈 春, 吴振森, 孙树计, 等. 利用神经网络提前24小时预报电离层foF2[J]. 电波科学学报, 2009, 24(1): 152-156.

CHEN Chun, WU Zhensen, SUN Shuji, et al. Forecasting the ionosphericfoF224 hour in advance by neural network techniques[J]. Chinese Journal of Radio Science, 2009, 24(1): 152-156. (in Chinese)

[8] 陈 春, 吴振森, 孙树计, 等. 利用F10.7短期预报电离层foF2[J]. 空间科学学报, 2009, 29(4): 383-388.

CHEN Chun , WU Zhensen, SUN Shunji, et al. Short-term Forecasting of the IonosphericfoF2With F10.7[J]. Chin J Space Sci, 2009, 29(4): 383-388. (in Chinese)

[9] 乐新安, 万卫星, 刘立波, 等. 基于Gauss-Markov 卡尔曼滤波的电离层数值同化现报预报系统的构件-以中国及周边地区为例德观测系统模拟实验[J]. 地球物理学报, 2010, 53(4): 787-795.

YUE Xinan, WAN Weixing, LIU Libo, et al. Development of an ionospheric numerical assimilation nowcast and forecast system based on Gauss-Markov Kalman filter-an observation system simulation experiment taking example for China and its surrounding area[J]. Chinese Journal of Geophysics, 2010, 53(4): 787-795. (in Chinese)

[10] 徐 彤, 吴振森, 吴 健, 等. Kalman滤波电离层短期预报[J]. 电波科学学报, 2007, 22(增刊): 146-149.

XU Tong, WU Zhensen, WU Jian, et al. Ionospheric short-time forecasting using the Kalman filter[J]. Chinese Journal of Radio Science, 2007, 22(Sup): 146-149. (in Chinese)

[11] 柳 文, 焦培南, 冯 静, 等. 电离层参数的相似日短期预测方法[J]. 电波科学学报, 2010, 25(2): 240-247.

LIU Wen, JIAO Peinan, FENG Jing, et al. Method of similar days for ionospheric parameter short-term forecasting[J]. Chinese Journal of Radio Science, 2010, 25(2): 240-247. (in Chinese)

[12] 冯 静, 柳 文, 凡俊梅, 等. 一种电离层短期预报的新方法[J]. 装备环境工程, 2009, 6(3): 15-20.

FENG Jing, LIU Wen, FAN Junmei, et al. A new method for short-term ionospheric forecast[J]. Equipment Environmental Engneering, 2009, 6(3): 15-20. (in Chinese)

[13] 孙树计, 丁宗华, 陈 春, 等. 一种爆时电离层foF2的经验预报方法[J]. 电波科学学报, 2009, 24(4): 637-644.

SUN Shuji, CHEN Chun, DING Zonghua, et al. A storm-time empirical ionosphericfoF2prediction method[J]. Chinese Journal of Radio Science, 2009, 24(4): 637-644. (in Chinese)

[14] 刘成思, 薛纪善. 关于集合Kalman滤波的理论和方法的发展[J]. 热带气象学报, 2005, 21(6): 628-633.

LIU Chengsi, XUE Jishan. The ensemble kalman filter theory and method development[J]. Journal of Tropical Meteorology, 2005, 21(6): 628-633.(in Chinese)

[15]BOUTTIER F, COURTIER P. Data Assimilation Concepts and Methods, Meteorological[R/OL]. 1999[2012-08-31]. http://www.ecmwf.int/newsevents/training/lecture_notes/pdf_files/ASSIM/Ass_cons.pdf

[16] BERGTHORSSON P, DONS B. Numerical weather map analysis. [J]. Tellus, 1955, 7(3): 329-240．

[17]GANDIN L S. Objective Analysis of Meteorology Fields[M]. Leningrad: Leningrad Hydromet Press, 1963.

[18] SASAKI Y. Some basic formalisms in numerical variational analysis[J]. Mon Wea Rev, 1970, 98: 875-883.

[19] KALMAN R E. A new approach to linear filter and prediction problems[J]. J Basic Eng, 1960, 82(3): 35-45.

[20] EVENSEN G. The ensemble Kalman filter: theoretical formulation and practical implementation[J]. Ocean Dymamics, 2003, 53: 343-367.

[21] 王 健, 惠守强, 付 炜, 等. 高频单站定位误差特性分析及精度优化构想[J]. 电波科学学报, 2010, 25(5): 925-933.

WANG Jian, HUI Shouqiang, FU Wei, et al. Error characteristic analysis and accuracy optimizing idea of HF single site location[J]. Chinese Journal of Radio Science, 2010, 25(5): 925-933. (in Chinese)