趙椏松，許輝群，王澤峰，聶 榮，楊夢瓊，李欣怡，魏文齋

(1.長江大學 地球物理與石油資源學院，湖北 武漢 430100;2.遼河油田公司 錦州采油廠,遼寧 凌海 121209)

1 引 言

地震數(shù)據(jù)在野外采集中會受到各種噪聲的影響，提高地震資料信噪比是地震資料處理中的基本問題，因此隨機噪聲處理在地震數(shù)據(jù)處理和解釋起著重要的作用[1,2]。目前常見的信號降噪方法有f-x域預測濾波[3]，小波變換[4,5]，曲波變換[6]，基于經(jīng)驗模態(tài)分解(Empirical Mode Decomposition，EMD)，集合經(jīng)驗模態(tài)分解(Ensemble Empirical Mode Decomposition，EEMD)，自適應噪聲完備集合經(jīng)驗模態(tài)分解(Complete Ensemble Empirical Mode Decomposition with Adaptive Noise，CEEMDAN)與基于小波和CEEMDAN的地震信號去噪等方法[7]，取得了較好的應用效果。EMD是黃鍔提出的一種針對非線性、非平穩(wěn)信號序列處理的一種分析方法，該方法是依據(jù)數(shù)據(jù)自身的時間尺度特征來進行信號分解，具有自適應性[8]。但是，EMD在分解過程中存在模態(tài)混疊的現(xiàn)象，不同時間尺度成分出現(xiàn)在同一特征模態(tài)函數(shù)中[9]；EEMD可以解決模態(tài)混疊，這是一種加入白噪聲輔助的數(shù)據(jù)分析方法，該方法在信號分解過程中白噪聲沒有完成消除[10,11]。為了解決這些問題，提出了CEEMDAN的方法[12]。在此方法上，提出了改進的帶有自適應白噪聲的完全集合經(jīng)驗模態(tài)分解(Improved Complete Ensemble Empirical Mode Decomposition with Adaptive Noise，ICEEMDAN)。針對地震信號噪聲處理問題，使用ICEEMDAN與曲波閾值相結(jié)合的方法，對信號進行降噪。曲波變換具有多尺度性、多方向性和各向異性特征[13]，本文利用 ICEEMDAN將非線性非平穩(wěn)的地震信號分解為不同頻率的IMF(IMF, Intrinsic Mode Function,本征模函數(shù))分量，并按照從高頻到低頻的順序依次排列，根據(jù)分解的IMF分量與原信號的相關(guān)性判斷噪聲，利用曲波閾值進行濾波，再對其各個IMF進行重構(gòu)并求出信噪比和均方差，判斷其去噪效果，通過與實際地震數(shù)據(jù)信號去噪進行對比分析，取得較好的效果。

2 方法原理

2.1 ICEEMADN原理

為了解決EEMD等方法出現(xiàn)的問題，Torres等提出了CEEMDAN的方法。該方法在分解的時候在每個階段都加入一種特定的噪聲，然后在每個階段計算并得到唯一的IMF和相應的余項。從而解決了EEMD分解加噪信號產(chǎn)生的不同數(shù)量的IMF，無法精確進行重構(gòu)、計算效率低等問題，由于CEEMDAN在實際分解的初期存在一些剩余的噪聲，存在虛假模式。進而在此方法上進一步提出了ICEEMDAN[14,15]，它利用在分解m層IMF的時候加入特殊的噪聲Em[w(i)]，第一次添加噪聲和分析信號之間所期望的信噪比(SNR，Signal to Noise Ratio)的對等。在分解的后期獲得振幅較小的噪聲，在其余模式中，使用EMD預處理產(chǎn)生的噪聲，即不通過標準偏差使其正?；?。其分解原理如下：

1)在原始信號y加入白噪聲E1[w(i)]，得到

y(i)=y+β0E1[w(i)]

(1)

其中，w(i)是添加的第i個白噪聲；

2)ICEEMDAN的第一階分量為：

(2)

其中，M(·)為局部均值函數(shù)，j是加入的白噪聲次數(shù)，x是原始信號，r1是余項。

3)第二階IMF分量為：

(3)

4)接下來計算第m個IMF分量：

(4)

其中，m=2,3,…,N。

2.2 曲波變換原理

1999年Candès和Donoho在脊波變換的基礎(chǔ)上提出了曲波變換(Curvelet)的方法，它繼承和發(fā)展了小波變換和脊波變換的理論[16,17]。小波變換是一種具有較強時、頻局部分析功能的非平穩(wěn)信號分析方法[18],然而曲波變換是一種多分辨、帶通、具有方向性的分析方法,因此在表達圖像中的曲線時明顯優(yōu)于小波變換。為了改善第一代Curvelet算法的速度減少冗余度，Candès等人于2005年又提出了實現(xiàn)更簡單、更便于理解的第二代Curvelet變換算法[19，20]，離散Curvelet變換公式如下：

