基于卡爾曼濾波的磁偶極子目標(biāo)張量定位

楊明明，連麗婷，高守勇

( 中國人民解放軍91388 部隊，廣東 湛江524022)

0 引 言

水中鐵磁物體，如鋼鐵航行體、水雷等，在地磁場的作用下會引起局部地磁異常，磁定位技術(shù)就是利用磁異常信號來確定這些磁性目標(biāo)的位置和磁矩等參數(shù)［1-3］?；诖排紭O子模型的磁性目標(biāo)定位常見于航空探潛、水下地磁反潛等領(lǐng)域。在航空探潛領(lǐng)域，通常先采用聲探測技術(shù)來確定可疑海域，再用裝載于飛機(jī)上的磁探儀接收到的潛艇磁場信號推算目標(biāo)潛艇的運(yùn)動要素及其與飛機(jī)的相對位置等參數(shù)，實現(xiàn)目標(biāo)精確定位。

在各種磁性目標(biāo)定位方法中，按采用的物理量可以分為基于總量磁場(梯度)［4］、磁場分量［5-7］和單點梯度張量［8-9］等方法，其中單點梯度張量定位方法基于磁偶極子目標(biāo)的線性探測方程，具有計算量小、實時性好的優(yōu)點，一直以來是研究的熱點［8-10］。文獻(xiàn)［8］對單點梯度張量定位方法進(jìn)行了較為深入的研究，指出目標(biāo)磁矩特征、測點與目標(biāo)距離、磁探儀分辨率和磁梯度測量系統(tǒng)的基線距離是影響定位結(jié)果的4個主要因素，但卻忽略了測量噪聲對定位精度的影響分析。在航空探潛的實際工作中，磁場測量噪聲不可避免地存在，其主要來源為測量系統(tǒng)噪聲和飛機(jī)磁場干擾［12］，這些測量噪聲對定位結(jié)果具有重要影響。為此，本文以航空探潛和水下小目標(biāo)磁場探測為應(yīng)用背景，在考慮磁場測量噪聲情況下，采用卡爾曼濾波理論對磁偶極子目標(biāo)梯度張量定位進(jìn)行研究，并設(shè)計了相應(yīng)的仿真實驗對其進(jìn)行驗證。

1 目標(biāo)定位數(shù)學(xué)模型

在卡爾曼濾波框架下，磁性目標(biāo)張量定位的數(shù)學(xué)模型稱為狀態(tài)空間模型，狀態(tài)空間模型由用于描述系統(tǒng)動態(tài)特性的狀態(tài)模型和描述系統(tǒng)測量信息的量測模型2 部分組成。

1.1 狀態(tài)模型

狀態(tài)模型的3個要素為狀態(tài)變量，狀態(tài)轉(zhuǎn)移矩陣和狀態(tài)噪聲。構(gòu)建狀態(tài)模型的基本原則是所建立的狀態(tài)模型既要符合系統(tǒng)動態(tài)特性實際，又要便于數(shù)據(jù)處理。在航空探潛過程中或小目標(biāo)磁場探測中，傳感器平臺的位置、速度和姿態(tài)一般由平臺的導(dǎo)航系統(tǒng)提供，為已知量，在進(jìn)行張量定位方法分析時可假定測量平臺靜止而目標(biāo)運(yùn)動，如圖1所示，其中坐標(biāo)系原點o 為探測系統(tǒng)的中心，x 軸，y 軸和z軸分別平行于地磁坐標(biāo)系相應(yīng)的3 軸。

圖1 目標(biāo)張量定位示意圖Fig.1 Schematic diagram of target gradient tensor localization

由于探測過程持續(xù)時間較短，一般假定在定位過程中潛艇目標(biāo)作深度不變的勻速直線運(yùn)動。基于此，可采用目標(biāo)在水平方面航速vx，vy以及目標(biāo)與探測系統(tǒng)的相對位置r 來描述其運(yùn)動過程的動態(tài)特性，即可定義狀態(tài)變量x=(x，y，z，vx，vy)T，則

式中T 為測量設(shè)備采樣周期;k=0，1，2，…N 為采樣時刻。

式(1)寫成矩陣的形式為Ak+1=Ak，k－1xk+ wk，其中Ak，k－1為狀態(tài)轉(zhuǎn)移矩陣，狀態(tài)噪聲wk=(wx，k，wy，k，wz，k，wvx，k，wvy，k)T為相互獨立的高斯白噪聲，其方差為R。

1.2 量測模型

磁偶極子是磁性目標(biāo)定位中使用最為廣泛且最有效的數(shù)學(xué)模型。對位于源點S的磁偶極子，其空間中任意場點P 處的磁感應(yīng)強(qiáng)度B (Bx，By，Bz)可按下式計算:

式中:r 為由源點S 指向場點P的位置矢量;m 為磁偶極子磁矩;μ 為介質(zhì)磁導(dǎo)率。

磁感應(yīng)強(qiáng)度的3個分量(Bx，By，Bz)在空間3個方向(x，y，z)的變化率構(gòu)成一個張量矩陣G=［Gij］3×3，Gij的表達(dá)式為:

式中:i，j=1，2，3;B1，B2和B3分別表示磁感應(yīng)強(qiáng)度的3個分量Bx，By和Bz;r1，r2和r3分別表示矢徑r的3個分量rx，ry和rz;δij為Kronecker-Delta 函數(shù)。

由文獻(xiàn)［8-9］可知，磁偶極子磁場梯度張量與磁場的關(guān)系為:

在磁偶極子目標(biāo)定位實際工程中，目標(biāo)磁場一般采用精度較高的超導(dǎo)量子磁力儀進(jìn)行測量，其測量精度可達(dá)到pT 級，磁場測量值中包含目標(biāo)磁場值、地磁噪聲和載體磁場干擾等。磁場梯度張量可采用5個三軸磁傳感器［9］或10個單軸磁傳感器［13］的輸出經(jīng)差分后得到，該測量方式極大地降低了同模測量噪聲的干擾，因此可認(rèn)為梯度張量測量值即為目標(biāo)梯度張量真值。為此，觀測方程可寫為:

式中:Yk和xk分別為k 時刻的觀測值(磁場測量值)和狀態(tài)向量;vk為k 時刻的測量噪聲，一般認(rèn)為其為方差為Q的高斯白噪聲;Hk為觀測矩陣，具體表達(dá)式為:

此外，在估計得到目標(biāo)與探測系統(tǒng)的相對位置后，可采用式(7)計算目標(biāo)磁矩參數(shù)，該參數(shù)可用于目標(biāo)的分類識別。

1.3 可觀測性分析

采用卡爾曼濾波求解式(1)和式(5)描述的磁性目標(biāo)定位模型的前提是定位模型可觀測。由線性系統(tǒng)理論，可采用如下的觀測性矩陣Qo來分析其可觀測性:

若Qo矩陣的秩大于等于待觀測向量維數(shù)，則系統(tǒng)可觀測。將H，A 代入上式后，取Qo陣的某一六階子陣如下:

明顯，若| | G| | ≠0，則TG 至少存在一個行列式不為0的2 階子陣;從而Qo1至少存在一個行列式不為0的5 階子陣，即若| | G| | ≠0，系統(tǒng)可觀測。

2 卡爾曼濾波基本原理

卡爾曼濾波是一種線性最小方差統(tǒng)計估算方法，它是通過處理一系列帶有誤差的實際測量數(shù)據(jù)而得到物理參數(shù)的最佳估算。其實質(zhì)是以預(yù)測方程和測量方程為基礎(chǔ)，用遞推的方法計算系統(tǒng)狀態(tài)變量xk的最小均方估計值。對于式(1)和式(5)描述的狀態(tài)空間模型，其遞推公式為:

狀態(tài)預(yù)測:

方差預(yù)測:

增益計算:

狀態(tài)估計:

方差更新:

給定初值x0和P0后，根據(jù)k 時刻(k=1，2，…)的觀測值Yk，就可以遞推計算得k 時刻的狀態(tài)估計的最小均方估計值。

3 數(shù)值算例分析

采用與文獻(xiàn)［9］中類似的數(shù)值算例來檢驗本定位方法的有效性，數(shù)值算例各參數(shù)具體取值如表1所示。

表1 數(shù)值算例參數(shù)Tab.1 Numerical example parameters

卡爾曼濾波迭代起始點x0由相對位置和水平方向的速度兩部分組成，其中相位位置可通過直接求解式(4)得到，而在迭代起始時刻的速度無法獲知，因此取速度的2個分量為0。P0隨著遞推的進(jìn)行，將較快收斂，因此實際使用時，P0可以取得稍大一些。在該算例中取P0=diag (10 m，10 m，10 m，5 m/s，5 m/s)。

