王蓉蓉,戴永壽,李 闖,張漫漫,張 鵬

(1.中國(guó)石油大學(xué)(華東)信息與控制工程學(xué)院,山東青島266580;2.中國(guó)石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院,山東青島266580)

基于奇異值分解的時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)方法

王蓉蓉1,戴永壽1,李 闖2,張漫漫1,張 鵬1

(1.中國(guó)石油大學(xué)(華東)信息與控制工程學(xué)院,山東青島266580;2.中國(guó)石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院,山東青島266580)

子波提取準(zhǔn)確性的評(píng)價(jià)在地震數(shù)據(jù)處理中占有重要地位，但是傳統(tǒng)的評(píng)價(jià)準(zhǔn)則受噪聲影響較大。為此，提出一種基于奇異值分解(SVD)的時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)方法，考慮非平穩(wěn)地震記錄的子波提取準(zhǔn)確性評(píng)價(jià)方法中Parsimony準(zhǔn)則、豐度準(zhǔn)則和絕對(duì)峰度準(zhǔn)則對(duì)噪聲環(huán)境的承受能力較強(qiáng)，選用Parsimony準(zhǔn)則與奇異值分解技術(shù)結(jié)合，構(gòu)造了一種抗噪/容噪能力更強(qiáng)的時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)準(zhǔn)則SVD_P。將SVD_P準(zhǔn)則、Parsimony準(zhǔn)則和豐度準(zhǔn)則同時(shí)應(yīng)用于仿真數(shù)據(jù)和實(shí)際資料處理，對(duì)比分析了時(shí)頻域時(shí)變子波提取方法與自適應(yīng)分段時(shí)變子波提取方法的準(zhǔn)確性，結(jié)果表明：SVD_P準(zhǔn)則、Parsimony準(zhǔn)則和豐度準(zhǔn)則都能對(duì)兩類子波提取方法進(jìn)行正確的評(píng)價(jià)，時(shí)頻域子波提取法提取子波的準(zhǔn)確性高于自適應(yīng)分段法提取的子波，但是SVD_P準(zhǔn)則評(píng)價(jià)的結(jié)果相對(duì)誤差最小，評(píng)價(jià)精度最高。

子波提取;準(zhǔn)確性;評(píng)價(jià)準(zhǔn)則;奇異值分解

地震子波提取是地震反演和地震解釋的基礎(chǔ)，其準(zhǔn)確性直接影響后續(xù)處理的結(jié)果[1-2]。由于子波在傳播過程中受地層吸收作用的影響導(dǎo)致高頻成分的缺失和相位特征的畸變，因此地震記錄具有非平穩(wěn)性[3-4]。針對(duì)非平穩(wěn)地震記錄中子波提取準(zhǔn)確性的判別，目前并無成熟有效的評(píng)價(jià)方法，因此，非平穩(wěn)地震記錄中時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)方法的研究對(duì)提高地震資料分辨率具有重要意義。

在非平穩(wěn)地震記錄中提取時(shí)變子波，現(xiàn)有的比較有效的兩種思路是：①將非平穩(wěn)地震記錄劃分為多個(gè)時(shí)窗后，每段記錄視為近似平穩(wěn)地震記錄提取子波；②將非平穩(wěn)地震記錄變換到時(shí)頻域，逐點(diǎn)提取時(shí)變子波。直接通過波形、能量或其它地震子波自身的屬性，都不能對(duì)子波提取精度有一個(gè)精確的評(píng)估[5]，目前較為有效的方法是通過非平穩(wěn)地震記錄反褶積得到反射系數(shù)序列，結(jié)合現(xiàn)有的評(píng)價(jià)準(zhǔn)則對(duì)反射系數(shù)序列進(jìn)行評(píng)價(jià)，進(jìn)而間接評(píng)價(jià)子波提取準(zhǔn)確性。

目前主要評(píng)價(jià)準(zhǔn)則有Parsimony準(zhǔn)則[6]、豐度準(zhǔn)則[7]、絕對(duì)峰度準(zhǔn)則[8]、相似系數(shù)準(zhǔn)則[9]和基于互信息率的判別準(zhǔn)則[10]。但這些評(píng)價(jià)準(zhǔn)則對(duì)改變其信號(hào)時(shí)域特性的噪聲較敏感，在有噪聲污染的信號(hào)中應(yīng)用這些準(zhǔn)則進(jìn)行子波提取準(zhǔn)確性評(píng)價(jià)效果不佳。為了解決這個(gè)問題，本文提出一種基于奇異值分解的子波提取準(zhǔn)確性評(píng)價(jià)方法，結(jié)合反褶積結(jié)果的主奇異值成分和Parsimony準(zhǔn)則，構(gòu)造一種抗噪/容噪能力較強(qiáng)的SVD_P準(zhǔn)則，利用模型分析和實(shí)際地震資料對(duì)方法進(jìn)行了驗(yàn)證。

1 基于奇異值分解的子波提取準(zhǔn)確性評(píng)價(jià)方法

1.1 非平穩(wěn)地震記錄中時(shí)變子波提取

1.1.1 時(shí)頻域時(shí)變子波提取

