未知姿態目標航磁異常探測的最優飛行方向分析

石 昱, 卞雷祥,*, 莊志洪, 付夢印, 王春娥, 秦 杰

(1. 南京理工大學機械工程學院, 江蘇 南京 210094; 2. 南京理工大學電子工程與光電技術學院,江蘇 南京 210094; 3. 南京理工大學自動化學院, 江蘇 南京 210094;4. 北京自動化控制設備研究所, 北京 100074)

0 引 言

隨著安靜型潛艇技術的發展,大多數潛艇航行噪聲已接近甚至低于海洋背景噪聲[1],聲探測等傳統探測方法面臨極大難度。磁異常探測(magneic anomaly detection, MAD)被認為是有效的非聲探潛手段之一[2]。此外,MAD在資源勘探[3-5]、未爆彈探測[6-8]、與地震預測[9-10]等領域也發揮著重要作用。

目前,MAD算法的定位精度與探測概率都依賴于信號的信噪比,基于磁張量縮并量正交基分解的目標探測定位算法,在信噪比低于-20 dB條件下,檢測效果受噪聲影響較大,無法確定檢測能量的峰值點[11];在噪聲均方根為0.01 nT的條件下,總場梯度(total field gradient, TFG)-張量縮并量梯度(tensor contraction gradient, TCG)定位方法中三軸的最大定位誤差接近0.8 m[12]。Wang等人[13]利用迭代算法對傳統三角定位測量(scalar triangulation and ranging, STAR)方法中的非球面誤差進行了補償,該方法在信噪比為40 dB、20 dB時,定位誤差分別在1%、8%左右,而傳統STAR方法此時的定位誤差分別在15%、20%附近。Yin等人[14]對基于歸一化磁源強度(normalized source strength, NSS)的STAR方法進行了仿真分析,仿真中該方法的探測范圍約為15 m,當信噪比高于19.4 dB時,相對位置誤差小于5%;當信噪比為-3 dB時,基于噪聲檢測的隨機共振法探測概率降到了70%[15];信噪比為-3 dB時,高階過零法(high order zero crossing, HOC)探測概率為0.74～0.82[16]。

對于隱匿目標的磁異常探測,需要選擇合適的探測路徑方向來獲得高信噪比的磁異常信號,提升目標的探測識別概率與定位精度。曲曉慧等人[17]利用蒙特卡羅法仿真了幾個特定飛行航向下磁探儀的搜潛概率,仿真結果表明飛機沿南北方向飛行搜索概率更大。吳芳等人[18]仿真了多個飛行航向下的目標磁異常信號,結果表明反潛機探測航向為0°(南北方向)、45°時,磁異常信號變化幅度皆大于其他探測航向。 Moser[19]研究了12個飛行航向下的目標磁異常信號,不同航向下磁異常信號峰峰值差異明顯,其中30°方向時信號變化幅度最大,但南北方向信號變化幅度與其極為接近。Schluckebier[20]分析了磁異常探測中飛機飛行航向與水下目標航向對探測距離的影響,無論水下目標的航向如何,飛行航向為南北方向時,磁異常探測距離都較大。

上述的研究[17-20]通過對多個特定探測路徑方向的磁異常信號進行比較,說明了沿著南北方向的磁感應強度變化較大。實際上,鐵磁性目標磁矩與其結構、姿態有關,運動目標的姿態是變化的,建立探測路徑上磁異常信號與目標姿態的定量關系,對于規劃路徑和目標特性反演定位均具有重要意義。

本文針對典型目標建立了磁矩精確計算模型,利用該模型建立了目標磁異常信號與其姿態之間的定量關系;計算并分析了航磁探測時飛行平面內磁異常變化最大路徑與背景磁場(地磁場)之間的關系,以此得到磁異常探測的最優飛行方向,最后借助多物理場仿真軟件進行了驗證。

1 理論分析

1.1 鐵磁性目標磁矩計算模型

1.1.1 實心鐵磁性目標的磁矩

在分析目標磁性時,只考慮固定磁性與感應磁性。固定磁性在短期內不會發生變化,感應磁性會受到磁性目標的姿態、形狀、體積、材料特性以及地磁場等因素的影響。因此,需要根據上述因素建立磁矩計算模型。

鐵磁性目標的磁矩M可以表示為

M=VJtotal

(1)

式中:Jtotal為目標的總磁化強度,由感應磁化強度與剩余磁化強度兩部分組成;V為磁性目標的體積。

鐵磁性目標的載體坐標系如圖1所示,l為其長軸方向;t為其短軸方向;v為其豎直方向。

