頻率域航空電磁異常的數字化識別篩選技術

李文杰

(中國地質科學院 地球物理地球化學勘查研究所， 廊坊 065000)

頻率域航空電磁法在進行礦產勘查時沿測線連續測量，平均每15 m一個數據，測線長度十幾千米至上百千米，一個測區測線總量從幾百測線千米到幾十萬測線千米不等。一般每測線千米數據中平均會出現2 處～3處明顯的異常，一個測區則會產生數千至數十萬個異常。在如此巨量的航空電磁異常中，排除飛行過載、工業電、飛行高度變化、地質噪聲等干擾引起的假異常，快速和高效地識別和篩選出有意義的基巖異常，是有效發揮航空電磁法快速高效勘查優勢的關鍵技術。

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

1 國內外研究情況

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

在國內，八十年代初，周鳳桐[7]總結了固定翼頻率域航空電磁資料的人工處理流程，李文杰[8]研究了固定翼頻率域航空電磁異常軟件識別技術。在航空電磁異常的軟件數字化自動快速拾取和識別方面，國內未見其他相關報導。

2 方法步驟

這里討論的航電異常數字化識別篩選，主要是指識別局部異常，篩選出其中的有意義的良導基巖異常。異常的計算機自動識別和交互判別，簡單地說就是將人工識別經驗和技術數字化為計算機標準識別參數和程序代碼。航電異常軟件自動識別的過程與人工異常識別大體一致，需要經過如下幾個步驟：

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

2)異常識別篩選：用約束條件來剔除待選異常中明顯無意義的異常。

3)交互輔助判別：進一步剔除地質噪聲、地形效應異常等難于自動識別的異常。

2.1 異常自動拾取

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

假設航電某頻道曲線為Eij(i為測線號，j為測點號)(圖1)，則它的梯度為

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

(1)

如果Δe(i,j)=0，則在Eij上出現一個拐點，當在拐點相鄰點的梯度滿足式(2)時

(2)

則拐點對應E(i,j)曲線上的一個谷點(如圖1中的點A、C)。同理，當拐點相鄰點的梯度滿足：

(3)

則拐點對應E(i,j)曲線上的一個峰點(如圖1中的點P)，一個波峰加相鄰兩個波谷就是一個異常。

異常拾取過程中，附帶地可以計算出異常的半極值寬度Pw、異常的峰值高度Pp、實虛分量比，以及峰值點處的視電阻率、視電阻率峰值等信息，以備后續步驟篩選。

2.2 自動異常識別

根據基巖異常的特征，設定異常峰值、半極值寬度、背景場值、飛行高度、多頻道交集、視電阻率峰值、覆蓋層判別等篩選條件，對待選異常進行自動篩選識別，給出異常分類。

雖然數據預處理時進行了噪聲濾波，但系統噪聲在異常篩選階段依然需要被考慮。由于噪聲的存在，異常峰值高度Pp一般會有一個最小限制即Ppmin≥3·vN，vN是系統的噪聲水平，即當異常的峰值高度Pp大于三倍噪聲水平時才被考慮。

圖1 異常的拾取原理

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

2)異常背景。有意義的礦體異常一般不會出現在大面積的低阻區，因此異常波谷的場值可以作為異常背景加以約束，波谷谷點處的場值大于某個限制后，異常可以不被考慮。

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

4)工業電干擾。航電接收的二次場響應異常有一部分是由于50 Hz工業電干擾及飛行過載引起的，這類異?？梢岳幂o助設備提供的輔助參數來判別和剔除。一般航電系統50 Hz監測輸出的信號小于100 mV；當測線上有工業電高壓線等干擾時，50 Hz監測輸出的信號會遠大于100 mV，并隨干擾強度的變化而改變。

5)過載干擾。系統自帶的過載計在飛行相對平穩時的輸出通常小于±0.15 g，當飛機遇到氣流或突然的機動動作會使過載計輸出產生正負變化，當過載輸出大于±0.15 g時所引起的異?？梢耘袛酁檫^載異常，并加以剔除。

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

在異常自動識別運算過程中，約束條件如果采用上限和下限雙邊限制，結合特別給定的約束條件組合，就可以對某類特定類型的異常進行識別，而不光是良導基巖異常。

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

2.3 交互輔助識別

數據解釋人員根據異常信息，參考地質、化探等相關資料，對自動識別拾取的異常逐個進行解釋推斷，剔除地形效應、負地形、覆蓋層等地質噪聲引起的假異常，以及過載、工業電等引起的干擾異常，對推斷的基巖異常進行標識。同時對未被自動識別拾取的異常進行補查，以免遺漏可能的基巖異常；最后將所有已標識的異常信息進行輸出和匯總。

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

2.4 異常自動識別的軟件流程

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

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

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

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

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

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

6)輸出異常篩選結果。將異常結果輸出為文本文件，其中包含異常的各種信息，及位置坐標、異常種類、評級等，供方法解釋人員進一步篩選。

7)輔助信息解釋。根據地形剖面曲線等輔助條件對異常篩選結果進行進一步分類，剔除負地形異常及地形效應異常。

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

圖3 異常自動識別流程Fig．3 Flow chart of the digital anomaly detection and classification

3 應用實例

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

圖4 內蒙某測區航電異常自動識別篩選結果Fig．4 Digital detection and classification results of a survey area in Inner Mongolia

圖5 內蒙某測區8040線上的已知礦上航電異常響應

4 結論

頻率域航空電磁異常的數字化自動識別篩選技術補充完善了國產HDY-402三頻航空電磁系統數據處理技術，開發的異常識別程序實用、高效，可根據解釋人員定義的不同約束條件檢出不同類型的航電異常，提高了數據處理的效率和質量，填補了國內該領域的技術空白。但是受到數據質量不佳及數據信息量有限等條件的限制，航電異常的自動識別篩選結果中包含有干擾異常，僅能實現異常的初步篩選，還不能完全替代解釋人員的人工判別。

參考文獻：

[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] 周鳳桐. 補償式航空電磁法[J]. 物化探研究報道，1982(11):1-78.

[8] 李文杰. 頻率域航空電磁數據處理技術研究 [D] .北京：中國地質大學, 2008.