唐曉紅, 趙魯陽, 李魯明,3, 何 為, 王 沛

(1.中國科學院 上海微系統(tǒng)與信息技術研究所 寬帶無線技術實驗室,上海 200050；2.上海師范大學 信息與機電工程學院,上海 200234；3.中國科學院大學,北京 100049)

0 引 言

微機電系統(tǒng)(micro-electro-mechanical system,MEMS)陀螺儀[1,2]是慣性導航系統(tǒng)中測量載體角速率的關鍵性器件,性能主要是受到常值誤差及隨機誤差的影響,抑制隨機誤差,通常需要建模及濾波等方式來進行補償。本文采用Allan方差分析方法對MEMS陀螺儀輸出數(shù)據(jù)進行隨機誤差辨識,建立自回歸滑動平均(auto-regressive moving average,ARMA)數(shù)學模型。針對模型中的有色噪聲,無法采用狀態(tài)擴充法進行白化[3],本文提出將其引入系統(tǒng)狀態(tài)方程中,作為方程的控制項,對解耦自適應Kalman濾波算法進行改進,實現(xiàn)對MEMS陀螺儀隨機噪聲的濾波補償。

1 Allan方差分析法

1.1 Allan方差的基本原理

對于時間序列{xt},其樣本數(shù)據(jù)的總長度為N,采樣時間為τ0,把樣本數(shù)據(jù)平均分為K組,每組樣本數(shù)據(jù)的長度為M<(N-1)/2。將τ=mτ0,(m=1,2,…,M)定義為相關時間,每一個相關時間下,每一組樣本數(shù)據(jù)的平均值為xj(M),則Allan方差計算方法為

(1)

1.2 Allan方差的應用

1)隨機誤差辨識

MEMS陀螺儀隨機噪聲主要包括:量化噪聲Q(μrad)、角度隨機游走N((°)/h1/2)、零偏不穩(wěn)定性B((°)/h)、角速率隨機游走K((°)/h3/2)及速率斜坡R((°)/h2)[4]。各隨機誤差項在Allan方差曲線中,直線斜率不同。

假設各項隨機噪聲相互獨立,根據(jù)Allan方差σ(τ)與MEMS陀螺儀隨機噪聲的對應關系,則有

(2)

將式(3)中的多項式系數(shù)采用最小二乘法進行擬合有

(3)

2)量測噪聲方差估計

Allan方差濾波器可直接濾除系統(tǒng)噪聲中部分低頻噪聲,量測噪聲的方差近似等于寬頻噪聲的Allan方差[6]。相關時間為采樣時間τ0時,Allan方差為

(4)

則量測噪聲方差為

(5)

式中k=2,3,…;初值1可任意選取。

2 數(shù)據(jù)采集、預處理與檢驗

2.1 數(shù)據(jù)采集

本實驗所測試的MEMS陀螺儀為MPU6050,使用STM32F030CCT6通過I2C總線讀取MPU6050的數(shù)據(jù)。實驗數(shù)據(jù)通過串口模塊、藍牙模塊進行接收、發(fā)送,由MATLAB讀取并記錄串口數(shù)據(jù)。在AeroSmith高精度單軸轉臺上,安裝MEMS陀螺儀,如圖1所示,靜止預熱20 min。設定采間隔時間為10 ms,以100 Hz的采樣頻率連續(xù)采集90 min數(shù)據(jù)。

圖1 MEMS陀螺儀轉臺測試

2.2 數(shù)據(jù)預處理

在建立MEMS陀螺儀樣本數(shù)據(jù)的時間序列模型前,對MEMS陀螺儀原始數(shù)據(jù)進行預處理以確保樣本數(shù)據(jù)的零均值、平穩(wěn)性以及正態(tài)性,通過提取樣本均值的方式進行消除常值誤差。由于周圍環(huán)境以及MEMS陀螺儀自身因素的干擾,會產(chǎn)生一定的趨勢性變化,可采用多項式擬合法[7]提取趨勢項。依據(jù)PauTa準則[8],若超過3倍測量值的標準偏差,則認為該數(shù)據(jù)是異常值,從測量數(shù)據(jù)中剔除。本文采用樣本數(shù)據(jù)均值代替異常值的方式來剔除異常值,預處理前后MEMS陀螺儀的角速度數(shù)據(jù)如圖2所示。

圖2 MEMS陀螺儀預處理前后數(shù)據(jù)

2.3 平穩(wěn)性及正態(tài)性檢驗

1)平穩(wěn)性檢驗

預處理后的陀螺儀數(shù)據(jù)采用輪次檢驗法[9]進行平穩(wěn)性檢驗。樣本數(shù)據(jù)的游程總數(shù)為r=16,查游程檢驗分布表,在顯著水平α=0.05下,其接收域為6≤r≤16。故樣本數(shù)據(jù)為平穩(wěn)的隨機信號。

2)正態(tài)性檢驗

將樣本數(shù)據(jù)分為K組等間隔區(qū)間,統(tǒng)計其在各組區(qū)間的個數(shù)。MEMS陀螺儀輸出角速度數(shù)據(jù)的頻數(shù)直方圖基本服從正態(tài)分布。

