多道聯(lián)合約束的匹配追蹤強(qiáng)反射軸壓制方法

楊子鵬 宋維琪* 劉 軍 陳俊安 劉 群 武 迪

(①中國(guó)石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院，山東青島 266580； ②中國(guó)石化西北油田勘探開發(fā)研究院，新疆烏魯木齊 830011)

0 引言

地震剖面中強(qiáng)反射軸的存在，對(duì)位于其內(nèi)部或附近的目的層有效信號(hào)造成了較嚴(yán)重的影響，尤其是強(qiáng)反射軸附近存在隱蔽油氣儲(chǔ)層時(shí)，儲(chǔ)層的反射信息會(huì)被強(qiáng)反射軸遮蓋或融合而難以識(shí)別和區(qū)分。但這些儲(chǔ)層往往又因?yàn)槁癫匚恢煤?、產(chǎn)量高，而成為重要勘探目標(biāo)。因此消除強(qiáng)反射軸影響、提高油氣儲(chǔ)層刻畫的精確度變得越來越重要。

經(jīng)過多年探究，多子波分解與重構(gòu)、匹配追蹤成為強(qiáng)反射軸壓制的兩種常用技術(shù)。劉炳楊等[1]認(rèn)為單一子波的假設(shè)會(huì)造成部分有用信息的丟失，由此應(yīng)用多子波分解與重構(gòu)技術(shù)成功剝離了目標(biāo)區(qū)的強(qiáng)屏蔽層。佘剛等[2]利用多子波分解與重構(gòu)技術(shù)去除了煤層強(qiáng)反射，提高了儲(chǔ)層預(yù)測(cè)精度。

匹配追蹤算法是由Mallat等[3]于1993年首次提出，并詳細(xì)推導(dǎo)了理論公式，自此該算法開始應(yīng)用于地震勘探多個(gè)領(lǐng)域，如薄層識(shí)別[4-5]、時(shí)頻分析[6-7]、地震道缺失重建[8]等。Wang[9-10]通過改進(jìn)提出了多道匹配追蹤算法，將其應(yīng)用于強(qiáng)反射軸剝離，取得較好效果。宋維琪等[11]將復(fù)數(shù)子波應(yīng)用于匹配追蹤算法，識(shí)別了薄層砂體。李海山等[12]對(duì)匹配追蹤算法中的過完備字典進(jìn)行優(yōu)化，通過匹配煤層強(qiáng)反射信息消除了其屏蔽作用。張?jiān)诮鸬萚13]提出層位與子波約束下的匹配追蹤去除強(qiáng)軸技術(shù)，剝離強(qiáng)軸后在儲(chǔ)層下方發(fā)現(xiàn)明顯的伴影異常。許璐等[14]針對(duì)瞬時(shí)頻率的“負(fù)頻率”問題，將局部頻率引入匹配追蹤算法，利用反射能量公式剝離了強(qiáng)反射軸。何峰等[15]從井旁地震道提取匹配子波，利用改進(jìn)的匹配追蹤算法去除了煤層強(qiáng)反射，凸現(xiàn)目標(biāo)砂體。

前人的研究大多卓有成效，但鮮見針對(duì)“迭代終止”方面的研究文獻(xiàn)，且迭代次數(shù)若選擇不當(dāng)，信號(hào)分解重構(gòu)時(shí)就會(huì)引入噪聲； 其次，上述研究主要聚焦于強(qiáng)反射軸附近被壓制的信號(hào)，卻忽略或繞開了強(qiáng)反射軸內(nèi)部含有的有效信號(hào)問題。

為此，本文利用殘差比閾值約束匹配追蹤算法中的迭代，避免了迭代次數(shù)選擇不當(dāng)引入的噪聲，并在一定程度上提高了信號(hào)分解的效率。為了保留強(qiáng)反射軸中的有效信號(hào)，利用從多道提取的強(qiáng)反射信號(hào)求取能量加權(quán)系數(shù)，增加了強(qiáng)反射軸中有效信號(hào)的能量占比，再引入強(qiáng)反射軸壓制參數(shù)λ，從而壓制了強(qiáng)反射軸。理論模型和實(shí)際地震資料的處理結(jié)果表明，與常規(guī)層位約束匹配追蹤強(qiáng)反射軸壓制方法相比，本文方法具有更好效果。

