基于水下運動平臺的磁性目標跟蹤與特征反演

馬劍飛， 丁凱， 顏冰， 林春生

(1.海軍工程大學 兵器工程學院， 湖北 武漢 430033；2.近地面探測技術重點實驗室， 江蘇 無錫 214035)

0 引言

相比于直接測量目標的三分量磁場，磁梯度張量信息包含5個獨立變量，具備不受地磁場波動的影響以及可以更加全面地反映磁場細節(jié)等突出的優(yōu)勢。在淺海環(huán)境下，聲納探測易受水文條件以及地形等因素的影響，磁梯度張量被動探測相對而言具有抗干擾能力強、定位精度高以及隱蔽性好的優(yōu)勢[1]。在軍事應用中，特別是在港口地區(qū)，可以有效探測敵方的潛艇、水雷目標或者未爆彈(UXO)。在民用領域，對于沉船打撈或者海洋礦物勘探等應用都有一定的應用價值[2-3]。

目前，基于磁梯度張量的探測技術逐步成為磁探測領域的研究熱點，德國、美國、澳大利亞等國家已開發(fā)出磁梯度張量探測系統(tǒng)，并開展了對應的航空航海試驗[4-6]。相比于航空磁探，基于水下運動平臺的磁探測能夠更加接近磁性目標，因此其測量的信噪比較高，特別是對于微弱的磁性目標檢測而言意義尤為重要，而且利用有限的幾個測量點能夠很快反演出磁性目標的位置，進一步可以反演出目標的形狀以區(qū)分目標種類。水下運動平臺受風浪的影響較小，比在海面上運行更加穩(wěn)定[1]，因此水下運動平臺的運動噪聲也會比較低。

現(xiàn)有的磁梯度定位方法大多都是將磁性目標視為偶極子模型[7-9]，在遠場情況下這種簡化模型是可靠的，但在近場情況下這種忽視目標尺度的模型導致定位結果出現(xiàn)較大的誤差。另外，遠場的磁梯度信號較微弱，定位結果受測量噪聲影響較大，而近場情況下信噪比較高，目標尺度因素(即模型不匹配)則是定位誤差的主要來源。Clark[10]與Beiki等[11]研究了常見幾何體的構造指數(shù)、幾何項以及其位置的反演方法，但是一般情況下目標的構造指數(shù)和幾何項是未知的，因此很難應用于實際。

在磁性目標模型選擇方面，吳志東等[12]研究了基于磁性橢球體模型的磁梯度張量水面艦船目標跟蹤問題，利用提出的高斯混合采樣粒子濾波算法取得了較好的跟蹤結果，在擬合目標磁梯度時假定橢球體的焦距是已知的，然而焦距參數(shù)是一個未知的參數(shù)，即使給定一個初值，由于焦距觀測性比較低[13]，濾波器更新過程中也很難趨近于真值，因此利用該模型去反演目標尺度特征是不合適的。戴忠華等[14]利用偶極子模型與橢球體的混合模型實現(xiàn)對艦船磁場模擬，已知目標與測量節(jié)點的距離，利用大量節(jié)點的磁場測量數(shù)據(jù)去反演混合模型的參數(shù)，實現(xiàn)了與艦船磁場很好的擬合，得到了相對較好的反演效果。張宏欣[15]以線列陣偶極子模型模擬艦船的磁場，并利用水下雙節(jié)點測量的三軸磁場對實現(xiàn)了艦船目標的準確定位跟蹤。事實上，利用多個偶極子完全可以對尺度目標的磁場特性進行精確逼近，磁偶極子的分布方式也可以有效反演目標的尺度特征。

在磁性目標分類的研究中，葛健等[16]研究了不同姿態(tài)航彈對應的垂直磁梯度空間分布特征，證明了利用磁梯度特征識別掩埋磁性目標的可行性。張敬東[17]利用小波多尺度分解數(shù)據(jù)處理方法有效提高了弱磁梯度異常信息的分辨率。鄭建擁等[18]以磁總量模值、磁梯度5個獨立變量和3個特征值不變量作為特征，利用改進的支持向量機方法提升了目標分類的準確率。Zhang等[19]基于多面體磁模型采用最陡下降法對目標的形狀進行了反演，證明了目標磁梯度能夠彌補磁場信息的不足，提升反演的精度。然而，上述兩種方式只局限于在特定距離條件下對特定目標的分類，在水下目標分類應用中不具有實際意義。

