韓穎
(中國空空導(dǎo)彈研究院,河南洛陽,471000)
將MEMS 器件應(yīng)用在飛航導(dǎo)彈系統(tǒng)中是未來微型化、小型化的發(fā)展趨勢。因技術(shù)受限,MEMS 陀螺儀目前多為商業(yè)級[1]。對于精度要求高的穩(wěn)像系統(tǒng),MEMS 陀螺儀噪聲在彈載飛行過程中引起運動載體輸出軸的抖動以及角度的漂移,會導(dǎo)致視頻圖像質(zhì)量的下降。為了使穩(wěn)像系統(tǒng)有高的跟蹤精度,在本文中設(shè)計使用卡爾曼濾波算法,對MEMS 陀螺儀的噪聲進行抑制或去除。
MEMS 陀螺儀隨機噪聲信號屬于非平穩(wěn)的時間序列,對于載體機動頻繁,不確定因素較多的系統(tǒng),為避免將真實數(shù)據(jù)當(dāng)作噪聲過濾掉,在卡爾曼濾波前,應(yīng)對陀螺儀噪聲信號進行預(yù)處理,保證其為零均值、平穩(wěn)、正態(tài)時間序列后,建立ARMA 模型[2,3]。
趨勢項是是MEMS 陀螺儀線性漂移的主要影響因素,一般常采用最小二乘法將其剔除。
在本研究中使用Matlab 運算工具中polyfit()函數(shù)和polyval()函數(shù)代替一般代數(shù)方法求得趨勢項系數(shù),在對陀螺儀數(shù)據(jù)進行3 階、4 階、5 階、6 階、7 階多項式趨勢項去除后發(fā)現(xiàn),去除5 階、6 階和7 階趨勢項后的數(shù)據(jù)方差很接近,見表1。

表1 去趨勢項后數(shù)據(jù)方差
考慮到在實際應(yīng)用中,模型階數(shù)大,計算量會隨之增大。因此,本文選擇采用剔除5 階趨勢項方法使陀螺的隨機漂移滿足平穩(wěn)性要求。
采集信號中偶然出現(xiàn)的異常數(shù)據(jù),被稱為野值(outl ier)。考慮到計算量的問題,設(shè)計中使用3σ 準(zhǔn)則循環(huán)3 次去除野值項,如圖1。

圖1 野值點與3 次去野值后的波形
利用Matlab 工具箱中dtrend()函數(shù)應(yīng)對去除野值之后的陀螺儀隨機噪聲信號進行平穩(wěn)化,零化處理,使噪聲信號為零均值序列,處理前后均值見表2。

表2 零化處理前后的均值
平穩(wěn)性檢驗是對時間序列進行ARMA 建模的前提條件。可以用來檢驗陀螺儀漂移噪聲序列是否不隨時間推移而變化。本設(shè)計中采用平穩(wěn)性的Daniel 檢驗方法對陀螺隨機噪聲序列進行檢驗。結(jié)果如表3:

表3 Daniel 平穩(wěn)性檢驗方法檢驗結(jié)果
總體分布的正態(tài)性檢驗一般采取Jarque-Bera 檢驗,觀察偏態(tài)系數(shù)g3與峰態(tài)系數(shù)g4是否滿足正太隨機變量的特性。對兩組去除趨勢項和野值點的陀螺隨機噪聲信號做正態(tài)性檢驗,其結(jié)果如表4。

表4 正態(tài)性檢驗結(jié)果
由上表可以看出,g3、g4都約等于0,因此可以認(rèn)此時MEMS 陀螺隨機信號序列均滿足正態(tài)性要求。
將MEMS 陀螺漂移信號轉(zhuǎn)化為零均值、正態(tài)、平穩(wěn)時間序列,符合建模要求后,選用AIC 準(zhǔn)則和BIC 準(zhǔn)則對陀螺隨機信號序列進行ARMA(p,q)模型階數(shù)的確定。
本文可以采用的時間序列模型一共有五種[4],即AR(1)、AR(2)、AR(3)、ARMA(1,1)、ARMA(2,1)。對這五種模型同時進行AIC 準(zhǔn)則和BIC 準(zhǔn)則的適用性檢驗,檢驗結(jié)果如表5。

