基于紋理分析的航磁異常邊界增強技術研究及應用

胡 聰

(江蘇省地質勘查技術院，南京 210049)

0 引言

位場數據邊界增強技術在地球物理界是經久不衰的話題。在地球物理資料的定性、半定量解釋中，引起地質學者特殊興趣的就有地球物理圖像中的線性異常，因為這類異?？赡軐叵碌臄嗔?、接觸帶或其他構造特征。一般來講，由重磁異常劃分斷裂的主要標志有線性梯度帶、異常特征的分界線，異常線性過渡帶、線性異常等[1]，但由于地質構造演化不同，以及地質與地球物理條件的差異，上述這些特征標志在圖像上并不一定都有明顯反映，這樣就需要應用突出斷裂信息的方法來輔助地質解釋。

目前航磁異常面積資料處理與轉換的方法很多，每種處理方法因其目的不同而各有其特點[2-3]。針對位場數據的特點，我們采取灰度圖像紋理分析技術，突出異常中的線性構造特點。該技術應用于航磁異常的線性構造特征增強方面具有實際應用價值。

早期地球物理數據的濾波在空間域中進行，人們用量板進行運算，效率很低。隨著計算機技術的發展，某些算法應運而生(如快速Fourier變換)，它促進了頻率域濾波技術的發展，并使該技術在地球物理學中得到廣泛應用[4]。在重力、地磁圖像處理中，常用的濾波方法如低通濾波和高通濾波，總梯度模、斜導數、歸一化標準差等方法應運而生[5]。

1 基本原理

1.1 紋理分析[6]

八十年代興起的圖像處理技術在勘查地球物理和勘查地球化學的信息處理中，己得到日益廣泛地應用。處理后的各種圖像，既能形象地表達數據中不同級次的強度信息，又能直觀地顯示出場的二維結構特征，特別對等值線圖中反映不明顯的微弱信息更具明顯效果。利用圖像處理方法對數據圖像實施增強和變換，增強和突出目標信息，壓制背景信息、干擾信息和噪聲，并通過對初始信息的選擇、演繹、分解和綜合，形成新的標志特征和組合特征，更有效地反映場的特征。同時圖像將場的信息轉化為生動的視覺形象，把場的強度、變化率、形態等轉化為亮度與色調、粗糙度和空間幾何形狀等直觀視覺特征，更便于分析解釋。圖像增強方法很多，如對比度增強、直方圖均衡化、邊緣銳化、平滑、中值濾波、高通低通濾波、帶通帶阻濾波、多譜增強和彩色增強等。紋理分析作為數字圖像處理技術中的一種，已經發展得非常成熟，目前還沒有統一和公認的有關紋理的確切定義，紋理圖像在局部區域內可能呈現不規則性，但整體上則表現出一定的規律性，其灰度分布往往表現出某種周期性，紋理圖像所表現出的這種特有的性質稱為紋理。紋理分析方法可大致分為兩大類：統計分析法和結構分析法。前者從圖像有關屬性的統計分析出發；后者則著力找出紋理基元，然后從結構組成上探索紋理的規律。

紋理分析對于異常突變具有較好地識別能力，有利于任何方位的線性特征的增強，可以用于提取航磁異常線性異常特征。在航磁異常圖中，線性異常特征的增強可以為斷裂體系分析、弱磁異常特征增強與目標提取起到很好的作用?；叶裙采仃囀且环N常用的數字紋理特征，灰度共生矩陣用灰度的空間分布來表征紋理。由于粗紋理的灰度分布隨距離的變化比細紋理緩慢得多，因此二者有完全不同的灰度共生矩陣。根據灰度共生矩陣在紋理分析中的這種特性，就可以將它用于航磁異常的異常特征分析。

1.2 灰度共生矩陣

灰度共生矩陣可以計算圖像紋理的一致性和反差性?；诨叶裙采仃囂崛〖y理特征的方法，是一個經典的統計分析方法，對它的研究也已經有了很久的歷史，是目前公認的一種紋理分析方法[7]。Haralick[9]于 1973 首先提出灰度共生矩陣(GLCM)，該方法優于灰度游程長度法和光譜方法，是一種得以廣泛應用的、常用的紋理統計分析方法和紋理測量技術。1992 年P.P.Ohanian[8]給出了對幾種紋理測量技術的比較結果。