圖1 鐵磁性目標的載體坐標系Fig.1 Carrier coordinate system for ferromagnetic targets

鐵磁性目標在l,t,v方向上的磁矩Ml,Mt,Mv可以表示[21-25]為

(2)

式中:Nl,Nt,Nv分別為鐵磁性目標在l,t,v方向上的退磁系數;HEl,HEt,HEv為地磁場在l,t,v方向上的分量;Jtl,Jtt,Jtv為總磁化強度在l,t,v方向上的分量;Jrl,Jrt,Jrv為剩余磁化強度在l,t,v方向上的分量;χ為材料的磁化率;χ=μr-1,μr為相對磁導率。

不規則形狀鐵磁體的退磁因子難以通過理論計算得到精確值。本文利用有限元軟件仿真計算出典型目標在地磁場中產生的磁矩值,將結果代入式(2),得到l,t,v3個方向上的退磁因子Nl,Nt,Nv,可以表示為

(3)

1.1.2 空心鐵磁性目標的磁矩

水下目標一般為空心鐵磁體,空心鐵磁體與實心鐵磁體產生的磁場成一定比值,可以根據實心鐵磁體的磁矩來計算空心鐵磁體的磁矩。假設空心球體內外徑之比為k,空心球體H空與實心球體H實之間產生的磁場比值R可以表示[26-27]為

(4)

當μr較大,k接近于1(球殼很薄)時R可以表示[21]為

(5)

空心鐵磁性目標各方向上的內外徑之比不統一,其在l,t,v方向上的磁場比值Rl,Rt,Rv可以表示為

(6)

式中:kl、kt、kv分別為鐵磁體在l,t,v方向上的內外徑之比。

可以借助多物理場仿真軟件計算出實心與空心鐵磁體在l,t,v方向上的磁矩之比,將結果作為Rl,Rt,Rv。

那么,空心鐵磁性目標在其載體坐標系上產生的磁矩Mls,Mts,Mvs可以表示為

(7)

1.1.3 鐵磁性目標在不同姿態下的磁矩

目標的姿態信息以及地磁場方向與探測坐標系的關系如圖2所示。假設探測坐標系x軸正方向為地理北,y軸正方向為地理西,z軸豎直向上。HE為地磁場,HExy為地磁場在xy平面(水平面)上的投影;I為HE與HExy之間的夾角,稱為磁傾角,I下傾時磁傾角為正,上仰時為負;D為HExy與地理北方向之間的夾角,稱為磁偏角,磁偏角東偏時為正,西偏時為負;l為目標的長軸方向,lxy為長軸在xy平面上的投影,α為lxy與地理北之間的夾角,稱為方位角,左偏為正;β為l與lxy之間的夾角,稱為俯仰角,下傾為正。

圖2 地磁場與目標姿態示意圖Fig.2 Schematic diagram of geomagnetic field and target attitude

目標姿態變化相當于載體坐標系在探測坐標系下發生旋轉。目標方位角與俯仰角分別為α、β時,載體坐標系相對于探測坐標系的運動如下:

(1) 繞z軸旋轉角為α;

(2) 繞y軸旋轉角為β。

繞z軸旋轉的變換矩陣Rot(z,α)可以表示為

繞y軸旋轉的變換矩陣Rot(y,β)可以表示為

對應的變換矩陣R(α,β)可以表示為

(8)

與R(α,β)對應的逆變換矩陣R(-β,-α)可以表示為

(9)

地磁場在目標載體坐標系(l,t,v)上的分量HEl,HEt,HEv可以表示為

(10)

式中:HE為地磁場值。

那么,剩磁在目標載體坐標系(l,t,v)上的分量Jl,Jt,Jv可以表示為

(11)

通過聯立式(2)、式(7)～式(11)可以得到目標在探測坐標系下的磁矩。

鐵磁性目標在探測坐標系與載體坐標系上的磁矩之間的關系可表示為

(12)

式中:MT=[Mx,My,Mz]T,MT為實心鐵磁性目標在探測坐標系上的磁矩;MC=[Ml,Mt,Mv]T,MC為實心鐵磁性目標在載體坐標系上的磁矩;MTS=[Mxs,Mys,Mzs]T,MTS為空心鐵磁性目標在探測坐標系上的磁矩;MCS=[Mls,Mts,Mvs]T,MCS為空心鐵磁性目標在載體坐標系上的磁矩。

1.2 鐵磁體磁場的分布規律

利用磁偶極子模型可以計算出鐵磁性目標在空間中的磁場分布情況,空間中某測點(x,y,z)的磁場B可以表示[28]為

