張冰， 徐嘉亮* ， 王維紅， 石穎， 王鵬

1 東北石油大學(xué)地球科學(xué)學(xué)院， 大慶 163318 2 油氣資源與勘探技術(shù)教育部重點實驗室(長江大學(xué))， 武漢 430100

0 引言

斷層及斷裂構(gòu)造對于改善儲層物性、聚集成藏等具有重要作用，因此在有利勘探區(qū)開發(fā)過程中，斷層的精細(xì)刻畫與有效識別占有重要地位(徐嘉亮等，2018a；余博文，2021).在復(fù)雜斷裂地區(qū)的勘探開發(fā)，由于受多期構(gòu)造運動影響，斷層較為發(fā)育，斷層類型較多，在解釋過程中對斷層的識別和刻畫的精細(xì)程度要求較高，就導(dǎo)致斷層精細(xì)解釋和斷層平面組合變得更加困難(馬玉歌等，2020；李海晨等，2020).隨著地球物理及計算機相關(guān)技術(shù)的發(fā)展，斷層識別技術(shù)取得較大進步，出現(xiàn)不同的斷層識別方法.現(xiàn)在經(jīng)過國內(nèi)外學(xué)者的大量研究，已經(jīng)形成較多斷層識別方法，如沿層相干屬性技術(shù)、相干體技術(shù)、方差體技術(shù)、曲率體屬性、優(yōu)勢頻帶相位分析技術(shù)、螞蟻追蹤算法、邊緣檢測技術(shù)、面包切片技術(shù)等(Marfurt et al.，1998；Gersztenkorn and Marfurt，1999).相干體技術(shù)是現(xiàn)在使用最廣泛的斷層識別技術(shù)，第一代相干體技術(shù)(C1)主要利用互相關(guān)；第二代相干體技術(shù)(C2)是基于地震道相似性；第三代相干體技術(shù)(C3)是提取協(xié)方差矩陣的特征向量(Bahorich and Farmer，1995)；在此基礎(chǔ)上又衍生出一些其他相干體技術(shù).螞蟻追蹤算法是根據(jù)螞蟻覓食行為而提出的斷層自動識別技術(shù)(陳桂和劉洋，2021).邊緣檢測技術(shù)來源于圖像處理，從最初應(yīng)用到地震勘探領(lǐng)域到現(xiàn)在已經(jīng)發(fā)展了較多的邊緣檢測技術(shù)，并在實際應(yīng)用中取得較好的成果.Hale(2013)在研究斷層掃描、提取斷面和估算斷距時提出最大似然屬性，并在識別斷層時效果較好.此外，在成像過程中斷面波會在斷層處成像，因此增強斷面波能量可以在一定程度上對斷層進行精細(xì)識別與刻畫(徐嘉亮等，2015).

目前對于復(fù)雜構(gòu)造斷層及小尺度斷層，利用識別度較高的斷面波進行斷層及斷面解釋是比較常用的方法(徐嘉亮等，2018b).如何提高斷面波識別度，降低繞射波等干擾成為提高斷面波識別度的關(guān)鍵(劉保金等，2012).本文利用斷層識別度更高的大似然屬性方法對斷層進行了有效識別，并利用拉普拉斯金字塔增強方法對識別的斷面波進行了增強，進而提高了斷面波的識別度.通過在X地區(qū)進行本文方法應(yīng)用，結(jié)果表明在地震剖面上斷面波能量有所提高，得到更加清楚、連續(xù)性更好的斷裂，斷層組合關(guān)系更加清晰，取得較好的應(yīng)用效果.

1 方法原理與關(guān)鍵步驟

1.1 最大似然屬性原理

最大似然法是用統(tǒng)計方法建立一組區(qū)分函數(shù)，在多類別分析中計算出各個類別樣本的歸屬概率，樣本屬于哪類由歸屬概率的大小決定(楊午陽，2017).

基于最大似然屬性識別斷層的原理源于對地震圖像的相似性分析.Hale提出面向斷面波的相似性理論，即應(yīng)用矩形求取視窗中傾角和傾向的相似度.相似度屬性由相似性通過構(gòu)造導(dǎo)向平滑得到，相似性計算的基礎(chǔ)是對地質(zhì)體上沿反射界面不同方向變化率的估計，因此應(yīng)先估計沿反射界面不同方向的反射斜率，然后進行最大似然屬性計算.

