石戰(zhàn)戰(zhàn) ,夏艷晴 ,周懷來 ,王元君
(1.成都理工大學(xué)工程技術(shù)學(xué)院,四川樂山614000;2.成都理工大學(xué)地球物理學(xué)院,成都610059)
隨機噪聲壓制是地震資料處理的關(guān)鍵環(huán)節(jié)之一[1]。Anvari等[2]將常用的地震信號去噪方法分為4類:①預(yù)測濾波。基于時空域或變換域中有效信號是可預(yù)測的,而隨機噪聲不可預(yù)測,如t-x域預(yù)測濾波[3]、f-x 域反褶積[4]。②中值濾波。該類方法以信號和隨機噪聲的統(tǒng)計規(guī)律不同為基礎(chǔ),包括時變中值濾波[5]和矢量中值濾波[6]等方法。③基于信號相干性的隨機噪聲壓制方法,如各向異性擴散濾波[7],字典學(xué)習(xí)[8]等。④基于信號分解的去噪方法,如經(jīng)驗?zāi)B(tài)分解[9]、奇異值分解[10]、稀疏表示[11-12]等。其中稀疏表示是一種現(xiàn)行有效的隨機噪聲壓制方法,該方法認為地震數(shù)據(jù)是可壓縮的,在某一字典下有效信號能夠獲得最佳的稀疏表示,而隨機噪聲與有效信號不同源,同一字典下不能實現(xiàn)有效的稀疏表示。因此,字典的有效性就決定了稀疏表示算法的去噪能力。根據(jù)構(gòu)建方法的不同,常用的字典分為:①分析型字典。由數(shù)學(xué)基函數(shù)構(gòu)建而成,如Fourier字典、Morlet小波字典、Curvelet字典和Ricker子波字典等,實際應(yīng)用中,這類字典常難以精確描述地震信號特征。②字典學(xué)習(xí)。這類方法是數(shù)據(jù)驅(qū)動下的地震去噪,優(yōu)點是具有自適應(yīng)性,但該類方法缺點是計算成本高。③地震物理子波構(gòu)建子波字典。該方法由地震資料提取統(tǒng)計子波或結(jié)合井震資料提取子波,構(gòu)建字典,能夠有效描述地震信號特征[12],物理意義明確,且計算效率高。
稀疏表示算法的核心思想是使用少量的基函數(shù)就能實現(xiàn)信號的有效表示,揭示信號的本質(zhì)特征[13]。被廣泛應(yīng)用于地震數(shù)據(jù)分析,如地震數(shù)據(jù)采集[14]、去噪[12]、重建[15]、最小二乘偏移[16]和反射系數(shù)反演[17]等領(lǐng)域。受壓縮感知[18-19]、線性表示及相關(guān)數(shù)學(xué)理論研究成果的促進,國內(nèi)外學(xué)者[20-22]提出了一系列稀疏表示算法。Zhang等[13]將其分為為4類:貪婪算法[20]、約束優(yōu)化算法[21-22]、同倫算法[23]和臨近算法[24]。這些算法均以單測量模型(Single Measurement Vector,SMV)為基礎(chǔ),多采用L1范數(shù)正則項,僅適用于一維信號稀疏表示[25]。受陣列信號處理需求推動,Cotter等[25]提出了基于多測量模型(Multiple Measurement Vector,MMV)的聯(lián)合稀疏表示(Joint Sparse Representation),采用 L2,1范數(shù)正則項,在同一支撐下,參與計算的各道信號同時獲得最優(yōu)稀疏表示,并指出采用MMV模型能夠放松解唯一性條件的上界。Eldar等[26]證明重建精度與測量矢量數(shù)成指數(shù)相關(guān)。Jin等[27]也證實了MMV模型能夠顯著提高信號重建精度。然而,聯(lián)合稀疏表示是一個非確定性多項式(Nondeterministic Polynomial-time,NP)問題,計算效率是制約該算法應(yīng)用于地球物理數(shù)據(jù)處理的瓶頸。交替方向乘子法(Alternating Direction Method of Multipliers,ADMM)是公認為計算效率最高的優(yōu)化方法,能夠適用于大尺度數(shù)據(jù)優(yōu)化處理[21]。Yang 等[21]證明了 ADMM 算法的收斂性,提出一種基于ADMM的聯(lián)合稀疏表示算法,并證實具有優(yōu)越的計算效率。
傳統(tǒng)基于稀疏表示的地震數(shù)據(jù)去噪方法均采用一維算法逐道處理,無法有效利用信號的空間相干性,道集內(nèi)各地震道的去噪處理是相互獨立的,處理結(jié)果一致性較差,去噪的同時會損害有效波,而疊前共偏移距道集中各波形均為水平同相軸,滿足聯(lián)合稀疏表示的共稀疏性(Common Sparsity)假設(shè)條件,可采用聯(lián)合稀疏表示對該道集中各地震道同時進行稀疏分解。針對上述問題,提出基于聯(lián)合稀疏表示的共偏移距道集隨機噪聲壓制方法,利用地震資料提取子波作為基函數(shù),構(gòu)造子波字典,再利用基于ADMM的聯(lián)合稀疏表示對共偏移距道集進行去噪處理,以期兼顧有效信號的道間相干性和空間變化,達到壓制噪聲的同時還能夠保持有效信號不被損害的目的。
隨著重點勘探領(lǐng)域由構(gòu)造油氣藏向深層巖性、地層等復(fù)雜油氣藏轉(zhuǎn)變,勘探開發(fā)研究更加依賴高保真地震數(shù)據(jù),對隨機噪聲壓制處理提出了更高的要求。傳統(tǒng)的一維稀疏表示去噪方法受多解性和單道處理方法局限性的制約,難以兼顧地震道間相干性和信號空間變化。共偏移距道集具有水平同相軸結(jié)構(gòu),表現(xiàn)出良好的空間相干性,滿足共稀疏性假設(shè)條件,而隨機噪聲空間分布不規(guī)律,將聯(lián)合稀疏表示用于疊前共偏移距道集就能提高去噪處理的有效性。
稀疏表示去噪方法以層狀介質(zhì)模型為基礎(chǔ),地震記錄滿足褶積模型
式中:s為地震道;r為地下介質(zhì)的單位沖擊響應(yīng);W為子波核矩陣(1.3節(jié)將詳述其構(gòu)建方法)。
有效地震信號和隨機噪聲由不同的場源產(chǎn)生,因此有效信號能夠被子波核矩陣W充分表示,而噪聲與反射信號不同源,在同一組基函數(shù)下不能獲得有效表示,成為稀疏表示殘差。如果以W為字典,利用稀疏表示算法對地震道s進行稀疏表示,再由稀疏表示結(jié)果r重構(gòu)出地震記錄s*,就能實現(xiàn)信號和噪聲分離。傳統(tǒng)算法常采用基追蹤(Basis Pursuit,BP)對地震道進行稀疏分解,目標(biāo)函數(shù)為
式中:λ >0為參數(shù);‖?‖1為L1范數(shù);‖?‖2為L2范數(shù)。
目標(biāo)函數(shù)式(2)采用L1范數(shù)正則項實現(xiàn)稀疏約束,反演出單位沖激序列r由少量的大反射脈沖和呈Gaussian分布的小反射脈沖組成,這也與實際的地層結(jié)構(gòu)相吻合。該方法為一維去噪算法,逐道處理地震數(shù)據(jù),沒有利用有效信號的空間相干性。同時算法受多解性影響,稀疏表示結(jié)果只是滿足條件的一個實例,相鄰地震道較小的變化可能會引起稀疏表示結(jié)果和重構(gòu)信號較大的區(qū)別,造成處理結(jié)果一致性較差,去噪同時會損害有效波。
針對一維算法的不足,受陣列信號處理啟發(fā),將聯(lián)合稀疏表示算法用于共偏移距道集地震數(shù)據(jù)去噪處理,其目標(biāo)函數(shù)為
式中:S為共偏移距道集,每一列表示一個地震道;R為單位沖擊響應(yīng)矩陣;‖? ‖2,1為 L2,1范數(shù);‖?‖F(xiàn)為Frobenius范數(shù)。
采用L2,1范數(shù)正則項實現(xiàn)共稀疏性約束,反演結(jié)果R僅有少量大反射脈沖行。
當(dāng)信號中混有異常值時,采用Frobenius范數(shù)擬合項會放大異常值的影響,容易造成算法不穩(wěn)定,從而導(dǎo)致算法不收斂或收斂到局部次優(yōu)解。采用L2,1范數(shù)擬合項代替式(3)中的Frobenius范數(shù),構(gòu)造出一種穩(wěn)健的 L2,1-L2,1范數(shù)稀疏表示方法,其目標(biāo)函數(shù)為
L2,1-L2,1范數(shù)稀疏表示的目標(biāo)函數(shù)由一個 L2,1范數(shù)正則項和一個L2,1范數(shù)擬合項構(gòu)成,因此,式(4)是一個非平滑凸函數(shù),難以直接迭代優(yōu)化。求解思路是:通過變量代換和組合,將目標(biāo)函數(shù)式(4)轉(zhuǎn)化為一個等價聯(lián)合基追蹤(Joint Basis Pursuit)問題,然后通過ADMM算法進行迭代優(yōu)化?;贏DMM的 L2,1-L2,1范數(shù)稀疏表示算法描述如下:
可以看出,式(6)是變量可分離的,可用ADMM進行優(yōu)化處理,其增廣Lagrange函數(shù)如下
式中:Λ1和 Λ2均為 Lagrange乘子;β1和 β2均為懲罰參數(shù)。
應(yīng)用ADMM算法優(yōu)化式(7),得到如下迭代公式
式中:Row_Shrink(·,·)為組迭代收縮閾值算子;γ1和γ2均為迭代步長。
字典的優(yōu)劣決定了去噪算法的性能[12],要求字典的原子能夠精確描述地震記錄。由地震資料提取的地震子波具有明確的物理意義,符合地震波場規(guī)律,能夠有效地描述地震信號,滿足稀疏表示方法需要。用地震子波構(gòu)造字典對地震信號進行稀疏分解。
基本的單子波字典構(gòu)建流程如下:
(1)由地震資料提取地震子波w;
(2)構(gòu)造子波字典為一個Toeplitz矩陣,其第(i,j)元素為Wij=wi-j+1(式中:wi-j+1為子波w的第i-j+1個元素);
(3)對字典的每一列做歸一化處理。
地震記錄不是穩(wěn)態(tài)信號,其振幅、相位隨著時間、空間發(fā)生變化,因此地震子波也是時變和空變的。Edgar等[28]、馮瑋等[29]的研究結(jié)果表明,子波的微小變化會對處理、解釋結(jié)果產(chǎn)生較大影響。隨著油氣藏勘探開發(fā)難度的增加,采用單子波構(gòu)造字典常難以滿足生產(chǎn)需求,可以采用時變子波[29]、空變子波[30]構(gòu)造字典,或由多子波字典串聯(lián)[12]構(gòu)造字典。當(dāng)研究薄層儲層預(yù)測時,可以利用雙極子分解[16,31-33]構(gòu)造過完備楔形子波字典對地震記錄進行稀疏分解。
圖1 地質(zhì)模型和地震正演Fig.1 Geological model and seismic forward modeling
通過數(shù)值模擬驗證所提基于聯(lián)合稀疏表示的共偏移距道集隨機噪聲壓制方法的有效性,設(shè)計一個5 層地質(zhì)模型[圖 1(a)],地層速度分別為:2 000 m/s,2 300 m/s,2 200 m/s,2 400 m/s和 2 200 m/s;界面埋深分別為:200 m,400 m,600 m和800 m。地震正演采用30 Hz雷克子波,采用中間放炮兩邊接收觀測系統(tǒng),道間距和最小炮檢距均為50 m,每炮48道接收。從正演地震記錄[圖1(b)]可以看出,受地質(zhì)構(gòu)造和觀測系統(tǒng)影響,疊前共炮點道集中反射波同相軸具有雙曲線特征,不同深度地層反射波同相軸傾角不同。因此,共炮點道集不滿足聯(lián)合稀疏表示的共稀疏性假設(shè),不能直接應(yīng)用聯(lián)合稀疏表示去噪方法。為了驗證算法的噪聲魯棒性和有效性,對正演地震數(shù)據(jù)混入隨機噪聲和異常值,加噪共炮點道集[圖1(c)],信噪比為10 dB。抽取50 m共偏移距道集[圖1(d)],圖中各反射波均為水平同相軸,各道中同一界面的反射波出現(xiàn)時間相同,滿足聯(lián)合稀疏表示的共稀疏性假設(shè)。因此,可以采用聯(lián)合稀疏表示在該道集中進行去噪處理。
圖2展示了基于聯(lián)合稀疏表示的共偏移距道集去噪結(jié)果。實驗中取λ=0.4,其中圖2(a)為處理后50 m共偏移距道集,與圖1(d)對比可以發(fā)現(xiàn),處理后隨機噪聲得到壓制,信噪比提高3.074 dB。圖2(b)為共偏移距道集殘差,濾波殘差呈隨機噪聲分布規(guī)律,說明去除隨機噪聲的同時沒有損害有效波。處理后,再由共偏移距道集抽取共炮點道集,圖 2(c),(d)分別為第 1 炮去噪結(jié)果和濾波殘差,處理后剖面信噪比提高到3.105 dB,隨機噪聲得到有效壓制,殘差剖面呈隨機噪聲分布規(guī)律,無有效波泄露。說明信號與噪聲不同源,以地震子波為基函數(shù),反射信號能夠被有效表示,而噪聲不能。同時,共偏移距道集中,聯(lián)合稀疏表示能夠有效利用信號空間的相干性,提高處理效果。
圖2 共偏移距道集聯(lián)合稀疏表示去噪結(jié)果Fig.2 De-noise based on joint sparse representation in common offset gathers
為了突出共偏移距道集去噪方法的優(yōu)越性,對比了傳統(tǒng)的基于一維信號稀疏表示的單道去噪方法(圖3)。與圖2實驗相同,取λ =0.4。圖3(a)為第1炮經(jīng)傳統(tǒng)的單道算法處理后剖面,雖然處理后隨機噪聲得到壓制,但信噪比提高有限,僅為0.849 dB。圖3(b)為第1炮濾波殘差,剖面上呈現(xiàn)明顯的雙曲線同相軸特征(圖中箭頭所示),說明單道處理算法不能有效利用信號的空間相干性,同時受多解性干擾,去噪同時會損害有效波。
為了進一步說明所提算法的優(yōu)越性,以第1炮第24道為例,對比2種算法處理結(jié)果(圖4)。圖4(a)和(b)分別為第24道正演地震記錄和加噪信號。圖4(c)為傳統(tǒng)算法處理結(jié)果,圖4(d)為共偏移距道集中聯(lián)合稀疏表示算法處理結(jié)果,對比2種算法處理結(jié)果可以看出,2種算法都能有效壓制隨機噪聲,處理后信噪比分別提升0.986 dB和2.781 dB,說明共偏移距道集算法具有更好的去噪效果。圖4(e)為傳統(tǒng)算法濾波殘差,與圖4(a)和(b)對比發(fā)現(xiàn),殘差中出現(xiàn)明顯的反射信號(圖中箭頭所示),與原始信號及加噪信號相似度高,說明有效信號被壓制。圖4(f)為共偏移距道集濾波殘差,呈現(xiàn)隨機噪聲分布特征,說明去噪同時沒有損害有效波。分析結(jié)果說明,所提的共偏移距道集聯(lián)合稀疏表示算法具有明顯優(yōu)勢。
圖3 第1炮一維算法單道處理結(jié)果Fig.3 Processing results of the first shot gather based on 1-D algorithm in trace-by-trace mode
圖4 第24道地震信號2種去噪方法對比分析Fig.4 Comparative analysis of two de-noising methods for the 24 th seismic trace
受采集因素制約,地震信號?;煊须S機噪聲,而稀疏表示方法是一個病態(tài)的逆問題,噪聲和異常值會影響算法的穩(wěn)定性。要求去噪方法具有較好的噪聲魯棒性,適應(yīng)信號的空間變化規(guī)律,去噪同時不損害有效波。為了驗證共偏移距道集中聯(lián)合稀疏表示去噪方法的有效性,以開放數(shù)據(jù)White Sea Line 5(數(shù)據(jù)由俄羅斯Deco Geophysical Software CO.提供)為例進行試算。由地震數(shù)據(jù)提取統(tǒng)計子波構(gòu)造子波字典,對傳統(tǒng)單道處理和共偏移距道集聯(lián)合稀疏表示2種去噪方法進行對比分析。圖5(a)為第2炮地震數(shù)據(jù),該數(shù)據(jù)受隨機噪聲和涌浪噪聲干擾嚴(yán)重,處理前必須對地震噪聲進行壓制。同時剖面上反射波為雙曲線同相軸,直達波時距曲線為直線,不同地層時距曲線斜率不同,互相干涉疊加,隨機噪聲壓制困難。圖5(b)為18 m 共偏移距道集,剖面各波形均為水平同相軸,各反射波雙程旅行時相同,滿足聯(lián)合稀疏表示共稀疏性假設(shè),信號與噪聲分布特征不同,有利于聯(lián)合稀疏表示算法進行去噪處理。因此,該道集中進行聯(lián)合稀疏表示壓制隨機噪聲,能夠有效挖掘地震數(shù)據(jù)的空間相干性,降低算法的多解性,提高算法的去噪性能。
圖6為共偏移距道集中聯(lián)合稀疏表示去噪結(jié)果。其中圖6(a)為共偏移距道集去噪結(jié)果,處理后信噪比提高4.364 dB。圖6(c)為去噪后第2炮,處理后信噪比提高 3.972 dB,與圖 5(a)和(b)比對可以看出,處理后道集同相軸連續(xù)性得到增強,隨機噪聲和涌浪噪聲被抑制。圖6(b)和(d)分別為共偏移距道集和炮集濾波殘差,剖面符合隨機噪聲分布規(guī)律,說明該方法能夠有效利用信號、噪聲空間分布特征的區(qū)別,多道數(shù)據(jù)同時參與計算,相當(dāng)于增加約束條件,縮小解空間,減小多解性,準(zhǔn)確反演出地下地層單位沖激響應(yīng),再由子波字典重構(gòu)出地震記錄,就能較好地實現(xiàn)信噪分離。
圖5 原始地震數(shù)據(jù)Fig.5 Raw seismic section
圖6 共偏移距道集聯(lián)合稀疏表示去噪方法實例分析Fig.6 Case study of seismic de-noise based on joint sparse representation in common offset gathers
圖7 傳統(tǒng)單道去噪算法實例分析Fig.7 Case study of seismic de-noise based on 1-D algorithm in trace-by-trace mode
圖7為傳統(tǒng)單道去噪算法處理結(jié)果。圖7(a)處理結(jié)果的信噪比提高2.397 dB,與圖5(a)對比可以看出,去噪后剖面上同相軸連續(xù)性得到一定程度的增強,噪聲水平降低。同時,噪聲也降低了單道算法的魯棒性,重建信號存在誤差,如圖中箭頭所示。圖7(b)為濾波殘差,剖面上存在明顯的連續(xù)同相軸(圖中箭頭所示),說明傳統(tǒng)算法去噪同時會損害有效波。其原因是,稀疏表示算法具有多解性,地震噪聲會降低算法魯棒性,使算法難以收斂到最優(yōu)解,子波字典表示地震信號的能力降低。同時,單道處理方法不能利用有效信號的空間相干性對稀疏表示算法進行約束,與聯(lián)合稀疏表示算法相比,其計算結(jié)果的準(zhǔn)確度必然要降低。
圖8以第8道地震信號為例進一步對2種算法進行定性分析。圖8(a)為原始地震道,信號同時受涌浪噪聲和隨機噪聲干擾。由于資料未振幅補償,深層地震信號能量弱,信噪分離困難。圖8(b)和(c)分別為傳統(tǒng)算法、共偏移距道集中聯(lián)合稀疏表示算法去噪結(jié)果,處理后資料信噪比分別提高2.263 dB 和 3.544 dB。圖 8(d)和(e)分別為 2種算法濾波殘差。對比2種算法處理結(jié)果和殘差可以發(fā)現(xiàn),受多解性和單道算法局限性制約,傳統(tǒng)算法去噪的同時會損害有效波,殘差與原始信號及處理后信號相似度較高[圖8(b)和(d)中紅色箭頭所示]。同時單道算法噪聲魯棒性不好,重建信號失真[圖8(b)中黑色箭頭所示],而共偏移距道集中聯(lián)合稀疏表示算法能夠有效利用信號空間相干性,實現(xiàn)有效信噪分離,濾波殘差滿足隨機噪聲分布規(guī)律,與原始信號、處理后信號不相干[圖8(c)和(e)中紅色箭頭所示]。
圖8 第8道信號2種去噪方法對比分析Fig.8 Comparative analysis of two de-noising methods for the 8 th seismic trace
(1)在共偏移距道集中采用聯(lián)合稀疏表示壓制隨機噪聲的方法,是基于共偏移距道集中各波形均為水平同相軸,滿足聯(lián)合稀疏表示的共稀疏性假設(shè)條件。該道集中聯(lián)合稀疏表示算法能夠兼顧地震道間的相干性和信號空間變化,提高去噪性能。
(2)地震信號常同時混有隨機噪聲和異常值,采用L2,1范數(shù)擬合項能夠有效提高算法的魯棒性,使隨機噪聲壓制方法穩(wěn)定可靠。L2,1-L2,1范數(shù)聯(lián)合稀疏表示需要采用變量代換和組合將其轉(zhuǎn)換為等價基追蹤問題,再通過ADMM進行優(yōu)化求解,大規(guī)模數(shù)據(jù)處理時需要考慮其計算效率。
(3)稀疏表示算法去噪能力取決于字典性能。地震資料提取地震子波物理意義明確,滿足地震波場規(guī)律。由地震子波構(gòu)建字典能夠有效描述地震信號,提高去噪能力。當(dāng)研究區(qū)有測井資料時,可以采用井震聯(lián)合提取地震子波。
(4)地震子波是空間變化的,當(dāng)研究區(qū)采集條件或地下地質(zhì)條件存在較大變化時,可以考慮采用時變子波或多子波構(gòu)建子波字典。針對薄儲層研究中反射信號弱、噪聲污染嚴(yán)重的問題,可以采用雙極子分解構(gòu)建楔形子波字典,提高去噪性能。