1 方法原理

1.1 匹配追蹤算法原理及實(shí)現(xiàn)方法

匹配追蹤算法在信號(hào)稀疏分解方面具有較好效果，常用的有單道匹配追蹤算法和多道匹配追蹤算法[16-17]。由于單道匹配追蹤算法結(jié)果的橫向連續(xù)性較差，因此本文選擇多道匹配追蹤算法。

假設(shè)一道地震信號(hào)x可被分解為一系列子波與隨機(jī)噪聲的組合，即可表達(dá)[18-19]為

x=〈x，ψγ1〉ψγ1+R1(x)

(1)

式中：〈x，ψγ1〉為第1次迭代的地震信號(hào)x與其最優(yōu)匹配子波ψγ1的內(nèi)積，γ1表示第1次迭代后最優(yōu)匹配子波的參數(shù)集合；R1(x)為第1次迭代后遺留的殘差項(xiàng)。由于R1(x)與ψγ1是正交的，因此滿足

(2)

K次迭代后，信號(hào)x可表示為

(3)

式中ψγk為第k次迭代出的匹配子波。由于匹配追蹤算法是一個(gè)迭代過程，此過程中為了使殘差項(xiàng)RK(x)足夠小，則必須使地震信號(hào)x與其最優(yōu)匹配子波ψγk∈Ψ的內(nèi)積盡量大，這就需要找到與信號(hào)x最匹配的子波，最終將x分解為K-1個(gè)子波與第K次迭代后殘差項(xiàng)的線性表示。因此，最優(yōu)匹配子波的參數(shù)求取關(guān)系到信號(hào)的分解和重構(gòu)精度。本文采用三步法[10]求取最優(yōu)匹配子波的控制參數(shù)，具體實(shí)現(xiàn)過程如下。

第一步，初始參數(shù)估計(jì)。以Morlet小波作為基本原子，其表達(dá)式[20]為

exp[iω(t-μ)+φ]

(4)

式中： μ表示位移參數(shù)；ω表示子波圓頻率；σ表示尺度參數(shù)；φ表示子波相位。

(5)

式中Ψ={ψγk}Γ，是子波組成的過完備字典，其中Γ表示γk的集合。

(6)

由式(6)得到的最優(yōu)子波參數(shù)適用于相鄰L道中的任意道。

第三步，匹配子波振幅估計(jì)。求取最優(yōu)子波參數(shù)后，利用

(7)

可得到相鄰L道中各單道對(duì)應(yīng)的振幅al，k，最后即得最優(yōu)匹配子波wl，k=al，kψγk。

1.2 迭代終止條件約束下的匹配追蹤法(TC-MMP)

一般來說，傳統(tǒng)的匹配追蹤算法迭代終止條件主要有兩類：一是給定一個(gè)硬閾值門限，如設(shè)定最大迭代次數(shù)K，那么迭代K次后的信號(hào)即可表示為K-1個(gè)匹配子波與其殘差項(xiàng)的線性組合。但不同地震道的最佳迭代次數(shù)難以確定。最大迭代次數(shù)設(shè)置得過小，會(huì)造成殘差較大，影響信號(hào)的重構(gòu)精度；設(shè)置得過大時(shí)，對(duì)于較低信噪比的信號(hào)，分解到一定迭代次數(shù)時(shí)噪聲方差可能比有效信號(hào)方差大，若繼續(xù)迭代分解，信號(hào)重構(gòu)時(shí)就會(huì)引入一定噪聲，且耗時(shí)較長(zhǎng)，因此硬閾值門限的迭代終止條件不一定可取。二是給定一個(gè)殘差閾值門限，即當(dāng)殘差小于給定閾值時(shí)終止迭代，這種方法對(duì)于較高信噪比信號(hào)較適用，最終重構(gòu)結(jié)果誤差較??；當(dāng)信號(hào)信噪比較低時(shí)，給定一個(gè)特定閾值，殘差會(huì)以指數(shù)形式減小，但降至一定階段殘差會(huì)隨迭代次數(shù)增加幾乎保持不變，此時(shí)子波能量遠(yuǎn)小于殘差，因此給定殘差閾值門限的終止條件也不一定可取。

