基于Kalman滤波技术的电离层电子密度反演研究

桂晓纯，洪振杰，赵运超

(温州大学数学与信息科学学院，浙江温州 325035)

基于Kalman滤波技术的电离层电子密度反演研究

桂晓纯，洪振杰†，赵运超

(温州大学数学与信息科学学院，浙江温州 325035)

利用Kalman滤波技术，实现了基于6个卫星接收站所获取三频信标数据的电离层电子密度反演，并对影响反演效果的因素进行了分析．结果表明：初始值、水平相关度以及垂直相关度等因素都会对电离层电子密度反演效果产生重要影响，其中水平相关设为14°、垂直相关设为200 km时反演效果较好．

Kalman滤波；电离层；卫星三频信标；总电子含量

三频信标是一种重要的空间环境探测星载设备，其工作机制是发射机在VHF、UHF和L频段上输出相对稳定相位的载频信号，并经全向天线向预定覆盖区域辐射，通过地面接收设备在单站或者多站进行信号幅度和相位测量，再进行反演可以获得电离层电子总含量、电离层闪烁和电子密度及电离层不规则体的信息[1-3]．

美国海军实验室的Bernhardt等人[4]创造性的工作使得三频信标技术有了重大突破．他们巧妙地利用数论知识，经过严格的数学推导，得到了TEC（总电子含量，Total Electron Content，简称TEC）的计算方法，使得其相位模糊系数值大幅度提高，大大降低了相位积分常数的求解难度．中国台湾和美国于 1997年联合制定了 COSMIC计划（气象、电离层和气候的卫星观测系统，Constellation Observing System for Meteorology, Ionosphere and Climate，简称COSMIC），并于2006年成功发射卫星．COSMIC携带有美国海军实验室研制的TBB（三频信标，Tri-Band Beacon，简称TBB）发射机，三个载波分别为150 MHz、400 MHz和1 067 MHz，它们发射的信号可以被地面或者其它低轨卫星接收机接收．通过不同频率之间测量的相位差就可以算出发射机到接收机路径上的总电子含量．将地面接收机得到的总电子含量数据作为输入参量利用CIT技术（电离层层析成像，Computerized Ionospheric Tomography，简称CIT）就可以反演电离层电子密度[5-9]．三频信标的精确探测的特点也逐渐引起我国地震电磁工作者的关注，经过几年的论证和实验，中国电波传播研究所吴健研究员、甄卫民研究员等人[10-12]在三频信标探测原理和 CT（层析成像，Computerized (or Computed) Tomography，简称CT）算法上有所突破，他们的工作为我国三频信标的研究奠定了基础．虽然三频信标TEC探测和CT成像技术在各方面都取得了重大的进展，但是我们需要看到的是不管是理论还是算法，三频信标CT技术还远没有达到理想的要求，都需要进一步的发展．

本文基于Kalman滤波技术，在介绍三频卫星信标测量TEC原理的基础上，实现了三频信标数据资料的反演，得到了2维的电离层电子密度反演图像，并对影响反演效果的一些因素进行了讨论，本研究对于三频信标数据获取站位布设及反演相关条件设定具有重要现实指导意义．

1 Kalman滤波反演算法及原理

Kalman滤波是Kalman在1960年提出的从与被提取信号有关的观测量中通过算法估计出所需要的信号的一种滤波算法[13]．Kalman滤波方法是一种时域方法，对于具有高斯分布噪声的线性系统可以得到系统状态的递推最小方差估计（Recursive Minimum Mean-Square Estimation, RMMSE）．

Kalman滤波是一种最优估计理论，与传统的估计方法相比，它具有以下几个特点：

第一，在时域内设计滤波器，并且算法是递推的，每一步只需要处理一个时刻的量测信息，既综合利用了t时刻以前的全部量测信息又不会使计算量随时间增大．

第二，设计Kalman滤波器采用了状态空间法，对于多维随机过程的估计也适用．采用动态方程来描述被估计量的动态变化规律，被估计量的动态统计信息由激励白噪声的统计信息和动态方程确定．由于激励白噪声是平稳过程，而且动态方程是已知的，所以被估计量可以是平稳的也可以是非平稳的，所以Kalman滤波也适用于非平衡过程．

假设Y为电子密度的观测值，对于电离层电子总含量，它和电子密度的关系为：

其中T(r1,r2)表示从卫星位置r1到接收器位置r2的电子总含量，Ne为电子密度．

把电子密度作为一个三维网格区域值可以得到一个离散近似为：

H中的元素Hi,k表示射线i经过网格点k处射线长度．

