李 倩,狄?guī)妥?魏建新
(1.中國地震局地球物理勘探中心,河南鄭州 450002;2.中國石油大學(北京)油氣資源與探測國家重點實驗室,北京102249;3.中國石油大學(北京)CNPC物探重點實驗室,北京102249)
基于稀疏約束反演譜分解的縫洞儲層疊后數據去噪應用效果分析
李 倩1,2,3,狄?guī)妥?,3,魏建新2,3
(1.中國地震局地球物理勘探中心,河南鄭州 450002;2.中國石油大學(北京)油氣資源與探測國家重點實驗室,北京102249;3.中國石油大學(北京)CNPC物探重點實驗室,北京102249)
高質量的地震剖面對碳酸鹽巖縫洞儲層的流體識別、儲量預測具有重要作用,因此對疊后地震數據進行去噪處理很有必要。為此,分析研究了基于稀疏約束反演譜分解的去噪方法在碳酸鹽巖縫洞儲層疊后地震數據中的應用效果。首先介紹了該方法的基本原理,然后分別將該方法用于數值模擬數據、物理模擬數據以及實際地震數據的處理,并分析去噪效果。在數值模擬測試中,將該方法與小波閾值去噪法的處理結果進行了比較,表明該方法具有較好的去噪性能和保幅性。物理模擬測試以及實際資料的處理結果均表明,該方法可以較好的壓制碳酸鹽巖縫洞儲層疊后地震數據中的隨機噪聲和偏移噪聲,使地震剖面中碳酸鹽巖孔洞的弱反射特征得到突顯,提高了地震數據處理成果的質量,為后續(xù)的定量解釋奠定了良好的基礎。
去噪;稀疏約束反演譜分解;碳酸鹽巖儲層;疊后地震數據;串珠狀特征;稀疏特性
高質量的地震剖面更有利于碳酸鹽巖縫洞儲層的解釋,因此有必要進行疊后地震數據的去噪處理。常用的地震數據去噪方法有中值濾波法、頻域濾波法、f-k濾波法、多項式擬合法、傅里葉變換法、小波變換法和波動方程法等[1-4]。2010年,BONAR等[5]提出了將時間域數據轉換到稀疏時頻域的一種譜分解方法進行噪聲壓制。當使用該方法來達到去噪目的時,它就相當于基追蹤去噪方法[6]。
目前對基于稀疏變換的去噪方法研究很多[7-11],其基本原理是地震信號在某一轉換域被認為是稀疏的,這就意味著可以用少量高能量的值來代表幾乎整個地震信號的能量,而噪聲在該轉換域被認為不聚焦,因此利用稀疏成分重構信號就可以有效地去除噪聲能量。ELBOTH等[7]等將時頻域去噪方法應用于海洋地震數據處理,取得了良好的效果。RODRIGUEZ等[9]提出了一種基于稀疏約束時頻變換的非相關噪聲衰減方法,并將之用于合成和實際的微地震數據中,去噪效果證明了該方法的優(yōu)越性。HAN等[10-11]在時頻域內地震數據重構方法的基礎上,在指定域內設置一個門檻值,使得地震同相軸比噪聲更集中且具有更強的能量,使用該時頻譜對地震信號進行重構,最終成功降低了隨機噪聲。
地震信號在稀疏約束譜域中是稀疏的,而隨機噪聲與子波庫無關,不會在稀疏譜中產生大的相關,因此使用稀疏約束反演譜對信號進行重構可以壓制地震信號中的隨機噪聲。本文將基于稀疏約束反演譜分解的去噪方法用于碳酸鹽巖縫洞儲層疊后地震數據處理,利用數學模型、物理模型以及實際資料進行驗證,結果表明該方法能夠有效壓制隨機噪聲,減少有效信號的損失。
1.1 基于稀疏約束反演譜分解的去噪方法
通常,地震褶積模型能夠被描述為一個地震記錄由地震子波w(t)和反射系數序列r(t)的褶積再加上隨機噪聲n(t)組成:
(1)
式中:“*”表示褶積運算。地震褶積模型在地震勘探中有重要意義,卻無法表示地震信號頻率隨時間的變化。因此,褶積模型被拓展為一個多子波褶積模型,即地震記錄由一系列不同頻率的地震子波和與之相對應的頻率依賴的偽反射系數序列褶積之后相加再加上隨機噪聲得到,此模型被稱之為非平穩(wěn)地震褶積模型。
(2)
式中:k可以認為與子波主頻線性相關;J表示參與計算的反射系數向量或者子波的總個數;wk(t)是主頻與角標k相對應的地震子波;rk(t)是頻率依賴的反射系數序列。
根據數學理論,方程(2)可以改寫為:
(3)
即:
(4)
其中,Wk和rk(k∈[1J]分別是中心頻率為fi的復雷克子波wk(t)形成的褶積矩陣以及與其有關的反射系數序列。矩陣D表示褶積矩陣庫,由一系列不同頻率的復雷克子波組成,也叫子波庫。m則被認為是包含了所有偽反射系數序列的矩陣。地震記錄由一系列中心頻率唯一確定的子波庫在特定時間被分解為不同的頻率成分,得到頻率依賴的反射系數矩陣m,這個反射系數矩陣就可以被認為是地震記錄的時間—頻率譜[5,11-14]。
如果不考慮隨機噪聲,頻率依賴的反射系數矩陣可以直接通過對D求逆得到,即:
(5)
然而,由于D可能不可逆,用這種方法求解不太可行。另外,由于矩陣m的元素個數遠大于地震記錄的元素個數,因此該方程式是欠定的。為了解這個欠定方程,得到稀疏時頻譜,本文采用了對反問題的解加額外約束項的方式將欠定方程轉化為一個適定方程[13-17]。本文采用L1范數作為約束項來產生稀疏解,最終的目標函數為:
(6)
方程右邊第一項是不適定項,起到了匹配分解結果與地震記錄的作用,而右邊第二項是約束項,確定了反射系數最終分布的形狀。λ(λ>0)是正則化算子,主要起調節(jié)偽反射系數序列稀疏度以及控制解的穩(wěn)定性的作用。
本文采用快速迭代軟閾值算法[16]求得L1范數約束反問題的解m之后,再將m代入方程(4)中,即可得去噪后的地震數據:
(7)
為了測試和分析基于稀疏約束反演譜分解的去噪方法在碳酸鹽巖儲層中的應用效果,首先將該方法用于數值模擬數據中,并與小波閾值去噪方法進行比較,分析該方法的優(yōu)越性。
2.1 數值模擬數據測試
設計一個用于數值模擬的碳酸鹽巖孔洞儲層的地質模型,如圖1a所示。圖中的孔洞寬為100m,高度依次是20,30,40,50,60,80,100,120m,圍巖的速度和密度分別是6000m/s和2.73g/cm3,孔洞的速度和密度分別為4000m/s和2.56g/cm3。經過波動方程正演模擬和數據處理,得到了疊后地震數據。圖1b 為圖1a黑框部分的疊后地震數據。模型各層和孔洞的詳細參數見表1。
圖1 碳酸鹽巖孔洞的地質模型和疊后地震剖面a 地質模型; b 對應圖1a黑框內數據的疊后地震剖面
層名P波速度/(m·s-1)S波速度/(m·s-1)密度/(g·cm-3)層1450025002.50層2600030002.77層3400023002.50孔洞400020002.56
首先從圖1b的疊后地震數據中提取第5個洞中心的單道數據(圖2)進行處理,圖3,圖4和圖5分別對比了不同信噪比情況下稀疏約束反演譜分解去噪與小波閾值去噪結果。由圖2可見,孔洞處的地震反射特征十分明顯。然后在單道數據上分別加上高斯白噪聲使數據信噪比為2.0,1.0和0.5,如圖3a,圖4a以及圖5a所示,較強的噪聲對地震信號產生了影響,容易對之后的解釋工作造成不良影響。對加噪后的單道地震數據進行稀疏約束反演譜分解,由稀疏約束反演譜分解時頻圖(圖3b,圖4b和圖5b)可以看出,地震信號的稀疏分布具有非常好的局部性特征,地震信號的有效部分具有較強的聚焦性,而隨機噪聲由于與子波庫不相關,在時頻圖上很少能夠得到呈現。最后利用公式(7)對稀疏約束反演譜分解結果進行信號重構,得到重構后的地震信號(見圖3c,圖4c和圖5c)。由于該方法只利用少量非零系數進行信號恢復,因此產生的估計信號相對干凈。為了對比,利用小波閾值去噪法分別對圖3a,圖4a和圖5a的含噪信號去噪,得到結果如圖3d,圖4d和圖5d所示。采用的小波是db6小波函數,使用啟發(fā)式閾值原則(即heursure準則)產生閾值,并進行軟閾值處理。對信噪比為2.0的含噪信號進行5層分解;對信噪比為1.0的含噪信號進行6層分解;對信噪比為0.5的含噪信號進行4層分解。當初始數據的信噪比較高時,從圖3和圖4中可以看出基于稀疏約束反演譜分解的去噪方法去噪后孔洞處的反射振幅保持較好,反射振幅得到較好的恢復,僅在部分點有微小的振蕩,而且比較光滑;而小波閾值去噪法去噪后孔洞處的反射振幅保持稍差,振蕩點也較多。當初始數據的信噪比較低時,從圖5可以看出,雖然基于稀疏約束反演譜分解的去噪方法去噪后振蕩點較多,但是對于目標區(qū)域也就是孔洞處的反射振幅恢復較好,信號形變較小;而小波閾值去噪法去噪后雖然振蕩點稍微少一些,但是孔洞處振幅恢復較差,信號形變也較大??偟膩碚f,稀疏約束反演譜分解方法較精確地恢復了單道信號的大部分頻率成分。盡管局部信號的能量稍有降低,孔洞處的反射仍然得到了較好的恢復,證明了該方法對于降低噪聲的有效性。
圖2 從圖1b的疊后地震數據中提取的第5個洞中心的單道數據
圖3 信噪比為2的單道數據的稀疏約束反演譜分解去噪與小波閾值去噪結果對比a 信噪比為2的單道含噪數據; b 稀疏約束反演譜分解的時頻圖; c 從稀疏譜域重構的地震信號; d 小波閾值去噪結果
圖4 信噪比為1的單道數據的稀疏約束反演譜分解去噪與小波閾值去噪結果對比a 信噪比為1單道含噪數據; b 稀疏約束反演譜分解的時頻圖; c 從稀疏譜域重構的地震信號; d 小波閾值去噪結果
信號的信噪比以及均方根誤差通常被作為信號去噪性能的評價標準。信噪比越高,原始信號和估計信號的均方根誤差越小,說明估計信號就越接近原始信號,去噪效果就越好。因此我們給出了基于稀疏約束反演譜分解去噪方法和小波閾值去噪法對不同含噪信號去噪后的信號信噪比(SNR)和均方根誤差(RMSE)(如表2所示)。對比2種方法的信噪比和均方根誤差,發(fā)現基于稀疏約束反演譜分解去噪方法的信噪比較小波閾值去噪法高,而均方根誤差比小波閾值去噪法小,說明基于稀疏約束反演譜分解去噪方法去噪性能較好。
圖5 信噪比為0.5的單道數據的稀疏約束反演譜分解去噪與小波閾值去噪結果對比a 信噪比為0.5的單道含噪數據; b 稀疏約束反演譜分解的時頻圖; c 從稀疏譜域重構的地震信號; d 小波閾值去噪結果
信號類型基于稀疏約束反演譜分解去噪小波閾值去噪SNRRMSESNRRMSE信噪比為2.0的原始信號15.98100.034413.02900.0483信噪比為1.0的原始信號14.57310.040510.03250.0683信噪比為0.5的原始信號9.29110.07438.19640.0843
為了測試該去噪方法在整個地震剖面上的應用效果,在圖1b的地震剖面中提取無噪地震數據,如圖6 所示,然后加上高斯白噪聲使其信噪比為2.0,1.0和0.5,分別利用基于稀疏約束反演譜分解去噪方法和小波閾值去噪法進行處理,結果如圖7,圖8和圖9所示。小波閾值去噪法選擇的小波、閾值、閾值處理方法以及分解層數與處理單道數據時選擇的參數相同。表3中給出了對整個地震數據分別利用2種去噪方法去噪后結果的信噪比和均方根誤差。結合圖7,圖8,圖9和表3進行綜合分析,發(fā)現基于稀疏約束反演譜分解去噪方法壓制噪聲的效果整體上比小波閾值去噪法效果好。在初始數據的信噪比較高時,基于稀疏約束反演譜分解的去噪方法和小波閾值去噪法的去噪效果都較好,2種方法去噪后數據的信噪比和均方根誤差十分接近;隨著初始數據的信噪比降低,基于稀疏約束反演譜分解的去噪方法表現相對穩(wěn)定,去噪后數據的信噪比和均方根誤差值變化較小,而小波閾值去噪法去噪后數據的信噪比和均方根誤差變化較大,去噪效果逐漸變差。在計算過程中,發(fā)現基于稀疏約束反演譜分解的去噪方法雖然壓制噪聲效果較好,但存在計算量大以及在求解稀疏反演譜的過程中正則化算子難以確定的缺點;而小波閾值去噪法計算速度快,但對低信噪比數據的去噪效果稍差,且閾值選擇對去噪效果至關重要,要根據具體情況來選擇合適的閾值??偟膩碚f,不同的去噪方法都有其優(yōu)點以及適用條件,沒有哪一種方法具有普遍適用性。
圖6 從圖1b的疊后地震剖面中提取的無噪地震數據
圖7 信噪比為2.0的數值模擬數據的稀疏約束反演譜分解去噪與小波閾值去噪結果對比a 信噪比為2.0地震數據; b 小波閾值去噪結果; c 基于稀疏約束反演譜分解去噪結果
圖8 信噪比為1.0的數值模擬數據的稀疏約束反演譜分解去噪與小波閾值去噪結果對比a 信噪比為1.0的地震數據; b 小波閾值去噪結果; c 基于稀疏約束反演譜分解去噪結果
圖9 信噪比為0.5的數值模擬數據的稀疏約束反演譜分解去噪與小波閾值去噪結果對比a 信噪比為0.5的地震數據; b 小波閾值去噪結果; c 基于稀疏約束反演譜分解去噪結果
信號類型基于稀疏約束反演譜分解去噪小波閾值去噪SNRRMSESNRRMSE信噪比為2.0的原始信號13.93120.011913.23210.0129信噪比為1.0的原始信號13.82730.012010.54070.0175信噪比為0.5的原始信號13.10890.01306.57240.0277
3.1 物理模型數據測試
將基于稀疏約束反演譜分解的去噪方法應用于物理模型數據,并分析其應用效果。該模型是塔里木油田哈拉哈塘地區(qū)碳酸鹽巖孔洞儲層地震物理模型[18],是根據相似比原理設計的與實際地層比較接近的復雜孔洞模型,其相似比分別是vm∶v=1∶2,lm∶l=1∶20000,fm∶f=10000∶1(其中,vm是模型材料的速度,v是實際地層速度,lm是模型的尺度,l是實際地層的尺度,fm是物理模擬實驗中數據采集時所使用的震源的主頻,f是實際地震數據的主頻)。模型重點考慮了哈拉哈塘地區(qū)主要含油氣層系中奧陶統(tǒng)一間房組至鷹山組一段上部地層,主要具有以下4個特征:①根據該地區(qū)的實際地層情況,目的層上部地層被簡化為水平層狀地層;②目的層之上3個相鄰的上覆層整體從南向北逐漸變薄,形成尖滅;③縫洞儲集體主要分布在深度為6300~7000m的地層中;④200多個不同類型的縫洞儲集體被布置在水平方向上的14個縫洞區(qū)。
圖10為物理模型Inline線的垂向剖面,圖兩側標示了各層的深度,圖中各層中的數字則表示該層的層速度和密度,并在圖中標出了相應的地質層位(斷層和河道未在圖中標出)。
圖11a是利用復雜三維地震物理模型得到的地震剖面,圖11b是去噪后的地震剖面。由模型制作引起孔洞反射之間的規(guī)則干擾以及在剖面最右側的偏移噪聲(如圖11中箭頭所示)在很大程度上被消除,使得孔洞反射更加清晰。
3.2 實際資料應用
將基于稀疏約束反演譜分解去噪方法應用于實際碳酸鹽巖縫洞儲層地震剖面,并分析其應用效果。圖12是目的層在3850ms到4300ms間的一個實際地震剖面,去噪后(如圖12b所示),在疊加剖面上的隨機噪聲和剖面最右側的偏移噪聲(如箭頭所示)得到了很大程度的衰減,目的層段內的構造變的更加清晰,層與層之間的信息更加豐富。剖面中間位置來自于某些孔洞的弱振幅響應和噪聲混雜在一起,去噪之后,“串珠狀”特征清晰可見。圖13為圖12中間部分的放大顯示,在圖13a中看到的是一片雜亂反射,無法判斷是孔洞的弱反射還是干擾,而在圖13b中可以清晰的看到一個個來自于孔洞的串珠狀的弱反射,對后續(xù)的解釋以及儲層預測都能提供更好的依據。
圖10 物理模型Inline線的垂向剖面
圖11 物理模型數據的稀疏譜分解去噪a 去噪前的地震剖面; b 去噪后的地震剖面
圖12 實際數據的稀疏譜分解去噪a 去噪前的地震剖面; b 去噪后的地震剖面
圖13 圖12中間部分的放大顯示a 去噪前的地震剖面; b 去噪后的地震剖面
本文在時頻變換域中對地震信號進行稀疏分解,得到有效信號的稀疏約束反演譜,再通過稀疏約束反演譜對信號進行重構,得到衰減了隨機噪聲的地震記錄。將基于稀疏約束反演譜分解的去噪方法用于碳酸鹽巖縫洞儲層疊后地震數據正演,通過數值模擬數據與小波閾值去噪法相比較,證明該方法具有較好的去噪性能。而物理模型數據及實際地震數據的去噪結果表明,該方法能夠有效地去除地震剖面上的隨機噪聲和偏移噪聲,使得孔洞反射更加清晰,而地震有效信號的能量損失較小,具有一定的保幅性。
致謝:感謝中國石油天然氣股份有限公司塔里木油田分公司研究院對縫洞地震物理模型實驗的支持;感謝中國石油天然氣股份有限公司石油勘探開發(fā)研究院西北分院對縫洞地震物理模型數據處理的幫助!
[1] 劉婷婷,陳陽康.f-x域經驗模式分解與多道奇異譜分析相結合去除隨機噪聲[J].石油物探,2016,55(1):67-75 LIU T T,CHEN Y K.Random noise attenuation based on EMD and MSSA in f-x domain[J].Geophysical Prospecting for Petroleum,2016,55(1):67-75
[2] 劉洋,王典,劉財,等.局部相關加權中值濾波技術及其在疊后隨機噪聲衰減中的應用[J].地球物理學報,2011,54(2):358-367 LIU Y,WANG D,LIU C,et al.Weighted median filter based on local correlation and its application to poststack random noise attenuation[J].Chinese Journal of Geophysics,2011,54(2):358-367
[3] 王姣,李振春,王德營.基于CEEMD的地震數據小波閾值去噪方法研究[J].石油物探,2014,53(2):164-172 WANG J,LI Z C,WANG D Y.A method for wavelet threshold denoising of seismic data based on CEEMD[J].Geophysical Prospecting for Petroleum,2014,53(2):164-172
[4] 胡天躍.地震資料疊前去噪技術的現狀與未來[J].地球物理學進展,2002,17(2):218-223 HU T Y.The current situation and future of seismic data prestack noise attenuation techniques[J].Progress in Geophysics,2002,17(2):218-223
[5] BONAR D C,SACCHI M D.Complex spectral decomposition via inversion strategies[J].Expanded Abstracts of 80thAnnual Internat SEG Mtg,2010:1408-1412
[6] CHEN S B,Donoho D L,SAUNDERS M A.Atomic decomposition by basis pursuit[J].SIAM Journal on Scientific Computing,1998,20(1):33-61
[7] ELBOTH T,PRESTERUD I V,HERMANSEN D.Time-frequency seismic data de-noising[J].Geophysical Prospecting,2010,58(3):441-453
[8] SEJDIC E,DJUROVIC I,JIANG J.Time-frequency feature representation using energy concentration:an overview of recent advance[J].Digital Signal Processing,2009,19(1):153-183
[9] RODRIGUEZ I V,BONAR D,SACCHI M.Microseismic data denoising using a 3C group sparsity constrained time-frequency transform[J].Geophysics,2012,77(2):V21-V29
[10] HAN L,SACCHI M D,HAN L G.Spectral decomposition and de-noising via time-frequency and space-wavenumber reassignment[J].Geophysical Prospecting,2014,62(2):244-257
[11] 韓利.高分辨率全譜分解方法研究[D].長春:吉林大學,2013 HAN L.Research on the methods of high-resolution full spectrum decomposition[D].Changchun:Jilin University,2013
[12] PORTNIAGUINE O,CASTAGNA J.Inverse spectral decomposition[J].Expanded Abstracts of 74thAnnual Internat SEG Mtg,2004:1786-1789
[13] HAN L,HAN L G,ZHAO L.Inverse spectral decomposition with the SPGL1 algorithm[J].Journal of Geophysics and Engineering,2012,9(4):423-427
[14] 劉炳楊,李勝軍,高建虎,等.最小平方約束反演譜分析方法的應用效果分析[J].石油物探,2014,53(5):562-569 LIU B Y,LI S J,GAO J H,et al.The application effects of constrained least-square spectrum analysis method[J].Geophysical Prospecting for Petroleum,2014,53(5):562-569
[15] PURYEAR C L,PORTNIAGUINE O N,COBOS C M,et al.Constrained least-squares spectral analysis:application to seismic data[J].Geophysics,2012,77(5):V143-V167
[16] BECK A,TEBOULLE M.A fast iterative shrinkage-thresholding algorithm for linear inverse problems[J].SIAM Journal on Imaging Sciences,2009,2(1):183-202
[17] GARDNER T J,MAGNASCO M O.Sparse time-frequency representations[J].Proceedings of the National Academy of Sciences of the United States of America,2006,103(16):6094-6099
[18] 李倩,狄?guī)妥?魏建新.碳酸鹽巖儲層孔洞體積的地震物理模擬估算[J].石油地球物理勘探,2014,49(6):1147-1156 LI Q,DI B R,WEI J X.Using seismic physical simulation estimate the pore volume of carbonate reservoir[J].Oil Geophysical Prospecting,2014,49(6):1147-1156
(編輯:朱文杰)
Applicationofdenoisingmethodbasedonsparseconstrainedinversespectraldecompositioninpoststackseismicdataofcave-fracturedreservoirs
LI Qian1,2,3,DI Bangrang2,3,WEI Jianxin2,3
(1.GeophysicalExplorationCenter,ChinaEarthquakeAdministration,Zhengzhou450002,China;2.StateKeyLaboratoryofPetroleumResourceandProspecting,ChinaUniversityofPetroleum,Beijing102249,China;3.CNPCKeyLabofGeophysicsExploration,ChinaUniversityofPetroleum,Beijing102249,China)
High-quality seismic sections are very important for fluid identification and reserve prediction of carbonate reservoirs.Therefore,it is necessary to suppress the noise for poststack seismic data from carbonate reservoirs.The application effect of a denoising method based on sparse constrained inverse spectral decomposition in poststack seismic data from carbonate reservoirs was analyzed.First,the basic principle of the above proposed method was introduced.Then,the proposed method was applied to numerical simulation data,physical modeling data,and field seismic data,and the application effect was analyzed.In the numerical simulation test,the proposed method was compared with the wavelet threshold denoising method,and the results proved that the proposed method has better denoising performance and amplitude preservation than the latter.The application effect of physical modeling data and field seismic data shows that this method can effectively suppress the random noise and migration noise in the stacked seismic section of a carbonate reservoir,highlight the weak reflections from some caves of a carbonate reservoir,improve the quality of the seismic data,and lay a good foundation for the subsequent quantitative interpretation.
denoisng,sparse constrained inverse spectral decomposition,carbonate reservoir,post-stack seismic data,string beads characteristic,sparse property
P631
:A
1000-1441(2017)05-0684-10DOI:10.3969/j.issn.1000-1441.2017.05.009
李倩,狄?guī)妥?魏建新.基于稀疏約束反演譜分解的縫洞儲層疊后數據去噪應用效果分析[J].石油物探,2017,56(5):693
LI Qian,DI Bangrang,WEI Jianxin,et al.Application of denoising method based on sparse constrained inverse spectral decomposition in poststack seismic data of cavefractured reservoirs
[J].Geophysical Prospecting for Petroleum,2017,56(5):693
2016-03-28;改回日期:2016-10-09。
李倩(1985—),女,工程師,博士,主要從事地震資料處理與解釋方面的研究工作。
國家科技重大專項(2011ZX05007-006)項目資助。
This research is financially supported by the National Science and Technology Major Project of China (Grant No.2011ZX05007-006).