基于上述分析，本文采用殘差比閾值作為迭代終止條件[21-22]，有效保證信號(hào)分解重構(gòu)的精度。

含噪聲地震信號(hào)可表示為

s=x+NΔf

(8)

式中：Δf為地震信號(hào)的頻寬；NΔf表示Δf內(nèi)的噪聲； 此處x為不含噪聲地震信號(hào)。

根據(jù)匹配追蹤算法的重構(gòu)理論[3]，殘差項(xiàng)Rk(x+NΔf)以指數(shù)收斂到0，即

(9)

式中：α為字典Ψ={ψγk}Γ中匹配子波的數(shù)量；β為過完備字典庫Ψ的相關(guān)系數(shù)，可表示為

進(jìn)一步可得

(10)

由式(10)可知，第k次與第k+1次的迭代殘差之差也會(huì)以指數(shù)形式減小，因此看似也可將殘差之差作為迭代終止條件，但分析式(1)～式(3)可知，第k次與第k+1次迭代的殘差之差還有如下關(guān)系

(11)

因此，為了進(jìn)一步增強(qiáng)迭代終止條件的魯棒性，擬將殘差比q[Rk(s)]作為迭代終止條件，即有

(12)

圖1 傳統(tǒng)迭代方式與殘差比閾值迭代信號(hào)重構(gòu)對(duì)比(a)不同迭代方式信號(hào)重構(gòu)結(jié)果； (b)殘差及殘差比迭代結(jié)果

以9個(gè)不同頻率、相位、尺度、時(shí)間的Morlet小波合成理論信號(hào)(圖2)，驗(yàn)證TC-MMP算法的分解效果。對(duì)不含噪信號(hào)(圖2a)用本文算法迭代9次，每次迭代都可得到最優(yōu)匹配子波，且最終重構(gòu)信號(hào)與理論信號(hào)毫無差別，殘差為零。對(duì)理論信號(hào)加入噪聲(信噪比為5∶1)，再對(duì)該含噪信號(hào)進(jìn)行迭代分解(圖2b)，可見用本文算法迭代9次分解出的最優(yōu)匹配子波都是構(gòu)成理論信號(hào)的子波，重構(gòu)信號(hào)與理論信號(hào)高度相似，且殘差為加入的隨機(jī)噪聲，說明TC-MMP算法在噪聲環(huán)境下依然較穩(wěn)定。

圖2 TC-MMP算法理論信號(hào)迭代分解過程(a)未加噪信號(hào)； (b)加噪信號(hào)

1.3 多道聯(lián)合壓制強(qiáng)反射軸方法

通常利用層位約束的匹配追蹤算法識(shí)別并剝離強(qiáng)反射軸，但這類方法同時(shí)也壓制了強(qiáng)反射軸中的有效信號(hào)。本文以目標(biāo)強(qiáng)反射軸為中心，上下截取時(shí)窗，并確保時(shí)窗內(nèi)強(qiáng)反射軸位置的振幅最強(qiáng)； 利用TC-MMP算法對(duì)地震信號(hào)進(jìn)行分解后，得到一系列匹配子波及其殘差項(xiàng)，其中能量最強(qiáng)的匹配子波可視為強(qiáng)反射信號(hào)。針對(duì)強(qiáng)反射軸中被淹沒的有效信號(hào)，采用多道聯(lián)合壓制強(qiáng)反射軸方法。首先通過多道提取的強(qiáng)反射信號(hào)求取能量加權(quán)系數(shù)ε，其作用是增加強(qiáng)反射軸中有效信號(hào)的能量占比，再調(diào)整壓制參數(shù)λ(0～1之間)，即可得到保留有效信號(hào)的強(qiáng)反射軸壓制結(jié)果dnew

