陳雷　張曉明　檀杰　邱江濤

摘要：針對(duì)彈體磁場(chǎng)嚴(yán)重影響旋轉(zhuǎn)彈捷聯(lián)地磁傳感器的測(cè)量精度這一問題，提出一種基于卡爾曼濾波算法（Kalmanfilter，KF）的彈體磁場(chǎng)校正方法。利用固定磁場(chǎng)和感應(yīng)磁場(chǎng)模型，將彈體磁場(chǎng)誤差系數(shù)轉(zhuǎn)換到橢球參數(shù)方程上，從而得到卡爾曼濾波的觀測(cè)方程。為提高算法的魯棒性，采用事先標(biāo)定法建立初始條件。根據(jù)卡爾曼濾波原理，給出辨識(shí)參數(shù)在線更新的實(shí)現(xiàn)步驟，推導(dǎo)彈體磁場(chǎng)的校正過程。仿真試驗(yàn)中，通過事先標(biāo)定法選取初值提高待估參數(shù)2倍的收斂速度。轉(zhuǎn)臺(tái)試驗(yàn)中，彈體磁場(chǎng)校正后的磁測(cè)誤差接近磁傳感器的測(cè)量噪聲，滾轉(zhuǎn)角解算精度優(yōu)于1°。試驗(yàn)驗(yàn)證該算法可在線更新彈體磁場(chǎng)誤差系數(shù)，實(shí)現(xiàn)彈體磁場(chǎng)的高精度補(bǔ)償。

關(guān)鍵詞：彈體磁場(chǎng);磁傳感器;誤差校正;濾波

中圖分類號(hào)：TN919.5

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1674–5124（2019）02–0048–06

0 引言

制導(dǎo)彈藥精確打擊的前提是精確獲取高精度的姿態(tài)導(dǎo)航信息，地磁姿態(tài)系統(tǒng)對(duì)常規(guī)彈藥的制導(dǎo)有著重要意義[1]。地磁姿態(tài)系統(tǒng)利用地磁場(chǎng)信息進(jìn)行載體姿態(tài)解算，與慣性系統(tǒng)組合導(dǎo)航，可補(bǔ)償陀螺漂移，提供全天候、全地域、全天時(shí)的導(dǎo)航姿態(tài)信息[2]。

捷聯(lián)地磁姿態(tài)系統(tǒng)利用地磁信息進(jìn)行姿態(tài)解算時(shí)，作用于磁傳感器上的磁場(chǎng)不僅有地磁場(chǎng)信息，還有干擾磁場(chǎng)信息[3]。盡管在彈藥設(shè)計(jì)制造階段可以通過消磁技術(shù)來減少彈體磁場(chǎng)的干擾，但彈體磁場(chǎng)仍可達(dá)到數(shù)百到數(shù)千納特[4]。而且由于彈體磁場(chǎng)長時(shí)間存儲(chǔ)磁化以及發(fā)射過載下彈體磁場(chǎng)的改變[5]，事先標(biāo)定的干擾磁場(chǎng)參數(shù)無法真實(shí)反映彈體磁場(chǎng)的變化，進(jìn)而嚴(yán)重影響姿態(tài)解算的精度。因此需要一種實(shí)時(shí)準(zhǔn)確的載體磁場(chǎng)校正方法，以提高地磁場(chǎng)的實(shí)時(shí)測(cè)量精度。

