基于平面波解構(gòu)的三維體曲率計(jì)算方法

李海山,楊午陽(yáng)

(1.中國(guó)石油天然氣股份有限公司勘探開(kāi)發(fā)研究院西北分院,甘肅蘭州730020;2.CNPC油藏描述重點(diǎn)實(shí)驗(yàn)室,甘肅蘭州730020)

表示地層彎曲程度的曲率屬性通常被用于構(gòu)造解釋和儲(chǔ)層分析,因其對(duì)各種復(fù)雜斷層、裂縫、河道及構(gòu)造彎曲的刻畫(huà)能力優(yōu)于相干屬性[1],近年來(lái)得到了廣泛的關(guān)注。曲率屬性自20世紀(jì)90年代中期被引入地震解釋,最初采用層面計(jì)算方式得到層位曲率[2]。ROBERTS[3]給出了詳細(xì)的層位曲率計(jì)算方法,明確了曲率屬性的物理意義并對(duì)其進(jìn)行了系統(tǒng)分類(lèi)。但層位曲率沒(méi)有直接利用地震資料的振幅信息,層位解釋的質(zhì)量對(duì)計(jì)算結(jié)果有著嚴(yán)重影響,容易產(chǎn)生構(gòu)造假象[4],同時(shí)也不能實(shí)現(xiàn)對(duì)整個(gè)構(gòu)造體及構(gòu)造內(nèi)部細(xì)節(jié)的精細(xì)刻畫(huà)[5]。

為克服層位曲率的局限性,AL-DOSSARY等[6]在ROBERTS[3]研究的基礎(chǔ)上,提出了三維多譜體曲率計(jì)算方法。由于體曲率利用了地震數(shù)據(jù)體的振幅信息,而振幅信息與構(gòu)造和巖性變化密切相關(guān),且計(jì)算不受人為因素的干擾,因此體曲率屬性相對(duì)于層面曲率屬性表現(xiàn)出諸多的優(yōu)點(diǎn)[4],很快得到了推廣應(yīng)用[7-11]。解釋人員可以從沿層曲率屬性上識(shí)別出小的擾曲、褶皺、凸起、差異壓實(shí)特征,而這些在常規(guī)解釋時(shí)無(wú)法追蹤,在相干上也呈現(xiàn)為連續(xù)的高相干特征[12]。體曲率的計(jì)算關(guān)鍵在于傾角數(shù)據(jù)體的計(jì)算,目前有多種傾角數(shù)據(jù)體的計(jì)算方法,主要有局部?jī)A斜疊加法[13]、復(fù)數(shù)道分析法[14]、相干傾角掃描法[15-16]、梯度結(jié)構(gòu)張量算法[17]以及平面波解構(gòu)法(PWD)[18]等。

平面波解構(gòu)法最初由CLAERBOUT[19]于1992年提出,其主要思想來(lái)源于描述地震數(shù)據(jù)的局部平面波模型,目前已經(jīng)廣泛應(yīng)用于波場(chǎng)分離與成像[20]、局部斜率或傾角計(jì)算[21-22]、構(gòu)造導(dǎo)向?yàn)V波[23]、反演參數(shù)建模[24]、地層自動(dòng)拾取與拉平[25]等。考慮到平面波解構(gòu)法具有簡(jiǎn)單快速、計(jì)算精度高、穩(wěn)定可靠等優(yōu)點(diǎn),本文將其引入三維體曲率計(jì)算。首先利用平面波解構(gòu)法將三維地震數(shù)據(jù)體轉(zhuǎn)化為視傾角數(shù)據(jù)體,再利用視傾角數(shù)據(jù)體計(jì)算三維空間中任意點(diǎn)的曲率,從而獲得三維體曲率屬性,最后將該方法用于實(shí)際資料的體曲率屬性的計(jì)算,結(jié)果表明本文方法有效且可靠,三維體曲率屬性能夠精細(xì)刻畫(huà)裂縫發(fā)育帶和河道等地質(zhì)特征。

1 方法原理

1.1 基于PWD的傾角體計(jì)算

