利用全球地磁場模型研究ΔT與Tap之間的差異性特征

陳 康 胡正旺 杜勁松,3

1 廣西壯族自治區第七地質隊，廣西壯族自治區柳州市柳堡路3號，545100 2 中國地質大學(武漢)地球物理與空間信息學院地球內部多尺度成像湖北省重點實驗室，武漢市魯磨路388號，430074 3 中國地質大學(武漢)地質過程與礦產資源國家重點實驗室，武漢市魯磨路388號，430074

將標量磁力儀測量的總磁場強度T的大小減去主磁場T0的大小，即可得到總磁場強度異常(ΔT)。在實際數據處理與正反演及定量解釋時，往往將ΔT近似為磁力異常矢量Ta在主磁場方向的投影分量Tap[1-3]，此時ΔT就和其他分量異常一樣，與場源的磁化強度大小呈簡單的線性關系，而與主磁場T0無直接關系，并且在場源外部空間滿足拉普拉斯方程，從而具備調和性質，進而可以簡化其數據處理、正反演與定量解釋過程。但當Ta幅值較大時，這種近似將會產生較大誤差。針對此誤差，部分學者[4-7]分別提出利用實測ΔT數據計算投影分量Tap的最優化方法、等效源方法、迭代方法，并且其理論模型的實驗和實際應用均表明，在高磁環境下考慮將ΔT近似為Tap的誤差具有必要性。

但上述關于ΔT與Tap之間差異性問題的相關研究均是針對局部磁力異常在笛卡爾直角坐標系中完成，而對于區域乃至全球性中、長波長ΔT與Tap磁力異常之間的差異程度、分布規律及影響因素研究尚未見到相關報道。因此，本文基于全球地磁場球諧模型，從不同波長范圍與不同觀測高度探討將ΔT近似為Tap的誤差的全球分布特征，旨在為全球巖石圈磁場模型的應用以及基于地面、航空與衛星ΔT數據構建區域乃至全球巖石圈磁力異常場模型提供參考資料。

1 ΔT與Tap差異性的理論分析

如圖1所示，總磁場強度異常(ΔT)為總磁場強度(T)的大小與主磁場強度(T0)的大小之差，Tap為磁力異常矢量(Ta)在主磁場方向的投影分量，即

圖1 ΔT、Ta與Tap之間的關系示意圖

ΔT=|T|-|T0|

(1)

Tap=|Ta|cosθ

(2)

式中，θ為T0與Ta之間的夾角，其計算公式為：

(3)

根據文獻[8]理論進行推導，ΔT與Tap之間的差異E為:

(4)

由于ΔT≤|Ta|，因此誤差E≥0。分析式(4)可知，當ΔT為0時，誤差E達到最大值，即

Emax=|Ta|2/(2|T0|)

(5)

此時，|T|與|T0|等長且與|Ta|形成等腰三角形，T0與Ta之間的夾角θ′為：

(6)

袁曉雨等[8]認為，誤差E與測點和場源之間的相對位置密切相關，誤差比例最大的測點位于磁力異常矢量Ta接近垂直于主磁場T0的區域，隨著Ta值增大，誤差按照近似指數關系迅速增大，且誤差最大值的極值主要分布在中緯度地區，對應的磁力異常矢量Ta與主磁場T0的夾角主要分布在90°～120°范圍內。但上述關于誤差E的幅值大小與空間分布特征以及影響因素的部分結論是基于局部平面直角坐標系的理論模型所獲得，可能與實際情況存在偏離。另外，對于不同的波長范圍與觀測高度，誤差E的幅值大小及其變化規律是構建與應用全球巖石圈磁場模型的重要因素。因此，需要從實際情況出發，采用全球巖石圈磁場模型研究誤差E的幅值大小與空間分布特征以及影響因素。

2 基于全球巖石圈磁場模型的ΔT與Tap及其差異計算方法

為分析實際情況下ΔT與Tap之間的差異性特征，本節給出基于全球巖石圈磁場模型計算相關物理量的方法。