(13)

(14)

式中：dorigin，Λ代表原始地震信號(hào)；dstrong，Λ表示強(qiáng)反射軸能量集合；Λ={1，2，…，N}表示地震道編號(hào)，N為總的地震道數(shù)；p為平衡參數(shù)(一般取為1)； M(·)表示取平均。

圖3a是截取的含有強(qiáng)反射軸的實(shí)際地震剖面，且強(qiáng)反射軸中存在有效信號(hào)(黃色箭頭標(biāo)記處)。針對(duì)該區(qū)域內(nèi)特征剖面(圖3b)，分別得到用常規(guī)方法(圖3c)和本文方法(圖3d)提取的強(qiáng)反射軸。從紅色標(biāo)記處看到，與常規(guī)方法相比，應(yīng)用本文能量加權(quán)法后，強(qiáng)反射軸中的有效信號(hào)并未被提取，說明強(qiáng)反射軸壓制以后，其中的有效信號(hào)可凸顯出來，證明能量加權(quán)法確實(shí)有效。

圖3 不同方法強(qiáng)反射軸提取對(duì)比(a)強(qiáng)反射軸中含有效信號(hào)的實(shí)際剖面； (b)局部剖面； (c)常規(guī)方法提取的強(qiáng)反射軸； (d)本文方法提取的強(qiáng)反射軸

1.4 整套方法實(shí)現(xiàn)步驟

基于上述研討，本文方法包含如下步驟(圖4)。

圖4 基于TC-MMP算法的強(qiáng)反射軸壓制方法流程

(1)選取一段含強(qiáng)反射軸地震剖面，以目標(biāo)強(qiáng)反射軸為中心截取上下時(shí)窗，得到一組地震信號(hào)X={x1，x2，x3，…，xL}，L一般取5。

(3)經(jīng)式(6)求取的子波參數(shù)適用于該組地震信號(hào)L道中任一道，故將從步驟(2)提取的強(qiáng)反射信號(hào)當(dāng)作這L道信號(hào)的強(qiáng)反射，并存儲(chǔ)這L道強(qiáng)反射。

(4)滑動(dòng)L道選擇另一組地震信號(hào)X′，重復(fù)步驟(1)～步驟(3)，直至遍歷所有地震道(N道)，可得N道強(qiáng)反射信號(hào)，即從多道提取的強(qiáng)反射信號(hào)。

(5)利用式(13)求取能量加權(quán)系數(shù)，并確定最優(yōu)壓制參數(shù)，聯(lián)合壓制強(qiáng)反射軸。

2 模型驗(yàn)證

構(gòu)建一個(gè)二維速度模型(圖5)，其具體參數(shù)如表1所示，密度是由Gardner公式[23]求得。模型中溶洞②和③距高阻層頂界面(強(qiáng)反射軸)很近，可近似代表強(qiáng)反射軸中的有效信號(hào)； 其他溶洞和速度異常體可近似地表征強(qiáng)反射軸附近的有效(反射)信息。為了確定最優(yōu)的強(qiáng)反射軸壓制參數(shù)，采用30Hz的Morlet小波作為地震子波，選取模型的第135道(強(qiáng)反射軸中含有溶洞②信息)進(jìn)行單道模擬(對(duì)應(yīng)圖6中黑色曲線)。

圖5 二維速度模型

表1 速度模型參數(shù)