針對(duì)非平穩(wěn)地震記錄,基于時(shí)頻域的時(shí)變子波提取方法如下:首先應(yīng)用改進(jìn)的廣義S變換將非平穩(wěn)地震記錄變換到時(shí)頻域,在時(shí)頻譜的點(diǎn)譜上擬合子波;然后在頻率域采用譜模擬法獲取每個(gè)時(shí)刻的子波振幅譜，采用基于高階累積量[11-12]的雙譜法預(yù)估計(jì)子波相位譜,結(jié)合局部相似度[13]逐點(diǎn)拾取時(shí)變子波相位譜,進(jìn)而實(shí)現(xiàn)非平穩(wěn)地震記錄混合相位時(shí)變子波的提取。

改進(jìn)的廣義S變換的窗函數(shù)表達(dá)式為[14]：

(1)

式中:q,p為大于0的調(diào)節(jié)因子。從公式(1)中可以看出,窗函數(shù)的寬度與頻率呈正比例變化,即在低頻處可獲得較高的時(shí)間分辨率,在高頻處可獲得較高的頻率分辨率,符合地震記錄動(dòng)態(tài)衰減的特性。利用該方法對(duì)地震記錄進(jìn)行時(shí)頻分析,可獲得更高的時(shí)頻分辨率和很好的時(shí)頻聚焦,能夠有效地分辨出地震資料中頻率成分的變化。利用該方法還可進(jìn)行無能量損失的反變換,準(zhǔn)確重構(gòu)時(shí)間域地震信號(hào)。

1.1.2 自適應(yīng)分段子波提取

對(duì)非平穩(wěn)地震記錄分段處理的核心思想是將地震記錄劃分為多個(gè)時(shí)窗并將每段記錄視為平穩(wěn)地震記錄，然后對(duì)每段記錄提取子波。假設(shè)在同一層段中反射界面上的地震子波是時(shí)不變的,則可認(rèn)為每個(gè)層段產(chǎn)生的地震記錄滿足Robinson模型[15],隨著地層深度的增加,不同層段中的子波是時(shí)變的。

本文采用高靜懷等[16]所提出的自適應(yīng)分子分解法對(duì)非平穩(wěn)地震記錄進(jìn)行分段,通過譜擬合得到每段地震記錄的子波振幅譜,結(jié)合基于高階累積量的雙譜重構(gòu)法提取的相位譜[17-18],完成混合相位子波提取。

1.2 含噪反褶積建模及反褶積處理

1.2.1 含噪反褶積建模

地震記錄含有噪聲的情況下,可以用Robinson褶積模型[19]表示為：

(2)

其中,r(n)為反射系數(shù)序列,w(n)為地震子波,v(n)為隨機(jī)干擾噪聲。上述信號(hào)應(yīng)當(dāng)滿足如下假設(shè):

1) 反射系數(shù)序列為獨(dú)立同分布過程;

2) 環(huán)境噪聲v(n)為隨機(jī)干擾噪聲,且與r(n)統(tǒng)計(jì)獨(dú)立,因此與x(n)也相互獨(dú)立。

(3)

由公式(3)可以看出,由于噪聲的影響,反褶積結(jié)果并不僅僅為原始反射系數(shù)序列,還含有噪聲與逆子波的成分，反褶積結(jié)果的性質(zhì)已經(jīng)偏離原始反射系數(shù)序列,若采用傳統(tǒng)方法對(duì)此結(jié)果進(jìn)行子波提取準(zhǔn)確性評(píng)價(jià),通常得不到準(zhǔn)確的結(jié)果。因此,尋找一種抗噪能力較強(qiáng)的評(píng)價(jià)準(zhǔn)則對(duì)子波提取準(zhǔn)確性評(píng)價(jià)至關(guān)重要。

1.2.2 反褶積處理

由于反褶積是檢驗(yàn)提取子波準(zhǔn)確性的重要工具,其結(jié)果直接影響子波評(píng)價(jià)的準(zhǔn)確性,因此,對(duì)適合于非平穩(wěn)地震記錄的反褶積方法的研究是子波提取準(zhǔn)確性評(píng)價(jià)方法的重要環(huán)節(jié)。

