磁偶極子梯度張量的幾何不變量及其應用

尹剛， 張英堂， 李志寧， 張光， 范紅波

軍械工程學院七系， 石家莊　050003

?

磁偶極子梯度張量的幾何不變量及其應用

尹剛， 張英堂*， 李志寧， 張光， 范紅波

軍械工程學院七系， 石家莊050003

摘要磁梯度張量系統姿態(tài)的變化將影響梯度場測量和數據解釋的精度，使得具有坐標變換不變性特點的張量不變量成為磁梯度張量數據解釋的研究熱點.本文在對磁偶極子產生的磁梯度張量進行特征值分析的基礎上得到了：測量點與磁偶極子位置形成的位置矢量、磁偶極子磁矩矢量與絕對值最小的特征值對應的特征向量垂直；位置矢量和磁矩矢量與最大及最小特征值對應的特征向量共面，且兩矢量間的夾角可由磁梯度張量矩陣的特征值表示.最后，將本文所得磁偶極子梯度張量的幾何不變量用于磁性目標的跟蹤中，取得了較好的實時跟蹤效果.

關鍵詞磁偶極子； 磁梯度張量； 幾何不變量； 特征值； 張量不變量

1引言

磁梯度張量測量及解釋應用被視為磁法勘探工作的一次突破(張昌達，2006； 2007)，其在資源勘探、軍事、環(huán)境等領域都有著重要的應用價值(Schmidt and Clark，2000； 2006).磁梯度張量系統測量的對象是磁場矢量的梯度場，測量的結果受磁化方向影響小，能夠反映目標體的矢量磁矩信息，且能更好地反演場源參數(方位、磁矩等)并對場源進行定位和追蹤(張鳳旭等，2007)，提高了磁源體的分辨率(馬國慶等，2012).

磁梯度張量具有良好的數學性質，但其搭載在運動載體上使用時，載體姿態(tài)的變化將影響其測量和數據解釋的精度(張光等，2013a)，由張量元素計算得到的張量不變量因其不受姿態(tài)變化的影響而得到了較多的關注，被廣泛應用于磁性目標的定位及磁異常源的邊界識別等多種場合(張光等，2013b).當探測距離大于2.5倍磁性物體長度時，可將磁性目標簡化為磁偶極子模型(張朝陽等，2010)，因此，本文在張量不變量的基礎上，以磁偶極子產生的磁梯度張量為對象，研究磁梯度張量矩陣、磁偶極子位置及磁矩間的本質關系，推導出了蘊含在其中的、不受姿態(tài)變化影響的固定幾何關系，本文稱之為磁梯度張量的幾何不變量.

2磁梯度張量測量理論及張量不變量

磁梯度張量是磁場矢量的3個分量在相互正交的3個方向的空間變化率.若B為磁場矢量場，則磁梯度張量可以表示為兩個三元素矢量的矩陣乘法運算：

(1)

其中,U為磁標勢，Bx，By和Bz是磁場矢量B在空間相互正交的三個方向上的分量.

磁法勘探中地磁場及鐵磁物質產生的異常場可看作無源的靜磁場，因此，磁感應強度的散度和旋度為0，即：

(2)

由式(1)和(2)可知，磁梯度張量矩陣G為對稱矩陣，其中9個元素中，只有5個是獨立的，即只需得到其中的5個元素就可測得磁梯度張量.在磁法探測中，磁性目標簡化為磁偶極子模型時可由6個未知量描述，即磁性目標的三維位置和三維磁矩.因此，磁矩矢量為m=(mx,my,mz)的磁偶極子，在距離其r=(x,y,z)處產生的磁矢量場(卞光浪等，2011)及磁梯度張量矩陣中的5個獨立分量(Rimetal.，2012)可分別由下式計算得到

(3)

(4)

易知，磁梯度張量矩陣的特征值與坐標的選擇無關，為張量的基本不變量(余天慶和毛為民，2006)，且張量矩陣的其余不變量可由特征值表示.設矩陣G的3個特征值為λ1，λ2和λ3，則由G的特征方程可得磁梯度張量的另外4個不變量，分別表示為

(5)

3幾何不變量的證明

3.1磁梯度張量矩陣的特征值分析

(6)

將式(4)及式(5)代入式(6)可得磁偶極子產生的磁梯度張量矩陣的特征值：

(7)

令：δ=ec2+bcd+abc+d2e+ade+e3，ε=a2c+adc+abe+c3+ce2+bde，φ=bc2+cde-be2-ace，則磁偶極子產生的梯度張量的三個特征值對應的特征向量V1，V2，V3分別為