(13)

式中:μ0=4π×10-7;r為目標與測點之間的相對位置矢量;M為目標的磁矩,可表示為

其中,[χHEl+Jrl,χHEt+Jrt,χHEv+Jrv]T可以表示為

(14)

該點產生的標量磁異常信號S可以表示[3]為

S=B·ie

(15)

式中:ie為地磁場方向的單位向量。

實際探測過程中載體平臺工作時所處的平面如圖3所示。假設載體平臺在z=600 m的xy平面內飛行,利用式(13)～式(15)可以計算出在探測平面內磁場分布情況,從而找到磁異常信號變化最大的路徑及其方向。假設探測平面內磁場強度的最大值、最小值位置分別為A(xmax,ymax),B(xmin,ymin),探測平面內磁異常信號變化最大的路徑方向為AB(xmax-xmin,ymax-ymin)。在進行磁異常探測時,測線一般只能布置在載體平臺所處的飛行平面內,由式(13)～式(15)可知飛行平面內磁場變化最大路徑主要與目標磁矩和地磁場有關,為此需要通過計算地磁場與磁矩在飛行平面內的投影以及平面內磁異常信號變化最大路徑之間的夾角來分析三者之間的關系。

圖3 飛行探測示意圖Fig.3 Flight detection diagram

飛行平面內磁異常信號變化最大的路徑所在直線與磁矩投影的夾角σxy可以表示為

(16)

其中,[χHEl+Jrl,χHEt+Jrt,χHEv+Jrv]T可以表示為

(17)

飛行平面內磁異常信號變化最大的路徑所在直線與地磁場投影的夾角γxy可以表示為

(18)

(19)

2 數值計算與仿真分析

2.1 水下目標磁矩的數值計算與仿真

假設水下目標為如圖4所示的圓柱形薄殼鐵磁體,將目標幾何參數設置為:半徑4.5 m,總長度84 m,兩端半球體半徑4.5 m,薄殼的厚度0.2 m;目標磁導率設置為200,剩余磁通密度設置為50 A/m[29];參考文獻[30]對地磁場參數進行設置:地磁場強度BE=50 000 nT,I=44°41′,D=-3°07′,假定地磁場值恒定;方位角α設置為0°～360°,俯仰角β設置為-90～90°,角度變化為1°。

圖4 水下目標模型Fig.4 Underwater target model

根據磁矩計算模型,得到水下目標磁矩隨其姿態的變化特性如圖5所示。目標磁矩的模值隨目標姿態變化范圍為3.55×105～2.42×106A·m2,最小值是最大值的14.67%,當目標長軸與地磁場方向一致時取最大值,當目標長軸與地磁場方向垂直時取最小值,目標長軸方向與地磁場方向偏離越大,目標的磁矩值越小。

圖5 水下目標磁矩在不同姿態下的數值計算結果Fig.5 Calculation results of magnetic moments of underwater targets in different attitudes

使用多物理場仿真軟件對水下目標磁矩隨其姿態的變化特性進行仿真,由于有限元仿真計算量過大,將角度變化設置為30°,其余參數與數值計算中的設置保持一致。仿真結果如圖6所示,磁矩模值變化范圍為3.57×105～2.34×106A·m2,目標在不同姿態的磁矩仿真結果與對應的磁矩模型計算結果差異為0.003%～0.04%,造成差異的原因可能是仿真軟件中網格劃分的不夠細致。

圖6 水下目標磁矩在不同姿態下的仿真結果Fig.6 Simulation results of magnetic moments of underwater targets in different attitudes

2.2 磁異常探測最優飛行方向分析

將z=600 m處的xy平面設置為載體平臺所處的飛行平面,其余參數設置與第2.1節中保持一致。利用磁矩計算模型與式(13)～式(19)可以計算出飛行平面內磁異常信號變化最大的路徑所在直線與磁矩投影的夾角σxy、以及與地磁場投影的夾角γxy隨目標姿態的變化特性,進而分析出飛行平面內磁矩投影方向、地磁場投影方向與磁場變化最大路徑之間的關系,得到磁異常探測最優飛行方向。

σxy隨目標姿態的變化特性如圖7所示。當目標處于不同姿態時,σxy的變化范圍為0°～89.996 4°。當目標處于某些姿態時,例如方位角為135°,俯仰角為70°時,σxy為87.66°,已接近兩個不同直線之間夾角所能達到的最大值(90°)。此時,沿著磁矩投影方向進行磁異常探測顯然是不合理的。當俯仰角為-30°～30°時,σxy最大值為44.158 8°,且35°以上區域面積很小,此時沿著磁矩方向探測可獲得變化幅度較大的磁異常信號曲線,但在實際磁異常探測中目標的磁矩信息是無法預知的。

