基于噪聲擾動下重力二階垂向導數與歸一化磁源強度的重磁相關性分析
——以相山鈾礦田為例

王彥國，田野，鄧居智，葛坤朋

（1. 東華理工大學核資源與環境國家重點實驗室，江西南昌，330013；2. 東華理工大學地球物理與測控技術學院，江西南昌，330013）

0 引言

基于巖（礦）石的密度、磁化率、電阻率及波阻抗等物性特點，地球物理方法可分為重、磁、電、震等不同勘探方法。由于觀測數據體在空間上的局限性，單憑某一種地球物理方法很難準確推測復雜的地下地質結構［1］。但不同類型的地球物理數據間具有互補性，可聯合多種地球物理數據進行綜合解釋，降低多解性，提高解釋的可靠性［2-4］。重力場和磁場均是地下所有場源產生異常的總和，具有嚴重的體積效應，多解性明顯。因而，利用重、磁異常進行綜合分析與解釋，有助于提高地質解釋的準確性。

重磁對應分析法是重磁聯合定性解釋的一種常用方法，該方法基于重磁泊松定理，通過研究重力垂向一階導數和化極磁異常的相關性開展巖性識別與構造劃分［5］。Garland［6］首先提出了重磁對應分析方法，應用于計算單一地質體的磁化強度與剩余密度的比值（泊松比），從而確定場源性質；Cordell 等［7］結合重力異常和磁源重力異常振幅譜，估計北極地區的柯尼希斯貝格比（即剩余磁化強度與感應磁化強度之比）及總磁化方向；Chandler等［8］利用滑動窗口實現了多源重磁場的相關性研究；Chandler 等［9］應用滑動窗口下的重磁對應分析法對美國明尼蘇達中東部至佩尼奧克統造山帶前寒武系開展地質填圖，取得了較好的應用效果。

劉心鑄［10］將重磁對應分析法引入中國，并進行了單一和疊加異常體的模型試驗及應用分析；黎益仕等［11］以重磁對應分析為約束，開展了重磁聯合反演，獲得了良好的應用效果；范正國等［12］分析了重磁對應分析法的局限性，提出了自適應重磁對應分析法，在一定程度上提高了方法的可靠性與實用性；宋景明［13］利用磁源重力異常剔除磁性體引起的重力異常，得到了反映沉積構造特征的剩余重力異常；在此基礎上，顏廷杰等［14］利用多尺度窗口滑動線性回歸修正方法將磁場換算為磁源重力異常，得到了純粹由無磁性地質體產生的剩余重力異常；王彥國等［15］基于重磁泊松定理推導了總磁化方向的估計方法，并在相山鈾礦田獲得了良好的應用效果。一些學者還在不同地區開展了重磁對應分析法的應用，取得了很好的應用效果［16-19］。然而，常規重磁對應分析法需要進行化極處理，且僅適用于重磁異常高度相關的情況。此外，該方法易產生虛假信息，嚴重限制了該方法的實際應用［12，18］。鑒于此，本文從重力二階垂向導數與歸一化磁源強度的表達式出發，推導出新形式下的重磁泊松公式。在此基礎上，提出了基于重力二階垂向導數與歸一化磁源強度互相關系數的重磁相關性分析方法，并通過噪聲擾動處理消除重磁同源體之外的虛假強相關信號。由于本方法使用了重力二階導數，很大程度上削弱了區域重力場的影響；同時，歸一化磁源強度基本不受磁化方向影響。因此，本文方法無需做化極處理，具有較強的適用性和實用性。通過模型試驗驗證了方法的可行性與優越性。最后，將該方法應用于相山鈾礦田重磁資料處理，獲得了良好的效果。

1 基本原理

典型地質體的重力異常可表示為［20-21］

式中：G是萬有引力常數；Δρ為地質體與圍巖的密度差；kg是與坐標位置無關但與地質體幾何形狀有關的參數；(x，y，z)和(x0，y0，z0)分別是計算點和場源點坐標；Ng是重力源構造指數，與地質體幾何形狀有關。

對重力異常g求取z方向的二階導數

