林坤, 馬朝選, 彭曉光, 李軍科

(中國船舶重工集團公司第七一八研究所, 河北 邯鄲 056027)

0 引 言

在自然伽馬測井中,由于核衰變的隨機性質(zhì),導(dǎo)致自然伽馬測井曲線上出現(xiàn)許多與地層性質(zhì)無關(guān)的統(tǒng)計漲落和毛刺干擾,因此在對自然伽馬測井曲線分析使用之前,需要對這些統(tǒng)計漲落和毛刺干擾進行濾波處理,只保留反映地層特性的信息[1]。本文通過幾種降噪常用方法分析,使用快速傅里葉變換FFT降噪方法和小波閾值降噪方法對自然伽馬測井曲線進行降噪分析,以信噪比、相關(guān)系數(shù)、降噪信號能量比和降噪信號與原信號均方差4個性能指標,對降噪效果進行對比,得到最佳的降噪方案。

1 小波變換算法的理論分析

小波變換在許多工程領(lǐng)域都得到了廣泛的應(yīng)用[2]。對于任意ψ(t)∈L2(R),即ψ(t)是平方可積函數(shù),如果ψ(t)的傅里葉變換滿足可容許條件

(1)

式中,R*=R-{0}表示非零實數(shù)全體,則稱ψ(t)是小波母函數(shù)。對于任意的實數(shù)對(a,b),其中參數(shù)a必須為非零實數(shù),稱式(2)形式的函數(shù)為小波基函數(shù)

(2)

對于任意能量有限信號f(t)其小波變換定義為

(3)

1.1 一維離散小波分解和重構(gòu)算法

一維離散小波變換給定一個長度為N的信號s,第1步從原始信號s開始,產(chǎn)生2組數(shù)據(jù),一組是作用于低通濾波器Lo_D得到的近似系數(shù)cA1;另一組是作用于高通濾波器Hi_D得到的細節(jié)系數(shù)cD1。這2個系數(shù)都是原信號在濾波器作用下以尺度為2的下采樣。一維小波分解算法見圖1。

圖1 一維小波分解算法

第2步采用同樣做法,把其中的低頻部分作為信號再次分解,以3層小波分解為例,圖形表示見圖2。

圖2 3層小波分解示意圖

首先對信號s進行3層小波分解,cA為近似系數(shù),cD為細節(jié)系數(shù)。噪聲信號通過包含在細節(jié)系數(shù)cD中,通過門限閾值對細節(jié)系數(shù)進行處理;對信號進行重構(gòu)可以達到降噪的目的。對信號s(t)降噪的目的就是抑制信號的噪聲,從而在s(t)中恢復(fù)真實信號f(t)。

重構(gòu)算法是小波變分解算法的逆變換,把分解得到的近似系數(shù)和細節(jié)系數(shù)通過低通濾波器和高通濾波器重構(gòu)得到原信號。

2 小波降噪算法判定準則和模型

2.1 小波降噪算法判定準則

測井信號的降噪效果可以通過信噪比、相關(guān)系數(shù)、能量比和均方差來衡量和評價[3-5]。

(1) 信噪比。信號信噪比RSN定義為

(4)

(2) 相關(guān)系數(shù)。用于計算小波函數(shù)與原測井信號波形的相似程度,相似程度越高,則對原數(shù)據(jù)信號進行逼近的時候產(chǎn)生的加權(quán)求和項就越小。ρxy為小波函數(shù)x與信號y的相似系數(shù)

(5)

(3) 能量比。降噪后的測井信號與原測井信號大部分情況下要有相同的光滑性,降噪后的測井信號與原測井信號的能量比越大,光滑性越好。

(4) 均方差。降噪后的測井信號和原測井信號的方差估計應(yīng)該是最壞情況下的方差最小,任意算法的降噪偏差均方差EMS為

(6)

2.2 小波降噪的基本模型

在實際工程應(yīng)用中,通過所分析的信號具有非線性、非平穩(wěn)性和奇異點較多的特點[6]。含噪聲的一維信號模型可表示為

s(t)=f(t)+σ·e(t)

(7)

式中,s(t)為含噪聲的信號;f(t)為有用信號;e(t)為噪聲;σ為噪聲強度。

小波變換降噪方法與步驟見文獻[7]。

3 小波基函數(shù)選取