由于存在不确定的模糊度，TEC的测量通常看作是相对测量，而不是绝对测量．通过将每个台站的测量值减去一个特定的测量值，转换成绝对测量量T1：

T0为给定接收器对指定卫星的参考测量，H0为参考测量的射线路径．

为简化推导过程，在接下来的讨论中默认H即为H-H0．在本文讨论中都是经度固定在120°时讨论的，所以H可以看成一个二维的矩阵．

在数据空间，分析电子密度X和初始电子密度X0具有如下的关系：

其中P为背景协方差矩阵，X0为初始电子密度，R为观测协方差．（5）式中PHT[R+HPHT]-1(Y-HX)在Kalman滤波中通常称为Kalman增益．

分析协方差与背景协方差具有如下关系：

P为背景协方差，aP为分析协方差，在迭代的过程中为下一时刻的背景协方差或者将其作为背景协方差的一个重要参考．

网格点k和l之间的背景协方差可以写为：

其中r表示经验误差估计，x0k和x0l分别表示k和l处的初始电子密度，Z表示高度，rkl表示点k和l之间纬度之差，Lz表示垂直相关度，Lα表示水平相关度．

2 实验过程及其结果分析

2.1 模拟反演实验准备

假设低轨卫星在800 km高空沿120°经度圈自北向南飞行，在中国境内沿卫星飞行轨道（即120° E）设置6个卫星接收站，站台纬度依次设为2°、5°、8°、11°、14°和17°，因此可以以高度和纬度进行网格划分，得到一个2维的区域．高度选取100 km至800 km，间隔20 km．纬度选为0°到20°，间隔1°．R矩阵我们通过经验给出，本实验取R=k·E，其中，E表示单位阵，其维数为所有卫星信号射线数量之和，k= 0.1× dim（E）．

对于系数矩阵的算法，此处假定的算法是：从第一条射线的 TEC开始，将每一条射线的TEC减去所有射线中最小的TEC，从而得到系数矩阵A．初始水平相关设为10°，垂直相关设为30 km，经验误差估计设为0.2．卫星的高度设为800 km，卫星纬度跨越为南纬3°至北纬23°，TEC采样间隔约为0.064°．在本文研究中，真实的电子密度以及初始的电子密度均由NeQuick模型[14]给出．

2.2 初始场对反演结果的影响分析

首先我们分析初始值对反演结果的影响．真实值为1月份UT6时，f 10.7（太阳10.7 cm辐射通量）为63；初始值为1月份UT6时，f 10.7分别设为68、73和78．得到电离层电子密度真实值与初始值以及它们的误差图像，如图1所示．从图1可以看出：初始值的设置对反演结果影响很大，初始值与真实值偏差越大，其反演结果的误差越大．

2.3 误差协方差矩阵的水平相关因子与垂直相关因子分析

误差协方差矩阵通常是根据经验给出，对反演结果的影响也比较大，其中水平和垂直相关因子与地理位置和昼夜变化相关．地理位置分低、中及高纬度讨论，假设低纬度的范围是0°到20°，中纬度是20°到40°，高纬度是40°到60°．对于不同的纬度，分为UT6时（下午）和UT18时（晚上）分别讨论．

讨论水平相关性时，取水平相关分别为4°、6°、8°、10°、12°、14°、16°、18°和20°；讨论垂直相关性时，取垂直相关分别为30 km、50 km、70 km、100 km、130 km、150 km、170 km、200 km、230 km、250 km、270 km和300 km．

图1 真实值图像、f 10.7分别为68、73及78时初始值的图像以及二者的误差图像

为了直观的看出水平相关和垂直相关对电子密度反演结果的影响，截取200 km到500 km这一个区域，计算每个网格点的反演误差，并统计总的平均误差，用来比对反演结果的好坏，分别给出低、中以及高纬度相关性误差分析的图像，如图2 – 图4所示：

图2 低纬度电子密度相关性的误差分析图像

图3 中纬度电子密度相关性的误差分析图像

图4 高纬度电子密度相关性的误差分析图像

对于低纬度，UT6时，随着水平相关和垂直相关的增大，误差整体上都呈现先增后减的趋势；UT18时，误差整体呈现减小的趋势．

对于中纬度，UT6时，随着水平相关和垂直相关的增大，误差整体呈现减小的趋势．UT18时，水平相关越大，误差会随之先减小后增大；垂直相关越大，误差会随之减小．

对于高纬度，UT6时，水平相关越大，误差整体呈现减小的趋势；垂直相关越大，误差先减小后增大．UT18时，随着水平相关和垂直相关的增大，误差整体呈现先减小后增大的趋势．

