小波閾值去噪方法在重力資料處理中的應用

高 鐵，徐英奎，顏 萍，吳燕岡

（沈陽地質調查中心方法室，沈陽 110034）

小波閾值去噪方法在重力資料處理中的應用

高 鐵，徐英奎*，顏 萍，吳燕岡

（沈陽地質調查中心方法室，沈陽 110034）

簡要介紹了小波變換的基本理論與小波閾值去噪實現過程。通過單球體與雙球體模型試驗，探討了閾值去噪時小波基以及分解層數的選擇問題；最后利用小波閾值去噪方法，對某實測布格重力資料進行處理，有效濾除了隨機噪聲。去噪后計算的水平方向導數明顯比去噪前直接計算的水平方向導數效果好。

小波閾值去噪；小波基；分解層數；重磁資料

0 引言

在重磁勘探中，由于受到各種因素的影響，使得重磁資料中常含有較強的噪聲，這些伴隨在有用信號中的噪聲是影響重磁資料后續處理精度的重要因素（例如在各種反演方法中常需要求取異常的水平和垂向高階導數，但高階導數方法對高頻干擾噪聲具有較強的放大作用，會大大降低反演的精度）。因此在對重磁資料進行各種特殊處理之前，在不降低信號分辨率的基礎上，對其進行去噪處理是十分必要的。

重磁資料去噪傳統的方法主要是基于傅里葉變換（例如匹配濾波、維納濾波等），但是傅里葉變換是把信號從時域完全變換到頻域進行分析，計算過程中丟失了時域的全部信息，在時域方面分辨性很差，而在頻域方面不能反映信號的細節［1－4］。

小波變換是20世紀80年代發展起來的一種新的方法技術，由于其具有時頻分析、多分辨率和去相關性等特點，在信號去噪方面得到了廣泛地應用。目前，已有許多學者做了大量的試驗，研究利用小波變換對地震信號、測井曲線、探地雷達信號等方面的去噪［2－10］，但是在重磁資料的應用，研究還比較少［11］。本研究首先介紹了小波變換的基本理論以及小波閾值去噪的原理與過程，然后探討了在重磁資料去噪處理中，閾值函數、小波基以及空間尺度的選取問題，最后將小波閾值去噪方法應用到某實測布格重力異常中。

1 方法原理

1．1 小波變換

設ψ（t）∈L2（R），若其傅里葉變換滿足條件為一個基本小波或母小波。將母函數ψ（t）經伸縮和平移后得ψa，b（t），稱其為一個小波序列。

其中：a為伸縮因子；b為平移因子。

對于任意函數f（t）∈L2（R）的連續小波變換為

小波變換的一個重要特征是有一個靈活可變的時間－頻率窗，它在高的“中心頻率”時自動變窄，而在低的“中心頻率”時自動變寬，這樣它對高頻信號有較高的分辨率，而對低頻信號又能給出完整的信息［］。

1．2 小波閾值去噪

閾值去噪方法的思想就是對小波分解后的各層系數中模大于和小于某閾值的系數分別處理，然后對處理完的小波系數再進行反變換，重構出經過去噪后的信號［12］。具體步驟如下［5］：

1）選擇一個小波基并確定分解層數，對原始含噪數據進行小波變換，得到含噪小波系數。

2）選取一種去噪準則，確定噪聲的閾值，給出小波系數估值。

3）由小波系數估值進行小波反變換，重構信號得到信號估值。

常用的閾值函數有硬閾值函數和軟閾值函數：

1）硬閾值函數。它可以很好地保留信號邊緣等局部特征，其表達式為

2）軟閾值函數。它的處理結果會相對圓滑，但是會造成邊緣模糊等失真現象，其表達式為

3）半軟閾值函數。它可兼顧軟閾值和半軟閾值方法的優點，其表達式為

2 模型試驗

重力異常是由地下淺部與深部所有地質體的綜合反映，所以等值線會比較圓滑，不會出現跳變點，因此選用軟閾值比硬閾值更適合重力資料的去噪。

小波去噪時需要考慮的另外兩個參數是小波基和分解層數，選擇不同的小波解決同一個問題可能會產生不同的結果，因此去噪處理時必須考慮選擇最優小波基，而隨著分解層數的增加，噪聲的幅值越來越小，有用信號的幅值逐漸增大，因此要選擇合適的分解層數以達到既保持原有信號的有用部分，又能濾去干擾部分的效果。