典型地質體的磁場可用歸一化磁源強度統一表示為［22-23］

式中：μ0是真空中的磁導率；M是磁源磁化強度；kM是與坐標無關但與場源形狀相關的參數；NM是磁源構造指數，與場源形狀有關；λ1、λ2和λ3是磁梯度張量Γ的特征值。

磁梯度張量Γ［24］可表示為

式中：U表示磁位；Bx、By、Bz分別表示磁場B的x、y、z三個方向分量，其值可利用波譜關系由磁異常T換算得到。

當觀測點位于場源上方時，重力二階垂向導數和歸一化磁源強度公式可分別簡化為

從常規重磁泊松公式可知，對于重磁同源體，其磁異常與重力異常的一階導數形式一致，因此NM=Ng+1，式（5）可轉換為

上式即為基于重力二階垂向導數與歸一化磁源強度的重磁泊松公式。 對于同一場源，記k=為一常數，此參數僅與場源的幾何形狀相關。

依據式（6），可引入下面的互相關系數判斷重磁的同源性

式中W表示滑動計算窗口的總點數。上式中包含重力二階導數和磁場一階導數，可大幅度降低重磁背景場的影響。

根據柯西不等式可得

因此，相關系數C取值范圍為［-1，1］。

同樣，依據式（6）可以計算包含場源幾何形狀影響因素的重磁泊松比

直接使用重力二階垂向導數和歸一化磁源強度進行互相關計算時，效果并不理想，主要原因是遠離場源時，歸一化磁源強度太小且窗口內數值變化不明顯，此時窗口內各點的Tμ趨于窗口內的平均值，則式（7）可寫作

根據上式可知，當遠離磁源區且重力二階垂向導數處于非零值區域時，系數C接近±1，即完全正相關或負相關，顯然這是不合理的。事實上，常規重磁分析方法同樣存在這樣的問題［12］，目前對此尚未有理想的解決方案。

為解決這一問題，對重力二階垂向導數和歸一化磁源強度添加兩個不相關的噪聲干擾，然后再進行互相關計算。雖然添加隨機噪聲會降低數據信噪比，但互相關系數計算是一個平均化處理過程，可有效削弱噪聲對有效信號的干擾。在接近場源區內，場源的有用信號占比更大，互相關系數主要取決于重磁異常的相關性；在遠離場源處，噪聲占主導地位，且隨機分布，此時互相關系數趨于零。由于歸一化磁源強度具有非負性，因此互相關系數趨于1 時，意味著高密度、強磁性場源的存在；接近-1 時，表明場源具有低密度、強磁性的特點。

2 模型試驗

2.1 單一異常體模型試驗

為了測試方法的可行性與有效性，首先建立一個單一異常體模型，分別模擬重磁完全同源、部分同源及完全不同源三種情況。

模型包含一個正方體異常體，異常體邊長為2 km，上頂面埋深為1 km，剩余密度為1.0 g/cm3，磁化強度為1 A/m，磁傾角及磁偏角均為45°。圖1為該模型的重力異常與磁異常。

圖1 單一異常體模型的重力異常（a）和磁異常（b）

圖2 是該模型的理論重力異常一階垂向導數和理論垂直磁化磁異常，及含10%隨機干擾的異常，以此數據計算常規重磁對應分析法的相關系數。需要指出的是，該模型是重磁完全同源，與重磁不同源的差別是僅異常位置發生偏移，異常形態是一致的，因此文中不對其進行展示。

圖2 單一異常體模型的重力一階垂向導數（a）及垂直磁化磁異常（b）

圖3（上）是基于圖2（上）數據，采用常規重磁分析法分別對重磁異常同源、部分同源、完全不同源三種情況下計算的系數C。可以看出，重磁異常完全同源時，C在整個測區都接近于1，這不利于確定重、磁源的真實位置；當重磁異常部分同源和完全不同源時，C在遠離場源區存在大面積的正相關和負相關，在重力源與磁源質心連線上，則存在明顯的負相關。這說明無論是重磁同源還是不同源，理論重力一階垂向導數和垂直磁化磁異常的重磁對應分析法效果并不理想，均不能較準確地反映重磁場的相關性。范正國等［12］對常規重磁對應分析法開展了詳細的模擬實驗工作，得到了同樣的結論，這也是該方法未能廣泛應用的主要原因。