利用本文方法對(duì)模擬的單道數(shù)據(jù)做強(qiáng)反射軸壓制處理，通過調(diào)整壓制參數(shù)得到不同的λ(0.2、0.4、0.6、0.8、1.0)時(shí)模型第135道的壓制結(jié)果(圖6)。對(duì)比圖6a～圖6e發(fā)現(xiàn)，當(dāng)λ=0.6時(shí)(圖6c)強(qiáng)反射壓制效果較理想；圖6f中λ=0.6時(shí)強(qiáng)反射軸壓制結(jié)果與有效信號(hào)(紅線)之差(綠線)幾乎為0，也充分說明本文方法在壓制強(qiáng)反射軸后還能突顯被強(qiáng)反射軸淹沒的有效信號(hào)。

圖6 模型單道不同λ壓制強(qiáng)反射軸結(jié)果(a)λ=0.2； (b)λ=0.4； (c)λ=0.6； (d)λ=0.8； (e)λ=1.0

對(duì)整個(gè)模型進(jìn)行模擬得到圖7a所示正演記錄，圖中強(qiáng)反射軸中的有效信號(hào)幾乎看不見，強(qiáng)反射軸下方的有效信號(hào)若隱若現(xiàn)。對(duì)比常規(guī)方法(圖7b)與本文方法(圖7c，λ=0.6)壓制強(qiáng)反射軸后結(jié)果，可見經(jīng)兩種方法處理后強(qiáng)反射軸附近(紅色箭頭處)有效信號(hào)雖都得到凸顯，但常規(guī)方法處理后強(qiáng)反射軸中的有效信號(hào)難以尋覓，而本文方法較充分地保留了強(qiáng)反射軸中的有效信號(hào)。

圖7 不含噪正演模型壓制強(qiáng)反射軸前、后剖面(a)未加噪正演剖面； (b)常規(guī)方法處理結(jié)果； (c)本文方法處理結(jié)果(λ=0.6)

再對(duì)模型數(shù)據(jù)加入噪聲(圖8a，信噪比為4∶1)，分析、對(duì)比常規(guī)方法(圖8b)與本文方法(圖8c)壓制結(jié)果，可看出本文方法在噪聲環(huán)境下依然效果明顯(黃色箭頭處)，這無疑有助于目標(biāo)區(qū)儲(chǔ)層精細(xì)描述及后續(xù)地質(zhì)綜合解釋。

圖8 加噪正演模型壓制強(qiáng)反射軸前、后剖面(a)加噪正演剖面； (b)常規(guī)方法處理結(jié)果； (c)本文方法處理結(jié)果(λ=0.6)

(a)～(e)黑線為模型第135道信號(hào)，藍(lán)線為本文方法提取的強(qiáng)反射軸信號(hào)，紅線為強(qiáng)反射軸壓制后結(jié)果； (f)紅線為設(shè)置于模型中的有效信號(hào)，綠線為λ=0.6時(shí)的強(qiáng)反射軸壓制結(jié)果與有效信號(hào)之差。

3 實(shí)際資料處理

中國(guó)新疆S區(qū)的奧陶系碳酸鹽巖地層發(fā)育典型的縫洞型儲(chǔ)集體。因碳酸鹽巖內(nèi)部地層界面波阻抗差異較小，奧陶系內(nèi)幕反射明顯較弱，發(fā)育的巖溶縫洞在地震剖面上表現(xiàn)為“串珠”狀反射，若存在強(qiáng)反射軸，則會(huì)嚴(yán)重影響有效信號(hào)的顯現(xiàn)，加大有效信號(hào)的識(shí)別難度。因此，亟待研發(fā)有效的強(qiáng)反射軸壓制技術(shù)。

將本文方法應(yīng)用于該區(qū)強(qiáng)反射軸壓制，并列舉兩個(gè)實(shí)例驗(yàn)證本文方法的有效性和適用性。

