頻率域航空電磁異常的數(shù)字化識(shí)別篩選技術(shù)

李文杰

(中國(guó)地質(zhì)科學(xué)院 地球物理地球化學(xué)勘查研究所， 廊坊 065000)

頻率域航空電磁法在進(jìn)行礦產(chǎn)勘查時(shí)沿測(cè)線連續(xù)測(cè)量，平均每15 m一個(gè)數(shù)據(jù)，測(cè)線長(zhǎng)度十幾千米至上百千米，一個(gè)測(cè)區(qū)測(cè)線總量從幾百測(cè)線千米到幾十萬(wàn)測(cè)線千米不等。一般每測(cè)線千米數(shù)據(jù)中平均會(huì)出現(xiàn)2 處～3處明顯的異常，一個(gè)測(cè)區(qū)則會(huì)產(chǎn)生數(shù)千至數(shù)十萬(wàn)個(gè)異常。在如此巨量的航空電磁異常中，排除飛行過(guò)載、工業(yè)電、飛行高度變化、地質(zhì)噪聲等干擾引起的假異常，快速和高效地識(shí)別和篩選出有意義的基巖異常，是有效發(fā)揮航空電磁法快速高效勘查優(yōu)勢(shì)的關(guān)鍵技術(shù)。

作者基于中國(guó)地質(zhì)科學(xué)院地球物理地球化學(xué)勘查研究所(簡(jiǎn)稱物化探所，IGGE)研制的國(guó)內(nèi)唯一自主技術(shù)固定翼頻率域航空電磁系統(tǒng)—HDY-402三頻航電系統(tǒng)的特有技術(shù)，參考國(guó)內(nèi)、外相關(guān)技術(shù)，借鑒傳統(tǒng)的人工異常識(shí)別篩選方法，開發(fā)了計(jì)算機(jī)軟件異常識(shí)別、篩選和解釋算法(Computer-aided Anomaly Detection, CAD)，并編制了應(yīng)用軟件，應(yīng)用在內(nèi)蒙古地區(qū)開展的約42.3萬(wàn)測(cè)線千米航空物探綜合勘查項(xiàng)目中，有效地提高了海量異常識(shí)別篩選的工作效率和質(zhì)量。

1 國(guó)內(nèi)外研究情況

國(guó)外在上世紀(jì)七十年代初，針對(duì)固定翼頻率航電數(shù)據(jù)，開發(fā)了基于計(jì)算機(jī)的交互式異常識(shí)別技術(shù)[1-2]，提高了異常解釋效率。九十年代末，澳大利亞礦產(chǎn)堪查工業(yè)聯(lián)合研究中心[3](CRC AMET)開發(fā)了EM FLOW數(shù)據(jù)處理解釋軟件，可以在時(shí)間常數(shù)(tau)域內(nèi)根據(jù)設(shè)定的參數(shù)實(shí)現(xiàn)航電異常的快速自動(dòng)識(shí)別和反演計(jì)算。Valleau[4]系統(tǒng)地介紹了DIGHEM直升機(jī)頻率域航空電磁數(shù)據(jù)的處理流程，實(shí)現(xiàn)了航電異常的計(jì)算機(jī)提取、篩選及人工輔助判別。Huang[5-6]提出用多個(gè)頻率電磁耦合譜的不同特征來(lái)識(shí)別頻率域航電異常，用于UXO探測(cè)。

在國(guó)內(nèi)，八十年代初，周鳳桐[7]總結(jié)了固定翼頻率域航空電磁資料的人工處理流程，李文杰[8]研究了固定翼頻率域航空電磁異常軟件識(shí)別技術(shù)。在航空電磁異常的軟件數(shù)字化自動(dòng)快速拾取和識(shí)別方面，國(guó)內(nèi)未見其他相關(guān)報(bào)導(dǎo)。

2 方法步驟