分析步驟為:1)假設(shè)目標(biāo)與探測系統(tǒng)的相對速度vx=3 m/s和vy=1 m/s，由采樣周期和探測線L 計算各采樣時刻目標(biāo)與探測系統(tǒng)的相對位置(記為理論值);由相對位置和磁偶極子磁矩按式(2)計算各采樣時刻的磁場感應(yīng)強(qiáng)度，并由文獻(xiàn)［9］中方法計算相應(yīng)的目標(biāo)磁場梯度張量;2)分別采用文獻(xiàn)［9］中的直接計算法(記為直接計算值)和本文的卡爾曼濾波方法計算目標(biāo)位置 (記為卡爾曼濾波值)后，按式(7)計算相應(yīng)的目標(biāo)磁矩;3)比較相對位置和磁矩的直接計算值、卡爾曼濾波值與理論值，結(jié)果如圖2和圖3所示。

從圖中可以看出，磁場測量噪聲的存在導(dǎo)致直接計算的位置參數(shù)具有較大的振蕩，相應(yīng)的磁矩參數(shù)計算結(jié)果也具有較大的振蕩，說明直接計算的穩(wěn)定性較差。而采用卡爾曼濾波得到的位置參數(shù)和磁矩參數(shù)在濾波器收斂后，能有效克服測量噪聲對位置參數(shù)計算的影響，具有較強(qiáng)的穩(wěn)定性，且設(shè)計的濾波器具有較快的收斂速度。

圖2 定位結(jié)果比較圖Fig.2 Comparison chart of localization result

圖3 磁矩參數(shù)計算結(jié)果比較圖Fig.3 Comparison chart of magnetic moment calculation result

4 結(jié) 語

本文在磁偶極子模型的基礎(chǔ)上，建立了適用于磁梯度張量定位的狀態(tài)空間模型。在目標(biāo)磁感應(yīng)強(qiáng)度測量不可避免地存在測量噪聲的實際情況下，采用卡爾曼濾波對模型進(jìn)行了求解。設(shè)計了相應(yīng)的數(shù)值算例，計算結(jié)果表明，基于卡爾曼濾波的定位方法具有較高的定位精度，且具有較強(qiáng)的穩(wěn)定性。該研究工作對航空磁探測和水下小目標(biāo)磁場探測，都具有一定的指導(dǎo)意義。

［1］TAREK E T，AHMED S.Position detect ion of unexploded ordnance from airborne magnetic anomaly data using 3-D self organized feature map［J］.International Symposium on Signal Processing and Information Technology，2005:322-327.

［2］徐杰.一種高效的艦船磁性定位方法［J］.艦船科學(xué)技術(shù)，2011，3(12):69-71.

［3］唐莉莉，宋勇，童官軍，等.一種新的運(yùn)動磁目標(biāo)定位算法［J］.傳感技術(shù)學(xué)報，2011，24(7):996-1000.

［4］XIAO Chang-han，LIU Sheng-dao，ZHOU Guo-hua.Rea-l time localization o f a magnetic object with total field data［M］.US:IEEE，2008.

［5］李華，李慶民，龔沈光，等.基于磁性目標(biāo)磁場模型的在線磁定位研究［J］.武漢理工大學(xué)學(xué)報(交通科學(xué)與工程版)，2009，33(1):165-168.

［6］黃玉，郝燕玲.水下地磁異常反演中位置磁矩聯(lián)合迭代算法［J］.華中科技大學(xué)學(xué)報(自然科學(xué)版)，2011，39(7):95-98.

［7］趙文春，劉大明，周國華.基于微分進(jìn)化算法的磁性目標(biāo)參數(shù)估計方法［J］.艦船電子工程，2011，31(167):167-170.

［8］肖昌漢，何華輝.磁性物體探測的線性求解方法［J］.華中理工大學(xué)學(xué)報，1997，25(4):83-85.

［9］張朝陽，肖昌漢，閆輝.磁性目標(biāo)單點梯度張量定位方法［J］.探測與控制學(xué)報，2009，31(4):44-48.

［10］王光源，馬海洋，章堯卿.航空磁探儀探潛目標(biāo)磁梯度定位方法［J］.兵工自動化，2011，30(1):32-34.

［11］李居偉，徐以成，孫明太.航空磁探儀目標(biāo)運(yùn)動分析仿真研究［J］.指揮控制與仿真，2011，33(4):72-74.

［12］MARIUS B.Recursive bayesianmethod for magnetic dipole tracking with a tensor gradiometer［J］.IEEE Transactions onmagnetics，2011，47(2):409-415.

［13］譚斌，林春生，張寧.地磁場梯度對飛機(jī)磁場求解精度的影響分析［J］.武漢大學(xué)學(xué)報(信息科學(xué)版)，2011，36(12):1482-1485.

［14］黃玉，郝燕玲.磁場梯度張量定位中的矢量磁力計最簡配置［J］.華中科技大學(xué)學(xué)報(自然科學(xué)版)，2010，38(12):52-55.