1.1.1 估計反射斜率

計算反射斜率方法有很多，其中梯度結(jié)構(gòu)張量(GST)是一種有效的方法在精確表征和分析地震數(shù)據(jù)中反射結(jié)構(gòu)的特征上.其計算原理如下：

假設(shè)原始地震數(shù)據(jù)為S(x,y,t)，其復(fù)數(shù)域的解析信號可表示為

S(x,y,t)=s(x,y,t)+jsH(x,y,t)

=A(x,y,t)ejPhs(x,y,t)，

(1)

(2)

其梯度結(jié)構(gòu)張量為：

(3)

其中w(x,y,t)為權(quán)重系數(shù)，表示對梯度向量的各個分量進行光滑平均，以達到對地震數(shù)據(jù)平滑濾波處理的目的.其表達式為：

(4)

式中：σ為噪聲尺度參數(shù)，通常取值在0.1～3.0之間，x、y、t分別為地震數(shù)據(jù)中沿線方向、道方向和時間方向的變量.

經(jīng)過平滑濾波后地震數(shù)據(jù)的梯度結(jié)構(gòu)張量GST根據(jù)矩陣特征分解，可表示為：

GST(x,y,t)=

(5)

式中：λ1、λ2、λ3為梯度結(jié)構(gòu)張量GST的特征值，且λ1>λ2>λ3>0，α1、α2、α3為梯度結(jié)構(gòu)張量GST的特征值λ1、λ2、λ3對應(yīng)的特征向量.即x、y方向的視傾角θx、θy可表示為：

(6)

(7)

式(6)和(7)的α1x(x,y,t)、α1y(x,y,t)、α1t(x,y,t)分別為式(5)中得到的最大特征值λ1對應(yīng)的特征向量α1(x,y,t)分別在x、y和t方向上的三個分量.

1.1.2 相似性屬性的計算

根據(jù)估計的反射斜率，由式(8)計算相似性：

(8)

式中：u為矩形方向窗口中的地震數(shù)據(jù)；K為矩形方向窗口中地震數(shù)據(jù)的道數(shù)；下標(biāo)k表示時窗中第k道信號；xk、yk分別為第k道信號在x、y方向上的距離；θx、θy分別為x、y方向上的視傾角；θxxk+θyyk為第k道信號在時間軸的時移量.

然而，應(yīng)用此方法計算的相似性，在分子和分母較小的地方會有很大變化.針對這種不穩(wěn)定性需要進行相應(yīng)的平滑處理，平滑處理后即得到相似屬性.假設(shè)時窗大小為(2ω+1)ms，Z=ω/Δt，Δt表示采樣率，此時式(8)可變?yōu)椋?/p>

(9)

式(9)即為進行構(gòu)造導(dǎo)向平滑處理后得到的計算穩(wěn)定性增強的相似屬性.

1.1.3 最大似然屬性的計算

通過矩形方向窗口對斷層傾角和傾向進行掃描得到的相似系數(shù)計算似然屬性，如圖1所示(以2D地震數(shù)據(jù)為例).假設(shè)斷層傾角θ掃描范圍為[θmaxmin]，掃描間隔為Δθ，斷層傾向φ掃描范圍為[φmaxmin]，掃描間隔為Δφ.在求取最大似然屬性過程中，對于某一個采樣點，分別對每一個掃描的傾角和傾向計算相似性，直到沿著相似度為最小值，即得到最大似然屬性，公式為：

Maximun likelihood=max[1-semblance8].

(10)

由公式(10)可知，最大似然屬性的數(shù)值范圍在0～1之間.這樣通過最大似然屬性對斷層及斷裂的有效識別分析，就可以為解釋斷層提供更好的支持.

1.2 建立拉普拉斯金字塔增強原理

拉普拉斯金字塔增強是空域濾波能量增強中銳化濾波增強的一種，其主要是在圖像的拉普拉斯金字塔分解過程中加入映射函數(shù)，得到由高斯金字塔演化而來的中間過程的拉普拉斯金字塔，最后得到需要的細(xì)節(jié)進行增強.

1.2.1 拉普拉斯金字塔分解