這里討論的航電異常數(shù)字化識(shí)別篩選，主要是指識(shí)別局部異常，篩選出其中的有意義的良導(dǎo)基巖異常。異常的計(jì)算機(jī)自動(dòng)識(shí)別和交互判別，簡(jiǎn)單地說(shuō)就是將人工識(shí)別經(jīng)驗(yàn)和技術(shù)數(shù)字化為計(jì)算機(jī)標(biāo)準(zhǔn)識(shí)別參數(shù)和程序代碼。航電異常軟件自動(dòng)識(shí)別的過(guò)程與人工異常識(shí)別大體一致，需要經(jīng)過(guò)如下幾個(gè)步驟：

1)異常自動(dòng)拾?。杭催x出一般意義上的待選異常。

2)異常識(shí)別篩選：用約束條件來(lái)剔除待選異常中明顯無(wú)意義的異常。

3)交互輔助判別：進(jìn)一步剔除地質(zhì)噪聲、地形效應(yīng)異常等難于自動(dòng)識(shí)別的異常。

2.1 異常自動(dòng)拾取

對(duì)經(jīng)過(guò)噪聲濾波、相位改正、零漂改正等處理的三頻航電剖面數(shù)據(jù)，取梯度后判別拾取出曲線波峰、波谷，再標(biāo)記出所有大于三倍噪聲水平的待選異常；然后提取異常信息，計(jì)算每個(gè)異常的半極值寬度、異常幅值、實(shí)虛分量比、視電阻率、視深度、飛行高度、地形及地形匹配系數(shù)、工業(yè)電干擾、過(guò)載等相關(guān)信息。這些數(shù)字化了的異常信息，是下一步自動(dòng)識(shí)別的基礎(chǔ)信息。

假設(shè)航電某頻道曲線為Eij(i為測(cè)線號(hào)，j為測(cè)點(diǎn)號(hào))(圖1)，則它的梯度為

Δe(i,j)=grad(Eij)=E(i,j+1)-E(i,j)

(1)

如果Δe(i,j)=0，則在Eij上出現(xiàn)一個(gè)拐點(diǎn)，當(dāng)在拐點(diǎn)相鄰點(diǎn)的梯度滿足式(2)時(shí)

(2)

則拐點(diǎn)對(duì)應(yīng)E(i,j)曲線上的一個(gè)谷點(diǎn)(如圖1中的點(diǎn)A、C)。同理，當(dāng)拐點(diǎn)相鄰點(diǎn)的梯度滿足：

(3)

則拐點(diǎn)對(duì)應(yīng)E(i,j)曲線上的一個(gè)峰點(diǎn)(如圖1中的點(diǎn)P)，一個(gè)波峰加相鄰兩個(gè)波谷就是一個(gè)異常。

異常拾取過(guò)程中，附帶地可以計(jì)算出異常的半極值寬度Pw、異常的峰值高度Pp、實(shí)虛分量比，以及峰值點(diǎn)處的視電阻率、視電阻率峰值等信息，以備后續(xù)步驟篩選。

2.2 自動(dòng)異常識(shí)別

根據(jù)基巖異常的特征，設(shè)定異常峰值、半極值寬度、背景場(chǎng)值、飛行高度、多頻道交集、視電阻率峰值、覆蓋層判別等篩選條件，對(duì)待選異常進(jìn)行自動(dòng)篩選識(shí)別，給出異常分類。

雖然數(shù)據(jù)預(yù)處理時(shí)進(jìn)行了噪聲濾波，但系統(tǒng)噪聲在異常篩選階段依然需要被考慮。由于噪聲的存在，異常峰值高度Pp一般會(huì)有一個(gè)最小限制即Ppmin≥3·vN，vN是系統(tǒng)的噪聲水平，即當(dāng)異常的峰值高度Pp大于三倍噪聲水平時(shí)才被考慮。

圖1 異常的拾取原理

1)半極值寬度。大多數(shù)礦體規(guī)模比較小，尤其硫化礦體，其走向長(zhǎng)度大約在200 m～300 m以下的居多。即使是規(guī)模比較大的礦體(如磁鐵礦體或大型硫化礦體)，其走向長(zhǎng)度可以大于500 m ～600 m，但往往由于地形或者礦體本身傾伏等原因，并不總是能使礦體整個(gè)地位于系統(tǒng)的探測(cè)深度范圍之內(nèi)，因而實(shí)際上反映在航電異常曲線的半極值寬度一般都小于500 m。

