孫 歡，楊賓峰，李 馳，王 潤(rùn)

(空軍工程大學(xué)信息與導(dǎo)航學(xué)院，陜西 西安 710077)

0 引言

地磁導(dǎo)航隱蔽性好，自主性高，抗干擾能力強(qiáng)，能提供全地域，全天候，全天時(shí)的導(dǎo)航信息，相比于傳統(tǒng)導(dǎo)航方式有著明顯優(yōu)勢(shì)和更廣闊的發(fā)展前景[1-6]。

精確的導(dǎo)航離不開(kāi)對(duì)地磁場(chǎng)的精確測(cè)量，地磁傳感器的觀測(cè)結(jié)果往往受到多種因素的干擾，儀表誤差、安裝誤差、載體內(nèi)部軟硬磁誤差及渦流磁場(chǎng)誤差都會(huì)對(duì)觀測(cè)結(jié)果造成不良影響，從而導(dǎo)致較大的定位偏差。因此對(duì)于地磁場(chǎng)觀測(cè)值進(jìn)行誤差補(bǔ)償，對(duì)提高地磁匹配導(dǎo)航的精度有著重要的意義[7-14]。

針對(duì)地磁傳感器的誤差補(bǔ)償問(wèn)題，國(guó)內(nèi)外學(xué)者進(jìn)行了諸多研究，V.Renaudin等人在磁場(chǎng)域?qū)Ω飨虍愋源抛鑲鞲衅?Anisotropic Magnetoresistive Sensor，AMR)的校正問(wèn)題進(jìn)行了研究，其詳細(xì)分析了傳感器自身儀表誤差、軟磁以及硬磁誤差的來(lái)源和產(chǎn)生機(jī)理，建模時(shí)將這些誤差都予以充分考慮[15]。文獻(xiàn)[16]采用十二姿態(tài)補(bǔ)償法對(duì)地磁傳感器進(jìn)行建模，補(bǔ)償算法采用了橢圓擬合校準(zhǔn)法和神經(jīng)網(wǎng)絡(luò)算法；文獻(xiàn)[17]對(duì)地磁傳感器儀表誤差、軟磁誤差進(jìn)行建模，并通過(guò)磁通量補(bǔ)償來(lái)消除磁矩干擾；文獻(xiàn)[18—19]對(duì)三軸磁力儀測(cè)量誤差和軟硬磁誤差進(jìn)行一體化建模，前者采用橢球擬合校準(zhǔn)法進(jìn)行補(bǔ)償，后者考慮到觀測(cè)向量和數(shù)據(jù)矩陣均存在誤差，因此用總體最小二乘進(jìn)行了補(bǔ)償；文獻(xiàn)[20—22]對(duì)三軸磁力儀的儀表誤差、安裝誤差和軟硬磁誤差進(jìn)行建模，采用ANSYS、COMSOL Multiphysics等軟件仿真建模，并采用最小二乘法進(jìn)行補(bǔ)償。

通過(guò)分析，以往研究仍存在以下問(wèn)題：1)模型中普遍沒(méi)有考慮渦流磁場(chǎng)的影響，對(duì)于高速、頻繁進(jìn)行姿態(tài)變換的載體，其渦流磁場(chǎng)的影響將會(huì)使補(bǔ)償精度下降；2)以地磁傳感器為核心的觀測(cè)系統(tǒng)性能不穩(wěn)定會(huì)形成系統(tǒng)噪聲，該噪聲與觀測(cè)噪聲相互耦合會(huì)降低觀測(cè)系統(tǒng)的性能，進(jìn)而影響補(bǔ)償精度。本文針對(duì)此問(wèn)題，提出了基于一體化模型的雙噪聲聯(lián)合估計(jì)補(bǔ)償方法。

1 誤差建模

在磁力儀標(biāo)定過(guò)程中，地磁傳感器的三軸x′,y′,z′的非嚴(yán)格正交，將會(huì)造成非正交誤差CNO，可以表示為：

(1)

在磁力儀安裝過(guò)程中，原則上傳感器坐標(biāo)系應(yīng)該與載體坐標(biāo)系嚴(yán)格一致，但是在安裝過(guò)程中，將不可避免地出現(xiàn)工差，從而造成安裝誤差，設(shè)磁力儀三軸與載體坐標(biāo)系三軸的誤差角為εx,εy,εz，則安裝誤差CM為：

(2)