小波基和分解層數的選擇對去噪的結果有很大的影響，前人已經做過許多理論分析和實驗：通過對地震信號進行小波閾值去噪后的信號與原始信號的誤差大小進行分析，得出sym族小波優于db族小波，而sym N小波族中，N取4～6的中等數，得到的誤差最小，而地震信號中的噪聲只存在于前三尺度的分解結果中，故采用3次分解去噪即可消除地震中的噪聲［6］；db小波比較適合探地雷達數據處理［8］，分解層數取“3”既能保證去噪的效果，又能盡量減少處理的工作量［9］。

本研究通過對疊加隨機噪聲的單球體模型和雙球體模型重力異常，選用不同的小波基以及分解層數進行去噪試驗，利用去噪后的信號與原始信號的標準偏差作為誤差分析的標準，來探討重磁資料小波閾值去噪處理中小波基和分解層數的選擇問題。

2．1 單球體模型

選取單個球體為模型，其參數為：球體中心坐標為（0，0），埋深50 m，半徑100 m，剩余密度為1．0 g／cm3。其重力異常等值線如圖1（a）所示，最大值為3．5 g．u．。疊加幅度為（－0．1，0．1）的隨機噪音后，重力異常等值線如圖1（b）所示。表1為不同小波基以及分解層數的誤差結果。當選取的小波基為db6，分解層數為4層時的小波閾值去噪后的重力異常等值線見圖1（c）。

2．2 雙球體模型

在單球體模型的基礎上再疊加一個球體，其參數為：球體中心坐標為（50，40），埋深20 m，半徑15 m。此球體埋深相對淺，規模也小，產生的重力異常相當于淺部的高頻有效信號。其重力異常等值線如圖2（a）所示，最大值為4．42 g．u．，同樣疊加幅度為（－0．1，0．1）的隨機噪音后，重力異常等值線如圖2（b）所示。表2為其不同小波基以及分解層數的誤差結果。當選取的小波基為db6，分解層數為3層時的小波閾值去噪后的重力異常等值線見圖2（c）。

2．3 參數選擇

由表1和表2可知，Haar方波的誤差很大，而當N和分解層數任意時，coif族小波的誤差變化范圍比sym族小波及db族小波大，所以sym族小波與db族小波優于haar小波及coif族小波。sym N和db N小波族中，N不能取“1”，其他取值對結果影響很小，最好取3～6；分解層數取3層最優。

圖1 單球體模型小波閾值去噪Fig．1 The application of wavelet threshold de－noising method to single sphere model

圖2 雙球體模型小波閾值去噪Fig．2 The application of wavelet threshold de－noising method to double sphere model

表1 不同小波基以及分解層數的誤差結果比較（單球體模型）Tab．1 The error by different wavelet basis function and decomposition order（single model）

3 實測數據處理

圖3是某區實測的布格重力異常等值線圖（單位毫伽：10－5m／s2），從圖3可以看出，等值線總體呈東北方向展布，西北高東南低。由于測量和各項改正過程中各種誤差的存在，導致等值線雜亂，甚至在某些地方出現明顯的沿南北測線方向的條帶狀虛假異常，如果在處理之前不對這些干擾去除，將嚴重影響后續各種反演與解釋的精度。采用二維小波閾值去噪對實測資料進行處理，采用的小波基為sym7，尺度為3，得到結果如圖4所示。從圖4可以看出，去噪后的等值線完好地保持了原有的整體形態，并且在沒有產生虛假異常的基礎上，等值線變得相對圓滑，去除了各種隨機干擾以及南北向條帶狀虛假異常。

表2 不同小波基以及分解層數的誤差結果比較（雙球體模型）Tab．2 The error by different wavelet basis function and decomposition order（double model）

水平方向導數可以突出某個特定方向的異常，是重磁數據處理中的一種經典方法，在此引入用以說明去噪的重要性。圖5是直接用實測的布格重力異常利用頻率域方法得到的水平方向導數，由于噪音的存在，求導過程中將高頻成分放大，導致結果非常雜亂，肉眼很難辨別出規律。圖6是先做小波閾值去噪處理然后再求導得到的水平方向導數，四個方向的導數等值線都比較光滑，很容易辨別出東北向的異常最強。

4 結論

圖3 實測布格重力異常等值線圖Fig．3 The practical gravity data contours