地震疊后數(shù)據(jù)的拉普拉斯金字塔構(gòu)建是從高斯金字塔演變而來，因此首先應(yīng)對地震信號進行高斯金字塔分解.假設(shè)原地震信號為F，并以原信號F作為高斯金字塔最底層地震信號G0；再對其進行低通濾波平滑和向下采樣，常采用水平方向和垂向均為1/2，即得到高斯金字塔第一層；不斷重復(fù)上述過程，就獲得縮小尺寸的地震信號；層級越高，地震信號越小，從而構(gòu)成高斯金字塔，其表達式為：

圖1 最大似然屬性計算示意圖(a) 斷面位置及采樣點； (b) 掃描斷面計算相似性； (c) 最小相似性對應(yīng)位置.Fig.1 Schematic diagram of maximum likelihood attribute calculation(a) Section position and sampling point; (b) Similarity of scanning section calculation; (c) Minimum similarity corresponding position.

(11)

式中：N為層數(shù)；Rl為高斯金字塔第l層的行數(shù)；Cl為高斯金字塔第l層的列數(shù)；ω(m,n)是5×5的二維高斯內(nèi)核，其表達式為：

(12)

ω(m,n)實際是高斯低通濾波器,經(jīng)過n次重復(fù)之后得到Gn，此時Gn只包含級少數(shù)的像素點.

1≤l≤N,0≤i≤Rl,0≤j≤Cl，

(13)

式中：N為層數(shù)；ω(m,n)是5×5的二維高斯內(nèi)核.拉普拉斯金字塔第l層地震信號Ll表達式為：

(14)

通過對地震信號的拉普拉斯金字塔分解，可以把原地震信號分解到不同的空間頻帶上，根據(jù)需要對相應(yīng)的細(xì)節(jié)進行增強，如圖2所示.通過拉普拉斯金字塔也可對增強后的地震信號進行重建，其表達式為：

圖2 構(gòu)造過程(a) 高斯金字塔； (b) 對高斯金字塔進行變換； (c) 拉普拉斯金字塔.Fig.2 Construction process(a) Gaussian pyramid; (b) Transform of Gaussian pyramid; (c) Laplace pyramid.

(15)

當(dāng)l=0時，得到重建后地震信號，假如在變換過程中沒有做相關(guān)參數(shù)修改，得到的重建信號與原始地震信號完全相同.

1.2.2 地震數(shù)據(jù)的重映射

在對地震剖面進行初始變化時，首先假設(shè)一個能量參數(shù)ε，如果地震剖面斷面部分的能量變化小于ε，就認(rèn)為是細(xì)節(jié)部分，如果大于ε，就認(rèn)為是邊緣部分.在重映射算法中，當(dāng)計算輸出金字塔系數(shù)時，某點強度(強度是指單通道像素值的大小)可表示為g0=Gl0(x0,y0).在構(gòu)建拉普拉斯金字塔過程中地震斷面波能量在ε附近的被認(rèn)為是細(xì)節(jié)處理，其表達式為：

φd(i)=g0+sign(i-g0)εfd(|i-g0|/ε)，

(16)

式中：g0為當(dāng)前要處理的像素點；sign(i-g0)為亮度差異(i-g0)的符號值；ε可通過限定細(xì)節(jié)部分像素點的個數(shù)≤nd來確定當(dāng)前層中ε的值；fd(x)=ex是映射函數(shù)，用于控制增強細(xì)節(jié)部分；|i-g0|為像素點i與數(shù)據(jù)周圍像素點的差異值；下標(biāo)d表示細(xì)節(jié)detail的首字母.

對重映射處理后的地震剖面構(gòu)造其拉普拉斯金字塔，然后再對增強后的地震剖面進行重建，即可得到斷面波能量增強的地震剖面，有利于后續(xù)的斷層解釋與有利圈閉的識別.

1.3 關(guān)鍵步驟

在實際應(yīng)用過程中，通過最大似然屬性進行斷裂識別和拉普拉斯金字塔斷面增強包括以下步驟：

1.3.1 斷裂特征分析與成像加強

通過對工區(qū)內(nèi)地震數(shù)據(jù)的分析，明確目的層內(nèi)斷裂構(gòu)造在地震上的反射特征.了解工區(qū)的構(gòu)造發(fā)育史，和斷裂形成的時期，以及斷裂特征類型，以便于后續(xù)相關(guān)參數(shù)的合理選擇.