(8)

(9)

(10)

3.2幾何不變量1

由式(7)中磁偶極子張量矩陣的特征值公式可知，假設3個特征值已由磁梯度張量系統測得的實際信號求解得到，則磁偶極子位置矢量與磁矩矢量間的夾角可由下式求解得到

(11)

由式(11)可知，磁梯度張量的特征值為張量不變量，因此，在磁梯度張量系統發(fā)生姿態(tài)變化時，由該公式計算得到的磁偶極子位置矢量與磁矩矢量間的夾角保持不變.本文將其定義為磁梯度張量的幾何不變量1，其幾何示意如圖1所示.

圖1　磁偶極子梯度張量場的幾何不變量示意圖Fig.1　Sketch map of magnetic gradient tensor′s geometric invariants

由磁偶極子產生的磁矢量、磁梯度張量場公式及幾何不變量1可得：

(12)

因此，基于本文推導的幾何不變量1，若已知磁偶極子產生的總磁場強度及磁梯度張量信息，可通過特征值分析及公式(12)求得磁偶極子相對于測量點的距離R和磁矩模M.

3.3幾何不變量2

將式(4)代入式(10)化簡可得如下等式：

(13)

由式(13)可知，磁梯度張量矩陣G的絕對值最小的特征值對應的特征向量垂直于磁矩矢量和位置矢量，且不隨著張量系統姿態(tài)的變化而變化.本文將此固定的垂直關系定義為磁梯度張量的幾何不變量2.3.4幾何不變量3

由對稱矩陣的特征向量之間的幾何關系可知，特征向量V3垂直于V1和V2；由幾何不變量2可知，特征向量V3垂直于磁矩矢量m和位置矢量r，因此，假設向量V1，V2，m和r共面，則4個向量之間可互相表示，令

(14)

(15)

又磁偶極子產生的梯度張量場與矢量場滿足歐拉反褶積公式(Reidetal.，1990；Naraetal.，2006；NaraandIto，2014)，且構造指數為3，即

(16)

(17)

則由式(14)及式(17)可得：

(18)

故磁偶極子位置矢量的單位向量可表示為

(19)

將式(15)及式(19)代入式(11)表示的幾何不變量1，可求得式(15)中的等式系數：

(20)

由式(18)及式(20)可知，位置矢量及磁矩矢量與最大及最小特征值對應的特征向量共面，且不隨著系統的姿態(tài)變化而變化，本文將其定義為幾何不變量3.因此，幾何不變量3可用于磁性目標的定位和磁矩反演.

在已知磁偶極子產生的磁梯度張量矩陣的情況下，可通過特征值分析，利用式(12)求解磁偶極子相對于測量點的距離和磁矩模，然后利用式(18)及式(20)求解磁偶極子的位置及磁矩矢量的單位向量，最終實現磁偶極子位置及磁矩矢量的正確反演，其幾何示意如圖2所示.

圖2　幾何不變量3表示的磁矩及位置矢量示意圖Fig.2　Inverse solution of source-dipole vector and magnetic moment based on geometric invariant 3

由圖2可知，該方法得到真實位置及磁矩矢量的同時，也得到了3個虛假解，因此，在實際應用中，可根據所測目標在磁梯度張量系統的上方或下方，直接去除兩個虛假解，然后根據磁梯度張量系統測得的多點總磁場強度及磁矢量場進行另一個虛假解的去除，進而實現真實磁偶極子位置及磁矩矢量的正確反演.

4仿真實驗分析

利用本文提出的幾何不變量進行磁偶極子位置及磁矩矢量反演的三組仿真實驗，以驗證所提幾何不變量的準確性.其中磁梯度張量系統位于坐標原點，磁偶極子運動軌跡的具體參數設置如表1所示.仿真實驗流程如下：

(1)根據磁偶極子數學模型，計算磁偶極子位于不同位置時在坐標原點產生的磁梯度張量值.

(2)對計算得到的磁梯度張量數據進行特征值分析，求解其對應的特征值和特征向量.

表1　仿真參數設置

(3)利用公式(12)求解磁偶極子相對于測量點的距離及磁矩模.

(4)將計算得到的特征值及特征向量代入式(19)及式(20)，求解磁偶極子位置及磁矩矢量.

利用本文所提幾何不變量反演得到的磁偶極子在x-y平面上的運動軌跡如圖3—5所示.由圖可知，幾何不變量估計得到磁偶極子真實運動軌跡的同時也求解得到3個虛假解，必須采用有效方法準確去除虛假解.虛假解的去除流程如圖6所示，去除虛假解后的真實軌跡在圖3—5中以深色表示.

