彭莉紅，張偉盟，程莎莎，孫棟華，駱燕，陳偉

1 核工業(yè)航測遙感中心，河北 石家莊 050002

2 中核集團(tuán)鈾資源地球物理勘查技術(shù)（重點(diǎn)實(shí)驗(yàn)室），河北 石家莊 050002

3 河北省航空探測與遙感技術(shù)重點(diǎn)實(shí)驗(yàn)室，河北 石家莊 050002

邊界識(shí)別方法在重力數(shù)據(jù)解釋中起著至關(guān)重要的作用，國內(nèi)外學(xué)者對它進(jìn)行了廣泛研究，提出多種邊界識(shí)別方法［1-16］。眾多邊界識(shí)別方法以導(dǎo)數(shù)分析類為主，而導(dǎo)數(shù)分析類中最基本的方法是總水平梯度（THG）［1-3］、垂向?qū)?shù)（VDR）［4-6］和解析信號(hào)（AS）［7-9］，但上述三種方法均存在異常幅度隨地質(zhì)體埋深的增加而迅速衰減的情況，因此對于埋深較深、異常幅度較小的地質(zhì)體，邊界識(shí)別能力稍顯不足。近年來，多種通過求取重力數(shù)據(jù)導(dǎo)數(shù)比值來均衡不同幅度異常的邊界識(shí)別方法相繼被提了出來。Miller 和Singh［17］提出了均衡不同幅度異常的邊界識(shí)別方法——傾斜角法（TDR），該方法對于場源體的深度不敏感，對埋深不同的多場源異常體的邊界識(shí)別效果較好；Verduzco 等［18］針對TDR 的處理結(jié)果峰值位于場源體中心，識(shí)別結(jié)果不直接的問題，提出了傾斜角的總水平導(dǎo)數(shù)法（THDR），該方法采用極大值識(shí)別場源邊界，結(jié)果相較于TDR 簡單、直接，且具有更高的分辨率，然而也存在隨場源深度的增加、異常幅度降低的情況［19］；Wijns 等［20］提出另一種均衡不同振幅的識(shí)別方法，利用解析信號(hào)與總水平導(dǎo)數(shù)的比值來進(jìn)行地質(zhì)體的邊界識(shí)別，即Theta圖法（TM），該方法相對于解析信號(hào)法的分辨能力增強(qiáng)，缺點(diǎn)是所識(shí)別出的邊界較模糊。Ferreira 等人［19］提出了總水平梯度的傾斜角法（TTHG），該方法能較好地識(shí)別出不同埋深的場源體邊界信息；Cooper［21］提出了另一種改進(jìn)的傾斜角法，稱為解析信號(hào)的傾斜角法（TAS），雖然TAS 可以更好地識(shí)別出埋深較深的場源體邊界信息，但存在與AS 法同樣的不足，即隨埋深的增加，其極大值位置會(huì)向?qū)嶋H邊界的內(nèi)部收斂。王萬銀等［22］提出了歸一化總水平導(dǎo)數(shù)垂向?qū)?shù)邊界識(shí)別方法（NVDR_THDR），該方法采用零閾值技術(shù)消除地質(zhì)體邊界，結(jié)果直觀，易于識(shí)別，是均衡不同埋深場源邊界識(shí)別效果的有效方法。馬國慶等提出增強(qiáng)型均衡濾波器，該濾波器利用不同階導(dǎo)數(shù)之間的組合識(shí)別地質(zhì)體邊界。該方法能有效識(shí)別深、淺部地質(zhì)體邊界，但高階導(dǎo)數(shù)的引入使噪聲不可避免地被放大，從而影響整體邊界識(shí)別的穩(wěn)定性［23-24］。Luan 等提出基于Logistic 函數(shù)及總水平導(dǎo)數(shù)的重力數(shù)據(jù)邊界增強(qiáng)方法（LTHG 法），認(rèn)為該方法處理結(jié)果簡單、直接且對于場源邊界具有更高的分辨能力。LTHG 方法是目前最新提出的可用于重力數(shù)據(jù)的邊界識(shí)別方法，但還未有其與主流邊界識(shí)別方法效果對比的相關(guān)研究成果，也缺少實(shí)際應(yīng)用方面的公開論述。

為此，本文基于理論模型及實(shí)測數(shù)據(jù)，對比分析了不同情況下LTHG 法與常用的THG、AS、TDR、THDR、TTHG、TAS、TM 和NVDR_THDR8 種方法對不同埋深場源體的均衡能力、分辨能力及抗噪能力，進(jìn)而驗(yàn)證LTHG 法的適用性及應(yīng)用效果。