2)異常背景。有意義的礦體異常一般不會(huì)出現(xiàn)在大面積的低阻區(qū)，因此異常波谷的場(chǎng)值可以作為異常背景加以約束，波谷谷點(diǎn)處的場(chǎng)值大于某個(gè)限制后，異常可以不被考慮。

3)測(cè)量飛行高度。航電接收的二次場(chǎng)隨著測(cè)量高度的增加以三次冪或更大的速率衰減。測(cè)量高度超高時(shí)，航電響應(yīng)中噪聲、零漂等干擾因素所占比重成倍增加，測(cè)得異常的可信度很低。根據(jù)實(shí)踐經(jīng)驗(yàn)，一般可以取70 m～80 m作為異常峰值點(diǎn)測(cè)量高度的上限。

4)工業(yè)電干擾。航電接收的二次場(chǎng)響應(yīng)異常有一部分是由于50 Hz工業(yè)電干擾及飛行過(guò)載引起的，這類異常可以利用輔助設(shè)備提供的輔助參數(shù)來(lái)判別和剔除。一般航電系統(tǒng)50 Hz監(jiān)測(cè)輸出的信號(hào)小于100 mV；當(dāng)測(cè)線上有工業(yè)電高壓線等干擾時(shí)，50 Hz監(jiān)測(cè)輸出的信號(hào)會(huì)遠(yuǎn)大于100 mV，并隨干擾強(qiáng)度的變化而改變。

5)過(guò)載干擾。系統(tǒng)自帶的過(guò)載計(jì)在飛行相對(duì)平穩(wěn)時(shí)的輸出通常小于±0.15 g，當(dāng)飛機(jī)遇到氣流或突然的機(jī)動(dòng)動(dòng)作會(huì)使過(guò)載計(jì)輸出產(chǎn)生正負(fù)變化，當(dāng)過(guò)載輸出大于±0.15 g時(shí)所引起的異?？梢耘袛酁檫^(guò)載異常，并加以剔除。

6)多頻道信息。HDY-402航電系統(tǒng)共有463 Hz、1 563 Hz、8 333 Hz三個(gè)工作頻率，相對(duì)而言中頻和低頻的探測(cè)深度較大，因此在判斷異常時(shí)主要以中低頻為主、高頻為輔。如果中低頻同時(shí)具有異常，則是良導(dǎo)基巖異常的可能性較大。結(jié)合低、中、高頻的視電阻率、視飛行高度數(shù)據(jù)，依據(jù)視電阻率與二層大地模型的關(guān)系，還可以推斷出二層大地的電阻率高低，對(duì)覆蓋層異常的判別幫助很大，也可以作為一類約束條件。

在異常自動(dòng)識(shí)別運(yùn)算過(guò)程中，約束條件如果采用上限和下限雙邊限制，結(jié)合特別給定的約束條件組合，就可以對(duì)某類特定類型的異常進(jìn)行識(shí)別，而不光是良導(dǎo)基巖異常。

舉例說(shuō)明，如果將半極值寬度、峰值高度、異常背景、測(cè)量高度等約束條件的上限和下限都放寬，即相當(dāng)于不限制這些條件，而只將輔助信息中的50 Hz工業(yè)干擾下限設(shè)得較小，如小于100 mV，上限放寬，則可以得到滿足大于100 mV下限標(biāo)準(zhǔn)的的50 Hz工業(yè)電干擾引起的干擾異常。

2.3 交互輔助識(shí)別

數(shù)據(jù)解釋人員根據(jù)異常信息，參考地質(zhì)、化探等相關(guān)資料，對(duì)自動(dòng)識(shí)別拾取的異常逐個(gè)進(jìn)行解釋推斷，剔除地形效應(yīng)、負(fù)地形、覆蓋層等地質(zhì)噪聲引起的假異常，以及過(guò)載、工業(yè)電等引起的干擾異常，對(duì)推斷的基巖異常進(jìn)行標(biāo)識(shí)。同時(shí)對(duì)未被自動(dòng)識(shí)別拾取的異常進(jìn)行補(bǔ)查，以免遺漏可能的基巖異常；最后將所有已標(biāo)識(shí)的異常信息進(jìn)行輸出和匯總。