圖3 常規重磁對應分析法對無噪數據（上）和含噪數據（下）計算的互相關系數C（窗口大小5×5）

圖3（下）是基于圖2（下）數據，分別對重磁異常同源、部分同源、完全不同源三種情況計算的C。可以看出，由于噪聲干擾的引入，C值在遠離場源區處于零值上下波動。在重磁完全同源的場源上方接近1，但重磁部分同源及完全不同源時，重磁源質心連線上仍存在明顯的負相關。這表明，對重、磁場進行噪聲擾動有助于提高重磁對應分析方法在重磁完全同源時的可靠性，但在其他情況（重磁部分同源及完全不同源）下仍不適用，易得到虛假信息，不利于重磁綜合解釋。

圖4 是基于圖1 中的重力異常和斜磁化磁異常換算得到的重力二階垂向導數和歸一化磁源強度，以及分別添加10%隨機干擾的結果。圖5（上）是對圖4（上）數據分別計算重磁同源、部分同源及完全不同源時的無噪聲重力二階垂向導數與歸一化磁源強度的互相關系數C。可以看出，無論重磁是否同源，在重力源內外均存在強相關性，這顯然是不合理的。圖5（下）是基于圖4（下）計算的重力二階垂向導數與歸一化磁源強度的互相關系數C。可以看出：添加噪聲干擾后，C在重磁同源位置上相關性較強，而在重磁場源外則相關性較弱；當重磁完全不同源時，整個測區內的相關系數基本在零值附近，即重磁完全不相關。顯然，通過添加噪聲擾動可顯著提高方法的有效性與可靠性。相對于常規重磁對應分析法，本文方法能夠更好地反映重磁同源體的真實位置。

圖4 單一異常體模型的重力二階垂向導數（a）及歸一化磁源強度（b）

圖5 本文算法基于圖4 無噪數據（上）和圖4 含噪數據（下）計算的互相關系數C（窗口大小5×5）

2.2 疊加模型試驗

為了驗證方法對復雜情況的應用效果，構建了一個多場源組成的疊加模型，模型參數見表1。其中正方體與巖脈為重磁同源，棱柱體為部分同源，球體是無磁性的。由于常規重磁對應分析法需要將斜磁化磁異常轉換成垂直磁化磁異常，而對于磁化方向不同的多場源，難以進行有效的化極處理，因此這里不再進行常規算法的試驗。

表1 組合模型參數

圖6a、圖6b 是組合模型產生的重力異常和磁異常。可以看出，受規模較大的棱柱體異常的影響，巖脈上方的重力異常表現為等值線同行扭曲；受磁化方向影響，磁異常較復雜，難以直接識別磁源位置。

圖6 組合模型重磁異常

圖6c、圖6d 是無噪聲時重力二階垂向導數及歸一化磁源強度。可以看出，重力二階垂向導數和歸一化磁源強的極值均與自身場源位置對應很好。圖6e、圖6f是添加10%隨機干擾后的重力二階垂向導數和歸一化磁源強度。

圖7a 是利用無噪聲時重力二階垂向導數與歸一化磁源強度（圖6c 和圖6d）計算的互相關系數。可以看出，在重力場源分布區內存在明顯的相關性，在場源外也呈現強相關性，尤其值得注意的是，球體本身是無磁性的，但也反映出明顯的負相關，即表現為負密度、強磁性的特征，這顯然會給解釋帶來不利影響。圖7b 則是利用噪聲擾動后重力二階垂向導數與歸一化磁源強度（圖6e 和圖6f）計算的互相關系數。可以看出，互相關系數在正方體上方為明顯的負相關（接近-1），在巖脈和棱柱體同源部位表現為強正相關（接近1），在球體上方則顯示為完全不相關（接近0），除此之外的區域，互相關系數同樣在零值附近波動。這也再次表明在進行重磁相關性計算時，添加噪聲擾動有助于準確判斷重磁異常是否同源。