圖3　估計得到的磁偶極子運動軌跡(仿真實驗1)Fig.3　Estimated track of magnetic dipole coming from the first simulation

由于在實際應用中所探測目標在測量平面的上方或下方為已知信息，且單點測量的數據反演得到的四組磁偶極子的三維坐標關于原點對稱，因此，根據已知信息中的Z坐標信息可直接刪除四組解中的兩個虛假解.將剩下的兩組解代入式(3)和式(4)中可計算得到兩組不同的解在測量點產生的磁梯度張量值.若為真實解，則該計算值與實際測量值相吻合；若為虛假解，則計算得到的張量數據與測量數據存在較大的偏差，因此，根據計算值與測量值之間的差異可進一步去除第3個虛假解，進而得到真實的磁偶極子位置及磁矩矢量.值得注意的是，虛假解的去除是實時的，即每個測量點反演得到的四組解都需要進行虛假解的去除.

圖4　估計得到的磁偶極子運動軌跡(仿真實驗2)Fig.4　Estimated track of magnetic dipole coming from the second simulation

圖5　估計得到的磁偶極子運動軌跡(仿真實驗3)Fig.5　Estimated track of magnetic dipole coming from the third simulation

圖6　虛假解去除流程Fig.6　Flow chart of elimination of ghost solutions

圖7　存在測量噪聲時估計得到的磁偶極子運動軌跡Fig.7　Estimated track of magnetic dipole existing measurement noises

由于在實際應用中，磁梯度張量系統存在測量誤差將在一定程度上影響磁偶極子探測和位置反演的距離，因此，將上述理論模型仿真中加入隨機噪聲進行研究.考慮真實磁梯度張量系統多用來測量5個獨立的磁梯度張量分量及張量矩陣對稱的原則，在5個獨立的磁梯度張量分量信號中加入均值為0 nT·m-1，方差為0.5 nT·m-1的高斯白噪聲.圖7為表1中第一組仿真實驗加入噪聲后的運動軌跡反演結果.由圖可知，加入隨機噪聲后，磁偶極子的可跟蹤范圍為y軸方向上的[-20 m,20 m].在可跟蹤范圍內，本文所提算法可較為準確地跟蹤磁偶極子，而超出此跟蹤范圍后，反演得到的磁偶極子的位置坐標存在較大的偏差和不確定性，反演結果無法為磁偶極子的運動軌跡提供參考.同樣地對表1中的第二及第三組仿真實驗加入隨機噪聲，則第二組仿真實驗的可跟蹤的坐標范圍為[-20,-5]至[8,23]，第三組仿真實驗的可跟蹤的坐標范圍為[0,-25]至[-22,-14].由式(3)及式(4)可知，磁矢量信號及磁梯度張量信號隨著距離的增加衰減極為迅速，因此，距離磁偶極子目標較遠處的磁場信號較弱，有效信號可能被測量噪聲所淹沒，進而反演的得到的磁偶極子位置及磁矩矢量存在較大的偏差.即在實際應用中，磁偶極子的跟蹤距離不僅僅受到其自身磁矩大小的影響，還受到測量系統中測量噪聲的影響.

5實驗驗證

為驗證本文所提幾何不變量在磁性目標定位及跟蹤中的可行性，采用4個磁通門傳感器構建了平面十字形磁梯度張量系統進行實驗研究.實驗時某鐵磁性運動載體(在測量范圍內可近似等效為磁偶極子)以60 km·h-1的速度在道路上直線運動，磁梯度張量系統放置于道路一側的地面上，同步采集各傳感器的輸出數據，采樣頻率為2048 Hz.選取兩條行駛路線進行目標軌跡的估計，假設運動載體沿著y軸正方向前進，在x方向上距離磁梯度張量系統的距離(Closest Proximity Approach, CPA)分別為1.4 m和2 m，z軸方向豎直向下構建右手坐標系進行數據計算和分析.

由于磁梯度張量系統靜止不動，因此可在無磁性目標的情況下采集背景地磁場，實時采集得到的磁傳感器數據減去背景地磁場即為磁性目標產生的磁場數據.在此，利用總磁場強度的變化判斷是否有磁性目標的出現，然后截取采集得到的實時信號進行運動軌跡的估計.由于運動載體等效的磁偶極子位于測量系統的上方，可直接去除兩個虛假解，則反演得到的載體在x-y平面的運動軌跡如圖8、圖9所示.利用圖6所示的方法去除剩余的虛假解，則圖8、圖9所示右側深色曲線為估計得到的運動軌跡的真實解.