傳統(tǒng)的彈體磁場(chǎng)校正方法將磁傳感器歷史測(cè)量值統(tǒng)一處理，事先求解出彈體磁場(chǎng)的校正參數(shù)。這類參數(shù)估計(jì)方法稱為事先標(biāo)定法，有橢圓擬合法[6]、橢球擬合法、超定方程求解最小二乘參數(shù)法[7]等。其求解固定磁干擾參數(shù)較為準(zhǔn)確，但由于需要內(nèi)存空間大，計(jì)算復(fù)雜度高，不適用于旋轉(zhuǎn)彈彈體磁場(chǎng)參數(shù)的在線校正或彈體磁場(chǎng)參數(shù)變化較大的環(huán)境。在線彈體磁場(chǎng)的校正方法包括遞推最小二乘法、八航向校正法[8]、二步校準(zhǔn)法[3]等。遞推最小二乘法估計(jì)參數(shù)方法較適合于機(jī)載解算，但算法固有特點(diǎn)導(dǎo)致在彈體磁場(chǎng)參數(shù)變化時(shí)增益參數(shù)會(huì)過度收斂[9-10]，會(huì)使得所估計(jì)參數(shù)出現(xiàn)誤差。而八航向校正法對(duì)彈體磁場(chǎng)參數(shù)校準(zhǔn)模型建模準(zhǔn)確性要求高，求解準(zhǔn)確磁場(chǎng)校正系數(shù)難度較大[11]。二步校準(zhǔn)法是通過最小二乘法求取固定磁場(chǎng)和比例因數(shù)誤差的迭代算法，該算法魯棒性好，但由于兩步估計(jì)均為非線性參數(shù)估計(jì)，估計(jì)過程復(fù)雜，計(jì)算量較大[12]。

為了實(shí)現(xiàn)旋轉(zhuǎn)彈藥彈體磁場(chǎng)的在線校正，首先分析了旋轉(zhuǎn)彈藥彈體干擾磁場(chǎng)誤差來源，建立了測(cè)量誤差模型，最終得到了基于三維磁場(chǎng)的線性觀測(cè)模型。利用事先標(biāo)定法為參數(shù)辨識(shí)提供初始條件，再采用卡爾曼濾波算法（KF）估計(jì)誤差系數(shù)，最終實(shí)現(xiàn)了磁測(cè)信號(hào)的高精度校正。最后進(jìn)行了一系列的仿真試驗(yàn)和轉(zhuǎn)臺(tái)試驗(yàn)，驗(yàn)證了彈體磁場(chǎng)校正算法的性能。

1 彈體磁場(chǎng)誤差模型

1.1 彈體磁場(chǎng)誤差機(jī)理

定義捷聯(lián)地磁傳感器的彈體坐標(biāo)系（b系），其坐標(biāo)原點(diǎn)位于傳感器安裝位點(diǎn)上，Xb、Yb、Zb三正交軸分別沿彈體的縱軸向前、橫軸向右、豎軸向下，如圖1所示。