(5)

(6)

2.3 基于ICEEMDAN的曲波閾值去噪方法

有效的地震數(shù)據(jù)與隨機噪聲之和可以表示為含隨機噪聲的地震數(shù)據(jù)：

x+n=y

(7)

其中，x表示有效的地震信號，y表示初始數(shù)據(jù)，n為隨機噪聲。根據(jù)稀疏去噪的原理：把y當成觀測數(shù)據(jù)，是可稀疏的，把n當成不可稀疏的，把觀測數(shù)據(jù)去除系數(shù)進行重構(gòu)，噪聲的處理為觀測數(shù)據(jù)與重構(gòu)數(shù)據(jù)的殘差，在重構(gòu)中去除，達到去噪的效果，則公式表達為：

(8)

這就是求最優(yōu)化因子問題：

(9)

由于

(10)

其中,λ>0是正則參數(shù)。公式(10)基于Ne下降的方向迭代將逐漸達到最優(yōu)，公式如下：

xk+1=Tε[xk-kg′(xk+ψΤ(ψxk-y)]

(11)

其中，xk+1是k+1次迭代結(jié)果，k是迭代步長，Tε為閾值，ψΤ為曲波變換因子。

本文方法基本流程：通過對數(shù)據(jù)進行ICEEMDAN分解，將所得的每個模態(tài)分量與原信號的相關(guān)性進行求解，然后將曲波閾值進行去噪處理，并將處理后的模態(tài)分量進行重構(gòu)，得到去噪后的數(shù)據(jù)。

3 模型測試

通過對一維數(shù)據(jù)測試，對EMD和ICEEMDAN這兩種方法進行對比分析，驗證ICEEMDAN的優(yōu)勢所在。對合成的地震記錄加入隨機噪聲，通過EMD、小波變換設(shè)置閾值等方法進行去噪，與本文基于ICEEMDAN的曲波閾值地震數(shù)據(jù)去噪方法進行對比分析，得出相關(guān)結(jié)論。

由圖1(a)和圖1(b)對比分析，EMD分解存在模態(tài)混疊現(xiàn)象，通過ICEEMDAN的分解可以解決模態(tài)混疊的問題，使每個分量能夠保留更好的局部特征。為了進一步驗證ICEEMDAN的優(yōu)勢，利用頻域的方法進行測試，由圖1(c)和圖1(d)對比分析，在圓圈標紅的位置可以明顯的發(fā)現(xiàn)，ICEEMDAN可以較好地改善模態(tài)混疊現(xiàn)象，使特征得到更好的刻畫。

圖1 正弦信號Fig.1 Sinusoidal signal

通過對合成的地震記錄，加入隨機噪聲進行去噪處理，圖2(a)是原始數(shù)據(jù)，圖2(b)是原始數(shù)據(jù)加入隨機噪聲合成的新的數(shù)據(jù)。利用EMD和小波變換設(shè)置閾值進行去噪處理，圖2(c)和圖2(d)分別是上述方法進行去噪得到的效果圖，可以看出這兩種方法都可以壓制隨機噪聲，但是依舊會保留一些隨機噪聲。圖2(f)是采用本論文方法得到的結(jié)果，可以看出隨機噪聲得到很好的壓制，原始數(shù)據(jù)也得到較好的保留。

圖2 正演數(shù)據(jù)Fig.2 Forward data

4 實際數(shù)據(jù)測試

為了進一步說明基于ICEEMDAN的曲波閾值去噪方法的有效性，對含隨機噪聲、采樣點為100、道數(shù)為151、采樣時間為2 ms的實際地震剖面進行去噪，通過對信噪比和均方根差的綜合對比，可以得到信噪比越高，均方根差越小，降噪效果越好。由圖3對比可知，圖3(d)比圖3(b)和圖3(c)的去噪效果更好，紅色圓圈標注的同相軸更連續(xù)。降噪效果如表1所示。

圖3 實際數(shù)據(jù)Fig.3 The actual data

表1 數(shù)據(jù)降噪效果對比

5 結(jié) 語

經(jīng)過正演模型和實際數(shù)據(jù)根據(jù)信噪比和均方根差的結(jié)果分析表明，相較于EMD方法和小波方法，本文提出的基于ICEEMDAN的曲波閾值地震數(shù)據(jù)去噪方法，效果更加明顯，對地震數(shù)據(jù)噪聲壓制更佳，處理后的數(shù)據(jù)信噪比更高，橫向分辨率更高，為地震資料處理和解釋提供了一種方法。