1 LTHG 方法簡介

LTHG 法是Luan 等（2019）新提出的一種基于Logistic 函數(shù)及總水平導(dǎo)數(shù)的位場數(shù)據(jù)均衡邊界識(shí)別方法，其核心原理是基于Logistic 函數(shù)，通過對總水平導(dǎo)數(shù)的一階垂向?qū)?shù)和水平導(dǎo)數(shù)的比值進(jìn)行Logistic 函數(shù)運(yùn)算來達(dá)到突出顯示不同深度的場源體邊界的效果。

Logistic 函數(shù)是機(jī)器學(xué)習(xí)領(lǐng)域中常用的一種激活函數(shù)，定義域（-∞，+∞），值域?yàn)椋?，1），形式如下：

Logistic 函數(shù)（函數(shù)曲線形式見圖1）在定義域內(nèi)單調(diào)連續(xù)，呈現(xiàn)出先緩慢增長，然后加速增長，最后逐漸穩(wěn)定的趨勢，與常用于識(shí)別重力數(shù)據(jù)邊界的arctan 函數(shù)相似。

圖1 Logistic 函數(shù)變化曲線Fig. 1 The Logistic function curve

基于Logistic 函數(shù)，LTHG 法的計(jì)算公式如下：

式（2）中： THG—重力數(shù)據(jù)的總水平導(dǎo)數(shù)；a—常數(shù)。LTHG 法求解后得到一個(gè)比值，是無量綱，其平面峰值表征了地質(zhì)體的邊界。

圖2 為二維情況下，采用LTHG 法a分別取值為0.5、1.0、5.0、10、20、50 和100 時(shí)直立長方體密度模型的邊界識(shí)別結(jié)果，可以看出a 取不同參數(shù)時(shí)，均可以較好的識(shí)別出模型的邊界位置，且邊界位置處異常幅值基本相同，但不同取值計(jì)算的結(jié)果還是存在一定的差異，當(dāng)a＞10 時(shí)，計(jì)算結(jié)果當(dāng)中，邊界位置處異常峰型會(huì)變的更加尖銳，但當(dāng)a取值大于10 時(shí)，處理結(jié)果中會(huì)因存在畸變點(diǎn)而出現(xiàn)空值的情況，當(dāng)a＜1 時(shí)，模型邊界位置處異常峰型平緩，1≤a＜10 時(shí)LTHG 法應(yīng)用效果較好。

圖2 不同參數(shù)LTHG 法邊界識(shí)別結(jié)果Fig. 2 The LTHG filter of gravity anomaly due to the 2D block with different a.

2 邊界識(shí)別方法對比

為說明LTHG 法識(shí)別邊界的有效性，選取八種常用的重力數(shù)據(jù)的邊界識(shí)別方法進(jìn)行對比，即總水平梯度（THG）、解析信號(hào)（AS）、傾斜角法（TDR）、傾斜角的總水平導(dǎo)數(shù)法（THDR）、Theta 圖法（TM）、總水平梯度傾斜角法（TTHG）、解析信號(hào)傾斜角法（TAS）和歸一化總水平導(dǎo)數(shù)垂向?qū)?shù)（NVDR_THDR）法，各方法定義如下：

THG 法［25］：

AS 法［25］：

TDR 法［25］：

THDR 法［25］：

TM 法［25］：

TAS 法［25］：

TTHG 法［25］：

NVDR_THDR 法［22］：

計(jì)算總水平導(dǎo)數(shù)（THG）的n階垂向?qū)?shù)（VDRn），n—垂向?qū)?shù)階數(shù)，n=1，2，3…。階數(shù)越大，橫向分辨力越高，通常對于重磁異常，階數(shù)取2 時(shí)較為合適；

使用閾值大于0 來計(jì)算總水平導(dǎo)數(shù)峰值（PTHDR）：

計(jì)算總水平導(dǎo)數(shù)峰值及總水平導(dǎo)數(shù)的比值：

計(jì)算總水平導(dǎo)數(shù)垂向?qū)?shù)的最大值（VDR - THDRmax），并使用最大值進(jìn)行總水平導(dǎo)數(shù)垂向?qū)?shù)歸一化，最終得到歸一化總水平導(dǎo)數(shù)垂向?qū)?shù)（NVDR-THDR）：