圖2是兩處在內(nèi)蒙二連浩特測(cè)區(qū)用軟件自動(dòng)識(shí)別出來(lái)的航電典型異常。圖2中曲線自上至下依次為異常響應(yīng)的視電阻率、高頻8 333 Hz(HI/HQ)、中頻1 563 Hz(MI/MQ)、低頻463 Hz(LI/LQ)實(shí)虛分量、飛行高度、地形等剖面曲線。圖2(a)為三頻航電在一處已知熱液型銅礦化點(diǎn)上的異常響應(yīng)曲線，對(duì)應(yīng)有明顯的中低頻視電阻率異常。圖2(b)為一處典型的殘坡積(低洼地)引起的地表覆蓋層異常，對(duì)應(yīng)有高頻視電阻率異常，中低頻視電阻率無(wú)明顯異常。

2.4 異常自動(dòng)識(shí)別的軟件流程

1)準(zhǔn)備數(shù)據(jù)。用于進(jìn)行異常識(shí)別的航電數(shù)據(jù)需經(jīng)過(guò)噪聲濾波處理、相位改正及零漂改正。噪聲濾波處理是必須的，其意不言自明。未進(jìn)行相位改正和零漂改正的航電數(shù)據(jù)，在進(jìn)行異常信息提取時(shí)，將得到錯(cuò)誤的視電阻率、視深度結(jié)果，無(wú)法保證異常背景、覆蓋層異常判定等約束條件的有效性。

2)設(shè)定異常篩選條件。設(shè)定異常峰值、半極值寬度、背景、測(cè)量高度、多頻道交集、視電阻率峰值、覆蓋層判別等篩選條件。上述異常篩選的約束條件因?yàn)槎加蒙?、下限?shù)值賦值，因此可以交集使用，或通過(guò)將某一項(xiàng)或幾項(xiàng)約束條件的上、下限數(shù)值設(shè)定為完全滿足而不參與異常篩選的條件約束。

3)異常拾取。計(jì)算航電各頻道實(shí)、虛分量曲線的梯度，根據(jù)正、負(fù)峰組合判斷拾取波峰、波谷，一個(gè)波峰及相鄰兩個(gè)波谷構(gòu)成一個(gè)被選異常。

4)異常信息提取。計(jì)算被選異常的峰值、半極值寬度、視電阻率、視電阻率背景、測(cè)量高度、多頻道交集、視電阻率峰值、視飛行高度等信息，并根據(jù)半極值寬度等條件識(shí)別雙峰異常。

圖2 已知礦化點(diǎn)異常和殘坡積引起的覆蓋層異常綜合剖面圖Fig．2 Composite profiles of an ore anomaly and a residual sediment anomaly (a)礦化點(diǎn)異常；(b)覆蓋層異常

5)應(yīng)用約束條件篩選異常。應(yīng)用異常篩選條件，篩選出符合各約束條件交集的異常。如上所述，約束條件可根據(jù)任務(wù)需要進(jìn)行設(shè)定，因而可以區(qū)分出50Hz工業(yè)電干擾異常、覆蓋層異常等特定類型異常。

6)輸出異常篩選結(jié)果。將異常結(jié)果輸出為文本文件，其中包含異常的各種信息，及位置坐標(biāo)、異常種類、評(píng)級(jí)等，供方法解釋人員進(jìn)一步篩選。

7)輔助信息解釋。根據(jù)地形剖面曲線等輔助條件對(duì)異常篩選結(jié)果進(jìn)行進(jìn)一步分類，剔除負(fù)地形異常及地形效應(yīng)異常。

將異常識(shí)別算法用軟件代碼實(shí)現(xiàn)，然后再依據(jù)以上流程就可以編織出航電異常識(shí)別程序(圖3)。

圖3 異常自動(dòng)識(shí)別流程Fig．3 Flow chart of the digital anomaly detection and classification

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