主磁場與巖石圈磁場分別起因于地球外核磁流體運動與巖石圈磁性結構，兩者共同組成地磁場的內源場，在地球表面及其以上空間，可以將內源場的磁位表示為球諧級數展開形式。根據展開式可以計算球面局部指北直角參考系(x軸指北、y軸指東、z軸指向地心)中的磁場三分量(X、Y與Z)，其無奇異性球諧計算表達式參見文獻[9]。為避免地球扁率的影響，所有解算點位坐標與參考坐標系分別為基于WGS84參考橢球的大地坐標(大地緯度B、大地經度L、大地高度H)與橢球面局部指北坐標系，坐標與矢量旋轉方法參見文獻[10]。

根據地磁場的球諧能量譜特征[11]，一般將地磁場球諧展開式中1～15階作為主磁場，而將16階及其以上部分作為巖石圈磁場。因此，將相應階數范圍的球諧系數代入球諧表達式中，可計算出主磁場三分量(X0、Y0與Z0)、巖石圈磁場三分量(ΔX、ΔY與ΔZ)與總磁場三分量(X、Y與Z)，進而可得到|T|、|T0|、ΔT、|Ta|、Tap：

(7a)

(7b)

(7c)

(7d)

(7e)

將式(7)代入式(4)與式(5)即可分別得到將ΔT近似為Tap的誤差E與最大誤差Emax，將式(7)代入式(3)即可得到T0與Ta之間的夾角θ為：

(8)

3 結果與分析

采用EMM2017地磁場模型計算與分析實際情況下ΔT與Tap之間的差異性特征。該模型由最新的EMAG2v3[12]全球總磁場強度異常2′×2′網格數據通過球諧展開獲得，最高展開階數達790階次，同時包含2000～2022年期間1～15階主磁場及其長期變化。因此，采用該模型1～15階、16～790階分別作為主磁場與巖石圈磁場。

3.1 ΔT與Tap之間差異的大小與分布特征

首先采用EMM2017(16～790階)計算得到參考橢球面上的總磁場強度異常(ΔT)、ΔT與Tap之差(圖2))。從圖中可以看出，ΔT與Tap之差均為非負值，說明ΔT≥Tap，這與理論分析結論一致。此外，該差異在全球大部分地區均低于2 nT，僅在部分區域存在較大值。

圖2 ΔT及其與Tap之間的差異分布

由表1可以看出，對于參考橢球面與5 km高度面上16～790階的巖石圈磁場，該差異的最大值分別可達248.5 nT與120.2 nT。由于目前地面與航空磁測的精度均遠高于該差異，因此在圖2(b)所示差異性較大地區的實際地面與航空ΔT數據處理、正演與反演時不能忽略該差異。此外，隨著全球巖石圈磁場模型展開階數的增高與計算高度的降低，ΔT與Tap之差的差異性會更加顯著，這說明在構建或應用高階全球巖石圈磁場模型時需要考慮該差異，但對于衛星高度的ΔT磁力異常數據則無需考慮。

表1 ΔT、Tap及兩者之間的差異

3.2 ΔT與Tap之間差異的影響因素

由式(4)可知，誤差E與|T0|、|Ta|和ΔT均存在關系。圖3(a)與3(b)分別為采用EMM2017解算參考橢球面上的磁力異常模量|Ta|(16～790階)與主磁場強度模量之倒數(1/|T0|)(1～15階)。對比圖2(b)與圖3可以明顯看出，誤差E主要受|Ta|大小影響，而受主磁場強度|T0|影響很弱，并且誤差E未表現出隨緯度而變化的分布特征。

圖4(a)為誤差E(≥2.0 nT)與|Ta|2之間的統計關系圖，從圖中可以看出，在整體趨勢上誤差E與|Ta|2近似滿足線性關系(藍色虛線)。圖4(b)為E(≥2.0 nT)與|Ta|2/2|T0|之間的統計關系圖，從圖中可以看出，采用式(5)計算的誤差E確實為誤差最大值Emax(紅色虛線)。對比圖4(a)與4(b)也可以看出，誤差E受主磁場強度T0影響很弱，而圖中偏離直線的點受ΔT大小影響。進一步基于圖2和圖3中ΔT、|Ta|與誤差E之間空間分布的相關性可知，誤差E受|Ta|影響最大、ΔT次之、|T0|最弱。

圖3 |Ta|與1/|T0|分布

圖4 誤差E與|Ta|2、|Ta|2/2|T0|的統計關系