平面波解構(gòu)濾波器是一種時(shí)空域預(yù)測(cè)誤差濾波器,它通過(guò)相鄰地震道預(yù)測(cè)當(dāng)前地震道來(lái)估算同相軸局部斜率,其預(yù)測(cè)準(zhǔn)則是原始地震道與預(yù)測(cè)出的地震道間的殘差能量最小。

1.1.1 平面波解構(gòu)

局部平面波微分方程可表示為[18]:

(1)

式中:u=u(t,x)為波場(chǎng);t和x分別為時(shí)間和距離;σ=σ(t,x)為同相軸局部斜率。

FOMEL[18]在Z變換域中給出了求解公式(1)的隱式有限差分公式:

(2)

(3)

式中:濾波器B(zt)的系數(shù)可以通過(guò)在零頻率附近進(jìn)行泰勒級(jí)數(shù)展開(kāi)確定,計(jì)算中通常采用三點(diǎn)中心濾波器B3(zt)[18]:

(4)

式中:

(5)

1.1.2 局部?jī)A角估算

用C(σ)表示平面波解構(gòu)算子C(zt,zx)在時(shí)間域與地震數(shù)據(jù)u(t,x)進(jìn)行二維褶積運(yùn)算的算子,在最小平方意義下,局部斜率σ的估算可轉(zhuǎn)化為求解最小二乘目標(biāo)函數(shù)[18]:

(6)

由公式(5)可知最優(yōu)化問(wèn)題(公式(6))是局部斜率σ的非線性?xún)?yōu)化問(wèn)題,利用高斯-牛頓法可將最優(yōu)化問(wèn)題((6)式)轉(zhuǎn)化為[18]:

(7)

式中:Δσ為斜率增量;σ0為初始斜率;C′(σ)為C(σ)對(duì)σ的偏導(dǎo)數(shù)。求解公式(7)后,初始斜率σ0通過(guò)加上斜率增量Δσ得到更新,然后循環(huán)迭代求解公式(7)。該方法需要進(jìn)行數(shù)次非線性迭代獲得期望的局部斜率估計(jì),收斂速率依賴(lài)于給定的起始斜率值。得到局部斜率σ后,就可以得到相應(yīng)的局部?jī)A角θ:

(8)

迭代求解公式(7)得到的局部斜率σ隨著時(shí)間和空間位置的改變而變化,同時(shí)避免如相干傾角掃描法那樣的局部時(shí)空窗模糊效應(yīng)。為了提高計(jì)算結(jié)果的連續(xù)性以及對(duì)噪聲的壓制能力,在局部斜率更新過(guò)程中,可引入預(yù)條件算子和正則化方法[20-21]:

(9)

式中:ε為常系數(shù)擾動(dòng)因子;D為正則化算子。

對(duì)于三維傾角體計(jì)算,首先沿縱測(cè)線方向利用公式(9)對(duì)每一個(gè)剖面u(t,x)進(jìn)行計(jì)算,獲得視傾角體px=p(t,x,y),然后沿聯(lián)絡(luò)測(cè)線方向利用公式(9)對(duì)每一個(gè)剖面u(t,y)進(jìn)行計(jì)算,獲得視傾角體qy=q(t,x,y),本文計(jì)算時(shí)起始斜率值均置為0。

1.2 體曲率計(jì)算

利用平面波解構(gòu)法獲得的傾角數(shù)據(jù)體,可計(jì)算出三維體曲率屬性,下面詳細(xì)介紹各種體曲率的計(jì)算方法。

一個(gè)反射面可以用二次曲面方程表示為[3]:

S(x,y)=ax2+by2+cxy+dx+ey+f

(10)

該方程中的系數(shù)與前文計(jì)算的視傾角px和qy具有如下關(guān)系[6]:

(11)

其中,算子Dx和Dy分別表示沿著x和y方向一階導(dǎo)數(shù)的數(shù)值近似。在獲得二次曲面方程系數(shù)后,利用下列公式計(jì)算曲率屬性[6]。

均值曲率為:

(12)

高斯曲率為:

(13)

最大和最小曲率為:

(14)

最正曲率和最負(fù)曲率為:

(15)

傾角曲率為:

(16)

走向曲率為:

(17)

2 局部?jī)A角估算方法性能測(cè)試

為說(shuō)明基于PWD的局部?jī)A角估算方法的性能,下面利用Sigmoid模型[18,26]數(shù)據(jù)進(jìn)行測(cè)試(圖1)。圖1a所示的Sigmoid模型包括一個(gè)褶皺中發(fā)育的斷層和褶皺被剝蝕后上覆沉積形成的角度不整合面。利用圖1b所示的初始傾角(所有樣點(diǎn)均為0)進(jìn)行迭代求解,50次迭代后得到的解構(gòu)殘差和傾角分別如圖1c和圖1d所示。由圖1d可見(jiàn)估算的傾角剖面中,不整合面上覆緩傾角地層的傾斜角度一致,不整合面下伏褶皺的傾角變化與起伏情況一致,褶皺內(nèi)的斷層也得到了清晰的刻畫(huà);由圖1c可見(jiàn)解構(gòu)殘差非常小,合成沉積模型中的同相軸信息基本被解構(gòu)出來(lái),進(jìn)一步驗(yàn)證了傾角估計(jì)的準(zhǔn)確性。圖2給出了平面波解構(gòu)殘差能量隨迭代次數(shù)的變化曲線,盡管初始傾角為0,解構(gòu)殘差能量仍然隨迭代次數(shù)的增加而單調(diào)遞減,可見(jiàn)基于PWD的傾角計(jì)算方法比較穩(wěn)定。

圖1 Sigmoid模型及相應(yīng)的解構(gòu)殘差及傾角剖面(1rad≈57.296°)a Sigmoid模型; b 初始傾角; c 平面波解構(gòu)殘差; d 估算的傾角剖面

為進(jìn)一步說(shuō)明基于PWD的局部?jī)A角估算方法計(jì)算精度高的優(yōu)點(diǎn),將其應(yīng)用于實(shí)際資料的傾角場(chǎng)計(jì)算。圖3a為某工區(qū)疊后地震剖面,圖3b為利用平面波解構(gòu)法得到的相應(yīng)傾角剖面,圖3c為利用相干傾角掃描法得到的相應(yīng)傾角剖面。由圖可見(jiàn),利用平面波解構(gòu)法得到的相應(yīng)傾角剖面對(duì)同相軸傾角和斷層等不連續(xù)性有較好的刻畫(huà)能力,不僅將如圖中紅色箭頭所指位置處的傾角場(chǎng)進(jìn)行了準(zhǔn)確刻畫(huà),而且將如圖中黑色箭頭所指位置處的斷層以?xún)A角場(chǎng)突變的形式也刻畫(huà)了出來(lái)。

圖2 平面波解構(gòu)殘差能量隨迭代次數(shù)的變化

圖3 疊后地震剖面及相應(yīng)的傾角剖面(1rad≈57.296°)a 疊后地震剖面; b 平面波解構(gòu)法得到的傾角剖面; c 相干傾角掃描法得到的傾角剖面

3 實(shí)際資料應(yīng)用

為驗(yàn)證本文體曲率屬性計(jì)算方法的有效性,利用西北某油田H工區(qū)的地震資料進(jìn)行了體曲率計(jì)算的應(yīng)用研究。采用的地震資料為三維疊后縱波資料,主要研究目的層為深層裂縫型儲(chǔ)層。圖4給出了該工區(qū)目的層解釋層位及沿層振幅切片,圖中W1,W2,W3指示了3口開(kāi)發(fā)井位置。圖5給出了平面波解構(gòu)法得到的視傾角沿層切片,為了對(duì)比分析,同時(shí)采用目前廣泛使用的相干傾角掃描法進(jìn)行了體傾角計(jì)算,圖6給出了相干傾角掃描法得到的視傾角沿層切片。對(duì)比可見(jiàn)平面波解構(gòu)法得到的視傾角在空間連續(xù)性及局部細(xì)節(jié)的刻畫(huà)(如圖中紅色箭頭所示的河道及微斷裂位置)方面都明顯優(yōu)于相干傾角掃描法得到的視傾角。此外,在同等計(jì)算條件下,計(jì)算兩個(gè)視傾角數(shù)據(jù)體時(shí)相干傾角掃描法所用計(jì)算時(shí)間是平面波解構(gòu)法所用計(jì)算時(shí)間的16.51倍,可見(jiàn)平面波解構(gòu)法的計(jì)算效率較高。