由圖8和圖9可知，虛假解表示的鐵磁體運動軌跡與真實運動軌跡相差較大，利用圖6所示的虛假解去除流程有效地識別出了由磁梯度張量數據計算得到的真實解.在可探測范圍內，利用幾何不變量估計得到的鐵磁性載體運動軌跡的真實解可較好地跟蹤運動載體的實際軌跡，可實現基于磁梯度張量數據的運動載體實時跟蹤.由于磁傳感器自身的噪聲以及實驗場地中地磁梯度的存在，導致估計得到的目標位置與真實值存在一定的偏差，可通過進一步提高磁傳感器的精度實現目標的遠距離跟蹤.

圖8　估計得到的載體運動軌跡(CPA=1.4 m)Fig.8　Estimated track of mobile vehicle when CPA equal to 1.4 m

圖9　估計得到的載體運動軌跡(CPA=2 m)Fig.9　Estimated track of mobile vehicle when CPA equal to 2 m

6結論

(1) 基于數學解析方法對磁偶極子產生的磁梯度張量矩陣的特征值及特征向量進行了理論推導，得出了磁偶極子位置矢量與磁矩矢量間的夾角關于特征值的數學表達式，并在此基礎上給出了磁矩模及磁偶極子與測量點間的距離的求解公式.

(2) 基于磁偶極子產生的矢量場、磁梯度張量場及其歐拉反褶積公式，驗證了磁偶極子位置矢量、磁矩矢量、最大特征值對應的特征向量及最小特征值對應的特征向量共面的假設，并給出了位置矢量及磁矩矢量關于特征值及特征向量的數學表達式.

(3) 對本文提出的幾何不變量在磁性目標跟蹤上的應用進行了仿真和實測分析，實驗結果證明了公式推導的正確性及其在目標跟蹤中的實用性.

致謝筆者向匿名評審專家及編輯部老師致以謝意！

References

Bian G L, Zhai G J, Fan L, et al. 2011. Applying two-step iterative least square approach to determine the geometry and physical parameters of magnetic sphere sources.ChineseJ.Geophys. (in Chinese), 54(5): 1375-1383, doi: 10.3969/j.issn. 0001-5733.2011. 05.027.

Ma G Q, Du X J, Li L L. 2012. Comparison of the tensor local wavenumber method with the conventional local wavenumber method for interpretation of total tensor data of potential fields.ChineseJ.Geophys. (in Chinese), 55(7): 2450-2461, doi: 10.6038/j.issn.0001-5733. 2012.07.029.

Nara T, Suzuki S, Ando S. 2006. A closed-form formula for magnetic dipole localization by measurement of its magnetic field and spatial gradients.IEEETransactionsonMagnetics, 42(10): 3291-3293.

Nara T, Ito W. 2014. Moore-Penrose generalized inverse of the gradient tensor in Euler's equation for locating a magnetic dipole.JournalofAppliedPhysics, 115: 17E504.

Reid A B, Allsop J M, Granser H, et al. 1990. Magnetic interpretation in three dimensions using Euler deconvolution.Geophysics, 55(1): 80-91.

Rim H, Park Y S, Jung H K. 2012. Interpretation of magnetic gradient tensor for automatic locating a dipole source.∥ 22ndInternational Geophysical Conference and Exhibition. Australia, 1-3. Schmidt P W, Clark D A. 2000. Advantages of measuring the magnetic gradient tensor. Preview: Australian Society of Exploration Geophysicists, 85: 26-30.Schmidt P W, Clark D A. 2006. The magnetic gradient tensor: its properties and uses in source characterization.TheLeadingEdge, 25(1): 75-78.

Yu T Q, Mao W M. 2006. Tensor Analysis and Application (in Chinese). Beijing: Tsinghua University Press.

Zhang C D. 2006. Airborne tensor magnetic gradiometry—the latest progress of airborne magnetometric technology.ChineseJournalofEngineeringGeophysics(in Chinese), 3(5): 354-361.Zhang C D. 2007. Some problems concerning the magnetic anomaly detection (MAD).ChineseJournalofEngineeringGeophysics(in Chinese), 4(6): 549-553.

Zhang C Y, Xiao C H, Gao J J, et al. 2010. Experiment research of magnetic dipole model applicability for a magnetic object.JournalofBasicScienceandEngineering(in Chinese), 18(5): 862-868.

Zhang F X, Zhang F Q, Meng L S, et al. 2007. Magnetic potential spectrum analysis and calculating method of magnetic anomaly derivatives based on discrete cosine transform.ChineseJ.Geophys. (in Chinese), 50(1): 297-304.

