基于Kalman滤波器的卫星钟差预报精度分析＊

刘建成，杨睿峰，徐 赟，王茂磊，桑怀胜

（北京环球信息应用开发中心，北京100094）

0 引 言

星地时间同步技术对卫星导航系统的导航、授时以及定位精度有着直接的影响。星地无线电双向时间比对的工作原理是地面和卫星同时进行星地伪码测距，根据两组测距信息解算出卫星钟相对于地面站基准钟的钟差。这种方法抵消了两组测距共同的误差，所以，钟差测量精度很高。大多数实时导航用户采用的是广播星历给出的钟差，其精度将影响导航定位的精度，因此，研究星地无线电双向时间比对体制下的卫星钟差预报精度有着重要意义。

常用的钟差预报方法包括最小二乘方法［1－2］、FIR 滤 波［1，3］、Kalman 滤 波［1，4－7］。 当 钟 差 状 态 模型是线性的，且状态噪声统计特性和测量噪声统计特性已知时，Kalman滤波性能最优，Kalman滤波在钟差预报方面获得了越来越多的应用。对于Kalman滤波器，估计误差方差随着时间的递推逐渐减小，直至达到稳态，此时估计误差方差最小，在稳态条件下得到的预报误差也最小，因此，稳态情况下的钟差预报精度可作为分析星地时间同步性能的重要依据。

Meditch研究了卫星钟差状态噪声为频率白噪声情况下的Kalman滤波器钟差预报精度，利用Kalman滤波器状态转移矩阵得到了最小预报误差的解析表达式［8］。Beisner研究了卫星钟状态噪声为频率白噪声、频率闪烁噪声和频率随机游走噪声且不考虑测量噪声情况下的Wiener预测器的最小钟差预报误差，利用Wiener方法直接从频域得到了比Meditch得到的更具有普适性的解析形式。目前，Kalman滤波器采用的钟差状态噪声逐渐扩展到3阶［9－11］，不但包括频率白噪声、频率随机游走噪声，还包括频率漂移随机游走噪声。研究钟差状态噪声包括3阶白噪声且测量噪声为白噪声情况下的基于Kalman滤波器的卫星钟差稳态精度。

根据Kalman滤波器达到稳态时的条件建立方程组，然后解该方程组，可以得到稳态解。该方法的优点是概念直观，但当目标状态变量增多时解方程组存在很大困难。通过消除Kalman递归形式中的增益项可得到预测协方差矩阵的递归公式，即Riccati方程，Vaughan解决了离散Riccati方程的非递归代数解问题［12］，利用这个结果可直接得到稳态解，绕开了解方程组的问题。在建立卫星钟差的状态方程和测量方程的基础上，利用这种方法得到Kalman滤波器的稳态解，进一步得到预报误差，建立起稳态情况下的钟差预报误差与卫星钟状态误差、测量误差及采样时间间隔之间的解析关系，最后做了典型参数情况下的数值分析。

1 卫星钟差Kalman滤波方程

采用3状态卫星钟差模型，包括相位、钟速和频率漂移值。考虑3类钟差噪声［10］，包括频率白噪声、频率随机游走噪声和频率漂移随机游走噪声。

若采样时间间隔为T，卫星钟差状态向量为xk＝ ［a0a1a2］T，其中a0为钟差，a1为钟速，a2为频率漂移量，那么建立多项式状态方程为［10］

式中状态转移矩阵误差矩阵为

式中q1、q2、q3分别为频率白噪声、频率随机游走噪声和频率漂移随机游走噪声的谱密度。

根据星地无线电双向时间传递的原理，钟差测量方程用xk表示为

从状态方程和测量方程可以看出，两者均为线性方程，则最优状态估计就是卡尔曼滤波，其误差协方差阵递推过程为

稳态协方差矩阵˜P定义为

式中:σ211为钟差预测误差方差；σ212为钟差、钟速预测误差相关系数；σ213为钟差、频率漂移量预测误差相关系数；σ222为钟速预测误差方差；σ223为钟速、频率漂移量预测误差相关系数；σ233为频率漂移量预测误差方差。

2 卫星钟差的稳态精度和预报精度

推导基于Kalman滤波器的卫星钟差稳态精度，得到卫星钟差预报误差协方差阵与稳态协方差矩阵和卫星钟状态误差矩阵的关系。

2．1 Kalman滤波器的稳态解

根据文献［12］提出的离散Riccati方程的非递归代数解算法，构造矩阵

式中F、H、R、Q由式（1）和式（3）确定。

由于星钟的系统模型为三阶，所以Kf为6阶特征值。若矩阵的特征值为λi（i＝1，2，3），则其对应的特征向量为

稳态解是S1＝λ1＋λ2＋λ3、S2＝λ1λ2＋λ2λ3＋λ3λ1和S3＝λ1λ2λ3的函数。下面说明如何确定S1、S2、S3．

矩阵Kf的特征方程为

由矩阵Kf的性质得到

式中λ1、λ2、λ3均为模大于1的特征值。

把式（11）展开，并与式（10）比较，得到

令S3＝x2，采用文献［13］的做法，解上式方程组得到

式中的x是式（13）代入式（12）而得到的一元四次方程的解的最大值

从式（14）可以看出，基于Kalman滤波器的卫星钟差稳态精度取决于参数q1、q2、q3、T、σ2．

2．2 卫星钟差预报精度

假设Kalman滤波器开始预测工作前已经达

等式右端第一部分反映了Kalman滤波稳态误差对预报精度的影响，取决于参数q1、q2、q3、T、N、σ2．第二部分反映了卫星钟的状态误差对预报精度的影响，取决于参数q1、q2、q3、T、N．