圖7 σxy隨目標姿態的變化特性Fig.7 Variation characteristics of σxy with target attitudes

當目標處于不同姿態時,γxy隨目標姿態的變化特性如圖8所示。當目標處于不同姿態時,γxy的變化范圍為0°～50.045 7°,且變化范圍主要為0°～40°,僅當方位角在90°/270°(目標長軸垂直于地磁場方向)左右,俯仰角在±30°左右時,γxy的變化范圍為40°～50.045 7°。地磁場在平面內的投影方向與磁異常變化最大路徑方向基本一致。

圖8 γxy隨目標姿態的變化特性Fig.8 Variation characteristics of γxy with target attitudes

如圖7和圖8所示,在目標姿態、磁矩信息未知的條件下,飛行平面內地磁場投影與磁異常變化最大路徑方向偏差較小。此外,地磁場在平面內的投影方向可事先通過測量地磁場獲取,因此沿著地磁場在平面內投影方向的探測路徑是較優的選擇。

目標上方飛行平面內沿地磁場投影的磁異常信號變化量對于飛行平面內最大變化量的占比REM可以表示為

(20)

利用磁矩計算模型與式(13)～式(15)可以計算出目標正上方600 m處載體飛行平面內各點的磁感應強度大小。將計算結果代入式(20)可得REM隨目標姿態的變化特性如圖9所示,當目標處于任意姿態時,REM的變化范圍為78%～100%;當γxy減小時,沿地磁場投影方向與探測平面內磁異常變化最大的方向偏離程度降低,這兩個方向的磁異常變化幅度相差就越小,由式(20)可得REM增大,反之亦然,這使得圖8與圖9相似。對比圖8可知,當γxy處于最大值附近區域時REM取到最小值,當γxy處于最小值區域時REM取到最大值;γxy越大REM越小,γxy越小REM越大,γxy小于35°時,REM基本都在90%以上。無論目標處于何種姿態,在目標正上方沿地磁場飛行獲取的磁異常變化量總是飛行平面最大變化量的78%以上。因此,在目標姿態和磁矩信息未知的情況下,沿地磁場在飛行平面內的投影方向飛行是最優的選擇。

圖9 REM隨目標姿態的變化特性Fig.9 Variation characteristics of REM with target attitudes

利用多物理場仿真軟件獲取不同方向測線的磁異常信號,對上一段的結果進行驗證,為減小計算量,只針對目標特定的姿態進行仿真,將目標姿態角設置為:α=80°、β=26°,此時σxy為14.29°、γxy為50.1°。其他參數設置與第2.1節中保持一致。將測線設置在水下目標正上方600 m處,只改變測線與地磁場投影的夾角來獲得磁異常信號曲線。沿著磁矩在平面內投影、地磁場在平面內投影設置兩條測線;其他測線分別與沿著地磁場投影方向測線的夾角為30°、60°、90°、120°、150°。仿真結果如圖10所示,圖中磁異常信號變化最大的測線是與地磁場投影夾角為30°的測線,該條測線與磁異常變化最大路徑方向的夾角最小;沿著地磁場投影測線的磁異常信號變化幅度為該平面內最大變化的88.50%,理論計算結果為84%,差異可能是仿真軟件在計算磁場分布時,空間中點與點之間的間隔設置過大使某些點的磁場數據丟失所引起的。

圖10 磁異常信號仿真結果Fig.10 Simulation of magnetic anomaly signal

3 結 論

為了優化航磁探測中的飛行探測路徑以達到提高磁異常信號強度的目的,本文綜合考慮目標的固定磁性與感應磁性,根據其結構、姿態等參數以及地磁場值建立了目標磁矩的計算模型,利用該模型可精確地計算出典型目標在地磁場中產生的磁矩;利用目標磁矩的計算模型與磁異常計算公式建立了探測路徑上磁異常信號與目標姿態的定量關系。推導公式數值計算結果與商用多物理場仿真軟件仿真結果吻合,結果表明:當目標處于任意姿態時,其正上分沿著地磁場在平面內投影的磁異常信號變化量約為平面內最大變化量的78%～100%,這也就意味著在目標姿態和磁矩信息未知的情況下,沿地磁場在飛行平面內的投影方向飛行是最優的選擇。后續將進一步開展縮比目標探測的驗證實驗。