磁性目標跟蹤濾波器的觀測信號一般為磁矢量信號，目標磁矩、尺寸以及速度等信息隱含于觀測信號之中，而且濾波跟蹤算法先驗信息缺失，設定的初值范圍不可避免地引入了人為經(jīng)驗，大量的未知參數(shù)也導致其跟蹤精度有限。相比較而言，以磁梯度張量作為估值濾波器觀測信號的優(yōu)勢在于：1)磁梯度張量擁有5個獨立變量，可以提高觀測矢量的維數(shù)；2)磁梯度張量相對于磁矢量可以更好地抑制環(huán)境噪聲的干擾；3)磁梯度定位算法能夠為跟蹤濾波器提供較為準確的濾波器磁矩初值以及位置初值。值得注意的是，隨著估值濾波器觀測值的不斷引入，偶極子的空間分布特性也逐漸顯現(xiàn)，其分布特性能夠在一定程度上反映磁性目標的尺度和結構特征。

針對磁性尺度目標模型不匹配問題，本文旨在結合磁梯度張量與濾波估值算法的優(yōu)勢，以多磁偶極子為模型研究磁性尺度目標的跟蹤方法，并根據(jù)多偶極子的空間分布狀態(tài)研究目標尺度和結構特征的反演算法。

1 濾波估值算法

目前以目標磁場作為信號源對運動目標進行跟蹤的研究已經(jīng)十分廣泛，其估值濾波算法主要分為粒子濾波框架與卡爾曼濾波框架，然而這兩類濾波器的估值精度對濾波器初值比較敏感，不能保證任意初值條件下的濾波收斂性[15]，但利用磁梯度定位算法可以為濾波器提供較準確的距離、速度以及磁矩初值，從而提升磁性目標跟蹤的準確性。

魯棒估計是指在出現(xiàn)異常量測情況下，利用適當方法對目標狀態(tài)進行估計，減小異常量測值對狀態(tài)估計帶來的誤差。當出現(xiàn)觀測野值時，傳統(tǒng)卡爾曼濾波算法在對狀態(tài)進行估計時會產(chǎn)生較大的誤差。為了對目標狀態(tài)進行穩(wěn)健的估計，人們將魯棒統(tǒng)計學應用于濾波算法中，提出了基于Huber方法的線性卡爾曼濾波算法[20]，實質是一種廣義極大似然估計的卡爾曼濾波算法。文獻[21]在線性魯棒卡爾曼濾波的基礎上，揭示了Huber方法在濾波算法中對卡爾曼增益和狀態(tài)估計的影響，并基于此推導基于無跡卡爾曼濾波(UKF)器的非線性魯棒濾波算法，以避免線性化帶來的額外誤差。

考慮如下運動目標：

xk+1=f(xk)+wk,

(1)

其量測方程為

zk+1=h(xk+1)+vk+1,

(2)

式中：xk為k時刻的運動狀態(tài)矢量；f為狀態(tài)轉移函數(shù)；wk為k時刻的運動噪聲；zk+1為k+1時刻的觀測矢量；h為觀測函數(shù)；vk+1為k+1時刻的測量矢量。

構造如下非線性回歸方程：

(3)

定義如下變量：

式中：Rk+1為觀測噪聲的協(xié)方差矩陣。則回歸變量

yk+1=g(xk+1)+ξk+1.

(4)

據(jù)此可定義Huber代價函數(shù)定義為

(5)

式中：ek+1,i為yk+1-g(xk+1)的第i個分量；m為回歸變量的維數(shù)；ρ函數(shù)為

(6)

u為門限值，e為估計誤差。

對ρ函數(shù)關于e求導，可得φ函數(shù)為

(7)

由(6)式可知，當誤差|e|

對(5)式進一步求導，可得

(8)

(9)

(10)

(11)

由于真實的目標狀態(tài)xk+1未知，令殘差δxk+1|k=0，則Ψx為單位矩陣I，即

(12)

可得Huber方法改進后的UKF算法等效量測噪聲為

(13)

由上述推導過程可知，非線性Huber無跡卡爾曼濾波(NHUKF)算法通過最小化目標代價函數(shù)來抑制異常觀測值對目標狀態(tài)估計帶來的影響，而且沒有對非線性方程的線性化過程，從而避免了傳統(tǒng)算法的線性化誤差[20]。

2 建立跟蹤估值模型