最大似然屬性斷裂識別和拉普拉斯金字塔增強斷面的理論基礎(chǔ)仍然是在地震剖面上反射波同相軸的不連續(xù)性，噪聲、較低的分辨率和地層巖性的突然變化也可能導(dǎo)致反射波同相軸不連續(xù)，因此在地震數(shù)據(jù)處理過程中要注重疊前相關(guān)噪聲的有效壓制和提高地震數(shù)據(jù)的分辨率，在Kirchhoff 疊前時間偏移過程中，要注重偏移孔徑和反假頻因子的合理選擇，在偏移成像后通過隨機噪聲壓制和疊后提高分辨率處理技術(shù)，來進一步提高斷面波成像效果.

1.3.2 最大似然屬性斷裂識別

最大似然屬性原理在上文中已經(jīng)敘述，從地震數(shù)據(jù)體中計算最大似然屬性時，有兩個關(guān)鍵參數(shù)：矩形窗口的水平計算因子、垂直計算因子.水平計算因子和垂直計算因子越大，計算周期越長，計算出最大似然屬性連續(xù)性越好.

1.3.3 拉普拉斯金字塔增強斷面

根據(jù)最大似然屬性對斷層的識別，通過拉普拉斯金字塔進行斷面增強，首先對地震信號進行高斯金字塔搭建，再通過從每一層減去高斯金字塔上一層信號且進行內(nèi)核卷積的地震信號，就得到拉普拉斯金字塔，在搭建拉普拉斯金字塔過程中加入一個映射函數(shù)來增強斷面部分，并用拉普拉斯金字塔對地震信號進行重建，即得到斷面增強后的地震數(shù)據(jù).

2 模型試算

本文以圖3所示的帶有兩條斷層的四層速度模型進行方法測試.該模型在水平和垂直方向上分別有500個網(wǎng)格點，網(wǎng)格間距為2 m，炮間距為20 m，每炮共有400個檢波器接收記錄.第一炮與最后一炮的位置分別在200 m和800 m的位置.該觀測系統(tǒng)能夠保證左側(cè)斷層滿覆蓋接收，右側(cè)斷裂不滿覆蓋接收，進而驗證在滿覆蓋與不滿覆蓋情況下本文方法的實用性.

圖3 帶有兩條斷層的正演速度模型Fig.3 Forward velocity model with two faults

由于正演模型斷層傾角較大，斷面波能量較低.圖4a為Kirchhoff疊前時間偏移原始剖面，圖4b為利用本文方法增強斷層后的疊前時間偏移剖面，兩條斷層能量均有所增強，連續(xù)性變好且識別度變高，更有助于利用斷面波進行斷層解釋.

3 實際應(yīng)用對比

本文將最大似然屬性和拉普拉斯金字塔增強在X地區(qū)進行應(yīng)用，并取得較好效果.X地區(qū)主控及分支斷層在剖面上組合成花狀構(gòu)造、階梯狀構(gòu)造和反Y字型構(gòu)造，同時還伴有一些小的走滑斷裂.在原疊加剖面上這些構(gòu)造位置關(guān)系不明確，造成斷層解釋和平面組合困難.此外，由于河道交錯疊置復(fù)雜，河道邊界的異常響應(yīng)也干擾斷層識別.為落實斷層解釋，采用最大似然屬性識別斷層和拉普拉斯金字塔增強斷面，并取得一定效果.圖5a為增強前測線A的偏移剖面，斷面處連續(xù)性較差且斷層搭接關(guān)系不明確，且存在假斷層響應(yīng)影響斷層解釋精度，圖5b是增強后測線A的疊加剖面，斷面波能量有所增強，斷層連續(xù)性較好且搭接關(guān)系明確.圖6a和圖6b為圖5斷層增強前后剖面的局部方法，由局部放大比較中能夠明顯看出利用本文提出的方法能夠增強斷面波的能量，進而提高斷層的識別度.

圖5 增強前后測線A的疊加剖面對比(a) 原始疊加剖面； (b) 增強后疊加剖面.Fig.5 Comparison of stacked profiles A before and after enhancement(a) Original section； (b) Enhanced section.