圖7給出了沿目的層相干切片及利用不同方法計(jì)算得到的最負(fù)曲率屬性切片。圖7a為基于本征值結(jié)構(gòu)的相干切片;圖7b為利用層位數(shù)據(jù)計(jì)算的最負(fù)曲率屬性,受層位解釋與地震同相軸追蹤質(zhì)量的影響,計(jì)算精度和可靠性較低;圖7c為利用平面波解構(gòu)法所得視傾角計(jì)算的最負(fù)曲率屬性;圖7d為利用相干傾角掃描法所得視傾角計(jì)算的最負(fù)曲率屬性,在錯(cuò)斷等不連續(xù)部位,由于在時(shí)窗內(nèi)不能掃描到精確的傾角而導(dǎo)致曲率計(jì)算結(jié)果模糊。比較可見(jiàn),圖7c的計(jì)算精度最高,細(xì)節(jié)刻畫(huà)能力優(yōu)于沿層相干屬性,特別是紅色箭頭所示的河道和微斷裂位置及黑色圓圈所示的裂縫發(fā)育帶位置得到了較細(xì)致的刻畫(huà)。平面波解構(gòu)法所得三維體曲率屬性的計(jì)算精度和對(duì)地質(zhì)特征的精細(xì)刻畫(huà)能力在圖8給出的沿目的層傾角曲率屬性切片上得到了更加突出的體現(xiàn)。

圖4 目的層解釋層位(a)及沿層振幅切片(b)

圖5 平面波解構(gòu)法得到的視傾角沿層切片(1rad≈57.296°)a 縱測(cè)線方向視傾角沿層切片; b 聯(lián)絡(luò)測(cè)線方向視傾角沿層切片

圖6 相干傾角掃描法得到的視傾角沿層切片(1rad≈57.296°)a 縱測(cè)線方向視傾角沿層切片; b 聯(lián)絡(luò)測(cè)線方向視傾角沿層切片

圖7 沿目的層相干切片及最負(fù)曲率屬性切片a 基于本征值結(jié)構(gòu)的相干切片; b利用圖4a所示解釋層位計(jì)算的最負(fù)曲率屬性; c 平面波解構(gòu)法得到的最負(fù)曲率屬性; d 相干傾角掃描法得到的最負(fù)曲率屬性

圖8 沿目的層傾角曲率屬性切片a 平面波解構(gòu)法得到的傾角曲率屬性; b 相干傾角掃描法得到的傾角曲率屬性

4 結(jié)論

1) 本文方法利用平面波解構(gòu)法由三維地震數(shù)據(jù)體計(jì)算得到視傾角數(shù)據(jù)體,再利用視傾角數(shù)據(jù)體計(jì)算三維空間中任意點(diǎn)的曲率,實(shí)現(xiàn)基于平面波解構(gòu)的三維體曲率計(jì)算。

2) 本文方法結(jié)合了體曲率計(jì)算時(shí)利用振幅信息、不受人為因素干擾的優(yōu)點(diǎn)和平面波解構(gòu)法簡(jiǎn)單快速、計(jì)算精度高、穩(wěn)定可靠的優(yōu)點(diǎn),使得計(jì)算效率高、穩(wěn)定可靠,具有較強(qiáng)的實(shí)用性。

3) 實(shí)際資料應(yīng)用結(jié)果驗(yàn)證了本文方法的有效性,利用平面波解構(gòu)法得到的體曲率屬性明顯優(yōu)于利用解釋層位數(shù)據(jù)得到的層位曲率屬性和利用相干掃描法得到的體曲率屬性,對(duì)裂縫發(fā)育帶及河道等非均質(zhì)地質(zhì)體有較強(qiáng)的刻畫(huà)能力。