朱劍兵
(中國石化 勝利油田分公司物探研究院,東營 257022)
眾所周知,地震勘探由于噪音來源多,一些弱信號往往會淹沒在噪音里,常規(guī)方法很難去識別弱信號中的構造、儲層信息。為了有效識別弱信號,前人引入了很多方法,其中包括:基于非線性動力系統(tǒng)的混沌理論[1]、基于獨立分量分析的盲源分離技術[2]、基于奇異值分解(SVD)的方法[3]、基于經驗模式分解(EMD)[4]和基于高階統(tǒng)計量[5]等識別技術。
Aharon[6]提出了KSVD算法;韓文功等[7]提出的SVD分解方法可以較好進行弱信號的識別和去噪。KSVD算法由于需要進行多次SVD分解的操作,因此,耗時過長, Rubinstein[8]在KSVD的基礎上提出了AKSVD算法(Approximate KSVD),但是其效果比較差;Dumitrescud等[9]在KSVD基礎上引入了正則化,得到RKSVD與RAKSVD(正則化近似KSVD)方法,并進行了對比,RAKSVD在效果幾乎與KSVD一致的情況下,大大提高了效率。在模態(tài)分解去噪方面,Huang等[10]在1998年首次提出了經驗模態(tài)分解理論(EMD)。然而EMD也并非完美,仍存有一定的缺陷性,諸如部分模態(tài)之間的混疊現象、端點效應以及篩分準則等;Wu等[11]引入了總體經驗模態(tài)分解(EEMD),對EMD算法進行了的改進,但又出現了新的問題,恢復信號中仍含有一定程度的殘余噪音成分,無法將其消除干凈;Torres[12]對EEMD算法進行了改進,引入的完備總體經驗模態(tài)分解(CEEMD)擁有高精度去噪的優(yōu)勢,但也不能盲目去除高頻模態(tài),其中還保存著很多有效信息;楊濤[13]結合KSVD和CEEMD的特點,利用KSVD對弱信號的識別優(yōu)勢及CEEMD的高效去噪,設計了一個基于KSVD的CEEMD去噪識別方法,去噪效果較好,但計算效率相對較低。筆者將改進的CEEMD(ICEEMD)與正則化近似KSVD(RAKSVD)有機結合起來,形成了基于ICEEMD-RAKVSD的去噪方法,該方法能夠對地震弱信號有效識別,取得了較好的弱信號去噪效果,同時計算效率也得到了提升。
由于KSVD中每次原子更新都要進行一次SVD分解,從而導致了計算效率較低。Rubinstein等[8]提出了一種近似KSVD算法(AKSVD),該算法能夠明顯提高計算效率,但由于引入了近似,SVD分解結果存在一定的誤差,因而影響了分解效果。針對這一問題,Irofti等[14]提出了正則化近似KSVD(RAKSVD)算法,建立嶺估計模型,解決稀疏分解中的病態(tài)問題,從而使SVD分解結果的誤差得到很大程度上的控制,提高了稀疏分解效果。
KSVD字典是學習型字典之一,具有更適合反映數據結構特點的優(yōu)勢,能獲得較好的去噪效果。KSVD算法通過對初始化字典進行逐列更新,用來尋找一個最佳字典的過程,假定訓練的信號為Y∈Rm×n,其中m為信號的行數,n為信號的列數,每列對應于一個訓練樣本。學習型字典模型可以描述為式(1)。
(1)
式中:D∈Rm×p為待學習字典,p為字典的列數;X∈Rp×n為稀疏表示系數;‖·‖0表示l0范數;‖·‖2表示l2范數;T為每列個數的最大值。
首先固定字典D,利用正交匹配追蹤算法(OMP)對已知樣本數據進行分解,獲取初始稀疏分解系數矩陣。
字典更新實施步驟如下:
1)記當前更新的第k列字典原子為dk,X中使用到dk原子的索引集合為Ik。
2)計算殘差:
(2)
式中:xj為x的第j列;Ek為去除dk外的殘差。
(3)
式中:U表示左奇異矩陣;V表示右奇異矩陣;Λ是奇異值合集;r表示非零奇異值的數量;ui表示左奇異矩陣的第i行;σi表示第i個奇異值;νi表示右奇異矩陣的第i列。
4)原子和稀疏系數更新為:
xk=σ1ν1,dk=u1
(4)
Rubinstein等[8]在KSVD的基礎上提出了AKSVD算法,主要是在上述步驟3)中采用近似計算更新原子,AKSVD算法更新后的原子和相關系數分別為式(5)和式(6)。
(5)
(6)
稀疏分解是一種不適定問題,而正則化方法能有效提高分解效果。正則化字典習模型可描述為式(7)。
(7)
式中,μ>0為正則化參數。同樣的,記當前更新的字典原子為dk,X中使用到原子的索引集合為Ik,則公式(2)可以變?yōu)槭?8)。
μ‖Xk,Ik‖2
(8)
由于相關系數里除了第k行外,都為固定項。
由公式(2)則可得到:
(9)
(10)
根據式(10),從KSVD算法的正則化[9]可以得出,正則化約束后的更新原子方法和原算法一致,而正則化參數是應用到了稀疏系數的更新,使系數變?yōu)?/(1+μ)倍,通過給稀疏系數一個權值,每次通過這個權值影響后續(xù)的原子更新,最終完成字典的更新。所以RAKSVD的原子和稀疏系數更新[14]為式(11)和式(12)。
(11)
(12)
圖1為KSVD、AKSVD、RAKSVD三種訓練字典的誤差對比圖,其中訓練集Y的稀疏度為5,字典維數為64*256,正則化參數μ為0.1,整個分解過程重復迭代120次,取平均值得到均方根誤差RMSE(單位無量綱)隨迭代次數的關系。
圖1 誤差對比圖Fig.1 Error comparison
由圖1可以看出,這三種KSVD訓練字典在迭代20次后RMSE基本趨于穩(wěn)定,其中AKSVD訓練效果里相對最差,RAKSVD與KSVD訓練效果接近,經前人研究AKSVD相比于KSVD效率提高,故而RAKSVD在提高效率的同時,也能達到較好的分解效果。
首先,EEMD算法在解決EMD存在的模態(tài)混疊、端點效應等缺陷的同時,也導致了分解的IMF分量存在噪音殘余,而CEEMD通過在不同階段加入不同高斯白噪音進行調整,比例大小根據待分解信號自適應計算得到,同時每個符合條件的模態(tài)分量對應的計算殘差值是固定唯一的。不僅解決了模態(tài)混疊問題,同時使原始信號復原的精度變得更高。
按照Torres[12]所設計的CEEMD方法步驟對信號做處理,原始信號中各個頻率成分的信號,基本上都能夠得到較好的分解,但該算法仍存在可能出現虛假模態(tài)的問題,同時由于高斯白噪聲的加入,分解的模態(tài)中會混入一定程度的噪聲。為了解決上述問題,改進的CEEMD(Improved CEEMD,ICEEMD)算法重新定義了各個模態(tài)局部均值的求取。
假設:Ej(·)為CEEEMD分解得到的第j階模態(tài),s(t)為原始地震信號,ωi為加入的第i個高斯白噪聲,白噪聲加入的比例為εk。則改進后的ICEEMD實現步驟如下:
1)將不同正負成對的白噪聲分別加入原始信號中,合成得到T個新信號(T為偶數),計算集合平均得到IMF1(t)。
(13)
2)令r0(t)=s(t),對k=1、…、K,計算k階殘差rk(t)為式(14)。
rk(t)=rk-1(t)-IMFk(t)
(14)
3)對k階殘差添加不同比例系數的白噪聲:rk(t)+εkEk+1(ωi(t)),經EMD分解得到IMF1(t),計算總體平均并將其作為IMF(k+1)(t)。
(15)
4)重復上述步驟,直到達到迭代條件,得到最終的殘差E(t)為式(16)。
(16)
CEEMD去噪往往是直接去除IMF1或IMF1和IMF2分量,但其中往往包含有效信息,所以,雖能有效去除高頻噪音,但對于弱信號去噪的問題,容易丟掉真正的有效信號,使同相軸不連續(xù)。
相比于CEEMD算法,改進后的CEEMD(ICEEMD)算法可以加入更少的白噪聲個數,就能使信號得到很好的模態(tài)分離,同時,在相同參數的條件下,ICEEMD分解的效果更好,計算效率更高。
KSVD字典稀疏去噪方法是屬于對壓縮感知原理的改進算法,通過訓練樣本求取冗余字典,較好地刻畫出信號特征,從而更加精確地重構原信號,實現噪聲去除和保留有效信號的目標。KSVD字典雖然能將噪音基本壓制和消除,但在信號的波峰波谷處仍會存有殘留一些毛刺。ICEEMD算法不需要提前設置基函數,它是以信號本身的特點作為基礎來自動分解和重構信號,打破了測不準原理的限制,達到了高精度識別弱信號的能力。
具體實現步驟為:
1)利用ICEEMD算法對地震信號進行分解得到固有模態(tài)函數集,根據分解效果,將其劃分為噪聲主導模態(tài)、信號主導模態(tài)和過渡模態(tài)。
2)對于噪聲主導模態(tài),則直接舍棄。
3)對于信號主導模態(tài),將其與原始信號疊加,形成主導信號集,通過RAKSVD方法進行主導信號集的分解和去噪。
4)對于過渡模態(tài),先將其疊加,再進行ICEEMD分解;舍棄噪音模態(tài)后,與原始信號疊加,重構過渡信號集,通過RAKSVD方法進行過渡信號集的分解和去噪。
5)將步驟3)和步驟4)中的兩個分解結果重構并平均,得到最終的去噪結果。
含有弱信號模型如圖2(a)所示,模型共有25道,每道300個采樣點,Ricker子波主頻25 Hz。在200 ms處有一弱信號。圖2(b)為含噪模型,加入噪聲后的峰值信噪比(PSNR)為1.2,圖2中強軸能夠明顯看到,但弱軸在噪聲的埋沒中已經不明顯了。對圖2(b)分別采用五種算法進行去噪處理,得到的結果如圖3所示,圖4為模型的15道的單道方法去噪結果,表1是各方法的去噪效果與時間對比。
表1 五種方法針對模型數據的去噪效果對比Tab.1 Comparison of the effects of five methods for denoising the model data
圖3 五種方法針對模型數據的去噪結果Fig.3 The result of five methods for denoising the model data(a)KSVD去噪;(b)RAKSVD去噪;(c)CEEMD去噪;(d)CEEMD-KSVD去噪;(e)ICEEMD-RAKSVD去噪
圖4 單道去噪方法結果對比Fig.4 Comparison of single channel denoising results
從圖3、圖4和表1可以看到,如果CEEMD單獨使用,則對弱信號的識別和去噪效果最差;KSVD和RAKSVD對模型中弱信號的識別和去噪效果較好,說明該類方法相比于CEEMD,對弱信號識別有著獨特的優(yōu)勢,但仍有一些毛刺的殘留,在實際應用中,當地下介質不均勻或由薄互層形成復波時,其缺點會被放大;CEEMD-KSVD和ICEEMD-RAKSVD是將模態(tài)分解和字典訓練聯合在一起進行去噪,將模態(tài)分解的信號分解為主導模態(tài)和過渡模態(tài),利用KSVD和RAKSVD對弱信號識別的優(yōu)勢,最終使含噪弱信號中的噪聲得到了較好的壓制,且保留了有效弱信號,ICEEMD-RAKSVD方法能夠達到最好的去噪效果。由表1可以看到,由于CEEMD算法需要多次迭代,每次迭代需要加入大量成對的高斯白噪聲,因此,其計算效率低。KSVD和RAKSVD計算時間少,RAKSVD較KSVD效率提高了60%。由于CEEMD的引入,CEEMD-KSVD和ICEEMD-RAKSVD算法明顯增加了計算時間,ICEEMD-RAKSVD 相對于CEEMD-KSVD效率提高6%,對于簡單的稀疏模型,CEEMD的計算時間占主導,因此,計算效率提升幅度不是很大。
針對實際資料,對本方法的去噪效果進行了進一步測試。圖5(a)為原始局部疊前地震道集,圖5中隨機噪聲十分發(fā)育,信噪比非常低,其中弱同相軸不清晰、不連續(xù)。上述五種方法的去噪結果如圖5(b)~圖5(f)所示,圖6為單道去噪結果對比展示,圖7為去噪后殘差剖面效果對比,去噪效果與時間對比見表2。
表2 五種方法針對實際資料的去噪效果對比Tab.2 Comparison of the effects of five methods for denoising the actual data
圖5 五種方法針對實際資料的去噪結果Fig.5 The result of five methods for denoising the actual data(a)原始疊前地震道集;(b)KSVD去噪;(c)RAKSVD去噪;(d)CEEMD去噪;(e)CEEMD-KSVD去噪;(f)ICEEMD-RAKSVD去噪
圖6 單道去噪方法結果對比Fig.6 Comparison of single channel denoising results
從圖5和表2可以看出,對于弱地震信號,實際資料的處理結果與模型測試得到的結論相一致:CEEMD去噪效果較差,去噪后PSNR最低,將部分弱同相軸的有效信號充當噪音去除,弱同相軸損失嚴重;RAKSVD與KSVD稀疏去噪較明顯,但仍存留了許多噪音。CEEMD-KSVD和ICEEMD-RAKSVD算法都能夠很好地壓制噪聲,紅方框標記可以看出ICEEMD-RAKSVD有效弱信號得以保留。從圖5(a)可以更為清晰看出五種方法對比,而圖6可以得出 CEEMD-KSVD和ICEEMD-RAKSVD算法去噪效果幾乎一致,能有效去除噪音,圖7藍圓圈中可以看出,CEEMD-KSVD去除的殘余還仍然保留些有效信號,而本文方法去除的噪音沒有保留有效信號,去噪效果更好,同時由于實際地震數據非常復雜,訓練字典的時間大大增加,從表2可以看出,RAKSVD較KSVD計算時間縮短了65%,下降幅度較為明顯,ICEEMD-RAKSVD相對于CEEMD-KSVD時間縮短了40%。對實際地震資料的去噪處理,更能突顯本文方法在提高計算效率方面的優(yōu)勢。
圖7 去噪后殘差剖面效果對比Fig.7 Comparison of residual profile effects after denoising(a)CEEMD_KSVD;(b)ICEEMD_RAKSVD
由于弱地震信息往往受噪音影響更大,通過去噪能更好地識別弱信號。CEEMD去噪方法直接舍棄高頻模態(tài),容易造成高頻中有效信息的丟失,而完整地保留了其他含有噪音模態(tài),KSVD去噪方法和RAKSVD去噪方法,對于弱信號的識別效果較好,但去噪后殘留毛刺,其中RAKSVD計算效率得到明顯提高,CEEMD-KSVD去噪方法能達到較好的去噪效果,但計算效率相對較低,這里指出的基于ICEEMD-RAKSVD去噪方法能在保證較好的去噪效果情況下,還能有效提升了計算效率,在實際地震處理中提升效率尤為重要,該方法具有一定的應用價值。