基于Allan方差的MEMS陀螺儀隨機誤差辨識與抑制*

馬 群， 王 慶， 陽 媛， 盛 浩

(1.東南大學 儀器科學與工程學院，江蘇 南京 210096；2.北京航空航天大學 計算機學院，北京 100191)

0 引 言

微機電系統(micro-electro-mechanical system,MEMS)陀螺儀具有體積小、成本低、高可靠性等優點，非常適合應用于機器人、車輛、小型無人機的姿態測量和導航中[1]。在實際應用中，不同MEMS陀螺儀往往受到制作工藝、環境因素制約，精度難以達到理想要求。針對這一問題，目前有效的解決方法是通過算法提高陀螺儀輸出精度。本文針對一MEMS陀螺儀靜態采集的原始數據，運用Allan方差法辨識其隨機誤差成分，建立隨機誤差的自回歸滑動平均(auto-regressive and moving average,ARMA)(1,1)模型，并基于該擬合模型參數設計卡爾曼(Kalman)濾波器，對原始數據進行預處理，處理后數據中包含的隨機誤差成分均得到有效衰減。

1 Allan方差分析法

1.1 Allan方差計算過程

通過靜態實驗對陀螺儀連續采集一段時間的原始數據，得到總數為的數據樣本，然后將平均角速率樣本劃分為包含不同數量采樣點的樣本組，分組方法如圖1所示。

圖1 Allan方差樣本分組示意

設第一組平均角速率樣本序列為

(1)

式中τ0為原始數據采樣間隔。

計算第一組數據的Allan方差

(2)

擴大采樣間隔，使τ1=2τ0，在相繼奇偶序號角速率之間取算數平均，即

(3)

此時樣本數N1=[N/2](取整)，組成新的平均角速率序列

(4)

計算其Allan方差，即

(5)

1.2 Allan方差與隨機誤差組成關系

Allan方差與MEMS陀螺儀隨機誤差的功率譜密度(power spectral density,PSD)有唯一的確定關系

(6)

式中SΩ(f)為隨機誤差的功率譜密度函數，sin4(πfτ)/(πfτ)2為誤差的Allan方差帶通濾波器。當τ發生變化時，不同頻率的信號可以在不同時間上依次通過濾波器[3]。假設MEMS陀螺儀的隨機噪聲包含統計獨立的量化噪聲、角度隨機游走、零偏不穩定性、角速率隨機游走和速率斜坡，則Allan方差可表示為

(7)

式中Q為量化噪聲系數，N為角度隨機游走系數，B為零偏不穩定系數,K為角速率隨機游走系數，R為速率斜坡系數[4]。將式兩邊同取對數，近似得

(8)

與冪律譜的f-S雙對數曲線類似，也可近似將Allan標準差的τ～σA(τ)雙對數曲線看作是由多段不同斜率線段首尾連接而成[5]。圖2為曲線斜率與誤差成分關系。

圖2 曲線斜率與誤差成分關系

(9)

(10)

利用此方法可得每一段函數的系數ak[6]。隨機誤差系數及其與ak轉換關系如表1所示。

表1 隨機誤差系數與ak對應關系[7]

1.3 數據采集與誤差辨識實驗

本文的實驗對象是國產MEMS-IMU—LPMS-USBAL2，該IMU包含3軸陀螺儀、3軸加速度計、3軸磁力計。將其固定在3軸捷聯慣導測試轉臺上，在25°C恒溫條件下，以100 Hz采樣頻率靜態采集2 h以上。

根據前述分段最小二乘法擬合Allan方差曲線，得出的誤差項系數值如表2所示，將表中參數與該型IMU用戶手冊上的參數對照，兩者差別不大，表明辨識結果準確。

表2 誤差成分系數擬合結果

2 MEMS陀螺儀隨機誤差建模

經過異常值剔除，分離趨勢項后，陀螺儀輸出是一組零均值、平穩、正態的時間序列[8]，可以通過建立ARMA模型來擬合陀螺儀數據的隨機誤差特性。陀螺儀輸出的時域遞推形式為