建立了由5 個(gè)不同規(guī)模、不同埋深的長方體組成的組合模型，模型三維立體圖見圖3，參數(shù)見表1，正演計(jì)算了3 種情況下組合模型的理論重力異常。情況1 為所有長方體的剩余密度均為正值，長方體物體1 和2 剩余密度為0.5×10-3kg·m-3，3、4、5 體的剩余密度為0.4×10-3kg·m-3，計(jì)算的重力異常如圖4A 所示。情況2 為長方體1、4 剩余密度為負(fù)值，體1 剩余密度為-0.5×10-3kg·m-3，體4 剩余密度為-0.4×10-3kg·m-3，體2、3、5 剩余密度與第一種情況保存一致，計(jì)算的理論重力異常如圖5A 所示。情況3 為在情況2 數(shù)據(jù)中加入了2 %的高斯噪聲，結(jié)果如圖6A 所示。選用包括THG、AS、TDR、THDR、TM、TAS、TTHG、NVDR_THDR 及LTHG在內(nèi)的九種方法對以上三種情況下的理論重力異常數(shù)據(jù)進(jìn)行處理，對比分析各方法的不同振幅重力異常的均衡能力、邊界分辨能力及抗干擾能力。

表1 組合模型參數(shù)Table 1 Parameters of the synthetic model

圖3 模型三維立體圖Fig.3 A 3D view of the synthetic model

圖4 模型理論重力異常（重力正異常）及邊界識(shí)別結(jié)果Fig. 4 Synthetic gravity anomaly of the first model and the edges detected by different methods

圖5 模型理論重力異常（重力正、負(fù)異常）及邊界識(shí)別結(jié)果Fig. 5 Synthetic gravity anomaly of the second model and the edges detected by different methods

圖6 加入噪聲的模型理論重力異常（重力正、負(fù)異常）及邊界識(shí)別結(jié)果Fig. 6 Synthetic gravity anomaly of the third model and the edges detected by different methods

圖4B、C 為THG、AS 的處理結(jié)果，可以看出當(dāng)存在多個(gè)埋深不同的場源體時(shí)，THG 和AS不能均衡不同埋深的場源體邊界信息，隨著深度增加，場源體邊界變得越來越不清晰，同時(shí)雖然薄板1 和2 埋深較淺，但AS 仍無法很好的顯示出它們的邊界信息。圖4D～I(xiàn) 分別為TDR、THDR、TM、TAS、TTHG 和NVDR_THDR的處理結(jié)果，雖然TDR 均衡了不同埋深的場源體的邊界信息，但處理結(jié)果的峰值位于場源體中心，處理結(jié)果不直接。THDR 法克服了TDR法的不足，利用極大值來識(shí)別場源體的邊界，但對于埋深較深的場源體，邊界分辨能力仍稍顯不足。TM 法識(shí)別的邊界位置不夠準(zhǔn)確，且存在假的邊界信息。TAS 法同樣利用極大值來識(shí)別場源體邊界，然而對于較深的場源，其峰值會(huì)從實(shí)際的邊界位置向內(nèi)偏移，使識(shí)別到的邊界位于真實(shí)邊界內(nèi)；對于較薄場源，TAS 的峰值直接位于場源體中心。TTHG 法較好地均衡了不同埋深的場源體邊界信息，識(shí)別出的邊界信息接近真實(shí)的場源體邊界，但隨著場源埋深的增加，識(shí)別的邊界信息也存在一定程度的擴(kuò)散。NVDR_THDR 法可直接用峰值識(shí)別出邊界，效果比上述各方法更好。圖4J 為 LTHG 法的處理結(jié)果，可以看出LTHG 很好地克服了導(dǎo)數(shù)計(jì)算隨深度增加而迅速衰減的問題，均衡了不同埋深場源體的邊界信息，結(jié)果與NVDR_THDR 相似，但LTHG 方法對邊界的分辨能力要更強(qiáng)。對比分析結(jié)果表明：LTHG 不僅能更清晰、準(zhǔn)確地識(shí)別出場源邊界，而且能給出較高的邊界分辨率。