圖4為內(nèi)蒙古某測(cè)區(qū)航電異常自動(dòng)拾取的結(jié)果。圖中彩色底圖為視電阻率平面圖(暖色代表視低電阻率，冷色代表高視電阻率)，紫色三角點(diǎn)為自動(dòng)拾取的異常位置，黑色箭頭所指藍(lán)色圓點(diǎn)為已知礦或礦化點(diǎn)?？梢钥吹皆趫D中間最上部①、④號(hào)已知礦化點(diǎn)與自動(dòng)異常拾取的位置是完全對(duì)應(yīng)的。①號(hào)異常為已知的熱液型銅礦化點(diǎn)，自動(dòng)識(shí)別在8 040線上有明顯航電響應(yīng)異常與之對(duì)應(yīng)(圖5)。該異?？刂崎L(zhǎng)度2 km，半極值寬度297 m，飛行高度48 m，中頻和低頻實(shí)虛分量比分別為1.88和1.76，異常清晰。④號(hào)已知熱液型銅礦化點(diǎn)在8 170、8 180線間穿過(guò)，在8 180線礦化點(diǎn)附近自動(dòng)識(shí)別找到了航電異常，但異常峰值較小。②號(hào)已知熱液型銅礦化點(diǎn)在6 870線上方通過(guò)，但測(cè)線上未見明顯異常，在相鄰的6 850線、6 860線航電異常明顯。③號(hào)已知熱液型錳礦化點(diǎn)位于8 240線、8 250線間，為熱液型銅礦化點(diǎn)，在測(cè)線上未發(fā)現(xiàn)明顯異常，但在相鄰的8 270、8 280、8 290線上自動(dòng)識(shí)別出了異常條帶。

圖4 內(nèi)蒙某測(cè)區(qū)航電異常自動(dòng)識(shí)別篩選結(jié)果Fig．4 Digital detection and classification results of a survey area in Inner Mongolia

圖5 內(nèi)蒙某測(cè)區(qū)8040線上的已知礦上航電異常響應(yīng)

4 結(jié)論

頻率域航空電磁異常的數(shù)字化自動(dòng)識(shí)別篩選技術(shù)補(bǔ)充完善了國(guó)產(chǎn)HDY-402三頻航空電磁系統(tǒng)數(shù)據(jù)處理技術(shù)，開發(fā)的異常識(shí)別程序?qū)嵱谩⒏咝?，可根?jù)解釋人員定義的不同約束條件檢出不同類型的航電異常，提高了數(shù)據(jù)處理的效率和質(zhì)量，填補(bǔ)了國(guó)內(nèi)該領(lǐng)域的技術(shù)空白。但是受到數(shù)據(jù)質(zhì)量不佳及數(shù)據(jù)信息量有限等條件的限制，航電異常的自動(dòng)識(shí)別篩選結(jié)果中包含有干擾異常，僅能實(shí)現(xiàn)異常的初步篩選，還不能完全替代解釋人員的人工判別。

參考文獻(xiàn)：

[1] DAVID FOUNTAIN, RICHARD SMITH. 55 Years of AEM; 50 Years of KEGS! Fugro Web Site, 2003, www.fugro.com.au.

[2] PALACKY G J, WEST G F. Computer Processing of Airborne EM Data[J]. Geophysical Prospecting, 1974, 22:490-509.

[3] EM FLOW. User Guide and Reference Manua[M]. Encom Technology Pty Ltd.

[4] NICHOLAS VALLEAU. HEM Data Processing-A Practical Overview[M]. Proc. 14th ASEG Conf and Exhibition, 2000, Perth, Australia.

[5] HAOPING HUANG, WON I J. Identification of Mineral Deposits using airborne electromagnetic spectra[M]. SEG Annual Meeting, 2002.

[6] HOAPING HUANG. Automated anomaly picking from broadband electromagnetic data in an unexploded ordance (UXO) survey [J] .Geophysics, 2003, 68(6), p1870-1876.

[7] 周鳳桐. 補(bǔ)償式航空電磁法[J]. 物化探研究報(bào)道，1982(11):1-78.

[8] 李文杰. 頻率域航空電磁數(shù)據(jù)處理技術(shù)研究 [D] .北京：中國(guó)地質(zhì)大學(xué), 2008.