圖6 增強前后測線A的疊加剖面對比(a) 原始疊加剖面； (b) 增強后疊加剖面.Fig.6 Comparison of stacked profile A before and after enhancement(a) Original section； (b) Enhanced section.

圖7和圖8是對測線B疊前時間偏移剖面的增強及其布局方法的比較，可看到增強前后斷面波能量明顯提高，斷層連續(xù)性較好且搭接關(guān)系更明確.

圖7 增強前后測線B的疊加剖面對比(a) 原始疊加剖面； (b) 增強后疊加剖面.Fig.7 Comparison of stacked profile B before and after enhancement(a) Original section； (b) Enhanced section.

圖8 增強前后測線B的疊加剖面對比(a) 原始疊加剖面； (b) 增強后疊加剖面.Fig.8 Comparison of stacked profile B before and after enhancement(a) Original section； (b) Enhanced section.

圖9a為該區(qū)T2(館陶組)層段常規(guī)相干體水平切片，可以發(fā)現(xiàn)X地區(qū)主斷裂為近東西向展布，但該相干體橫向分辨能力較低，斷層邊界模糊影響其真實性.圖9b表示斷面波照明度增強后最大似然水平切片，最大似然體橫向分辨能力較高，次級斷裂與主斷裂的接觸關(guān)系較清楚，斷層細(xì)節(jié)更加清楚、連續(xù)，且刻畫斷層邊界更加清晰，斷層組合關(guān)系明確，有利于斷層識別，更加符合斷裂的實際地質(zhì)特征.

圖9 X地區(qū)地震資料平面分析(a) 常規(guī)相干體水平切片； (b) 斷面波照明度增強后相干體水平切片.Fig.9 X area seismic data plane analysis(a) Horizontal slice of conventional coherent body； (b) Horizontal slice of coherent body after enhancing fault-surface wave illumination.

4 討論及結(jié)論

通過分析最大似然屬性識別斷層和基于拉普拉斯金字塔算法增強斷面識別度的原理，以及在X地區(qū)的應(yīng)用，與直接應(yīng)用常規(guī)地震資料對斷層精細(xì)解釋和識別相比，應(yīng)用此方法，可以顯著提高斷面波在疊加剖面上的振幅能量，有利于后續(xù)分析斷層位置以及斷層展布特征、斷層平面組合關(guān)系，較好地實現(xiàn)斷層精細(xì)解釋和斷裂系統(tǒng)的有效識別.然而在最大似然屬性識別斷裂時，要注意矩形窗口的水平計算因子、垂直計算因子的合理選擇，保證運算時間的合理性，且計算的最大似然屬性具有較好的連續(xù)性，同時，在研究中發(fā)現(xiàn)最大似然屬性在某些高振幅的地質(zhì)體邊界處仍具有一定異常響應(yīng)，影響斷層識別，需要進一步改善提高.此外在基于拉普拉斯金字塔算法增強斷面時，由于地震圖像中存在與斷面部分強度變化相似的異常響應(yīng)，在拉普拉斯金字塔算法增強斷面時異常相應(yīng)也會得到增強，對于該算法的不足之處，還需要針對原地震剖面中斷層的特點進行改進，以獲得更好的斷面增強效果.

本文分析了最大似然屬性識別斷層和基于拉普拉斯金字塔算法增強斷面的原理，并在X地區(qū)進行應(yīng)用，即通過最大似然屬性對工區(qū)數(shù)據(jù)體進行掃描，并計算數(shù)據(jù)體采樣點之間的相似性，以獲得工區(qū)內(nèi)最可能發(fā)育斷層的位置及概率，然后應(yīng)用拉普拉斯金字塔增強算法進行斷面增強，即在構(gòu)造拉普拉斯金字塔的過程中加入一個映射函數(shù)，使斷面得到增強.在X地區(qū)的實際應(yīng)用表明，最大似然屬性在疊加剖面上識別斷層位置比較準(zhǔn)確，再通過拉普拉斯金字塔算法增強斷面，可得到斷面波能量較高的疊加剖面，在剖面上斷層平面組合關(guān)系更明確，提高了斷層解釋精度，為后續(xù)有利圈閉識別和井位優(yōu)選提供了依據(jù).