小波基函數(shù)選取一般原則[8-10]。①正交性,便于數(shù)學(xué)分析和工程應(yīng)用中理解操作;②緊支集,保證有優(yōu)良的時域局部特征,也利于算法的實現(xiàn);③對稱性,關(guān)系到小波的濾波特性是否具有線性相位,這與失真問題密切相關(guān);④平滑性,關(guān)系到頻率分辨率的高低。

4 自然伽馬測井曲線降噪

4.1 信號降噪常用方法分析

4.1.1 通過快速傅里葉變換FFT實現(xiàn)信號降噪

首先對原信號進行快速傅里葉變化FFT,求出原信號的頻譜;根據(jù)得到的信號頻譜,對所需要的頻譜成分保留,對不需要的信號頻譜成分進行抑制;對FFT變換后的頻譜做傅里葉逆變換,得到降噪后的信號。這個過程是對原信號在一定范圍濾波,還原到時域則相當于對信號進行卷積運算。設(shè)原信號為f(t),降噪后的信號為g(t),其傅里葉變換形式為F(ω)和G(ω)。則G(ω)=H(ω)F(ω),其中H(ω)為頻域的濾波器,用于抑制噪聲信號的頻譜。

4.1.2 小波閾值實現(xiàn)信號降噪

小波閾值降噪方法通過對小波分解后的各層系數(shù)中模大于閾值和小于某閾值的系數(shù)分別處理;對處理完的小波系數(shù)再進行小波重構(gòu)變換,重構(gòu)出降噪后的信號[11-12]。

(1) 使用penalty軟閾值實現(xiàn)信號閾值降噪。使用函數(shù)wnoisest獲取噪聲方差;使用函數(shù)wbmpen獲取小波去噪閾值;使用函數(shù)wdencmp實現(xiàn)信號去噪。

sigma=wnoisest(c,l,1)

thr3=wbmpen(c,l,sigma,2)

xd3=wdencmp(′gbl′,c,l,wname,lev,thr3,′s′,1);

(2) 使用Birge-Massart分層硬閾值實現(xiàn)信號閾值降噪。使用函數(shù)wdcbm獲取小波去噪閾值;使用函數(shù)wdencmp實現(xiàn)信號去噪。

[thr4,nkeep]=wdcbm(c,l,3,)

[xd4,cxd,lxd,perf0,perfl2]=wdencmp(′lvd′,c,l,wname,lev,thr4,′h′);

(3) 使用全局軟閾值實現(xiàn)信號閾值降噪。在Matlab中使用函數(shù)ddencmp獲取信號去噪閾值;采用函數(shù)wdencmp實現(xiàn)信號去噪。

[thr5,sorh,keepapp]=ddencmp(′den′,′wv′,x)

xd5=wdencmp(′gbl′,x,′db4′,4,thr5,′s′,keepapp);

4.1.3 基于樣本估計的自適應(yīng)閾值實現(xiàn)信號降噪

對信號做無偏似然估計,根據(jù)最壞情況下降噪信號與原信號方差最小的原則確定一個統(tǒng)一的閾值,然后截去超出這個閾值的系數(shù)。

(1) 使用rigrsure自適應(yīng)軟閾值降噪

thr6=thselect(x,′rigrsure′)

xd6=wdencmp(′gbl′,c,l,′db4′,4,thr6,′s′,1);

(2) 使用sqtwolog自適應(yīng)軟閾值降噪

thr7=thselect(x,′sqtwolog′)

xd7=wdencmp(′gbl′,c,l,′db4′,4,thr7,′s′,1);

(3) 使用heursure自適應(yīng)軟閾值降噪

thr8=thselect(x,′heursure′)

xd8=wdencmp(′gbl′,c,l,′db4′,4,thr8,′s′,1);

(4) 使用minimaxi自適應(yīng)軟閾值降噪

thr9=thselect(x,′minimaxi′)

xd9=wdencmp(′gbl′,c,l,′db4′,4,thr9,′s′,1);

4.2 自然伽馬測井曲線降噪效果分析

分別使用快速傅里葉變換降噪算法1、penalty軟閾值降噪算法2、Birge-Massart分層硬閾值降噪算法3、全局軟閾值降噪算法4、rigrsure自適應(yīng)軟閾值降噪算法5、sqtwolog自適應(yīng)軟閾值降噪算法6、heursure自適應(yīng)軟閾值降噪算法7和minimaxi自適應(yīng)軟閾值降噪算法8對自然伽馬測井曲線進行了降噪,降噪后的信噪比、相關(guān)性、能量比和均方差指標對比見圖3。

圖3 降噪算法指標對比