(11)

式(11)為ARMA(p,q)模型，其中,w(n)為時刻高斯白噪聲，x(n)為n時刻觀測到的陀螺儀輸出。當p=0時，有

x(n)=w(n)+b1w(n-1)+b2w(n-2)+…+

bqw(n-q)

(12)

此為MA(q)模型，表示x(n)為由w(n)與其既往連續q個時移序列加權平均而得。

而當q=0時，有

x(n)=a1x(n-1)+a2x(n-2)+…+apx(n-p)+

w(n)

(13)

此為AR(p)模型，表示x(n)是序列自身前p個觀測值的線性回歸再加上時刻白噪聲w(n)的結果[9]。

2.1 模型形式

平穩時間序列的自相關函數、偏自相關函數特性與ARMA模型形式、階次(p,q)之間有著明確的關系[10，11]，因此，可依據圖3的流程來判斷樣本所適用的模型形式。

圖3 ARMA模型形式判斷流程

圖4 z軸樣本自相關和偏自相關函數

2.2 階次判定

當樣本數量過多時，BIC的懲罰項比AIC的多，可有效防止模型精度過高造成的模型復雜度過高[12]。因此本文選用貝葉斯準則來判定模型階次。BIC定義為

BIC=kln(n)-lnL(X)

(14)

式中k為模型參數個數，n為樣本數量，kln(n)為對模型復雜度的懲罰項，lnL(X)為擬合模型在樣本集上的極大后驗似然度，當BIC取得最小值時的擬合模型即為最優模型[13]。實際應用中的階次p，q不超過4階，當p，q分別取1～4時，從z軸樣本的BIC值中數據可以看出:ARMA(1,1)模型的BIC值最小(-52524)，因此選取該模型擬合陀螺儀Z軸的隨機誤差特性最優。

2.3 參數估計

由于ARMA(1,1)模型的階次較小，所以,采用計算量較小的矩估計方法估計模型參數比較方便[14]。經計算得出陀螺儀z軸輸出模型參數為

代入式(11)得z軸輸出模型為

X(n)=0.2864X(n-1)+w(n)-0.1869w(n-1)

(15)

3 卡爾曼濾波

卡爾曼濾波器是一種基于信號和噪聲的時域統計特性進行最優估計的方法，可以實現最小均方估計誤差意義下隨機信號的最優線性濾波，結合建立的陀螺儀數據隨機誤差模型，將陀螺儀隨機誤差作為系統輸入，建立離散隨機信號的狀態空間模型

Xk=ΦXk-1+BWk,Zk=HXk+Vk

(16)

式中Xk為k時刻系統狀態向量，Z為量測向量，Φ為狀態一步轉移矩陣，B為噪聲分配矩陣，H為量測矩陣，Vk為量測噪聲向量，Wk為系統噪聲向量，Vk,Wk均為零均值的高斯白噪聲向量。卡爾曼濾波器遞推方程為

(17)

式中Kk為濾波增益矩陣，Pk,k-1為一步預測誤差方差陣 ，Pk為估計誤差方差陣，R,Q分別為過程噪聲和測量噪聲協方差[15]。

代入得到的陀螺儀z軸ARMA模型及其參數，可得濾波器輸入參數：Xk=[xkxk-1]T，Wk=[wkwk-1]T，Φ=[0.286 4 0;1 0]，B=[1-0.186 9;0 0],R=0.006 3。令初值X=[0 0]T,H[1 0]。圖5為濾波前后z軸靜態輸出對比。

圖5 卡爾曼濾波前后對比

4 結 論

經對濾波后的z軸數據進行Allan方差分析，對比了濾波前后z軸數據的隨機誤差項系數值，結果表明：濾波后的噪聲系數均有所減小，其中量化噪聲系數減小了2.8 %，角度隨機游走系數減小了19.8 %，零偏不穩定性系數減小了8.1 %。驗證了使用卡爾曼濾波器能夠有效抑制陀螺儀z軸原始輸出的隨機誤差，尤其是能夠減小角度隨機游走的影響。