由于制作工藝精度的問(wèn)題，三軸磁力儀三軸存在細(xì)微差異，由此產(chǎn)生刻度因子誤差，設(shè)三軸的刻度因子誤差分別為sx,sy,sz，則刻度因子誤差矩陣CSF為：

CSF=diag(sx,sy,sz)

(3)

制作工藝精度和安裝工差同樣也會(huì)造成零刻度漂移誤差w和剩磁誤差bs，矩陣表達(dá)式分別為：

w=[δ1,δ2,δ3]T，bs=[bx1,by1,bz1]T

(4)

硬磁誤差是載體內(nèi)部硬磁材料產(chǎn)生的磁場(chǎng)所造成的誤差，記硬磁誤差為：

Bh=[bh1,bh2,bh3]T

(5)

軟磁誤差是軟磁性材料被外部磁場(chǎng)激勵(lì)所產(chǎn)生的誤差，外部磁場(chǎng)包括地磁場(chǎng)和載體內(nèi)部所產(chǎn)生的磁場(chǎng)，可以表示為：

BSI=CSI(H+Bh)

(6)

載體做姿態(tài)變換的過(guò)程中，將會(huì)產(chǎn)生渦流磁場(chǎng)，將其考慮在內(nèi)，得模型：

B=CH+b+P(ΔB/Δt)+ε

(7)

式(7)中，C=CMCNOCSFCSI，b=CMCNOCSFCSIBh+w+bs，ε為零均值高斯地磁觀測(cè)噪聲，其方差為Σ。

由于觀測(cè)系統(tǒng)并非絕對(duì)穩(wěn)定，在地磁場(chǎng)觀測(cè)過(guò)程中，將不可避免地產(chǎn)生系統(tǒng)誤差，體現(xiàn)在各個(gè)系統(tǒng)參量的微小變化，因此地磁場(chǎng)觀測(cè)模型可以進(jìn)一步表示為：

B=(C+ΔC)H+(b+Δb)+(P+ΔP)(ΔB/Δt)+ε

(8)

式(8)中ΔC，Δb，ΔP是觀測(cè)系統(tǒng)不穩(wěn)定造成的系統(tǒng)參數(shù)波動(dòng)。

2 雙噪聲聯(lián)合估計(jì)補(bǔ)償算法

實(shí)際觀測(cè)過(guò)程中，除了觀測(cè)噪聲之外，觀測(cè)系統(tǒng)不可避免地受到各種擾動(dòng)因素影響，因此不能保證系統(tǒng)參量的絕對(duì)精確可信，根據(jù)上文式(8)所提模型：

z=B=(C+ΔC)H+(b+Δb)+(P+ΔP)(ΔB/Δt)+ε=f(x,v)+ε

(9)

式(9)中，z為實(shí)際觀測(cè)值，x為真值，v為含擾動(dòng)的系統(tǒng)參量，其值為v=w+m，w為原始系統(tǒng)參量，m為系統(tǒng)擾動(dòng)誤差，z,x均為3×1維向量，v,w,m為21×1維向量，f(?,?)為非線性觀測(cè)方程，ε為觀測(cè)誤差，各個(gè)量如下所示：

z=[Bx,By,Bz]T，x=[Hx,Hy,Hz]T，

w=[c11,c12,c13,c21,c22,c23,c31,c32,c33,b1,b2,b3,p11,p12,p13,p21,p22,p23,p31,p32,p33]T，

m=[Δc11,Δc12,Δc13,Δc21,Δc22,Δc23,Δc31,

Δc32,Δc33,Δb1,Δb2,Δb3,Δp11,Δp12,Δp13,

Δp21,Δp22,Δp23,Δp31,Δp32,Δp33]T。

以往研究表明，觀測(cè)誤差、系統(tǒng)誤差服從零均值高斯分布，且統(tǒng)計(jì)獨(dú)立，令其方差矩陣分別為Q1=E[εεT]，Q2=E[mmT]，其中，Q1為3×3維方差矩陣，Q2為21×21維方差矩陣。

由最小二乘算法原理， 當(dāng)(z-f(x,w))T·(z-f(x,w))取得最小值時(shí)，x取得最優(yōu)估計(jì)值，但當(dāng)觀測(cè)方程中含有高斯觀測(cè)噪聲ε時(shí)，進(jìn)行加權(quán)處理來(lái)抑制噪聲的影響，此時(shí)，x的最優(yōu)估計(jì)值為：

(9)