圖7 本文算法基于疊加模型重磁異常計算的互相關系數C（窗口大小5×5）

3 應用實例

為了檢驗算法對實際資料的應用效果，選取了江西相山鈾礦田地區的地面重、磁數據（比例尺1：50000）進行測試。江西相山火山盆地是中國第一大、世界第三大火山巖型鈾礦田［25］。

相山火山盆地位于華夏板塊湘桂贛地塊北東緣樂安—撫州斷隆帶上（圖8）。該區遭受了揚子—加里東期、海西—印支期造山作用，燕山期NEE 向贛杭構造火山巖帶西南端與近SN 向贛中南花崗巖帶的交接地帶發生了強烈的構造—巖漿—成礦作用［26-27］。相山鈾礦田中心區域出露大面積的下白堊統鵝湖嶺組碎斑熔巖，向外依次分布著條帶狀的下白堊統打鼓頂組流紋英安巖、下白堊統花崗斑巖及青白口系變質巖，西北側還發育上白堊統的沉積巖。相山地區鈾礦化主要受斷裂構造、火山構造、變質基底及巖石層間接觸面控制，多控礦因素交匯更有利于富礦的形成［28-31］。因此，推斷有利控礦因素的分布位置對相山地區鈾礦勘查具有重要的研究意義。

為方便重、磁異常解釋及其相關性分析，圖9給出了相山火山盆地主要巖石的密度、磁化率散點分布，其中碎斑熔巖、流紋英安巖、花崗斑巖及變質巖樣品分別為822、353、101、75 塊。從圖9 可以看出，碎斑熔巖和花崗斑巖密度均主要為2.60～2.65 g/cm3，均值分別為2.63、2.64 g/cm3，同屬于低密度體；流紋英安巖密度主要為2.65～2.75 g/cm3，均值為2.69 g/cm3，為中—高密度體；變質巖密度主要為2.7～2.8 g/cm3，均值為2.76 g/cm3，為高密度體。部分碎斑熔巖巖樣的磁化率為（10～100）×10-6SI，屬于低磁性，但大部分巖樣的磁化率分布在（100～1000）×10-6SI，相當于中強磁性；流紋英安巖巖樣磁化率則多數為（10～100）×10-6SI，部分為（100～1000）×10-6SI，還有少量樣本磁化率大于1000×10-6SI，即流紋英安巖磁性變化較大，弱、中、強磁性樣品均有；花崗斑巖約一半樣品磁化率為（10～100）×10-6SI，另一半為（100～1000）×10-6SI，說明花崗巖表現為弱、中磁性；變質巖的磁化率分布較集中，多數為（20～30）×10-6SI，表現為弱磁性。上述巖石密度、磁化率統計分析結果表明，碎斑熔巖、中—強度磁性花崗斑巖及中—強磁性流紋英安巖這三類巖石都可看作是重、磁同源體。

圖9 相山地區主要巖石密度（上）和磁化率（下）散點分布圖

圖10 是相山鈾礦田的布格重力異常和總磁力異常圖。布格重力異常整體表現為中心區域異常低、外圍異常高，這樣的分布特征主要與大面積的低密度碎斑熔巖和高密度變質巖有關。磁異常主要以局部異常為主，正磁異常主要分布在研究區的南側與中東側，以SN、NE、NEE 向為主，這些正磁異常在近正北側都可以尋找到相對應的負磁異常分布。然而，在鄒家山—鄒家一帶，存在大面積的負磁異常，且這些負磁異常幅值明顯大于周圍的正異常幅值，表明該區域的巖石有可能具有與當代地磁場方向不一致的剩磁［15］。

圖10 相山地區布格重力異常（左）及地面磁異常（右）

圖11（上）是研究區觀測面上的重力二階垂向導數和歸一化磁源強度。可以看出，地面上的重力二階垂向導數分布凌亂，以小規模條帶狀異常為主，各個方向異常相互切割嚴重；地面歸一化磁源強度主要在中西部、中南部及西北部存在明顯的異常，多數鈾礦點位于異常高—低分界位置。