在實(shí)例1原始地震剖面(圖9a)上，受T74強(qiáng)反射層影響，其附近的有效信號(hào)被掩蓋，若無測(cè)井資料，很難識(shí)別出有效信號(hào)。分別利用常規(guī)方法和本文方法對(duì)原始剖面進(jìn)行強(qiáng)反射軸壓制處理，發(fā)現(xiàn)兩種方法處理結(jié)果(圖9b，圖9c)都使強(qiáng)反射軸附近的信號(hào)得到了一定程度凸顯。進(jìn)一步將圖9a中紅色方框區(qū)域做局部放大顯示(圖10)，可見S5井漏失位置處有效信號(hào)弱，且強(qiáng)反射軸中也存在有效信號(hào)。對(duì)比常規(guī)方法(圖10b)與本文方法(圖10c)的強(qiáng)反射軸壓制結(jié)果，可見本文方法使強(qiáng)反射軸中的有效信號(hào)得到最大程度的保留。

圖9 實(shí)例1壓制強(qiáng)反射軸前、后剖面(a)原始地震剖面； (b)常規(guī)方法處理； (c)本文方法處理

圖10 實(shí)例1局部放大顯示(a)對(duì)應(yīng)圖9a； (b)對(duì)應(yīng)圖9b； (c)對(duì)應(yīng)圖9c

實(shí)例2的強(qiáng)反射軸壓制效果如圖11、圖12所示，可見S52井漏失位置處的有效信號(hào)原本被強(qiáng)反射軸淹沒，但經(jīng)本文方法處理后得以很好地凸顯，再次表明本文的強(qiáng)反射軸壓制方法確實(shí)有效且實(shí)用，可為該區(qū)縫洞儲(chǔ)層識(shí)別研究提供技術(shù)支撐。

圖11 實(shí)例2壓制強(qiáng)反射軸前、后剖面(a)原始地震剖面； (b)常規(guī)方法處理； (c)本文方法處理

圖12 實(shí)例2局部放大顯示(a)對(duì)應(yīng)圖11a； (b)對(duì)應(yīng)圖11b； (c)對(duì)應(yīng)圖11c

4 結(jié)論與認(rèn)識(shí)

通過本項(xiàng)研究，得到如下認(rèn)識(shí)和結(jié)論：

(1)當(dāng)目的層附近存在高阻層時(shí)，在地震剖面上形成強(qiáng)反射軸且會(huì)影響其附近有效信號(hào)的顯現(xiàn)，更難精確識(shí)別強(qiáng)反射軸內(nèi)部的有效信號(hào)，給儲(chǔ)層預(yù)測(cè)帶來巨大困難。

(2)將殘差比閾值作為迭代終止條件引入多道匹配追蹤，即形成TC-MMP算法。通過與兩種傳統(tǒng)迭代方式進(jìn)行對(duì)比，發(fā)現(xiàn)殘差比閾值迭代終止條件可較好地實(shí)現(xiàn)自適應(yīng)終止迭代，在一定程度上提高了信號(hào)分解的效率。

(3)以目標(biāo)強(qiáng)反射軸為中心建立時(shí)窗，用TC-MMP算法分解信號(hào)，將能量最強(qiáng)匹配子波作為強(qiáng)反射信號(hào)，通過多道提取的強(qiáng)反射信號(hào)求取能量加權(quán)系數(shù)，增大強(qiáng)反射軸中有效信號(hào)的能量占比，再通過試算方法選取理想強(qiáng)反射軸壓制參數(shù)，使最終強(qiáng)反射軸壓制結(jié)果最大程度上保留強(qiáng)反射軸中有效信號(hào)，同時(shí)也凸顯該強(qiáng)反射軸附近有效信號(hào)。

(4)模型數(shù)據(jù)和實(shí)際資料的處理結(jié)果表明，與常規(guī)層位約束匹配追蹤強(qiáng)反射軸壓制方法相比，本文方法壓制強(qiáng)反射軸效果更佳，實(shí)用性更強(qiáng)。

需要指出的是：應(yīng)用本文方法時(shí)需首先確保時(shí)窗內(nèi)強(qiáng)反射位置的振幅能量最大； 其次時(shí)窗的選取需根據(jù)實(shí)際資料特征確定，若取得過大就會(huì)延緩信號(hào)分解速度，進(jìn)而影響其使用效率。