鑒于f(?,?)的非線性特征，在解算過(guò)程中將進(jìn)行迭代求解。在認(rèn)為系統(tǒng)參量w精確可信(觀測(cè)量誤差僅由觀測(cè)噪聲產(chǎn)生)的情況下，經(jīng)常采用基于一階Taylor展開(kāi)的數(shù)值迭代算法，也稱為Gauss-Newton迭代法，其原理是對(duì)f(?,?)作線性近似，每一次的迭代值作為下一次迭代的更新公式，重復(fù)此過(guò)程直至收斂。

(10)

則k-1次迭代結(jié)果可通過(guò)如下線性加權(quán)最小二乘器獲得，如式(11)所示。

由于該問(wèn)題是關(guān)于x的二次優(yōu)化問(wèn)題，因此存在閉合最優(yōu)解，其值如式(12)所示。

這同樣是關(guān)于x和w的二次優(yōu)化問(wèn)題，因此存在閉合解，如式(16)所示。

(11)

(12)

(13)

(14)

(15)

(16)

(17)

3 模擬仿真與實(shí)驗(yàn)驗(yàn)證

3.1 模擬仿真

利用仿真數(shù)據(jù)對(duì)上述算法進(jìn)行驗(yàn)證，根據(jù)上文所述地磁場(chǎng)觀測(cè)值模型，設(shè)地磁場(chǎng)真實(shí)值為50 000 nT，相關(guān)參量取值如下：

產(chǎn)生仿真數(shù)據(jù)時(shí)，觀測(cè)噪聲、系統(tǒng)噪聲服從均值為0，方差分別為50和30的高斯分布。

首先，不加噪聲產(chǎn)生地磁場(chǎng)觀測(cè)數(shù)據(jù)，分別采用含渦流磁場(chǎng)模型和不含渦流磁場(chǎng)模型進(jìn)行補(bǔ)償，補(bǔ)償結(jié)果如圖1所示。

圖1 不同模型補(bǔ)償效果Fig.1 Different model compensation effects

由圖可知，含渦流磁場(chǎng)模型在地磁場(chǎng)補(bǔ)償中，不僅能夠顯著降低誤差峰值，而且可以有效濾除載體姿態(tài)轉(zhuǎn)動(dòng)所引起的波動(dòng)干擾，具有更高的補(bǔ)償精度和穩(wěn)定性。由此確定了含渦流磁場(chǎng)模型的優(yōu)越性。

其次驗(yàn)證Gauss-Newton迭代法在無(wú)系統(tǒng)噪聲時(shí)的可用性，此時(shí)仿真數(shù)據(jù)不含系統(tǒng)噪聲(即系統(tǒng)參量真實(shí)可信)。在補(bǔ)償過(guò)程中，相應(yīng)的觀測(cè)噪聲標(biāo)準(zhǔn)差設(shè)為σ=100δ nT，其中δ為觀測(cè)量擾動(dòng)參數(shù)(數(shù)值從1到29變化)，地磁場(chǎng)解算值的均方根誤差隨觀測(cè)量擾動(dòng)參數(shù)變化如圖2所示，整體補(bǔ)償效果如圖3所示。

圖2 地磁場(chǎng)解算值的均方根誤差隨δ變化Fig.2 The root mean square error of the geomagnetic field solution value changes with δ

圖3 地磁場(chǎng)補(bǔ)償效果(不含系統(tǒng)噪聲)Fig.3 Geomagnetic field compensation effect (without system noise)

由圖2和圖3可知，解算值的均方根誤差可以達(dá)到相應(yīng)的Gramer-Rao界(CRB)，在補(bǔ)償中將峰值誤差從400 nT降至30 nT，證明Gauss-Newton迭代法在系統(tǒng)噪聲不存在的地磁場(chǎng)補(bǔ)償中完全可用。

當(dāng)仿真數(shù)據(jù)中加入系統(tǒng)噪聲以后，分別以基于Gauss-Newton迭代法(A算法)和Taylor-B迭代法(B算法)的加權(quán)最小二乘算法進(jìn)行補(bǔ)償，設(shè)相應(yīng)的觀測(cè)噪聲標(biāo)準(zhǔn)差為σ=100δnT，系統(tǒng)噪聲標(biāo)準(zhǔn)差為υ=30εnT,其中，δ和ε為觀測(cè)量擾動(dòng)參數(shù)(數(shù)值從1到29變化)，則地磁場(chǎng)解算值的均方根誤差隨δ和ε變化如圖4和5所示，整體補(bǔ)償效果如圖6所示。