Zhang G, Zhang Y T, Fan H B, et al. 2013a. The effect of vehicle attitude variation on magnetic field gradient tensor localization precision for invisible object.JournalofDetection&Control(in Chinese), 35(1): 20-24.

Zhang G, Zhang Y T, Li Z N, et al. 2013b. Localizing method of magnetic field gradient tensor under carriers moving parallelly.JournalofHuazhongUniversityofScienceandTechnology(NatutalScienceEdition) (in Chinese), 41(1): 21-24.

附中文參考文獻

卞光浪, 翟國君, 范龍等. 2011. 用最小二乘兩步迭代法求解磁性球體幾何與磁性參數. 地球物理學報, 54(5): 1375-1383, doi: 10.3969/j.issn.0001-5733.2011.05.027.

馬國慶, 杜曉娟, 李麗麗. 2012. 解釋位場全張量數據的張量局部波數法及其與常規(guī)局部波數法的比較. 地球物理學報, 55(7): 2450-2461, doi: 10.6038/j.issn.0001-5733.2012.07.029.

余天慶, 毛為民. 2006. 張量分析及應用. 北京: 清華大學出版社.

張昌達. 2006. 航空磁力梯度張量測量-航空磁測技術的最新進展. 工程地球物理學報, 3(5): 354-361.

張昌達. 2007. 關于磁異常探測的若干問題. 工程地球物理學報, 4(6): 549-553.

張朝陽, 肖昌漢, 高俊吉等. 2010. 磁性物體磁偶極子模型適用性的試驗研究. 應用基礎與工程科學學報, 18(5): 862-868.

張鳳旭, 張鳳琴, 孟令順等. 2007. 基于離散余弦變換的磁位譜分析及磁異常導數計算方法. 地球物理學報, 50(1): 297-304.

張光, 張英堂, 范紅波等. 2013a. 載體姿態(tài)變化對磁張量定位精度的影響. 探測與控制學報, 35(1): 20-24.

張光, 張英堂, 李志寧等. 2013b. 載體平動條件下的磁梯度張量定位方法. 華中科技大學學報(自然科學版), 41(1): 21-24.

(本文編輯汪海英)

基金項目軍內科研資助.

作者簡介尹剛，男，1988年生，博士生，主要研究方向為磁異常探測、磁梯度張量測量及地磁導航.E-mail: gang.gang88@163.com *通訊作者張英堂，男，1960年生，教授，博士生導師，主要研究方向為測試技術與信號處理.E-mail: zyt01@malis.tsinghua.edu.cn

doi:10.6038/cjg20160232 中圖分類號P631

收稿日期2014-06-16，2015-12-23收修定稿

Research on geometric invariant of magnetic gradient tensors for a magnetic dipole source and its application

YIN Gang, ZHANG Ying-Tang*, LI Zhing-Ning, ZHANG Guang, FAN Hong-Bo

Department7th，MechanicalEngineeringCollege，Shijiazhuang050003，China

AbstractAs attitude change of the magnetic gradient tensor system may affect precision of gradient measurements and data interpretation, the rotational invariants which are invariant under rotation of the reference frame have been extensively studied for magnetic gradient tensor data. Based on the eigenvalue analysis of magnetic gradient tensors produced by the magnetic dipole, some invariant geometrical relationships, which are defined as geometric invariants, are deduced in this paper. Firstly, the eigenvector corresponding to the eigenvalue which has the smallest absolute value is perpendicular to both the dipole moment and source-dipole displacement vector. Secondly, source-dipole displacement, dipole moment and eigenvectors corresponding to the maximal and minimal eigenvalues are coplanar. Thirdly, the angle between dipole moment and source-dipole displacement vector can be calculated by eigenvalues of magnetic gradient tensors. At last, the proposed geometric invariants are used to track a magnetic target, and simulation and experimental results confirm the correctness and practicability of the proposed geometric invariants.

KeywordsMagnetic dipole; Magnetic gradient tensor; Geometric invariant; Eigenvalue; Tensor invariant

尹剛， 張英堂， 李志寧等. 2016. 磁偶極子梯度張量的幾何不變量及其應用.地球物理學報，59(2):749-756,doi:10.6038/cjg20160232.

Yin G, Zhang Y T, Li Z N, et al. 2016. Research on geometric invariant of magnetic gradient tensors for a magnetic dipole source and its application.ChineseJ.Geophys. (in Chinese),59(2):749-756,doi:10.6038/cjg20160232.