圖11（下）是向上延拓500 m 的重力二階垂向導數和歸一化磁源強度。可見此時重力導數及歸一化磁源強度的局部異常規模較大，走向明顯，均以NS、NE、NEE 向為主。鈾礦點主要位于重力導數梯度帶（或零值線）附近，及歸一化磁源強度異常高—低分界線位置這些位置往往與斷裂帶或巖性接觸帶關系密切。

圖12 是基于添加了10%噪聲擾動下相山地區重力二階垂向導數與歸一化磁源強度計算的互相關系數。可以看出，地面上互相關系數C（圖12a）僅在牛頭嶺北側存在一個呈近東西向展布的條帶狀高值區，該區對應大面積出露的流紋英安巖，其他位置基本都在-0.5～0.5，即大部分地區重磁相關性并不明顯，這是由于對重力場進行了二階導數計算，因而異常較為凌亂。向上延拓500 m 后的互相關系數C（圖12b）的分布反映出多個重磁強相關區域：①牛頭嶺北側、淮頭南側這兩個正強相關區域反映了高密度、強磁性的流紋英安巖具有一定規模；②云際—堯崗西沿線、芙蓉山西側、居隆庵西側存在強正相關區，推測為碎斑熔巖之下存在一定規模的流紋英安巖，其中居隆庵地區的大量鉆孔資料已證實了這一推測［32-33］，這為相山地區在層間接觸面尋找鈾礦提供了參考，尤其相山鈾礦田西部的云際—堯崗西沿線一帶；③牛頭嶺、居隆庵、相山、堯崗北等地區存在大面積負強相關區域，應是低密度、強磁性碎斑熔巖的反映，何家—戴家村一帶存在近東西向的條帶狀負相關區域，應與低密度的花崗斑巖有關，不過何家附近的負相關性更加顯著，區域更大，表明這里的花崗斑巖應具有較大規模，而位于相山鈾礦田北部的沙洲礦床產于花崗巖內外接觸帶［34-35］，因此重磁相關性分析結果為沙洲礦床開展深部鈾礦勘查也提供了依據；④云際東、鳳崗南及芙蓉山等區域的強負相關性區域，可能反映了低密度、強磁性的花崗斑巖分布；⑤幾乎所有的鈾礦點都位于強正相關與強負相關的過渡帶，反映了相山地區鈾礦點與巖性接觸面關系密切；⑥牛頭嶺、芙蓉山、相山、堯崗北分布的古火山口處于大面積的強負相關區域，表明這些火山口噴發的主要物質是碎斑熔巖［36］，而牛頭嶺北、淮頭南兩地的火山口則以流紋英安巖噴出為主，另外居隆庵東側大面積的負強相關則表明該區可能還存在一個碎斑熔巖的古火山口。顯然，這些豐富的重磁相關性信息為相山鈾礦田地質研究及鈾礦勘查提供了一定參考依據。

圖12 基于10%噪聲擾動下相山地區重力二階垂向導數與歸一化磁源強度的互相關系數C 分布（窗口大小3×3）

4 結論

基于重力二階垂向導數與歸一化磁源強度表達式關系，本文提出了一種新的重磁相關分析方法，削弱了背景場與磁化方向的影響。為提高相關分析結果的可靠性，本文對重力二階垂向導數和歸一化磁源強度進行較強的隨機噪聲擾動。模型試驗表明，此方法可以有效地識別重磁同源體的位置，且不存在虛假信息；同時，該方法適用于重磁場源位置不同、磁源磁化方向不一致的復雜環境。將本文方法應用于相山鈾礦田，重磁相關性分析結果較好地揭示了鵝湖嶺組碎斑熔巖、打鼓頂組流紋英安巖及花崗斑巖等重磁同源巖石的分布位置及規模，也反映了不同位置古火山噴出的主要巖石類型，同時揭示了鈾礦與巖性接觸面的密切關系。這些結論為研究相山地區深部鈾礦成礦環境提供了數據基礎。