3.3 ΔT與Tap之間的差異估計與壓制方法

根據式(5)，采用EMM2017地磁場模型(16～790階)計算得到Emax分布，結果如圖5(a)所示。對比圖5(a)與圖2(b)可以看出，兩者在空間分布上具有較好的一致性，兩者差異見圖5(b)。從圖中可以看出，Emax幅度大于實際誤差，因此采用式(5)計算E會過高估計誤差，但計算得到的結果一方面可以作為誤差估計的上限，另一方面可以用于判斷誤差較大的區域。

圖5 Emax及其與實際誤差E之間的差異分布

從式(4)也可以看出，在將Tap近似為實際ΔT數據進行磁力異常三分量與模量轉換處理時，為減小誤差E，可以引入Ta趨近于ΔT的約束條件。通過表1也可以看出，在統計意義上，|Ta|與ΔT的幅值比較接近，但這僅對磁力異常矢量的大小進行約束，并未約束其矢量方向，因此后文將分析誤差E與磁力異常矢量Ta方向之間的關系。

根據式(8)，采用EMM2017地磁場模型計算得到θ分布，結果如圖6(a)所示。對比圖6(a)與圖2(b)可以看出，誤差E與夾角θ無相關性。根據式(6)，采用EMM2017地磁場模型計算得到θ′分布，結果如圖6(b)所示，從圖中可以看出，誤差E與夾角θ′具有明顯的對應關系，即當E較小時，θ′接近90°；而當E較大時，θ′逐漸增大(本次計算中夾角θ′未超過95°)。由式(6)可知，θ′受控于磁力異常模量|Ta|與主磁場強度模量|T0|之比值，這說明在對總磁場強度異常數據ΔT進行轉換處理或磁場建模時，為減小將Tap近似為ΔT的誤差，可以對磁力異常矢量的方向進行約束，即要求磁力異常矢量與主磁場垂直的分量的幅度較小。這點從圖1也可以看出，即欲使誤差E最小，則要求磁力異常矢量Ta與主磁場T0方向平行。

圖6 T0與Ta之間的夾角θ分布以及當ΔT為0時T0與Ta之間的夾角θ′分布

Lesur等[13]在利用全球總磁場強度異常匯編ΔT數據反演巖石圈磁場球諧模型系數時采用了上述磁力異常矢量的方向約束，并且認為引入該約束是為壓制Backus效應[14]。但根據上述分析，本文認為該約束應該是為減小將Tap近似為ΔT的誤差，而非壓制Backus效應。實際上，所謂Backus效應，即由標量磁場值反演磁場矢量具有多解性；而Tap方向即主磁場方向，為磁力異常矢量在主磁場方向的投影。因此，Backus效應和Tap與ΔT的差異性為兩個不同的問題。現今根據ΔT可以精確計算Tap[4-7]，因此在利用總磁場強度異常ΔT數據反演構建巖石圈磁場模型時，無需引入磁力異常矢量的大小或方向約束。

4 結 語

本文采用全球地磁場模型EMM2017分別計算得到不同波長范圍和不同高度的ΔT與Tap及兩者差異E的全球分布，進而分析該差異的幅值大小與空間分布特征以及影響差異性的因素，得到以下結論：

1)對不同波長范圍與不同高度情況下誤差E的幅值大小與空間分布特征進行分析，可為構建與應用巖石圈磁場模型時是否需要考慮ΔT與Tap之間的差異性提供判斷依據。

2)在實際情況下，誤差E受|Ta|影響最大、ΔT次之、|T0|最弱，并且誤差E未表現出隨緯度而變化的分布特征，與T0和Ta之間的夾角也無相關性。

3)為壓制ΔT與Tap之間的差異性，可以引入|Ta|趨近于ΔT的約束條件，也可以將磁力異常矢量與主磁場垂直分量的幅度最小化作為約束，但目前研究表明，無論是等效源方法、最優化方法還是迭代方法均能夠基于ΔT精確計算Tap，因此在利用總磁場強度異常ΔT數據反演構建巖石圈磁場模型時，無需對磁力異常矢量的大小或方向進行約束。

致謝：EMM2017全球地磁場模型來源于https:∥www.ngdc.noaa.gov/geomag/EMM/，在此表示衷心感謝！