圖5 原始布格重力異常水平方向導數Fig．5 Directional derivatives calculated from the original gravity data

由于小波變換具有時頻分析、多分辨率和去相關性等特點，在信號去噪方面得到了廣泛地應用。通過建立簡單的球體模型，用小波閾值方法對其疊加了隨機噪聲后的數據進行去噪處理，探討了在重磁資料去噪處理中，軟硬閾值函數、小波基以及分解層數的選取問題。結論如下：

1）軟閾值函數比較適合于重磁數據去噪。

2）通過誤差分析比較，各類小波都可以很好地濾除噪聲，但db N族小波和sym N族小波優于其他小波基，N不能取1，其他取值對結果影響很小，最好取3～6。

3）分解層數取3層去噪效果最佳。

最后將小波閾值去噪方法應用到某實測布格重力異常中，有效地濾除了隨機噪聲，去噪后計算的水平方向導數明顯比去噪前直接計算的水平方向導數效果好。經理論和實測數據的處理結果表明，小波閾值去噪方法能夠有效地去除重磁資料中的隨機噪聲，是一種切實可行的重磁信號去噪新方法。

圖4 小波閾值去噪后的重力等值線圖Fig．4 The contours after de－noising

圖6 去噪后的布格重力異常水平方向導數Fig．6 Directional derivatives calculated from the denoised gravity data

［1］ 方興．小波分析及其去噪性能研究 ［J］．浙江師范大學學報：自然科學版，2003，26（4）：345－349．

［2］ 中國生，徐國元，江文武．基于小波變換的爆破地震信號去噪的應用［J］．中南大學學報：自然科學版，2006，37（1）：155－159．

［3］ 李英，張淑貞，許康生．小波降噪方法在地震信號處理中的應用［J］．西北地震學報，2006，28（2）：159－162．

［4］ 張旭東，詹毅，馬永琴．小波變換與傅立葉變換在地震資料去噪中的對比［J］．內蒙古石油化工，2007，7：29－44．

［5］ 姜弢，劉慶普，胡留軍．地震信號去噪的小波方法研究［J］．哈爾濱工程大學學報，2002，23（4）：86－90．

［6］ 柳建新，韓世禮，馬捷．小波分析在地震資料去噪中的應用［J］．地球物理學進展，2006，21（2）：541－545．

［7］ 鄒長春，楊欣德，潘令枝．一種基于小波變換的測井曲線去噪新方法［J］．物探與化探，1999，23（6）：462－466．

［8］ 鄒海林，寧書年，林捷．小波理論在探地雷達信號處理中的應用［J］．地球物理學進展，2004，19（2）：268－275．

［9］ 劉秀娟，鄧世坤，徐保林．探地雷達信號小波變換去噪［J］．工程勘察，2006，10：66－71．

［10］何蘭芳，王緒本，何展翔，等．MT時間序列的小波去噪分析［J］．地震地質，2001，23（2）：222－226．

［11］謝仕民，王子靜，鐘鈞宇．小波降噪技術在地磁數據庫預處理中的應用研究［J］．戰術導彈控制技術，2008，30（2）：15－22．

［12］周偉．MATLAB小波分析高級技術［M］．西安：西安電子科技大學出版社，2006．

Application of wavelet threshold de－noising method in gravity data processing

GAO Tie，XU Ying－kui，YAN Ping，WU Yan-gang
（Shenyang Center China Geological Survey，Shenyang 110034，China）

In this paper，the basic theory of wavelet transform and the process of wavelet threshold denoising are introduced．Then，we discussed the criterion for wavelet bases and scales to be selected by establishing single sphere and double spheres models and doing various experiments on the gravity of the models．Finally，the wavelet threshold de－noising method is applied to practical gravity data，and the result is completely denoised．Directional derivatives calculated from the denoised data are more effective than that calculated from the original data

wavelet threshold de－noising；wavelet bases；wavelet scales；gravity and magnetic data

P 631．1

A

10．3969／j．issn．1001-1749．2014．05．09

1001-1749（2014）05-0566-05

2014-03-03 改回日期：2014-07-23

沈陽地質礦產研究所自設項目（2012007）

高鐵（1983-），男，博士，主要研究方向為地球物理綜合方法，E-mail：39085496＠qq．com。

*通訊作者：徐英奎（1963-），男，工程師，研究方向為構造地質，E-mail：997182861＠qq．com。