表5 AIC 準(zhǔn)則和BIC 準(zhǔn)則的檢驗結(jié)果
通過比較后,選取AIC 和BIC 值相對較小的ARMA(2,1)模型作為陀螺儀的隨機漂移誤差模型。其表達式為:
在Matlab 中使用最小二乘法即可解得到相應(yīng)的ARMA(2,1)的模型參數(shù),如表6 所示,將此數(shù)據(jù)作為卡爾曼濾波算法的初值。

表6 陀螺儀參數(shù)模型
卡爾曼濾波是根據(jù)前一個估計值和最近一個觀測值,按一定的遞推方式求得新的狀態(tài)估計值,然后不斷循環(huán)的過程。基礎(chǔ)公式(1)之外,其狀態(tài)方程和測量方程分別如下[5,6]:


其中,R k是測量噪聲協(xié)方差矩陣,為正定矩陣;Q k是過程噪聲協(xié)方差矩陣,為非負(fù)正定矩陣,δ kj為克羅尼克(Kroneker)δ函數(shù),其特性是:

如果被估計狀態(tài)和觀測量是滿足式(2)、(3),則k 時刻的觀測的估計可按下述方程求解。
根據(jù)以上五條基本公式,給定初值X0和P0,根據(jù)k時刻的觀測值Zk,就可以遞推計算得k 時刻的狀態(tài)估計(k=1,2…,N)。
卡爾曼濾波算法框圖如圖2 所示。

圖2 一步預(yù)測卡爾曼濾波算法
在試驗中ARMA(2,1)模型參數(shù)如表6 所示。R k是觀測噪聲序列方差陣,取濾波數(shù)據(jù)的方差varX;Qk是系統(tǒng)噪聲序列方差陣,為正定矩陣,取其中σα為ARMA(2,1)模型系統(tǒng)激勵白噪聲方差。
在本設(shè)計中假定濾波估計從開始就是無偏的,且估計的誤差陣最小,初始值另有,Φk,k?1和Γk,k?1為tk,k?1時刻至tk時刻的狀態(tài)轉(zhuǎn)移矩陣,取
理想情況下,當(dāng)MEMS 陀螺輸入為零時,輸出也應(yīng)該為零,但是MEMS 陀螺存其它形式的隨機干擾和漂移[7]。在本研究中,選用IEEE 公認(rèn)的陀螺儀參數(shù)分析評價的標(biāo)準(zhǔn)方法——Allan 方差,對MEMS 陀螺儀濾波結(jié)果進行評價。
在25 ℃室溫靜態(tài)條件下,以500Hz 的采樣頻率對MEMS 陀螺樣品進行靜態(tài)數(shù)據(jù)采集,敏感軸指天,陀螺標(biāo)度因數(shù)為7.0mV/°/sec,并假設(shè)其不變。采樣時間1min。利用Matlab2007a 對采集到的數(shù)據(jù)進行預(yù)處理后,方可對陀螺儀輸出的隨機噪聲信號序列進行濾。濾波前后對比,如圖3 所示。

圖3 卡爾曼濾波前后數(shù)據(jù)以及濾除噪聲
卡爾曼濾波前后Allan 方差對比如圖4 所示。

圖4 經(jīng)典卡爾曼濾波前后Allan 方差對比圖
卡爾曼濾波前后Allan方差分析各噪聲分量的結(jié)果如7。
卡爾曼濾波前后陀螺隨機誤差標(biāo)準(zhǔn)差如表8 所示。

表8 卡爾曼濾波前后陀螺隨機誤差標(biāo)準(zhǔn)差
從圖3 及表7 中可以看出,在經(jīng)過卡爾曼濾波處理后,量化噪聲約為濾波前的50%,角速度隨機游走、零偏不穩(wěn)定性、速率隨機游走、速率斜坡約為濾波前的60%。在表8 中可以看出,通過卡爾曼濾波之后,陀螺儀隨機誤差標(biāo)準(zhǔn)差降為濾波前的50%。說明卡爾曼濾波對該MEMS 陀螺儀有明顯的濾波效果。

表7 卡爾曼濾波前后Allan 方差結(jié)果
將卡爾曼濾波算法與MEMS 陀螺儀相結(jié)合,運用在彈載系統(tǒng)中,既節(jié)省了存儲空間,又易于工程實現(xiàn),使濾波過程更加直接有效,非常適合應(yīng)用在穩(wěn)像系統(tǒng)中對MEMS 陀螺儀的噪聲處理中。