如果定义间隔为T时的状态噪声矩阵函数为Q［T］，由式（2）确定，那么NT 预报误差协方差阵还可表示为

式中NT卫星钟差的预报误差方差为σ211（N）＝P˜N（1，1）．

当不考虑频率随机游走噪声和频率漂移随机游走噪声时，即q2＝q3＝0，此时β＝γ＝0，x＝到稳态，则钟差参数预报初始误差协方差阵为˜P1＝˜P．根据Kalman滤波器预测方程，预报误差协方差阵的递推关系为

根据此递推关系，得到卫星钟差状态向量的N步预报误差协方差阵，即卫星钟差状态向量的NT预报误差协方差阵为度，再根据式（17）或（18）得到钟差预报误差方差，结果与文献［8］一致。

3 数值分析

根据上面得出的理论结果，采样时间间隔、伪码距离测量误差均方根对基于Kalman滤波器的钟差稳态误差均方根和钟差预报误差均方根的影响做了数值分析。数值分析采用的卫星钟状态误差参数为［11］:q1＝1．11×10－22s2／s、q2＝2．22×10－32s2／s3、q3＝6．66×10－45s2／s5．

3．1 钟差稳态精度与伪距测量精度的关系

分析了采样时间间隔分别为5s、10s、15s时，基于Kalman滤波器的钟差稳态误差均方根随伪距测量误差均方根的变化关系，如图1所示。从图1可以看出，伪距测量误差均方根越大，钟差稳态误差均方根就越大。采样时间间隔越大，钟差稳态误差均方根就越大。因此，提高钟差稳态精度的方法包括提高伪距测量精度和提高数据率。

3．2 钟差预报精度与伪距测量精度的关系

图1 钟差稳态精度与伪距测量精度的关系曲线

分析了采样时间间隔分别为5s、10s、15s时，基于Kalman滤波器的钟差预报误差均方根随伪距测量误差均方根的变化关系。1h和2h的钟差预报误差均方根与伪码测量误差均方根的关系曲线分别如图2（a）和图2（b）所示。

从图2可以看出，伪距测量误差均方根越大，钟差预报误差均方根就越大。采样时间间隔越大，钟差预报误差均方根就越大。

图2 钟差预报精度与伪距测量精度的关系曲线

4 结 论

研究了卫星钟差包括3阶白噪声情况下的基于Kalman滤波器的卫星钟差稳态精度。在建立卫星钟差的状态方程和测量方程的基础上，利用离散Riccati方程的非递归代数解，得到了Kalman滤波器的稳态解，进一步得到了卫星钟差预报误差，建立起钟差预报误差与卫星钟状态误差、测量误差及采样时间间隔之间的解析关系。做了典型参数情况下的数值分析。下一步的工作是分析卫星钟的闪烁噪声对基于Kalman滤波器的卫星钟差预报精度的影响程度。

［1］ 刘 利．相对论时间比对理论与高精度时间同步技术［D］．郑州:解放军信息工程大学，2004．

［2］ DAVIS J A，HARRIS P M，et al．Least－squares analysis of two－way Satellite time and frequency Transfer Measurements［C］∥Proceedings of the 33thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，Long Beach，California，USA，2001（11）:121－132．

［3］ SHMAL I Y，IBARRA－MANZHNO O．An optimal FIR filtering algorithm for a time error model of a precise clock［C］∥Proceedings of the 34thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，Reston，Virginia，USA，2002（12）:53－68．

［4］ DAVIS J A，STACEY P W，et al．Combining time transfer measurements using a Kalman filter［C］∥Proceedings of the 34thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，Reston，Virginia，USA，2002（12）:53－68．

［5］ DAVIS J A，GREENGHALL C A，STACEY P W．A Kalman filter clock algorithm for use in the presence of flicker frequency modulation noise［C］∥Proceedings of the 35thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，Reston，Virginia，USA，2003（12）:53－68．

［6］ 陈小敏，李孝辉．基于自适应Kalman滤波器的原子钟信号预测方法［J］．吉林大学学报·理学版，2009，47（3）:591－593．

［7］ GALLEANI L，TAVELLA P．Time and the Kalman filter－applications of optimal estimation to atomic timing［J］．IEEE control systems magazine，2010，30（4）:44－65．

［8］ BEUSMER H M．Clock error with a Wiener predic－tor and by numerical Calculation［J］．IEEE Transactions on Instrumentation and measurement，1980，29（2）:105－113．

［9］ HUTSELL S T，REID W G，GRUM J D，et al．Operational use of the Hadamard variance in GPS［C］∥Proceedings of the 28thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，Reston，Virginia，USA，1996（12）:201－214．

［10］ HUTSELL S T．Relating the Hadamard variance to MCS Kalman filter clock estimation［C］∥Proceedings of the 27thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，San Diego，California，USA，1995（12）:291－302．

［11］ HUTSELL S T．Fine tuning GPS clock estimation in the MCS［C］∥Proceedings of the 26thAnnual Precise Time and Time Interval（PTTI）Applications and Planning Meeting，Reston，Virginia，USA，1994（12）:63－74．

［12］ VAUGHAN D R．A non－recursive algebraic solution for the discrete riccati equation［J］．IEEE Transactions on Automatic Control，1970，15（1）:597－599．

［13］ RAMACHANDRA K V，MOHAN B R，GEETHA B R．A three－state Kalman tracker using position and rate measurements［J］．IEEE Transactions on Aerospace and Electronic Systems，1993，29（1）:215－222．