圖4 自然伽馬測井曲線濾波前后曲線

從圖3可見,在使用小波閾值降噪方法2～8進行噪聲抑制的時候,降噪后的信噪比、相關(guān)系數(shù)和能量比都比通過傅里葉變換降噪方法1的大,降噪后的均方差比通過傅里葉變換降噪方法1的小,說明通過小波閾值降噪方法2～8降噪效果好,降噪后的信號與原信號保持了很好的相似性,光滑性好。在小波閾值降噪方法2～8的7種降噪方法中,降噪方法5和7信噪比高,均方差小,說明降噪方法5和7的降噪結(jié)果最優(yōu)。

自然伽馬測井曲線使用rigrsure自適應(yīng)軟閾值小波降噪算法前后的曲線見圖4。從圖4可見,經(jīng)過該方法降噪后的測井曲線比較光滑,無尖棱毛刺。該方法在對自然伽馬測井曲線降噪中有較好的濾波效果,能有效去除自然伽馬測井曲線中與地層信息無關(guān)的統(tǒng)計漲落和毛刺干擾。

5 結(jié) 論

(1) 小波變換具有熵值低、多分辨率特性、去相關(guān)性和小波基函數(shù)選擇靈活的特性,可根據(jù)測井信號特點和降噪要求選擇適合的小波基函數(shù)。本文選用了db4小波基函數(shù)。

(2) 通過對小波變換算法理論分析,建立了小波降噪算法的判定準則:信噪比越大,相關(guān)系數(shù)越接近于1,降噪后的測井信號與原測井信號的能量比越大,降噪后的測井信號與原測井信號的均方差越小,相似性越好,降噪結(jié)果越接近原測井信號,光滑性越好,降噪效果越好。

(3) 通過快速傅里葉變換FFT降噪方法和小波閾值降噪方法對測井信號進行了降噪進行了分析,采用rigrsure自適應(yīng)軟閾值對自然伽馬曲線降噪,可以得到最優(yōu)的降噪性能。

參考文獻:

[1] 鐘雷文. 基于Matlab的煤田核測井信號小波濾波方法研究 [J]. 科學(xué)技術(shù)與工程, 2009, 9(2): 398-400.

[2] 冉啟文. 小波變換與分數(shù)傅里葉變換理論及應(yīng)用 [M]. 哈爾濱: 哈爾濱工業(yè)大學(xué)出版社, 2001.

[3] 鄭一, 孫曉峰, 陳健, 等. 基于集合經(jīng)驗?zāi)B(tài)的隨鉆脈沖信號優(yōu)良降噪算法 [J]. 石油勘探與開發(fā), 2012, 39(6): 750-753.

[4] 劉志松. 基于小波分析的信號去噪方法 [J]. 浙江海洋學(xué)院學(xué)報(自然科學(xué)版), 2011, 30(2): 150-153.

[5] 吳雅娟, 高興, 王輝, 等. 改進的小波閾值法在測井曲線去噪中的應(yīng)用 [J]. 計算機系統(tǒng)應(yīng)用, 2013, 22(3): 182-185.

[6] 李宏, 方世輝, 李族, 等. 泥漿脈沖信號的小波檢測方法研究及應(yīng)用 [J]. 長江大學(xué)學(xué)報(自然科學(xué)版), 2010, 7(1): 68-71.

[7] 殷雷. 基于小波變換的隨鉆數(shù)據(jù)去噪算法研究 [J]. 科學(xué)技術(shù)與工程, 2012, 12(32): 8663-8665.

[8] 汪新凡. 小波基選擇及其優(yōu)化 [J]. 株洲工學(xué)院學(xué)報, 2013, 17(5): 33-35.

[9] 林木, 劉明哲, 度先國, 等. 基于MATLAB GUI的γ能譜小波去噪方法對比研究 [J]. 核電子學(xué)與探測技術(shù), 2013, 33(10): 1245-1248.

[10] 呂瑞蘭, 吳鐵軍, 于玲. 采用不同小波母函數(shù)的閾值去噪方法性能分析 [J]. 光譜學(xué)與光譜分析, 2004, 24(7): 826-829.

[11] 周偉. 基于Matlab小波分析應(yīng)用 [M]. 西安: 西安電子科技大學(xué)出版社, 2011.

[12] 董長虹, 高志, 余嘯海. Matlab小波分析工具箱原理與應(yīng)用 [M]. 北京: 國防工業(yè)出版社, 2006.