3 模型建立

對于平穩(wěn)時間序列樣本{xt},t=1,2,…,N,其ARMA(p,q)為

(6)

式中N為數(shù)據(jù)樣本長度；εt～WN(0,σ2)；p,q為模型階數(shù)。通常,MEMS陀螺儀誤差建模不會超過3階,且p≥q。則,可建立模型ARMA(1,1),ARMA(2,1),ARMA(2,2)。

赤池信息(Akaike information criterion,AIC)準則為

(7)

則3種ARMA模型的AIC值分別為-6.623 309,-6.623 351,-6.623 347,且ARMA(2,1)模型的AIC值最小,可作為MEMS陀螺儀的隨機誤差模型。

本文采用自回歸逼近法[10]進行ARMA(2,1)模型的參數(shù)估計:

1)建立高階的AR(6)模型,逼近ARMA(2,1)模型

xt=0.022 7xt-1+0.009 1xt-2+0.004 9xt-3+

0.007 3xt-4+0.007 7xt-5+0.008xt-6+εt

(8)

0.007 3xt-4-0.007 7xt-5-0.008xt-6

(9)

3)建立的ARMA(2,1)模型為

xt=0.917 1xt-1-0.011 2xt-2+εt-0.894 5εt-1,

εt～WN(0,0.001 33)

(10)

4 Kalman濾波補償

設A為狀態(tài)轉移矩陣,B為系統(tǒng)噪聲矩陣,H為量測矩陣,Wk和Vk分別為系統(tǒng)過程噪聲和量測噪聲,根據(jù)所建的ARMA(2,1)模型,濾波器的參數(shù)設置如下

(11)

基于ARMA模型的系統(tǒng)狀態(tài)方程亦可寫為

Xk=A·Xk-1+Sk,Sk=[εk+b1·εk-10]T

(12)

可知系統(tǒng)噪聲Sk為有色噪聲。通常采用狀態(tài)擴充法,將系統(tǒng)有色噪聲擴充到狀態(tài)量中,進行白化有色噪聲。而該方法對有色噪聲進行白化的成形濾波器,其固定形式為

Sk+Π·Sk-1=εk

(13)

(14)

(15)

根據(jù)ARMA模型可得系統(tǒng)過程協(xié)方差Q,本文將固定Q陣,只對量測噪聲協(xié)方差R陣進行自適應,則改進的解耦自適應Kalman濾波方程為

(16)

5 實驗與結果分析

為了驗證本文所提的改進的解耦自適應Kalman濾波算法的正確性與有效性,根據(jù)建立的ARMA(2,1)模型及濾波方程,分別對Kalman濾波算法,解耦自適應Kalman濾波算法[11,12]和本文算法進行仿真,如圖3所示。

(a) 標準Kalman濾波前后結果

(c) 改進的解耦自適應Kalman濾波前后結果圖3 三種算法濾波前后仿真結果

可見本文所提改進算法、解耦自適應Kalman濾波的效果均好于標準Kalman濾波。MEMS陀螺儀濾波前后數(shù)據(jù)的統(tǒng)計特性均方差:表明3種算法濾波后的噪聲均方差為:Kalmax濾波算法為0.000 637,解耦Kalman濾波算法為0.000 291,本文算法為0.000 198,相比濾波前均方差(0.001 33)均減小了一個數(shù)量級。與濾波前相比,Kalman濾波后的噪聲均方差減小52.1 %,解耦自適應Kalman濾波后的噪聲均方差減小78.1 %,改進的解耦自適應Kalman濾波的噪聲均方差減小85.1 %。說明本文所提改進算法較其他2種Kalman濾波具有更好的濾波效果。

MEMS陀螺儀的短期誤差表現(xiàn)為角度隨機游走,長期誤差表現(xiàn)為角速率隨機游走[11]。為了更進一步驗證本文所提改進算法的正確性與有效性,對MEMS陀螺儀濾波前后的隨機誤差項進行對比分析,如表1所示。

根據(jù)表1,經(jīng)改進的解耦自適應Kalman濾波算法濾波后,2種誤差項皆可有效減小41.9 %。相較于Kalman濾波算法、解耦自適應Kalman濾波算法,本文所提改進算法的隨機誤差項分別減小30.5 %,10.1 %。表明本文所提改進算法可更好地提高MEMS陀螺儀的輸出精度,減小隨機噪聲對MEMS陀螺儀性能的影響。

6 結 論

本文對解耦自適應Kalman濾波算法改進,對改進算法、Kalman濾波算法、解耦自適應Kalman濾波算法進行仿真,從數(shù)據(jù)的統(tǒng)計特性均方差、MEMS陀螺儀的隨機誤差項,對濾波前后的數(shù)據(jù)進行對比分析。實驗結果表明,本文所改進的算法可有效減小MEMS陀螺儀隨機誤差,其濾波效果優(yōu)于Kalman濾波算法、解耦自適應Kalman濾波算法,對提高MEMS陀螺儀的精度具有一定的工程實用性。