圖1所示為一個水下運動平臺磁梯度跟蹤定位的典型示意圖，其中平臺的運動參數(shù)可以通過平臺的慣性導航準確獲取，因此平臺的速度參數(shù)不需要作為狀態(tài)估計量?？紤]到需要通過偶極子的分布特征來判斷目標形狀，在三維空間中最少需要4個點才能表達立體結構特征，而過多的偶極子會導致狀態(tài)估計量的維數(shù)激增從而大幅增加運算量，故選用4偶極子模型作為NHUKF算法的狀態(tài)估計模型。

圖1 載體平臺運動示意圖Fig.1 Schematic diagram of moving carrier platform

目標的磁性來源主要是鐵磁性材料，而且同一物體經(jīng)過的外磁場磁化過程也大致相同，同時也是出于降低觀測量的維數(shù)和利用磁偶極子分布特征反映物體形狀特征的考慮，假設4個偶極子的磁矩是相同的。綜上所述，設定濾波器狀態(tài)估計量為

(14)

由于磁梯度張量僅包含5個獨立變量，狀態(tài)變量是15維的，而觀測量維數(shù)過低導致卡爾曼增益產(chǎn)生的狀態(tài)量更新方向不一定是真實解的方向，容易陷入歧義解的范圍。在水下運動平臺磁定位的應用中，平臺速度是已知的，可以利用前面時刻的測量磁梯度擴充觀測量的維數(shù)，相當于變相地進行陣列跟蹤定位，因此設置觀測量為

zk=[Gk-2n,Gk-n,Gk],

(15)

式中：Gk、Gk-n、Gk-2n分別表示k時刻、k-n時刻和k-2n時刻5個獨立磁梯度變量組成的向量。不難看出，此時觀測量和狀態(tài)量都是15維的向量，可以進行準確地磁梯度跟蹤定位。

3 目標尺度反演方法

NHUKF估值濾波算法的偶極子位置估計值是4個離散的三維空間點坐標，對于三維離散點要提取其分布特征，直觀的想法是通過線性擬合的方式提取其主要的分布尺度與方向。然而使用多元線性回歸在三維空間擬合的結果是一條曲線，很難反映目標的主尺度特征，而且現(xiàn)有擬合方式大多通過計算預測值和觀測值的均方誤差來解算擬合函數(shù)，對于求解主尺度特征而言，更希望以各點到直線的距離為代價函數(shù)。

主成分分析(PCA)法是一種通過正交變換將一組可能存在線性相關性的變量轉換為一組線性不相關變量的統(tǒng)計方法，常用于特征降維以更好地發(fā)掘數(shù)據(jù)模式。下面結合PCA法的計算過程和目標主尺度分析的需要進行簡單說明。假設偶極子的坐標分別為(xi,yi,zi)，首先計算這些偶極子坐標的協(xié)方差矩陣，即

(16)

對C進行矩陣正交分解，可以得到3個特征值λ1>λ2>λ3，分別對應3個正交的單位特征向量u1、u2、u3.顯然，特征值越大表示這些離散點在該特征向量方向投影的方差就越大，最大特征值對應的特征向量為第1主成分。對于目標尺度分析的應用而言，第1主成分表示目標的最大尺度在此方向上通過各離散點在第1主成分上的投影就可以估計出目標的主尺度，各離散點在第1主成分方向的投影為

pi=[xi,yi,zi]u1.

(17)

常見的掩埋目標形狀主要為圓柱體以及圓盤兩種典型結構，這里約定當圓柱體的高H大于底面半徑時稱之為圓柱，當圓柱體的高H小于底面半徑時稱之為圓盤。對于圓柱體而言，第一主成分方向與圓柱體的軸線一致，其主尺度對應圓柱體的高，偶極子分布接近于一條直線，此時λ1應該遠大于λ2和λ3，因此目標的主尺度可估計為

(18)

式中：N為狀態(tài)估計模型的偶極子數(shù)。

對于圓盤而言，λ1和λ2的大小相差不多，且遠大于λ3，其第1主尺度和第2主尺度反映的是圓盤的直徑特征。鑒于本文利用4個偶極子模擬目標的磁場特征，在分析主尺度和圓盤直徑的代數(shù)關系時采用4個小圓來模擬大圓的磁場。根據(jù)等面積法，小圓盤的半徑r等于大圓盤半徑的一半，顯然，若小圓覆蓋的區(qū)域與大圓覆蓋的區(qū)域重合部分越多，則擬合的效果越好。圖2所示為兩種典型的覆蓋策略。其中第1種策略是小圓之間相交且與大圓相切，此時所估計的小圓半徑為