為對比分析上述各方法同時(shí)識(shí)別正負(fù)異常邊界的能力，選用包括THG、AS、TDR、THDR、TM、TAS、TTHG、NVDR_THDR 及LTHG在內(nèi)的9 種方法對第2 種情況下的理論重力異常數(shù)據(jù)進(jìn)行了處理。圖5B、C 為THG、AS 的處理結(jié)果，同樣THG 法和AS 法不能同時(shí)均衡振幅不同的異常。圖5D～I(xiàn) 分別為TDR、THDR、TM、TAS、TTHG 和NVDR_THDR 的處理結(jié)果，由圖5D 可知，當(dāng)同時(shí)存在正負(fù)異常時(shí)，TDR 法無法準(zhǔn)確識(shí)別場源的邊界，且結(jié)果中出現(xiàn)了假邊界信息。THDR 法也存在假邊界信息。TM法除不能很好地識(shí)別出埋深較深或薄的場源體邊界外假邊界信息更加顯著。TAS 法能避免假邊界信息的產(chǎn)生，但識(shí)別到的邊界位于真實(shí)邊界內(nèi)。TTHG 能避免假邊界信息的產(chǎn)生，但識(shí)別的邊界信息隨著場源埋深的增加有一定程度的擴(kuò)散。NVDR_THDR 法識(shí)別出的邊界信息輪廓清晰、位置準(zhǔn)確，無虛假邊界信息。圖5J 為 LTHG 法的處理結(jié)果，可以看出LTHG 法識(shí)別的邊界位置接近真實(shí)場源邊界，無假邊界信息，相較于NVDR_THDR，LTHG 方法對邊界的分辨能力要更強(qiáng)。

為對比分析各方法的抗噪能力，在第二種情況的數(shù)據(jù)中加入幅度為異常振幅2 %的高斯噪聲。圖6B～J 分別為THG、AS、TDR、THDR、TM、TAS、TTHG、NVDR_THDR 及LTHG 的處理結(jié)果。THG、AS 法對淺源體3 和4 的邊界識(shí)別效果相對較好，但對深源體5 的邊界識(shí)別效果卻相對模糊，同時(shí)對于薄板1 和2 ，AS 法也不能很好地發(fā)揮作用。TDR 和TM 依然無法對場源體邊界進(jìn)行準(zhǔn)確的識(shí)別，受噪聲的影響THDR、TAS 法識(shí)別的場源邊界變的模糊，分辨率降低，TTHG、NVDR_THDR 和LTHG 三種方法依然可以較準(zhǔn)確的識(shí)別出場源體邊界位置，且NVDR_THDR 和LTHG 受深度變化的影響小，未出現(xiàn)明顯的邊界向內(nèi)或外側(cè)偏移的情況，相對而言LTHG 法識(shí)別的位置更清晰，分辨率更高。

3 應(yīng)用實(shí)例

埃塞俄比亞歐加登盆地是非洲東海岸油氣勘探前景最好的盆地之一，是發(fā)育在古生代結(jié)晶基底之上的一個(gè)古內(nèi)陸裂谷—被動(dòng)大陸邊緣盆地，根據(jù)基底構(gòu)造特征，可進(jìn)一步劃分為6 個(gè)次一級構(gòu)造單元：北部隆起、北部斜坡、西部斜坡、中部坳陷、南部斜坡和南部隆起（圖7），目前在中部坳陷帶發(fā)現(xiàn)了C、H 和 D 三個(gè)氣田，展示了盆地良好的勘探潛力［26-30］。但盆地的勘探還處于早期評價(jià)階段，對斷裂的整體空間延伸及展布規(guī)律還認(rèn)識(shí)不清。對重力數(shù)據(jù)及其轉(zhuǎn)換參數(shù)進(jìn)行分析，可以為研究盆地?cái)嗔褬?gòu)造發(fā)育特征及其與油氣成藏的關(guān)系提供重要依據(jù)［31-33］。

圖7 歐加登盆地構(gòu)造單元?jiǎng)澐旨把芯繀^(qū)位置（據(jù)參考文獻(xiàn)［30］）Fig. 7 Location and simplified geological maps of the study area