事先標(biāo)定的校正參數(shù)，可以補(bǔ)償?shù)舻卮艂鞲衅鞯闹圃煺`差和安裝誤差，因此可以將其視為理想的捷聯(lián)地磁傳感器。捷聯(lián)地磁傳感器在彈體坐標(biāo)系下的量測(cè)方程可表示為：

式中：ζ為傳感器的測(cè)量噪聲誤差矢量;Hmb=HxHyHzT為地磁場(chǎng)測(cè)量矢量，Heb=HexHeyHezT為地磁場(chǎng)矢量，上標(biāo)“b”表示在彈體坐標(biāo)系下。固定磁場(chǎng)矢量B0=bxbybz是由彈體上硬磁材料所產(chǎn)生的。彈體在制造或存放期間，較長時(shí)間地停放在固定地點(diǎn)上，致使彈體上的硬磁材料受到不同程度的磁化而具有較強(qiáng)的磁性，而且這些硬磁材料本身具有一定的固定剩余磁場(chǎng)。這些固定磁場(chǎng)疊加在地磁場(chǎng)上，使得地磁傳感器的測(cè)量產(chǎn)生輸出偏置。感應(yīng)磁場(chǎng)系數(shù)Cs主要是由于彈體中軟磁材料在地磁場(chǎng)中磁化而產(chǎn)生的。軟磁材料具有較低的矯頑力和較窄的磁滯回線，可在較弱磁場(chǎng)中被磁化而產(chǎn)生感應(yīng)磁場(chǎng)。感應(yīng)磁場(chǎng)的大小及方向隨載體姿態(tài)和載體在地磁場(chǎng)中的位置變化而變化[13]。根據(jù)泊松方程，可以得到感應(yīng)磁場(chǎng)系數(shù)：

此外，彈體內(nèi)部的電子儀器還可能產(chǎn)生時(shí)變的干擾磁場(chǎng)，不過單根通電導(dǎo)線產(chǎn)生的電磁場(chǎng)相對(duì)較小，可以通過合理布線的方式減小其影響，因此此類干擾磁場(chǎng)基本可以忽略不計(jì)。

1.2 彈體捷聯(lián)磁場(chǎng)測(cè)量模型

若地磁姿態(tài)系統(tǒng)相對(duì)位置保持靜止，只改變姿態(tài)，地磁場(chǎng)矢量的模值Heb為常量，其矢量圖為標(biāo)準(zhǔn)球體;由于彈體磁場(chǎng)誤差的影響，磁場(chǎng)測(cè)量矢量Hmb的矢量圖為橢球體，可表述為二次曲面方程：

點(diǎn)的標(biāo)準(zhǔn)球面，磁場(chǎng)測(cè)量分量HxHyHzT的軌跡為偏離球心并帶有一定旋向的橢球，測(cè)量點(diǎn)（Hx，Hy，Hz）分布于理想橢球面兩側(cè)。彈體磁場(chǎng)環(huán)境的標(biāo)定就轉(zhuǎn)化為擬合橢球參數(shù)問題，關(guān)鍵在于橢球參數(shù)的精確辨識(shí)。將式（3）改寫成二次曲線方程的一般形式

2彈載磁場(chǎng)的校正

2.1 基于卡爾曼濾波的參數(shù)在線更新算法

考慮到噪聲誤差ζ是隨機(jī)變量，式（1）表示磁測(cè)系統(tǒng)誤差模型中沒有對(duì)噪聲誤差ζ標(biāo)識(shí)。通過硬件濾波和軟件濾波的方法可以在一定程度上消除噪聲誤差的影響，為了實(shí)現(xiàn)對(duì)橢球參數(shù)的高精度估計(jì)，傳統(tǒng)的方法采用直接最小二乘法進(jìn)行橢球擬合。直接橢球擬合是利用最小二乘法批量處理磁測(cè)數(shù)據(jù)，存在著求解結(jié)果受測(cè)量噪聲影響較大、不能在線實(shí)時(shí)更新、存儲(chǔ)數(shù)據(jù)量大、計(jì)算量大等缺點(diǎn)。利用卡爾曼濾波算法（KF）可實(shí)現(xiàn)系統(tǒng)的在線參數(shù)辨識(shí)，充分利用采樣數(shù)據(jù)減小隨機(jī)誤差的影響，并且由于計(jì)算量小，可有效節(jié)省存儲(chǔ)空間。

首先，濾波器的狀態(tài)向量取值定義為θ=abcdeglmnp，其元素為橢球方程各次項(xiàng)系數(shù)?？紤]到離散時(shí)間系統(tǒng)狀態(tài)向量恒定的一般形式從而確定了卡爾曼濾波器的狀態(tài)方程，由橢球方程建立了卡爾曼濾波器的觀測(cè)方程，如下式所示：

其中：Hk=[Hk2xHk2yHk2zHkxHkyHkxHkzHkyHkzHkxHkyHkz1]為觀測(cè)方程的量測(cè)矩陣，Zk=0為量測(cè)向量。隨機(jī)測(cè)量噪聲Vk滿足：

根據(jù)卡爾曼濾波原理，在線彈體磁場(chǎng)誤差系參數(shù)辨識(shí)的實(shí)現(xiàn)步驟為：

1）利用事先標(biāo)定法確定狀態(tài)向量和協(xié)方差矩陣的初始狀態(tài)θ0、P0;

2）狀態(tài)預(yù)測(cè)：

3）狀態(tài)估計(jì)：

4）返回到第2個(gè)步驟，直到狀態(tài)向量θk滿足：

時(shí)，即當(dāng)各待估參數(shù)變化不大時(shí)，可認(rèn)為參數(shù)辨識(shí)完畢。

2.2 彈體磁環(huán)境的校正

上述過程實(shí)現(xiàn)了橢球模型狀態(tài)向量θ的在線辨識(shí)，由式（4）可以得到固定磁場(chǎng)矢量：

由于M是實(shí)對(duì)稱矩陣，對(duì)該矩陣進(jìn)行奇異值分解（singular value decomposition，SVD）可得到M=UTSMU。其中U為對(duì)稱矩陣，SM是由M的特征值組成的對(duì)角陣。為簡(jiǎn)化計(jì)算，假設(shè)感應(yīng)磁場(chǎng)系數(shù)矩陣是實(shí)對(duì)稱矩陣，則該矩陣的逆矩陣也是一個(gè)實(shí)對(duì)稱矩陣。因此，矩陣A也可以被SVD分解[14]，再結(jié)合式（3）可以得到：

標(biāo)定的目的在于盡可能消除噪聲誤差ζ的影響，通過對(duì)誤差模型各參數(shù)矩陣進(jìn)行辨識(shí)，進(jìn)一步根據(jù)下式的誤差補(bǔ)償原理，實(shí)現(xiàn)誤差的校正。因此，利用地磁傳感器的輸出值計(jì)算得到真實(shí)的磁分量：

3 仿真分析及轉(zhuǎn)臺(tái)試驗(yàn)驗(yàn)證

利用Matlab環(huán)境仿真飛行彈體下的磁測(cè)數(shù)據(jù)，通過一組仿真試驗(yàn)檢驗(yàn)了初值選取對(duì)彈體磁場(chǎng)校正算法魯棒性和收斂速度的作用效果，然后通過轉(zhuǎn)臺(tái)試驗(yàn)驗(yàn)證算法的實(shí)時(shí)性和有效性。

3.1 仿真分析

根據(jù)旋轉(zhuǎn)彈藥在外彈道無控飛行狀態(tài)下的運(yùn)動(dòng)規(guī)律，彈體偏航角近似為零，理想無干擾磁場(chǎng)時(shí)捷聯(lián)三軸磁傳感器測(cè)量數(shù)據(jù)的特征為：軸向磁場(chǎng)分量隨彈體俯仰變化規(guī)律而變化;兩個(gè)徑向磁場(chǎng)呈類正弦變化：頻率為彈體滾轉(zhuǎn)角速度，幅值隨彈體俯仰變化規(guī)律而變化，兩者彼此正交。地磁場(chǎng)在載體坐標(biāo)下的模值Heb設(shè)定為0.5267Gauss。根據(jù)無控彈的飛行規(guī)律，假設(shè)彈體在該磁場(chǎng)下以5r/s的轉(zhuǎn)速飛行，彈體俯仰角由45°線性變化到–45°，整個(gè)過程持續(xù)5s，仿真出干擾環(huán)境下的磁場(chǎng)數(shù)據(jù)。

表1給出了仿真試驗(yàn)的試驗(yàn)條件，試驗(yàn)通過常規(guī)法和事先標(biāo)定法分別選取初始條件來評(píng)估初值選取對(duì)算法性能的影響。根據(jù)以往磁測(cè)數(shù)據(jù)，設(shè)定干擾磁場(chǎng)下的誤差矩陣系數(shù)，地磁場(chǎng)測(cè)量噪聲服從精度為0.0005Gauss、均值為零的高斯分布。常規(guī)法選取的參數(shù)初值均為1，事先標(biāo)定法選取的初值如表1的試驗(yàn)II所示。試驗(yàn)I的仿真結(jié)果如圖3（a）、圖3（b）所示，圖3（a）顯示了參數(shù)逐步收斂的全過程，可以看出在初始段參數(shù)估計(jì)出現(xiàn)了多次較大的震蕩，甚至發(fā)散。這是由于初始階段信息量的不足，估計(jì)參數(shù)無法滿足橢球方程的約束條件，從而導(dǎo)致了參數(shù)估計(jì)過程的大幅震蕩。估計(jì)參數(shù)震蕩后的收斂過程如圖3（b）所示，可以看出一旦估計(jì)參數(shù)滿足約束條件后，該算法便能使其始終維持在約束條件內(nèi)，從而在3s左右逐步完成了參數(shù)收斂。圖4顯示了試驗(yàn)II的仿真結(jié)果，可以看出通過事先標(biāo)定法得到的初始條件在參數(shù)估計(jì)的一開始就滿足約束條件，從而避免了參數(shù)估計(jì)初始階段的大幅震蕩，而之后除cxx外的各參數(shù)更是在1s的時(shí)間內(nèi)完成了收斂過程。

3.2 轉(zhuǎn)臺(tái)試驗(yàn)

為了進(jìn)一步驗(yàn)證在線校正算法的補(bǔ)償效果，使用如圖5所示的飛行轉(zhuǎn)臺(tái)模擬彈體飛行過程的姿態(tài)變化，進(jìn)行校正試驗(yàn)。試驗(yàn)轉(zhuǎn)臺(tái)設(shè)置為：轉(zhuǎn)臺(tái)內(nèi)框以5r/s的轉(zhuǎn)速轉(zhuǎn)動(dòng)，中框由45°線性變化到–45°，外框靜止，整個(gè)過程持續(xù)5s。為了模擬彈載下的磁干擾環(huán)境，將系統(tǒng)裝配在含有鐵磁材料的彈體尾翼中，再通過圓筒形工裝將系統(tǒng)固連到飛行轉(zhuǎn)臺(tái)的內(nèi)框上。本次試驗(yàn)中，在線校正算法得到的彈體磁環(huán)境誤差系數(shù)在線更新過程如圖6所示，各系數(shù)最終收斂值如表2所示，磁測(cè)誤差的補(bǔ)償效果如圖7所示?？梢钥闯鲂Ｕ皒、y、z各軸的測(cè)量誤差均值接近固定磁場(chǎng)誤差系數(shù)，校正后的測(cè)量誤差更接近測(cè)量噪聲。采用12位置點(diǎn)記錄滾轉(zhuǎn)角數(shù)據(jù)，利用表1標(biāo)定的誤差系數(shù)對(duì)數(shù)據(jù)補(bǔ)償，分別利用補(bǔ)償前后的磁場(chǎng)數(shù)據(jù)進(jìn)行滾轉(zhuǎn)角解算，對(duì)比滾轉(zhuǎn)角解算誤差如圖8所示。校正前滾轉(zhuǎn)角解算誤差達(dá)到了8°，殘差均值為0.547°，標(biāo)準(zhǔn)差為5.59°;校正后滾轉(zhuǎn)角解算精度優(yōu)于1°，殘差均值為–0.251°，標(biāo)準(zhǔn)差為0.367°。試驗(yàn)表明通過在線校正算法的補(bǔ)償，可以滿足高精度的磁場(chǎng)測(cè)量需求，驗(yàn)證了算法的有效性。

4 結(jié)束語

本文提出了一種用于校正旋轉(zhuǎn)彈捷聯(lián)地磁傳感器的在線算法。不同于傳統(tǒng)的直接橢球擬合標(biāo)定法和“二步校準(zhǔn)法”，該算法使用線性模型實(shí)時(shí)校正磁測(cè)誤差，能實(shí)現(xiàn)彈體層次化誤差系數(shù)的在線更新，且計(jì)算量小。該算法根據(jù)磁傳感器測(cè)量的誤差模型，將固定磁場(chǎng)和感應(yīng)磁場(chǎng)模型系數(shù)轉(zhuǎn)換到橢球方程的參數(shù)上，再通過卡爾曼濾波算法估計(jì)方程參數(shù)，從而實(shí)現(xiàn)誤差系數(shù)的在線更新。該算法提出了事先標(biāo)定法建立初始條件，并經(jīng)過了仿真試驗(yàn)和轉(zhuǎn)臺(tái)試驗(yàn)。試驗(yàn)結(jié)果表明，該在線算法具有非常好的估計(jì)精度，可有效補(bǔ)償彈體磁場(chǎng)干擾誤差，實(shí)現(xiàn)磁傳感器的高精度測(cè)量。該算法具有適用性，不僅可以實(shí)現(xiàn)彈體磁場(chǎng)誤差的校正，還實(shí)現(xiàn)磁傳感器誤差的校正。

參考文獻(xiàn)

[1]陳國彬.彈載三軸磁測(cè)系統(tǒng)設(shè)計(jì)與誤差標(biāo)定補(bǔ)償技術(shù)研究[D].太原：中北大學(xué)，2012.

[2]龍禮，張合，劉建敬.地磁陀螺復(fù)合測(cè)姿系統(tǒng)誤差補(bǔ)償方法[J].火力與指揮控制，2014（7）：102-105.

[3]GEBREE.Anon-linear， two-stepestimationalgorithmforcalibratingsolid-statestrapdownmagnetometers[J].InternationalStPetersburgConferenceonNavigationSystems， 2001： 28-30.

[4]王婕，郭子祺，劉建英.固定翼無人機(jī)航磁探測(cè)系統(tǒng)的磁補(bǔ)償模型分析[J].航空學(xué)報(bào)，2016，37（11）：3435-3443.

[5]GOPALBV， SARMAVN， RAMBABUHV.Realtimecompensationforaircraftinducednoiseduringhighresolutionairbornemagneticsurveys[J].JIndGeophysUnion， 2004， 46：

48-52.

[6]張曉明，趙剡.基于橢圓約束的新型載體磁場(chǎng)標(biāo)定及補(bǔ)償技術(shù)[J].儀器儀表學(xué)報(bào)，2009，30（11）：2438-2443.

[7]趙鑫爐，張曉明，白渚銓，等.基于磁阻傳感器的航姿測(cè)量系統(tǒng)羅差補(bǔ)償技術(shù)研究[J].傳感技術(shù)學(xué)報(bào)，2013（11）：1504-1507.

[8]Processforcompensatingthemagneticdisturbancesinthedeterminationofamagneticheading， anddevicesforcarryingoutthisprocess： US4414753A[P].[2018-06-08]

[9]郭鵬飛，華春紅，任章，等.基于遞推最小二乘的航姿系統(tǒng)羅差校正[J].中國慣性技術(shù)學(xué)報(bào)，2008，16（1）：24-27.

[10]孫宏偉，房建成，李艷.橢圓擬合方法在磁羅盤羅差校準(zhǔn)中的應(yīng)用[J].光學(xué)精密工程，2009，17（12）：3034-3039.

[11]韓蘭懿，米松林，潘建軍，等.一種無人機(jī)羅差修正方法[J].航空計(jì)算技術(shù)，2013，43（3）：69-72.

[12]GEBREEGZ.Calibrationofstrapdownmagnetometersinmagneticfielddomain[J].JournalofAerospaceEngineering，2014， 19（2）： 87-102.

[13]張曉明，趙剡.一種快速載體磁場(chǎng)補(bǔ)償方法研究[J].中北大學(xué)學(xué)報(bào)（自然科學(xué)版），2009，30（3）：286-291.

[14]GRANDVALLETB， ZEMOUCHEA， BOUTAYEBM， etal.Real-timeattitude-independentthree-axismagnetometercalibrationforspinningprojectiles： aslidingwindowapproach[J].IEEETransactionsonControlSystemsTechnology， 2013， 22（1）： 255-264.