圖2 偶極子擬合策略Fig.2 Fitting strategy of multiple dipoles

r1=L/2,

(19)

式中：L為4個小圓圓心組成的正方形的對角線長度。第2種策略是小圓之間相切且與大圓相交，此時所估計的小圓半徑為

(20)

對于正方形頂點組成的離散點集合，利用PCA法可以推導出其協(xié)方差矩陣的最大特征值λ1與其對角線長度L的關系為

(21)

綜合以上兩種策略，可以推導得到圓盤目標的主尺度為

(22)

4 仿真實驗

為了驗證算法的可行性，開展利用搭載磁梯度張量系統(tǒng)的運動載體平臺探測靜止磁性目標的仿真實驗，實驗內(nèi)容主要分為兩部分：第一部分檢驗NHUKF算法的定位性能；第二部分檢驗PCA法在目標尺度反演時的性能。

4.1 基于NHUKF算法的目標定位

設磁性目標的中心為空間直角坐標系的原點，x軸與水下運動平臺的軸線重合并指向自主水下機器人(AUV)平臺的艏部，y軸指向右舷，z軸豎直向下。仿真條件設置為：初始時刻目標的位置為(-40 m， 5 m， -5 m)；磁梯度張量測量系統(tǒng)基線Ls為0.5 m；采樣頻率fs為5 Hz；AUV平臺運動速度為(5 m/s，0 m/s，0 m/s)；磁梯度系統(tǒng)的測量精度為0.01 nT. 如圖3所示，目標是高為3 m、半徑為0.5 m的圓柱體，目標橫滾角、俯仰角以及方位角分別為π/6 rad、π/4 rad和π/3 rad，磁梯度正演基于密集的偶極子網(wǎng)格，網(wǎng)格邊長為0.1 m，每個磁偶極子的磁矩M為(-3 A·m2， 1 A·m2， 5 A·m2)。跟蹤濾波器的初值由3階張量法確定[8]。

圖3 磁性目標正演模型Fig.3 Forward modeling of magnetic target

如圖4所示，在無異常觀測的情況下，NHUKF算法與UKF算法的定位性能相近，隨著觀測值的不斷引入，其定位性能要優(yōu)于3階張量法，定位誤差小于0.1 m. 在接近目標時，多數(shù)時刻3階張量法都有較好的定位結果，但有的時刻其定位矩陣條件數(shù)很大(與磁偶極子的特征平面有關)，稍微的擾動或者觀測精度不足就會導致定位結果發(fā)生較大的偏差。如圖5所示，異常觀測的引入在距目標-20 m、-10 m以及0 m的時刻，不難看出NHUKF算法相對于UKF算法而言，能夠有效消除異常觀測對定位性能帶來的影響。另外可以發(fā)現(xiàn)3階張量法的定位性能受異常觀測影響十分大，在距目標-20 m、-10 m以及0 m的時刻，其定位誤差都出現(xiàn)了比較大的峰值。

圖4 無異常觀測時的定位性能比較Fig.4 Comparison of localization performances without abnormal observation

圖5 存在異常觀測時的定位性能比較Fig.5 Comparison of localization performances with abnormal observation

由上述分析可知，利用磁梯度張量的NHUKF算法夠有效定位尺度目標，其定位性能要優(yōu)于傳統(tǒng)磁梯度張量定位算法，而且能夠有效消除異常觀測的影響。

4.2 基于PCA算法的目標尺度反演

對于掩埋或者半掩埋的磁性目標而言，大多數(shù)目標都具有規(guī)則的幾何形狀，通過反演的偶極子分布特征、磁矩特征和目標尺度特征，在很大程度上對目標的種類進行區(qū)分，也可以排除環(huán)境的干擾。本節(jié)基本仿真條件與4.1節(jié)一致，分別檢驗PCA法在不同磁矩、不同姿態(tài)、不同尺度以及不同形狀情況下對目標尺度的反演性能。

4.2.1 不同磁矩的目標尺度反演

表1所示為不同目標磁矩[Mx,My,Mz]對應的反演結果，其中Mx、My和Mz分別表示單個偶極子在x軸、y軸和z軸方向的磁矩，主尺度方向表示偶極子分布協(xié)方差矩陣的最大特征值λ1對應的單位特征向量。從表1中可以看出，反演主尺度方向誤差在1°的范圍內(nèi)，反演主尺度模值的誤差小于0.3 m，正演偶極子的磁矩方向不會影響算法反演的性能，而且算法可以準確估計出目標的整體磁矩特征。