核工業(yè)航測遙感中心和中國石油集團(tuán)東方地球物理有限責(zé)任公司于2018 年在歐加登盆地中部共同開展了1：10 萬航空重力測量，獲得一套高精度航空重力數(shù)據(jù)。本文為驗(yàn)證LTHG 法實(shí)際應(yīng)用效果，選取歐加登盆地中部坳陷帶作為研究區(qū)（圖7），對研究區(qū)的航空布格重力數(shù)據(jù)進(jìn)行LTHG 處理，并與THG、AS、TDR、THDR、TM、TAS、TTHG 和NVDR_THDR在內(nèi)的幾種方法進(jìn)行對比分析，為降低數(shù)據(jù)噪聲的影響，對布格重力異常數(shù)據(jù)進(jìn)行了向上延拓，延拓高度為2 km，經(jīng)延拓處理后識(shí)別的場源體邊界會(huì)更加光滑，但場源體的形態(tài)不會(huì)發(fā)生改變；LTHG 法a取值為5（a為5 時(shí)，識(shí)別結(jié)果足夠尖銳且無空值現(xiàn)象）。從各方法的識(shí)別結(jié)果可以看出，研究區(qū)斷裂走向整體呈NE 和NEE 向，與地質(zhì)資料相吻合，但各方法的識(shí)別結(jié)果及識(shí)別精度上存在一定差異。由圖8A、B可見，THG 法和AS 法對幅值較大的異常識(shí)別效果較好，對幅值較低的識(shí)別效果稍顯不足；圖8C 顯示，TDR 雖然均衡了不同埋深的場源體的邊界信息，但無法產(chǎn)生銳利的邊界信息；圖8D、E 和F 顯示THDR、TM 和TAS 識(shí)別效果不甚理想，假邊界信息較多。圖8G、H 和I 顯示TTHG、NVDR_THDR 和LTHG 三種方法的識(shí)別效果較好，處理得到的重力異常線性特征更加明顯，同時(shí)又顯示更多的地質(zhì)細(xì)節(jié)，但相對而言LTHG 識(shí)別的結(jié)果更加直觀，分辨率更高。

圖8 研究區(qū)航空布格重力異常邊界識(shí)別結(jié)果Fig.8 The edge recognition results of airborne gravity data in the study area

結(jié)合研究區(qū)的區(qū)域地質(zhì)資料和其他物探資料，以LTHG 方法為主，TTHG、NVDR_THDR等方法為輔助，推斷了研究區(qū)內(nèi)的斷裂共10 條（表2），以近NE 向和近NW 向?yàn)橹鳎渲?，NE 向斷裂5 條，分別為 F1、F2、F3、F4、F5，NW 向斷裂5 條，為 F6、F7、F8、F9、F10。 將本次的邊界識(shí)別結(jié)果與前人文獻(xiàn)［30］進(jìn)行對比分析，如圖9 所示。本次推斷的NE 向斷裂與前人研究成果基本一致，同時(shí)新推斷了一組NW 向斷裂，從布格重力數(shù)據(jù)及邊界識(shí)別結(jié)果來看，NW 向斷裂基本錯(cuò)切了NE 向斷裂；前人研究成果中，中部凹陷、中央凸起及西部凹陷北界斷裂與本次推斷的F2 斷裂一致，南界斷裂與本次推斷F3 斷裂一致。

表2 研究區(qū)斷裂推斷成果表Table 2 The results of inference fault in the study area

圖9 研究區(qū)推斷斷裂與前人解釋成果對比圖Fig.9 Comparison of interpretation faults and previous interpretation results in the study area

4 結(jié) 論

1）理論模型的合成重力數(shù)據(jù)處理結(jié)果顯示，LTHG 法能很好的識(shí)別出不同埋深的場源體邊界，提供更豐富的場源體信息。與其他常用方法進(jìn)行比，LTHG 法的識(shí)別的結(jié)果準(zhǔn)確、清晰，分辨率更高，且具有一定的抗干擾能力；

2）理論模型試驗(yàn)表明，采用LTHG 法進(jìn)行邊界識(shí)別時(shí)，邊界位置不受a取值的影響，但當(dāng)a＞10 時(shí)，識(shí)別結(jié)果會(huì)出現(xiàn)數(shù)據(jù)缺失的情況，a＜1 時(shí)，模型邊界位置處異常峰型平緩，因此1≤a＜10 時(shí)，LTHG 法應(yīng)用效果較好；

3）對實(shí)測航空重力數(shù)據(jù)處理得到的重力異?？傮w上與研究區(qū)的斷裂有較好的相關(guān)性，同時(shí)顯示了更多的地質(zhì)細(xì)節(jié)，最終以LTHG 法處理的結(jié)果為主，并綜合其他方法對研究區(qū)的斷裂進(jìn)行了推斷解釋，推斷NE 向斷裂與前人研究成果基本一致，新推斷一組NW 向斷裂，且NW 向斷裂錯(cuò)切了NE 向斷裂。