唐曉紅, 趙魯陽(yáng), 李魯明,3, 何 為, 王 沛

(1.中國(guó)科學(xué)院 上海微系統(tǒng)與信息技術(shù)研究所 寬帶無(wú)線技術(shù)實(shí)驗(yàn)室,上海 200050；2.上海師范大學(xué) 信息與機(jī)電工程學(xué)院,上海 200234；3.中國(guó)科學(xué)院大學(xué),北京 100049)

0 引 言

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

1 Allan方差分析法

1.1 Allan方差的基本原理

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

(1)

1.2 Allan方差的應(yīng)用

1)隨機(jī)誤差辨識(shí)

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

假設(shè)各項(xiàng)隨機(jī)噪聲相互獨(dú)立,根據(jù)Allan方差σ(τ)與MEMS陀螺儀隨機(jī)噪聲的對(duì)應(yīng)關(guān)系,則有

(2)

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

(3)

2)量測(cè)噪聲方差估計(jì)

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

(4)

則量測(cè)噪聲方差為

(5)

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

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

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

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

圖1 MEMS陀螺儀轉(zhuǎn)臺(tái)測(cè)試

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

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

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

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

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

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

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

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

3 模型建立

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

(6)

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

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

(7)

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

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

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濾波補(bǔ)償

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

(11)

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

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

(12)

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

Sk+Π·Sk-1=εk

(13)

(14)

(15)

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

(16)

5 實(shí)驗(yàn)與結(jié)果分析

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

(a) 標(biāo)準(zhǔn)Kalman濾波前后結(jié)果

(c) 改進(jìn)的解耦自適應(yīng)Kalman濾波前后結(jié)果圖3 三種算法濾波前后仿真結(jié)果

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

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

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

6 結(jié) 論

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