表1 不同磁矩的目標反演

圖6所示為不同偶極子磁矩對應的反演結果。由圖6可以看出，反演的偶極子位置雖然略有差別，但基本都是沿圓柱體的軸線分布，反演的主尺度參數(shù)能夠很好地反映出其軸線的尺度和方向的特征。在濾波器狀更新過程中，有限的幾個偶極子通過不斷調(diào)整位置更好地擬合目標磁場特征，顯然當偶極子沿主尺度分布時能夠更好地減小濾波器的預測值和觀測值的誤差。

圖6 不同磁矩的目標反演Fig.6 Inversion of targets with different magnetic moments

4.2.2 不同姿態(tài)的目標尺度反演

表2所示為在不同目標姿態(tài)條件下對目標尺度的反演結果，其中[θ,β,φ]中θ為橫滾角、β為俯仰角、φ為方位角。從表2中可以看出，不同姿態(tài)目標的反演結果都比較準確，反演的主尺度方向誤差在2°范圍內(nèi)，反演的主尺度模值誤差小于0.3 m.

表2 不同姿態(tài)的目標反演

圖7所示為不同姿態(tài)的目標反演結果，反演的偶極子沿圓柱體軸線分布，反演的主尺度也對應于圓柱體軸線，反演的姿態(tài)角也與目標真實姿態(tài)角十分接近，表明本文的主尺度估計算法不受姿態(tài)角影響，而且該算法可以準確估計目標的姿態(tài)。

圖7 不同姿態(tài)的目標反演Fig.7 Inversion of targets with different attitudes

4.2.3 不同尺度的目標尺度反演

表3所示為不同目標尺度對應的反演結果，其中R、H分別表示圓柱的底面半徑和高。從表3中可以看出，不同尺度目標的反演結果都比較準確，反演的主尺度方向誤差在2°的范圍內(nèi)，反演的主尺度模值誤差小于0.3 m.

表3 不同尺度的目標反演

圖8所示為不同尺度目標的反演結果，其中反演的偶極子沿圓柱體軸線分布，在不同目標尺度情況下，反演的主尺度模值和方向與圓柱體軸線的模值和方向相近，表明反演算法可以準確估計不同大小磁性圓柱體的尺度特征。

圖8 不同尺度的目標反演Fig.8 Inversion of targets with different scales

表4 不同的圓柱目標反演

4.2.4 不同形狀的目標尺度反演

表4對比了實心以及空心圓柱體的反演結果，表5對比了實心以及空心圓盤的反演結果，其中實心表示無空腔的目標；外殼表示空心目標的外部殼體；空腔表示目標內(nèi)部的中空部分。由表4和表5可以看出，空腔對目標反演的結果影響不大，算法仍然可以準確估計目標主尺度的方向以及模值，圓柱體偶極子分布協(xié)方差矩陣的特征值滿足λ1?λ2>λ3的關系，圓盤偶極子分布協(xié)方差矩陣的特征值滿足λ1≈λ2≥λ3的關系，以此特征值關系可以區(qū)分目標為圓柱或者圓盤，聯(lián)合目標的磁矩特征以及尺度特征就可以推斷目標的類別。

表5 不同的圓盤目標反演

如圖9(a)、圖9(b)、圖9(c)、圖9(d)所示，圓柱狀結構目標的反演偶極子基本沿軸線分布，圓盤狀結構目標的反演偶極子沿中心切面對稱分布。對于圓柱狀結構目標，第1主尺度反映的是其軸線的尺度；對于圓盤狀結構目標，第1主尺度和第2主尺度反映的其橫截面的尺度。

圖9 不同形狀的目標反演Fig.9 Inversion of targets with different structures

5 結論

本文利用多個偶極子組合模型，結合具有魯棒性的NHUKF算法研究一種針對尺度磁性目標的磁梯度定位方式，該方法能夠有效提高磁性尺度目標的定位精度，特別是當出現(xiàn)異常觀測時仍然可以準確定位目標，同時也可以準確反演偶極子的分布特征和磁矩特征，具有良好的魯棒性。進一步提出一種基于偶極子分布模式反演目標尺度特征的PCA算法，該方法對不同磁矩、不同姿態(tài)、不同尺度以及不同形狀的目標，主尺度反演誤差小于0.3 m，對目標主尺度方向的反演誤差小于2°，亦可識別目標主尺度為線結構或者面結構，為水下目標種類判別提供了一定的參考。