設圖像水平和垂直方向上各有Nc×Nr個像元，將每個像元出現的灰度量化為Ng層，設Lx={1,2,…,Nc}為水平空間域，Ly={1,2,…,Nr}為垂直空間域，G={1,2,…,Ng}為量化灰度層集。集Lx×Ly為行列編序的圖像像元集，則圖像函數f可表示為一個函數:指定每一個像元具有Ng個灰度層中的一個值G，即f:Lx×Ly→G。灰度共生矩陣定義為在圖像域Lx×Ly范圍內，兩個相距為d，方向為θ的像元在圖像中出現的概率，即:

P(i,j|d,θ)=#{[(k,l),(m,n)]∈(Lx×Ly)×(Lx×Ly)|d,θ,f(k,l)=i,

f(m,n)=j}

(1)

式中 #{x}表示大括號中成立的像元對數；i、j分別為點(k，l)、(m，n)的灰度值。如距離為d，水平方向的灰度共生矩陣計算公式為式(2)：

P(i,j|d,0)=#{[(k,l),(m,n)]∈(Lx×Ly)× (Lx×Ly)|k-m=0,|l-n|=d,f(k,l)=i,f(m,n)=j}

(2)

同理，距離為d，對角方向的灰度共生矩陣為式(3)與式(4)：

P(i,j|d,45)=#{[(k,l),(m,n)]∈(Lx×Ly)×(Lx×Ly)|k-m=d,l-n=-d,

f(k,l)=i,f(m,n)=j}

(3)

P(i,j|d,135)=#{[(k,l),(m,n)]∈(Lx×Ly)×(Lx×Ly)|k-m=-d,l-n=d,

f(k,l)=i,f(m,n)=j}

(4)

灰度共生矩陣提供了圖像灰度方向、間隔和變化幅度的信息，可根據共生矩陣來計算一些相應的特征值，以相應特征值來表征圖像的紋理信息?；叶裙采仃嚨募y理分析方法需要選擇一定大小的移動窗口，這里對每一個特征都以此滑動窗口(如 3×3，5×5)進行計算，在數據處理過程中，移動窗口的選擇和灰度量化級的選擇都會對特征值產生影響，對于較粗糙的紋理，窗口可取大一些，而細膩的紋理窗口可取較小，合適的灰度級別并不容易權衡，且考慮到實時性的要求，所以在不影響紋理特征的前提下，可以盡量縮小灰度值的范圍。文獻[9]提出了十四種基于灰度共生矩陣的紋理特征參數。作者利用局部平穩度作為特征參數，進行了航磁異常邊界增強的實驗。局部平穩度的公式如式(5)：

(5)

局部平穩度反映出圖像紋理的同質性、度量圖像紋理局部變化的程度。其值大，則說明圖像紋理的不同區域間缺少變化，局部非常均勻。

2 模型試驗

球體是很常見的三度模型，如礦巢、巖株和穹隆構造等地質體在地面所產生的磁異常，可以近似看作球體的異常[10]。設一球體地質模型，磁化強度模值為M= 12 A/m，磁化方向的傾角為45°，磁化方向水平投影與x方向的夾角為45°，與y方向的夾角為45°。以m為單位，x軸坐標范圍是(-300，300)，y軸坐標范圍是(-300，300)，z方向埋深是100；球體中心坐標為(0，0，100)，球體半徑為55 m，計算面坐標為z= 0[11]。球體理論異常灰度圖和紋理分析結果如圖1所示。

圖1 球體模型理論異常圖與紋理分析結果圖Fig．1 Theoretical magnetic anomaly map of sphere model and texture analysis result(a)理論異常圖;(b)紋理分析結果圖件

圖1(a)為此球體模型磁異常灰度圖像，白色代表正值，黑色代表負值，航磁異常邊界信息較模糊，當選取窗口為3×3，灰度量化級為“20”時，圖像的局部平穩度達到臨界值，磁異常邊界信息得到最大的增強，如圖1(b)所示。在北東向方向，正異常、負異常交替變化，圍繞極值點，會出現水波紋型異常特征。

3 應用