現(xiàn)有的確定性反褶積方法已經(jīng)十分成熟,如最小二乘反褶積、最小熵反褶積[20]、稀疏脈沖反褶積[21]、譜稀疏貝葉斯學(xué)習(xí)反射系數(shù)反演[22]、多分辨率地震信號(hào)反褶積[23]、基追蹤反演[24]等。由于ARMA模型描述下的子波在時(shí)域表現(xiàn)為無限長(zhǎng)脈沖響應(yīng)序列,時(shí)域反褶積處理不可避免地要對(duì)子波進(jìn)行截?cái)?，為了減少截?cái)嗾`差對(duì)反演過程的影響,在提取時(shí)變子波后,我們選擇譜除法在頻率域做反褶積處理。

將Robinson褶積模型進(jìn)行傅里葉變換可得：

(4)

其中,X(ejω),W(ejω)和R(ejω)分別為地震記錄、地震子波和反射系數(shù)序列的頻域表示。則反射系數(shù)序列的估計(jì)可表示為：

(5)

1.3 基于奇異值分解的評(píng)價(jià)方法

SVD濾波[25-26]是通過利用地震記錄的相關(guān)性優(yōu)化提取信息來達(dá)到波場(chǎng)分離與去噪目的的[27-28]。沈鴻雁等[29]指出,將地震記錄特征值或奇異值作為正交基,通過矩陣分解與變換將其轉(zhuǎn)換到新的坐標(biāo)系,利用信號(hào)能量與噪聲的可分性將較大奇異值成分重構(gòu),即可實(shí)現(xiàn)對(duì)隨機(jī)噪聲的消除,增強(qiáng)相干能量。

設(shè)地震記錄為x,采樣點(diǎn)數(shù)為N,x=[x1,x2,…,xn],由此構(gòu)成的分解矩陣Dm可表示為:

(6)

式中:m為矩陣重構(gòu)的維數(shù)；n為時(shí)間延遲量,為每維記錄個(gè)數(shù)。各變量之間的關(guān)系可由(m-1)×τ+n=N確定。時(shí)間延遲量的計(jì)算采用單道地震記錄自相關(guān)歸一化方法求取,可表示為[30]：

(7)

對(duì)分解矩陣Dm進(jìn)行奇異值分解,其奇異值由兩部分組成,即：

(8)

式中:SDM是分解矩陣Dm對(duì)應(yīng)的奇異值;SD為有效信號(hào)對(duì)應(yīng)的奇異值;SV為隨機(jī)干擾對(duì)應(yīng)的奇異值。

采取保留地震信號(hào)對(duì)應(yīng)的奇異值，其它奇異值置零的方法對(duì)矩陣進(jìn)行重構(gòu),從而得到去除隨機(jī)噪聲后的地震信號(hào)。將該方法與傳統(tǒng)的子波提取評(píng)價(jià)準(zhǔn)則結(jié)合，可構(gòu)造有效去除隨機(jī)噪聲的評(píng)價(jià)準(zhǔn)則。對(duì)非平穩(wěn)地震記錄中子波提取準(zhǔn)確性評(píng)價(jià)方法的研究表明,Parsimony準(zhǔn)則、豐度準(zhǔn)則和絕對(duì)峰度準(zhǔn)則對(duì)噪聲環(huán)境的承受能力強(qiáng)。因此本文將Parsimony準(zhǔn)則與奇異值分解技術(shù)結(jié)合,構(gòu)造一種抗噪/容噪能力強(qiáng)的時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)準(zhǔn)則。

王友新等[31]指出,Parsimony的意思是“節(jié)約”或“吝嗇”,表示地震記錄經(jīng)過該準(zhǔn)則處理后可以使結(jié)構(gòu)更簡(jiǎn)單[32],即：

(9)

其中,xi為反褶積結(jié)果,滿足n>2。

新構(gòu)造的評(píng)價(jià)準(zhǔn)則稱之為SVD_P準(zhǔn)則,公式如下:

(10)

式中:xSVD表示對(duì)反褶積結(jié)果進(jìn)行奇異值分解并取主奇異值重構(gòu)反射系數(shù)序列,滿足n>2。當(dāng)只有一個(gè)尖脈沖時(shí),PSVD(x)=0。由于該準(zhǔn)則采用了xi的對(duì)數(shù)形式,因此可以使強(qiáng)反射和弱反射都能夠?qū)χ貥?gòu)結(jié)果產(chǎn)生一定的影響。

將SVD_P準(zhǔn)則與Parsimony準(zhǔn)則和豐度準(zhǔn)則進(jìn)行對(duì)比以驗(yàn)證SVD_P準(zhǔn)則的優(yōu)越性。豐度準(zhǔn)則公式為:

(11)

式中:xi為反褶積結(jié)果。該方差模準(zhǔn)則的最大值代表一個(gè)滿足非高斯性質(zhì)的概率密度分布函數(shù)的峰值,其值越大,概率密度分布形狀越尖銳。

2 實(shí)驗(yàn)分析

為了驗(yàn)證本文提出的基于奇異值分解的時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)準(zhǔn)則SVD_P的有效性,本文采用自回歸滑動(dòng)平均模型(ARMA)描述地震子波并合成地震記錄,通過對(duì)稱映射ARMA模型零極點(diǎn)的方式構(gòu)造因果性及相位性不同的地震子波并進(jìn)行反褶積處理。采用基于奇異值分解的時(shí)變子波提取準(zhǔn)確性評(píng)價(jià)準(zhǔn)則對(duì)反褶積結(jié)果進(jìn)行評(píng)價(jià),辨識(shí)出真實(shí)的反射系數(shù)序列,從而判斷出真實(shí)的地震子波。

采用Margrave等[33]提出的具有明確吸收效應(yīng)的非平穩(wěn)褶積模型構(gòu)造非平穩(wěn)地震記錄,公式為：

s(t)=w(f)αQ(u,f)r(u)exp[2πif(t-u)]dfdu

(12)

式中:s(t)表示地震記錄;w(f)為地震子波;r(u)

為反射系數(shù)序列;αQ(u,f)為衰減因子,可以表示為:

(13)

式中：Q為大地濾波因子。本實(shí)驗(yàn)中采用混合因果、混合相位子波作為真實(shí)子波,如圖1a所示。圖1b 為本文選用的反射系數(shù)序列,該反射系數(shù)序列滿足獨(dú)立同分布IID及Bernoulli-Gaussian分布,采樣間隔為1ms,長(zhǎng)度為1000ms,其與地震子波褶積合成的地震記錄如圖1c所示,構(gòu)造的非平穩(wěn)地震記錄如圖1d所示。

2.1 方法的有效性驗(yàn)證

2.1.1 未加入噪聲時(shí)評(píng)價(jià)準(zhǔn)則的準(zhǔn)確性

以圖1a所示混合相位子波為真實(shí)子波,采用時(shí)域?qū)ΨQ映射ARMA模型零極點(diǎn)的方式構(gòu)造出一系列因果性及相位性不同的地震子波[34],其時(shí)域波形如圖2a所示(第12道為真實(shí)子波)。應(yīng)用構(gòu)造出的不同子波對(duì)合成的非平穩(wěn)地震記錄(圖1d)進(jìn)行反射系數(shù)序列反演,通過譜除法在頻域進(jìn)行反褶積處理,最后將頻域反褶積結(jié)果映射到時(shí)域進(jìn)行比較。16個(gè)不同相位譜子波反褶積處理后的結(jié)果如圖2b所示,其中第12道為真實(shí)子波反演結(jié)果。

圖1 仿真實(shí)驗(yàn)分析

圖2 構(gòu)造的16道子波(a)及反射系數(shù)序列反演結(jié)果(b)

由圖2b可以看出,僅真實(shí)子波能夠完全恢復(fù)真實(shí)的反射系數(shù)序列,其余子波反演結(jié)果均不滿足真實(shí)反射系數(shù)序列的稀疏性質(zhì)。應(yīng)用SVD_P準(zhǔn)則、Parsimony準(zhǔn)則及豐度準(zhǔn)則對(duì)圖2b所示反射系數(shù)序列進(jìn)行評(píng)價(jià),確定出真實(shí)或準(zhǔn)確的反射系數(shù)序列,結(jié)果如圖3所示。

由圖3可知,SVD_P準(zhǔn)則及Parsimony準(zhǔn)則的最小值對(duì)應(yīng)真實(shí)反射系數(shù),而豐度準(zhǔn)則最大值對(duì)應(yīng)真實(shí)反射系數(shù)?？梢娫谖醇尤朐肼暤那闆r下,各種評(píng)價(jià)準(zhǔn)則均能識(shí)別出真實(shí)反射系數(shù)序列,進(jìn)而識(shí)別出真實(shí)地震子波(第12道)。

2.1.2 噪聲對(duì)評(píng)價(jià)準(zhǔn)則的影響

為了得到適用于實(shí)際地震資料,對(duì)噪聲環(huán)境承受能力較強(qiáng)的評(píng)價(jià)準(zhǔn)則,我們將合成的非平穩(wěn)地震記錄(圖1d)加入不同信噪比的高斯白噪聲,采用構(gòu)造的16道地震子波(圖2a)進(jìn)行反褶積處理,并用SVD_P準(zhǔn)則、Parsimony準(zhǔn)則及豐度準(zhǔn)則對(duì)16道反褶積結(jié)果進(jìn)行對(duì)比判別。3種評(píng)價(jià)準(zhǔn)則對(duì)噪聲環(huán)境的承受能力如表1所示。

圖3 不同評(píng)價(jià)準(zhǔn)則評(píng)價(jià)結(jié)果

表1 噪聲對(duì)評(píng)價(jià)準(zhǔn)則的影響

注:“√”表示可準(zhǔn)確判別;“×”表示無法準(zhǔn)確判別。

由表1實(shí)驗(yàn)結(jié)果可以看出,本文所提出的基于奇異值分解的子波提取準(zhǔn)確性評(píng)價(jià)準(zhǔn)則對(duì)所有噪聲強(qiáng)度是有效的,具有良好的抗噪/容噪能力,優(yōu)于Parsimony準(zhǔn)則和豐度準(zhǔn)則。

2.2 時(shí)變子波的提取

2.2.1 時(shí)頻域時(shí)變子波提取

在時(shí)頻域點(diǎn)譜上采用譜模擬法提取子波振幅譜,結(jié)合雙譜法及局部相似度重構(gòu)子波相位譜,提取圖1d中非平穩(wěn)地震記錄的時(shí)變子波。圖4分別為100,500,520,900ms時(shí)原始子波與提取的時(shí)變子波時(shí)域波形對(duì)比圖。

由圖4可知:①在t=100ms時(shí),時(shí)頻域提取的子波時(shí)域波形與理論值比較吻合,隨著時(shí)間的增加,擬合出的子波振幅逐漸減少,能夠反映出子波在傳播過程中的動(dòng)態(tài)衰減特性,說明了子波提取的正確性;②在500ms和520ms處,由于層位相距較近,子波相互疊加難以識(shí)別,但這兩個(gè)點(diǎn)仍有較好的子波提取結(jié)果,說明該方法能夠較準(zhǔn)確地提取出每一反射點(diǎn)處的子波。

2.2.2 基于自適應(yīng)分段提取時(shí)變子波

采用自適應(yīng)分子分解法[16]對(duì)圖1d所示非平穩(wěn)地震記錄自適應(yīng)分段,如圖5所示。圖5a表示非平穩(wěn)地震記錄,紅色曲線表示希爾伯特變換得到的振幅包絡(luò),綠色“*”表示包絡(luò)峰值。取圖5a中兩個(gè)相鄰峰值的中點(diǎn)位置作為每個(gè)初始分子窗的終點(diǎn)位置;將位于相鄰終點(diǎn)內(nèi)的原子窗累加,得到初始分子窗,如圖5b所示;采用自適應(yīng)分子分解法得到自適應(yīng)分子窗(圖5c),將地震記錄分為8段,從左至右分別為1,2,…,8。

將圖5c劃分的每段地震記錄進(jìn)行譜擬合,再進(jìn)行傅里葉反變換后得到每段的子波振幅，采用基于高階累積量的雙譜重構(gòu)法提取子波相位,每段子波提取結(jié)果如圖6所示。在圖6中,藍(lán)色實(shí)線為每段地震記錄估計(jì)得到的子波,紅色虛線為實(shí)際子波。

圖4 原始子波與提取的時(shí)變子波波形對(duì)比

圖5 自適應(yīng)分子分解法分段

由圖6可知:①第1～2段估計(jì)得到的子波振幅和相位與實(shí)際子波基本一致,說明地震波在淺層傳播時(shí)受到的地層吸收作用較小;②隨著深度的增加,地層吸收作用越來越明顯,第3～8段估計(jì)得到的子波時(shí)域波形變寬、頻帶變窄、振幅衰減,但仍能保持子波的基本形態(tài)。子波估計(jì)的結(jié)果真實(shí)反映了地震波在地下的傳播規(guī)律,從側(cè)面證實(shí)了基于自適應(yīng)分段處理的子波提取的正確性。

圖6 自適應(yīng)分子分解分段法提取的每段子波

2.3 反射系數(shù)序列反演

根據(jù)以上兩種方法估計(jì)得到子波后,通過公式(5) 對(duì)非平穩(wěn)地震記錄做反褶積處理。時(shí)頻域時(shí)變子波提取后通過逐點(diǎn)譜除得到反褶積結(jié)果;分段處理方法對(duì)各段地震記錄提取子波后,對(duì)每段記錄分別做反褶積處理。圖7a 為原始反射系數(shù)序列,圖7b為時(shí)頻域反褶積結(jié)果;圖7c為時(shí)頻域反褶積結(jié)果的誤差;圖7d為分段后反褶積結(jié)果;圖7e 為分段后反褶積結(jié)果的誤差。

由圖7可知:①分段后反褶積結(jié)果中噪聲較弱,可以分辨大部分層位,部分薄層及反射系數(shù)較小的層位(箭頭所示位置)難以分辨;②時(shí)頻域反褶積結(jié)果中基本沒有噪聲,且層位清晰,更接近真實(shí)反射系數(shù);③兩種反褶積結(jié)果與真實(shí)反射系數(shù)的差均較小,說明兩種方法估計(jì)得到的子波都能反演出較準(zhǔn)確的反射系數(shù)。

2.4 應(yīng)用評(píng)價(jià)準(zhǔn)則對(duì)反褶積結(jié)果的評(píng)價(jià)

對(duì)圖7中兩種反褶積結(jié)果及真實(shí)反射系數(shù)分別應(yīng)用SVD_P準(zhǔn)則、Parsimony準(zhǔn)則及豐度準(zhǔn)則進(jìn)行評(píng)價(jià),結(jié)果如表2所示。從2.1節(jié)的討論中可知,反褶積結(jié)果越接近真實(shí)反射系數(shù),SVD_P準(zhǔn)則和Parsimony準(zhǔn)則的值越小,豐度準(zhǔn)則的值越大。據(jù)此判斷,表2中帶“*”標(biāo)志的數(shù)據(jù)為最佳結(jié)果。

由表2可以看出,時(shí)頻域反褶積得到的反射系數(shù)序列更加接近真實(shí)值,這與2.3節(jié)反褶積對(duì)比的結(jié)果是一致的。因此,時(shí)頻域時(shí)變子波提取法提取子波的準(zhǔn)確性高于自適應(yīng)分段提取的子波。由表2還可以看出,雖然SVD_P準(zhǔn)則、Parsimony準(zhǔn)則及豐度準(zhǔn)則都能對(duì)兩類子波提取方法進(jìn)行正確的評(píng)價(jià),但是SVD_P準(zhǔn)則評(píng)價(jià)的結(jié)果相對(duì)誤差最小。

圖7 反褶積結(jié)果分析

表2 不同評(píng)價(jià)準(zhǔn)則對(duì)兩類反演結(jié)果的評(píng)價(jià)及其相對(duì)誤差

分類時(shí)頻域反褶積評(píng)價(jià)結(jié)果分段后反褶積評(píng)價(jià)結(jié)果原始反射系數(shù)序列時(shí)頻域反褶積評(píng)價(jià)結(jié)果的誤差分段后反褶積評(píng)價(jià)結(jié)果的誤差SVD＿P準(zhǔn)則3．3729?3．56483．32281．50%?7．28%Parsimony準(zhǔn)則3．5213?3．89833．46261．69%?12．58%豐度準(zhǔn)則0．0281?0．02350．02861．75%?17．83%

3 實(shí)際地震數(shù)據(jù)應(yīng)用效果分析

圖8a為某探區(qū)實(shí)際二維地震剖面,圖8b為第189道地震記錄。分別應(yīng)用時(shí)頻域時(shí)變子波提取法和自適應(yīng)分段子波提取方法提取地震子波,結(jié)果如圖9所示。

由圖9可見:自適應(yīng)分段子波估計(jì)和時(shí)頻域時(shí)變子波估計(jì)兩種方法提取的子波形態(tài)相似,都能反映出子波在傳播過程中的動(dòng)態(tài)衰減特性。為了判斷兩種子波提取方法的準(zhǔn)確性,根據(jù)圖9中的子波分別對(duì)地震記錄做反褶積,得到估計(jì)的反射系數(shù)序列,如圖10所示。分別應(yīng)用SVD_P準(zhǔn)則、Parsimony準(zhǔn)則和豐度準(zhǔn)則對(duì)兩種反褶積結(jié)果和測(cè)井得到的反射系數(shù)序列(圖10c)進(jìn)行了評(píng)價(jià),結(jié)果如表3所示。

根據(jù)2.1節(jié)討論的結(jié)果判斷,表3中帶“*”標(biāo)志的評(píng)價(jià)結(jié)果為最佳結(jié)果，因此時(shí)頻域時(shí)變子波提取法提取的子波準(zhǔn)確性高于自適應(yīng)分段提取的子波。由表3還可以判斷,雖然SVD_P準(zhǔn)則、Parsimony準(zhǔn)則及豐度準(zhǔn)則都能對(duì)兩類子波提取方法進(jìn)行正確的評(píng)價(jià),但是SVD_P準(zhǔn)則評(píng)價(jià)的結(jié)果相對(duì)誤差最小,評(píng)價(jià)精度最高。

圖8 實(shí)際地震剖面(a)與第189道地震記錄(b)

圖9 實(shí)際數(shù)據(jù)子波提取結(jié)果

圖10 實(shí)際數(shù)據(jù)反褶積結(jié)果分析

表3 不同評(píng)價(jià)準(zhǔn)則對(duì)兩類反演結(jié)果的評(píng)價(jià)及其相對(duì)誤差

分類時(shí)頻域反褶積評(píng)價(jià)結(jié)果分段后反褶積評(píng)價(jià)結(jié)果原始反射系數(shù)序列時(shí)頻域反褶積評(píng)價(jià)結(jié)果的誤差分段后反褶積評(píng)價(jià)結(jié)果的誤差SVD＿P準(zhǔn)則4．2756?4．43174．02166．31%?10．20%Parsimony準(zhǔn)則4．5352?4．59274．13269．74%?11．13%豐度準(zhǔn)則0．0171?0．01690．019311．40%?12．43%

4 結(jié)論

本文研究表明，通過譜除法在頻域進(jìn)行反褶積處理、提取反褶積結(jié)果的主奇異值成分并結(jié)合Parsimony準(zhǔn)則進(jìn)行子波提取準(zhǔn)確性評(píng)價(jià)的方法(SVD_P準(zhǔn)則)是有效且可行的。與Parsimony準(zhǔn)則和豐度準(zhǔn)則相比，SVD_P準(zhǔn)則對(duì)隨機(jī)干擾噪聲的承受能力更強(qiáng)，更適用于實(shí)際地震資料子波估計(jì)的準(zhǔn)確性評(píng)價(jià)。

[1] 袁三一,王尚旭.子波相位不準(zhǔn)對(duì)反演結(jié)果的影響[J].應(yīng)用地球物理,2011,8(1):48-59 Yuan S Y,Wang S X.Influence of inaccurate wavelet phase estimation on seismic inversion[J].Applied Geophysics,2011,8(1):48-59

[2] 王蓉蓉,戴永壽,張亞南,等.非平穩(wěn)地震記錄中時(shí)變子波提取方法研究[J].地球物理學(xué)進(jìn)展,2015,30(2):700-708 Wang R R,Dai Y S,Zhang Y N,et al.Time-varying wavelet extraction methods in non-stationary seismogram[J].Progress in Geophysics,2015,30(2):700-708

[3] 戴永壽,王俊嶺,王偉偉,等.基于高階累積量ARMA模型線性非線性結(jié)合的地震子波提取方法研究[J].地球物理學(xué)報(bào),2008,51(6):1851-1859 Dai Y S,Wang J L,Wang W W,et al.Seismic wavelet extraction via cumulant based ARMA model approach with linear and nonlinear combination[J].Chinese Journal of Geophysics,2008,51(6):1851-1859

[4] 劉喜武,劉洪,李幼銘,等.局域波分解及其在地震信號(hào)時(shí)頻分析中的應(yīng)用[J].地球物理學(xué)進(jìn)展,2007,22(2):365-375 Liu X W,Liu H,Li Y M,et al.Local waves decomposition and its applications to seismic signal time-frequency analysis[J].Progress in Geophysics,2007,22(2):365-375

[5] 陳健,戴永壽,張亞南,等.基于高階統(tǒng)計(jì)量的地震子波提取方法評(píng)價(jià)[J].石油地球物理勘探,2013,48(3):497-503 Chen J,Dai Y S,Zhang Y N,et al.Summary of the evaluation approaches for seismic wavelet pick-up based on higher order statistics [J].Oil Geophysical Prospecting,2013,48(3):497-503

[6] 李振春,王希萍,韓文功.地震數(shù)據(jù)處理中的相位校正技術(shù)綜述[J].地球物理學(xué)進(jìn)展,2008,23(3):768-774 Li Z C,Wang X P,Han W G.Review of phase correction in seismic data processing[J].Progress in Geophysics,2008,23(3):768-774

[7] White R E.Maximum kurtosis phase correction[J].Geophysical Journal,1988,95(2):371-389

[8] 何同林,鄭鵬,劉郁林,等.基于峰度準(zhǔn)則與判決引導(dǎo)的非線性盲解卷積[J].電子科技大學(xué)學(xué)報(bào),2007,36(2):186-189 He T L,Zheng P,Liu Y L,et al.Nonlinear blind deconvolution based on kurtosis criterion and decision directed algorithm[J].Journal of University of Electronic Science and Technology of China,2007,36(2):186-189

[9] 姚逢昌.振幅譜補(bǔ)償和相位校正[J].石油物探,1990,29(3):46-58 Yao F C.Amplitude spectra compensation and phase correction[J].Geophysical Prospecting for Petroleum,1990,29(3):46-58

[10] Pham D T.Generalized mutual information approach to multichannel blind deconvolution[J].Signal Processing,2007,87(9):2045-2060

[11] 戴永壽,鄭德玲,魏磊,等.高階統(tǒng)計(jì)量地震子波估計(jì)建模[J].石油地球物理勘探,2006,41(5):514-518 Dai Y S,Zheng D L,Wei L,et al.Primary research on modeling for wavelet estimation via high-order cumulants[J].Oil Geophysical Prospecting,2006,41(5):514-518

[12] 唐斌,尹成.基于高階統(tǒng)計(jì)的非最小相位地震子波恢復(fù)[J].地球物理學(xué)報(bào),2001,44(3):404-410 Tang B,Yin C.Non-minimum phase seismic wavelet reconstruction based on higher order statistics[J].Chinese Journal of Geophysics,2001,44(3):404-410

[13] 劉玉金,李振春,郭愷.基于局部相似度的非穩(wěn)態(tài)相位校正方法[J].石油地球物理勘探,2012,47(6):887-893 Liu Y J,Li Z C,Guo K.Non-stationary phase correction based on local similarity[J].Oil Geophysical Prospecting,2012,47(6):887-893

[14] 齊春艷,李彥鵬,彭繼新,等.一種改進(jìn)的廣義S變換[J].石油地球物理勘探,2010,45(2):215-218 Qi C Y,Li Y P,Peng J X,et al.An improved generalized S-transform[J].Oil Geophysical Prospecting,2010,45(2):215-218

[15] Ziolkowski A.Why don’t we measure seismic signatures? [J].Geophysics,1991,56(2):190-210

[16] 高靜懷,汪玲玲,趙偉.基于反射地震記錄變子波模型提高地震記錄分辨率[J].地球物理學(xué)報(bào),2009,52(5):1289-1300 Gao J H,Wang L L,Zhao W.Enhancing resolution of seismic traces based on the changing wavelet model of the seismogram[J].Chinese Journal of Geophysics,2009,52(5):1289-1300

[17] 張亞南,戴永壽,王少水,等.高效 ARMA 模型高分辨率地震子波提取方法[J].石油地球物理勘探,2011,46(5):686-694 Zhang Y N,Dai Y S,Wang S S,et al.High resolution wavelet estimation by ARMA modeling[J].Oil Geophysical Prospecting,2011,46(5):686-694

[18] Wang R R,Dai Y S,Li C,et al.A time-varying mixed-phase wavelet extraction and evaluation method based on adaptive segmentation[J].Electronic Journal of Geotechnical Engineering,2015,20(6):1403-1418

[19] Robinson E A.Predictive decomposition of time series with application to seismic exploration[J].Geophysics,1967,32(3):418-484

[20] Wiggins R.Minimum entropy deconvolution[J].Geoexploration,1978,16(5):21-35

[21] Sacchi M D.Re-weighting strategies in seismic deconvolution[J].Geophysical Journal International,1997,12(9):651-656

[22] Yuan S Y,Wang S X.Spectral sparse Bayesian learning reflectivity inversion[J].Geophysical Prospecting,2013,61(4):735-746

[23] 章珂,李衍達(dá),劉貴忠,等.多分辨率地震信號(hào)反褶積[J].地球物理學(xué)報(bào),1999,42(4):529-535 Zhang K,Li Y D,Liu G Z,et al.Multiresolution seismic signal deconvolution[J].Chinese Journal of Geophysics,1999,42(4):529-535

[24] Zhang R,Castagna J.Seismic sparse-layer reflectivity inversion using basis pursuit decomposition[J].Geophysics,2011,76(6):147-158

[25] Oropeza V E,Sacchi M D.Multifrequency singular spectrum analysis[J].Expanded Abstracts of 79thAnnual Internat SEG Mtg,2009,3193-3197

[26] Li C,Huang J P,Li Z C,et al.Simultaneous de-noising and regularization method for 3D field data based on the MSSA[J].CPS/SEG International Geophysical Conference,2014,296-299

[27] Yuan S Y,Wang S X.A local f-x Cadzow method for noise reduction of seismic data obtained in complex formations[J].Petroleum Science,2011,8(3):269-277

[28] Oropeza V E,Sacchi M D.A randomized SVD for multichannel singular spectrum analysis (MSSA) noise attenuation[J].Expanded Abstracts of 80thAnnual Internat SEG Mtg,2010,1989-1992

[29] 沈鴻雁,李慶春.奇異值分解 (SVD) 實(shí)現(xiàn)地震波場(chǎng)分離與去噪新思路[J].地球物理學(xué)進(jìn)展,2010,25(1):225-230 Shen H Y,Li Q C.New idea for seismic wave filed separation and denoising by singular value decomposition(SVD)[J].Progress in Geophysics,2010,25(1):225-230

[30] 胡永泉,尹成,潘樹林,等.基于單道奇異值分解的微地震資料去噪方法[J].石油天然氣學(xué)報(bào),2013,35(4):64-69 Hu Y Q,Yin C,Pan S L,et al.Denoising method for microseismic data based on single-channel SVD[J].Journal of Oil and Gas Technology,2013,35(4):64-69

[31] 王有新,周興元.最小熵反褶積和集約化準(zhǔn)則[J].石油地球物理勘探,1993,28(6):678-684 Wang Y X,Zhou X Y.Minimum entropy deconvelution and parsimonious criterion[J].Oil Geophysical Prospecting,1993,28(6):678-684

[32] 李大衛(wèi),尹成,趙維坤,等.基于相位掃描的地震子波提取方法研究[J].西南石油大學(xué)學(xué)報(bào),2007,29(3):17-19 Li D W,Yin C,Zhao W K,et al.Seismic wavelet extraction based on phase scan[J].Journal of Southwest Petroleum University,2007,29(3):17-19

[33] Margrave G F,Lamoureux M P,Henley D C.Gabor deconvolution:estimating reflectivity by nonstationary deconvolution of seismic data[J].Geophysics,2011,76(3):W15-W30

[34] 張亞南,戴永壽,陳健,等.用對(duì)稱映射ARMA模型的零極點(diǎn)研究子波相位對(duì)反射系數(shù)序列反演的影響[J].地球物理學(xué)報(bào),2013,56(6):2043-2054 Zhang Y N,Dai Y S,Chen J,et al.The research on the influence of wavelet phase on the inversion results of reflection coefficient sequences by using the ARMA model of symmetrical mapping pole-zeros[J].Chinese Journal of Geophysics,2013,56(6):2043-2054

(編輯：戴春秋)

An evaluation criterion on the accuracy of time-varying wavelet extraction based on singular value decomposition

Wang Rongrong1,Dai Yongshou1,Li Chuang2,Zhang Manman1,Zhang Peng1

(1.CollegeofInformationandControlEngineering,ChinaUniversityofPetroleum,Qingdao266580,China;2.SchoolofGeosciences,ChinaUniversityofPetroleum,Qingdao266580,China)

The accuracy evaluation of time-varying wavelet extraction plays an important role in seismic data processing.However,the conditional evaluation criterion is influenced seriously by noise.Therefore,we propose a time-varying wavelet accuracy criterion based on singular value decomposition (SVD).Since the Parsimony criterion,Kurtosis criterion and Absolute kurtosis criterion have good tolerability to noisy environment among the existing evaluation criteria for the non-stationary seismic wavelet extraction accuracy,the Parsimony criterion and SVD technology are combined to construct a SVD_P criterion which has better noise-tolerant ability; and the spectrum division is employed as the deconvolution method.The Parsimony criterion,Kurtosis criterion and SVD_P criterion are applied to the simulation experiment and field data processing to compare the precision of time-frequency domain time-varying wavelet extraction method and adaptive segmentation time-varying wavelet extraction method.The results show that all three criteria could provide valid evaluation of these two wavelet extraction method while the time-frequency domain wavelet extraction method is more accurate than the adaptive segmentation method.Additionally,the evaluation result of SVD_P criterion owns smallest error and highest evaluation precision.

wavelet extraction,accuracy,evaluation criterion,singular value decomposition

2015-01-13;改回日期：2015-04-24。

王蓉蓉(1990—),女,碩士在讀,主要從事地震信號(hào)處理方法研究工作。

戴永壽(1963—),男,教授,博士生導(dǎo)師,主要從事地震信號(hào)處理方法研究工作。

國(guó)家自然科學(xué)基金項(xiàng)目(40974072)和中國(guó)石油大學(xué)(華東)研究生創(chuàng)新工程資助項(xiàng)目(YCX2015050)聯(lián)合資助。

P631

A

1000-1441(2015)05-0531-10

10.3969/j.issn.1000-1441.2015.05.006