选取高纬度UT6时，将水平相关设为14°，垂直相关设为200 km，由公式（5）可以计算出每个网格点的电子密度，将其与真实电子密度先做差，再取绝对值，最后算出平均值为2.38%,反演效果图如图5所示．由图5不难看出，初始值与真实值的误差较明显，但反演值与真实值的差别较小，说明将水平相关设为 14°，垂直相关设为200 km时反演效果比较好．

图5 高纬度UT6时, 水平相关为14°、垂直相关为200 km的电子密度反演图像

3 结 论

本文介绍了三频信标技术的国内外研究现状、Kalman滤波原理与算法，利用模拟实验对电离层电子密度进行了反演，分析了初始值等因子对反演效应的影响．研究发现初始值、水平相关度以及垂直相关度等因素都会对电离层电子密度反演效果产生重要影响，其中水平相关设为14°、垂直相关设为200 km时反演效果较好．

[1] Na H, Shen J, Lee H. A fourier domain technique for ionospheric tomography [J]. Radio science, 1995, 30(3):747-754.

[2] Garriott O K, Nichol A W. Ionospheric information deduced from the doppler shifts of harmonic frequencies from earth satellites [J]. Journal of Atmospheric and Terrestrial Physics, 1961, 22(1): 50-63.

[3] Titheridge J E. Determination of ionospheric electron content from the faraday rotation of geostationary satellite signals [J]. Planetary and Space Science, 1972, 20(3): 353-369.

[4] Bernhardt P A, Siefring C L. New satellite-based systems for ionospheric tomography and scintillation region imaging [J]. Radio science, 2006, 41(5): 3-14.

[5] 徐继生, 叶公节. 东亚赤道异常区电离层CT诊断: 实验及初步结果[J]. 地球物理学报, 1995, 38(5): 553-563.

[6] 许正文, 吴健. Mexico草帽小波与电离层CT反演[J]. 电波科学学报, 1999, 14(3): 308-312.

[7] 吴雄斌, 徐继生. 一种改进的电离层层析成像算法[J]. 地球物理学报, 2000, 43(1): 19-28.

[8] 万卫星, 宁百齐, 袁洪, 等. 电离层扰动的GPS探测[J]. 空间科学学报, 1998, 18(3): 247-252.

[9] 徐继生, 邹玉华, 马淑英. GPS地面台网和掩星观测结合的时变三维电离层层析[J]. 地球物理学报, 2005, 48(4): 759-767.

[10] 吴健. 用三频卫星信标测量电离层天气新方法[J]. 中国科学: A辑, 2000, 30(S1): 111-114.

[11] 赵海生, 许正文, 吴健, 等. 三频卫星信标测量TEC方法探讨及数值模拟[J]. 装备环境工程, 2010, 7(4): 9-13.

[12] 许正文, 吴健, 曹冲. 电离层CT变换重建方法中的小波外推技术[J]. 空间科学学报, 1999, 19(2): 114-121.

[13] 秦永元, 张洪钺, 汪叔华. 卡尔曼滤波与组合导航原理[M]. 西安: 西北工业大学出版社, 1998: 33-48.

[14] Di Giovanni G, Radicella S M. An analytical model of the electron density profile in the ionosphere [J]. Advances in Space Research, 1990, 10(11): 27-30.

Research of Ionospheric Electron Density Inversion Based on Kalman Filter Technology

GUI Xiaochun, HONG Zhenjie, ZHAO Yunchao
(College of Mathematics and Information Science, Wenzhou University, Wenzhou, China 325035)

This paper probes into the experiment that a tri-band beacon data obtained from six satellites and earth stations has been developed to realize ionospheric electron density inversion based on Kalman filter technology and the factors that affect the inversion effect is also analyzed. It turns out that the factors like initial value, horizontal correlation and vertical correlation have a significant effect on ionospheric electron density inversion. When horizontal correlation is set as 14 degree and vertical correlation is set as 200 kilometer, the inversion effect is the best.

Kalman Filter; Ionosphere; Satellite Tri-band Beacon; TEC (Total Electron Content)

P353

A

1674-3563(2015)04-0056-06

10.3875/j.issn.1674-3563.2015.04.010 本文的PDF文件可以从xuebao.wzu.edu.cn获得

(编辑：封毅)

2014-12-22

公益性行业(地震)科研专项基金(201108004)

桂晓纯(1990- )，男，安徽宣城人，硕士研究生，研究方向：应用分析与最优化理论．† 通讯作者，hong@wzu.edu.cn