顾及有色噪声的光纤陀螺随机噪声自适应滤波方法
靳凯迪
信息工程大学地理空间信息学院, 河南 郑州 450001
基金项目：国家自然科学基金(42074014;41574010)
关键词：光纤陀螺 随机误差 ARIMA模型 有色噪声 自适应滤波
靳凯迪, 柴洪洲, 宿楚涵, 等. 顾及有色噪声的光纤陀螺随机噪声自适应滤波方法. 测绘学 ，2022，51(1)：80-86. DOI:
10.11947/j.AGCS.2022.20200614
JIN Kaidi, CHAI Hongzhou, SU Chuhan, et al. Adaptive Kalman filter method with colored noise for fiber optic gyroscope random drift. Acta Geodaetica et Cartographica Sinica, 2022, 51(1): 80-86. DOI: 10.11947/j.AGCS.2022.20200614
阅读全文：
http://xb.sinomaps.com/article/2022/1001-1595/2022-1-80.htm
引 言
由于元件构造和环境影响，光纤陀螺(fiber optic gyroscope，FOG)会产生复杂的噪声源，直接影响惯性导航系统初始对准和导航解算整个过程的精度。光纤陀螺噪声主要由系统性噪声与随机噪声组成，系统性噪声可通过试验标校有效补偿，随机噪声由于性质复杂、随机性强、补偿难度大，成为衡量光纤陀螺精度的重要指标[1-3]。因此，对光纤陀螺随机噪声进行数学建模并通过滤波等手段加以补偿，对准确分析光纤惯性导航系统误差和提高导航解算精度等方面具有重要意义。目前，用于光纤陀螺随机误差建模的方法主要包括Allan方差分析、功率谱密度分析和时间序列分析等[4-8]。文献[9]使用自回归滑动平均模型(auto-regressive moving average，ARMA)建模光纤陀螺随机噪声，并给出原始观测序列的Kalman滤波方程，但并未考虑MA模型存在时状态噪声为有色噪声的情况。文献[10]使用高阶AR模型计算ARMA模型的噪声估值，达到白化有色噪声的目的，但其作为一种近似方法并不严密。后续研究中或是将随机噪声建模为AR模型[11]，或是沿用文献[9-10]的思想[12-13]。针对随机噪声滤波模型先验信息不准确的问题，文献[11]使用Sage-Husa自适应滤波(Sage-Husa adaptive Kalman filter，SHAKF)在线估计系统噪声和量测噪声。文献[12]使用动态Allan方差分析法(dynamic Allan variance，DAVAR)估计量测噪声。文献[7]针对自适应滤波参数耦合问题，使用DAVAR和SHAKF分别估计状态噪声和量测噪声。然而，SHAKF易丧失矩阵正定性或非正定性[14-16]，尤其在ARIMA模型收敛初期，此现象更加显著。
1 光纤陀螺随机噪声建模
ARIMA模型实质上是差分运算与ARMA模型的融合。ARIMA建模光纤陀螺随机噪声序列的主要内容包括噪声随机序列平稳性检验、ARIMA模型阶数确定、模型参数估计及适用性检验等[17-20]。ARMA模型要求平稳的时间序列，因此，首先要对采样光纤陀螺原始随机噪声进行平稳性和随机性检验，如检验不通过，则需进行差分直至满足平稳随机过程。使用中海达iPos光纤惯导以100 Hz采样频率采集静态数据如图 1所示。采用单位根检验(包括ADF检验与KPSS检验)对信号检测通过，因此，ARIMA模型差分阶数为0。设ARIMA(p, 0, q)模型如下[21]
(1)
图 1 FOG随机噪声原始输出信号Fig. 1 Original random noise signals of FOG
式中，xk为历元k时随机噪声值；ai为自回归项系数；bj为滑动平均项系数；εk为未知方差的白噪声；c为常值。
使用前2000个历元数据建立ARIMA模型，计算过程分为模型阶数(p, q)确定和模型系数(ai，bj)计算。(p, q)可以通过时间序列的偏自相关函数图确定，但由偏自相关函数图得到的ARIMA模型阶数往往不准确。因此采用赤池信息准则(akaike information criterion，AIC)与贝叶斯信息准则(Bayesian information criterionm，BIC)确定随机噪声序列ARIMA阶数[6]，计算结果为p=3，q=3。模型阶数确定后，通过最小二乘法计算模型系数值，所得ARIMA模型见表 1。
表 1 光纤陀螺随机噪声ARIMA(p, 0, q)模型参数
Tab. 1 Parameters for ARIMA(p, 0, q) of FOG random nosie
图 2为对所建模型采用标准化残差检验结果。由图 2可知，所得残差序列符合正态分布。Durbin-Watson统计量是常用的自相关统计量，设残差序列为vk(k=1, 2, …, n)，构造DW统计量为[22]
(2)
图 2 标准化残差检验Fig. 2 Standardized residuals test
式中，ρ为自相关因子。ρ越接近0，等价于DW统计量越接近2，序列越不存在一阶相关性，本次残差序列相关性检验所得结果为1.953 9，认为残差序列符合高斯正态分布。综上，认为所建模型较为准确，通过适用性检验。
2 ARIMA模型有色噪声Kalman滤波
在光纤陀螺随机噪声建模为ARIMA模型后，使用Kalman滤波是进行随机噪声降噪的有效方法。以往处理存在MA系数的ARIMA模型时，往往忽视噪声的相关性或使用噪声估值进行白化。为此，建立基于Harvey法的Kalman滤波模型白化ARIMA有色噪声。
2.1 ARIMA模型滤波方程
获得如式(1)所示的随机噪声ARIMA(p, 0, q)模型后，令状态向量Xk=[xk xk-1 … xk-p+1 c]T，建立Kalman滤波状态方程和量测方程为[9]
(3)
式中，Vk为均值为0的高斯分布量测噪声，且与状态噪声Wk= [εk εk－1 … εk－q]T不相关。
2.2 Harvey法有色噪声滤波方程
当存在MA系数时，式(3)中的状态噪声为有色噪声，此时，Kalman滤波将会导致结果发散或失真[14, 21]。传统上，通常将有色噪声建模为一阶AR模型[23-24]，但ARIMA的有色噪声数学模型是明确的，对噪声再进行时间序列建模并不完全准确。在标准Harvey法滤波方程中[25]，针对如下所示ARIMA模型
(4)
令y1, k=xk，有
(5)
考虑y2, k
(6)
递推可得如下关系
(7)
式中，bi=0(i=q+1, q+2, …, p－1)；m=1, 2，…, p－1。构建状态方程为
(8)
(9)
继续进行递推
(10)
令Yk=[y1, k y2, k … yp, k c]T，得到的滤波方程为
(11)
式中
(12)
式中，εk和Vk均为白噪声且两者不相关，满足随机模型
(13)
扩展后的Harvey方程多出了一维由“常值+白噪声”构成的随机常值状态参数。实现了附带常数项，且q=p的ARIMA模型状态方程的构建。通过Harvey法构造的状态噪声为高斯白噪声，满足标准Kalman滤波条件。
3 自适应Kalman滤波
3.1 ARIMA模型及状态噪声在线更新
光纤陀螺随机噪声复杂且时变，实际应用中使用有限数量或不同时期的数据所得到的ARIMA模型随时间的推移，误差会逐渐增大。因此，有必要依据当前观测数据进行ARIMA参数在线更新。将式(1)改写为
(14)
式中
(15)
使用序贯平差实现模型参数在线更新的主要流程如下
(16)
实际应用中，
(17)
3.2 变分贝叶斯自适应滤波
ARIMA模型在线更新中，参数估值
基于Allan方差的滤波器是一种带通滤波器，可直接滤除部分低频系统噪声的特性，传统上，将Allan方差值近似为宽频量测噪声的方差[1]，选择取样间隔为最短采样时间的Allan方差为
(18)
对于上文所述的量测模型，量测值Zk：N(xk, R)，其中N(μ, ∑)表示均值为μ和协方差阵为∑的高斯概率密度函数，由式(18)可得
(19)
式中，Δk为当历元k时的真误差；
为动态Allan方差值；R为真实量测噪声方差。由于白噪声序列前后历元不相关，式(19)中第3行的
。由式(19)可知，Allan方差的量测噪声估值与真实量测噪声间存在偏差项
(20)
VBAKF是一种近似的贝叶斯法，它利用已知的模型信息、观测信息和先验信息来获得状态向量和未知参数联合后验密度的近似解。可同时修正Qk与Rk的VBAKF的主要流程[15]如下。
(1) 状态更新过程
(21)
式中，
为k－1历元参数平差值协方差阵，
为参数预测值协方差阵。
选择
(22)
式中，IW分布的两参数分别为自由度参数和逆尺度矩阵。
(2) 变分量测迭代更新过程。
①参数初始化：
② 更新
③ 更新Yk后验概率分布q1(Yk)=N(
返回②继续迭代。
(3) 输出
VBAKF通过估计参数预测值协方差阵
4 试验与分析
对图 1所示的光纤陀螺随机噪声进行试验分析，原始陀螺随机噪声Allan方差曲线如图 3所示，可见该陀螺随机漂移主要由角度随机游走、零偏不稳定性组成。图 4为通过序贯平差建立的ARIMA系数在线更新，在试验条件下，光纤陀螺随机噪声序列的ARIMA模型较为稳定，估值很快收敛，因此，在实际测量中，可仅于初始阶段进行模型参数更新，以提升计算效率。
图 3 光纤陀螺随机噪声Allan方差Fig. 3 Allan-variance of FOG random noise
图 4 光纤陀螺随机噪声ARIMA系数更新
Fig. 4 ARIMA coefficients update for FOG random noise
将动态Allan方差应用于量测噪声自适应估计是一种近似方法。为对比DAVAR与VBAKF对量测噪声的估计能力，分别在原始训练数据中，增加均方差为1×10-3、1.5×10-3、2×10-3(°)/s的高斯白噪声进行解算，量测噪声估计结果如图 5所示，应认为原始训练数据与所建ARIMA模型相符性较好，因此，量测噪声与添加噪声大小近似相同，在图 5中第1000历元，VBAKF量测噪声估值分别收敛至1.05×10-3、1.44×10-3、1.95×10-3(°)/s，DAVAR量测噪声估值收敛至1.46×10-3、1.84×10-3、2.36×10-3(°)/s。因此，可认为VBAKF较DAVAR对量测噪声有更好的估计性能。使用Harvey法构造Kalman滤波方程，结合VBAKF(以Harvey VBAKF代替)对光纤陀螺随机噪声进行实时滤波处理。同时，分别使用不考虑有色噪声的滤波方程(以KF代替)及不使用VBAKF的Harvey法滤波方程(以Harvey代替)对同一数据进行处理，原始数据与滤波结果如图 6所示，数据统计特性见表 2。
图 5 DAVAR与VBAKF噪声自适应对比
Fig. 5 Adaptive comparison of DAVAR and VBAKF measurement noise
图 6 滤波输出数据比较Fig. 6 Comparison of filtering output data
表 2 滤波结果统计特性对比
Tab. 2 Comparison of statistical results for filtering data
声明：本站部分文章内容及图片转载于互联 、内容不代表本站观点，如有内容涉及侵权，请您立即联系本站处理，非常感谢！