由圖4和圖5可知，當(dāng)系統(tǒng)噪聲存在時(shí)，基于Gauss-Newton迭代法的加權(quán)最小二乘解算值的均方根誤差遠(yuǎn)遠(yuǎn)高于相應(yīng)的CRB，已經(jīng)不適用于誤差補(bǔ)償；而基于Taylor-B迭代法的加權(quán)最小二乘解算值均方根誤差達(dá)到了相應(yīng)的CRB，理論上可以進(jìn)行補(bǔ)償；由圖6可知，在對(duì)仿真數(shù)據(jù)的補(bǔ)償中，基于Gauss-Newton迭代法的加權(quán)最小二乘將峰值誤差從700 nT降至400 nT，仍達(dá)不到期望精度，而基于Taylor-B迭代法的加權(quán)最小二乘由于去除了系統(tǒng)噪聲帶來(lái)的誤差，將峰值誤差降至20 nT，進(jìn)一步提高了補(bǔ)償精度，使得地磁場(chǎng)觀測(cè)值補(bǔ)償有了更好的適用性。

圖4 地磁場(chǎng)解算值的均方根誤差隨δ變化Fig.4 The root mean square error of the geomagnetic field solution value changes withδ

圖5 地磁場(chǎng)解算值的均方根誤差隨ε變化Fig.5 The root mean square error of the geomagnetic field solution value changes with ε

圖6 地磁場(chǎng)補(bǔ)償效果(含系統(tǒng)噪聲)Fig.6 Geomagnetic field compensation effect (with system noise)

3.2 實(shí)驗(yàn)設(shè)計(jì)

在設(shè)計(jì)實(shí)驗(yàn)中，采用一款Bartington公司的Mag-03地磁傳感器進(jìn)行補(bǔ)償，為模擬飛行載體姿態(tài)變化，選用3FHT30C無(wú)磁轉(zhuǎn)臺(tái)提供姿態(tài)信息，數(shù)據(jù)采集軟件用National Instruments，如圖7所示，后期用Matlab進(jìn)行數(shù)據(jù)處理。地磁場(chǎng)真實(shí)值由標(biāo)準(zhǔn)磁圖提供。

選取遠(yuǎn)離強(qiáng)磁干擾的實(shí)驗(yàn)場(chǎng)所，將地磁傳感器置于無(wú)磁轉(zhuǎn)臺(tái)，快速旋轉(zhuǎn)進(jìn)行全方位姿態(tài)信息測(cè)量。為保證測(cè)量噪聲擬合高斯分布，測(cè)量10組數(shù)據(jù)并選取相同姿態(tài)的數(shù)據(jù)作為實(shí)驗(yàn)數(shù)據(jù)，系統(tǒng)參量初值選用經(jīng)驗(yàn)值(先驗(yàn)值的不準(zhǔn)確度即為系統(tǒng)噪聲)。將軟件采集數(shù)據(jù)用本文所提兩種算法進(jìn)行處理，所得結(jié)果如圖8所示，未經(jīng)補(bǔ)償?shù)牡卮艌?chǎng)測(cè)量值峰值誤差可達(dá)3 000 nT，基于Gauss-Newton迭代法的加權(quán)最小二乘算法有著一定的補(bǔ)償效果，將誤差峰值降至1 600 nT，但仍然達(dá)不到地磁導(dǎo)航精度要求，基于Taylor-B迭代法的加權(quán)最小二乘算法由于濾除了系統(tǒng)噪聲帶來(lái)的干擾，將峰值誤差降至90 nT，因此實(shí)現(xiàn)了更高精度的補(bǔ)償，證明了所提算法的正確性和實(shí)用性。

圖7 實(shí)驗(yàn)裝置Fig. 7 Experimental device

圖8 實(shí)際補(bǔ)償效果Fig. 8 Actual compensation effect

4 結(jié)論

本文提出了地磁導(dǎo)航中采用一體化模型的地磁傳感器雙噪聲聯(lián)合估計(jì)補(bǔ)償方法。該方法有效抑制了渦流磁場(chǎng)和雙噪聲帶來(lái)的不良影響，提高了地磁傳感器補(bǔ)償穩(wěn)定度和補(bǔ)償精度。仿真和實(shí)驗(yàn)結(jié)果表明，采用一體化模型的雙噪聲聯(lián)合估計(jì)補(bǔ)償方法補(bǔ)償結(jié)果波動(dòng)明顯減小，且補(bǔ)償精度提高一個(gè)數(shù)量級(jí)，具有一定的工程應(yīng)用價(jià)值。