為了進一步驗證方法的可行性，選擇了非州某區域的航磁異常平面數據進行了處理。該區域地層大體分為下部基巖和上部第四系覆蓋層，區內植被十分茂密，地層走向總體呈NEE向，圖2為此區域的航磁異?；叶葓D，航磁異常范圍在 -500 nT到300 nT，正、負異常伴生，大致呈東高西低的趨勢，此區域的大部分位置磁異常變化較為平緩。將灰度共生矩的窗口固化為3×3，灰度量化級選作5、10和20，經過本文的方法對這塊區域的磁異常進行了計算，如圖3所示。

圖2 非洲某區域航磁異?；叶葓DFig．2 Aeromagnetic anomalies grayscale in a region of Africa

圖3 非洲某區域航磁異常數據紋理分析對比結果圖Fig．3 Texture analysis comparison of aeromagnetic anomalies in a region of Africa窗口都為3×3；(a)灰度量化級為5；(b)灰度量化級為10；(c)灰度量化級為20

圖4 非洲某區域斷裂位置圖Fig．4 Fracture location map in a region of Africa

從圖2中可以看出，磁異常主體走向是NEE，與地層走向基本吻合?；叶攘考壔x為“5”時，圖2中出現串珠狀磁異常帶，隨著灰度量化級的增大，細節得到逐步的增強，當灰度量級化選為“20”時，航磁異常細節信息達到最大。在磁異常梯級帶及磁異常極值點周圍出現水波紋狀紋理特征，紋理越密，則表明此處的磁異常梯度變化越大。根據串珠狀磁異常帶，異常等值線扭曲等，即以劃分斷裂的主要標志為原則，作者劃分出五條斷裂位置(圖4中以紅實線表示)，三條為NEE向的，一條為NW向的，一條為近EW向的，根據實際資料的顯示，該地區的斷裂位置是基本準確的。

4 結論

(1)理論實驗和實際運用表明，基于紋理分析的航磁異常邊界增強方法，是處理航磁異常數據的一種有效方法。

(2)將紋理分析運用到航磁數據中，可以有效地增強航磁異常邊界信息。模型實驗表明，紋理分析方法能夠在精度上、清晰度上，較穩定地提取航磁邊界信息。

(3)在實際運用中，也顯示了該方法的可行性和實用性，根據航磁異常判斷與推測礦區的斷裂分布與延伸，可以為下一步礦體靶區的圈定提供有利的證據，也可以解決地質構造其他相關問題。

(4)隨著航磁數據采集精度越來越高，我們可以獲取更加豐富的磁異常信息，顯然對于這些更加細致的信息進行分析具有很好的應用前景。在計算灰度共生矩陣和紋理特征值的過程中，灰度量化級和窗口的選擇需要根據經驗和不斷地嘗試。相信這種方法也可應用于其他大數據量高精度的數據信息增強和提取中。

參考文獻：

[1] 王家林．石油重磁解釋[M]．北京：石油工業出版社，1991．

[2] 柏冠軍，吳漢寧，趙希剛，等．重力資料識別鄂爾多斯盆地線性構造方法研究[J]．地球物理學進展，2007，22(5):1386-1392．

[3] 嚴良俊，胡文寶，姚長利．重磁資料面積處理中的濾波增強技術與應用[J]．勘探地球物理進展，2006， 29(2):102-103．

[4] 穆石敏，申寧華，孫運生．區域地球物理數據處理方法及其應用[M]．長春：吉林科學技術出版社，1990．

[5] 肖鋒，吳燕岡，孟令順．重力異常圖中的邊界增強和提取技術[J]．吉林大學學報：地球科學版，2011，41(4):1197-1203．

[6] 高成．Matlab圖像處理與應用[M]．北京：國防工業出版社，2007．

[7] COOPER G R J. The Antialiased Textural Analysis of Aeromagnetic Data[J]. Computers &Geosciences.2009，35: 586-591．

[8] OHANIAN P P,DUBES R C. Performance Evaluation for Four Class of Texture Features[J]. Pattern Recognition，1992，25(8):819-833．

[9] HARALICK R M, SHANMUGAM K. Texture features for image classification[J], 1973,3(6):610-621．

[10] 胡聰．利用重力資料進行邊界和質心提取的方法研究及應用[D]．長春：吉林大學，2012．

[11] 侯重初，劉奎俊，徐瑾．任意形體磁場及其各階導數的正演[J]．物探化探計算技術，1988， 10(1):7-20．