多尺度邊緣檢測技術(shù)在斷層識(shí)別及裂縫發(fā)育帶預(yù)測中的應(yīng)用
——以車排子地區(qū)排691井區(qū)為例

馬承杰

（中國石化勝利油田分公司信息化管理中心，山東東營 257000）

常規(guī)邊緣檢測技術(shù)可以識(shí)別砂體、斷裂、砂巖透鏡體、火成巖巖體等地質(zhì)體的邊界特征［1-2］，但是地下地質(zhì)條件的復(fù)雜性決定了這些地質(zhì)體既具有邊緣特征，又具有多級(jí)多尺度的特征，常規(guī)的邊緣檢查方法對(duì)于斷層識(shí)別與裂縫發(fā)育帶預(yù)測，特別是地震資料品質(zhì)較差條件下的識(shí)別表現(xiàn)出極大的不適應(yīng)性［3］。為此，需要根據(jù)不同的地質(zhì)體選擇特定的方法來檢測其邊緣特征，最大限度利用地震數(shù)據(jù)及其衍生出的多種信息識(shí)別其不同的地質(zhì)現(xiàn)象，提高對(duì)復(fù)雜地質(zhì)體的精準(zhǔn)描述能力。前期的邊緣檢測算法僅根據(jù)地震數(shù)據(jù)，提取相關(guān)屬性參數(shù)，如索伯算子、Canny 算子、Prewitt 算子等，基本沒有考慮地質(zhì)體的多尺度的性質(zhì)。近年來又有學(xué)者根據(jù)不同地質(zhì)體多尺度邊緣響應(yīng)特征提出了小波變換多尺度邊緣檢測方法［4］，在實(shí)際應(yīng)用中取得了較好的效果，但是這些方法僅局限于二維地震剖面的識(shí)別，而且需要進(jìn)行層位解釋的約束，應(yīng)用范圍具有一定的局限性。針對(duì)上述難點(diǎn)，研發(fā)三維檢測算子，建立三維地震數(shù)據(jù)邊緣檢測方法，即三參數(shù)小波變換與結(jié)構(gòu)導(dǎo)向梯度聯(lián)合約束的多尺度邊緣檢測算法，該算法同時(shí)考慮了傾角和方位角導(dǎo)向下沿層檢測。選取準(zhǔn)噶爾盆地西緣車排子地區(qū)排691井區(qū)作為方法應(yīng)用的靶區(qū)，研究區(qū)的三維全覆蓋面積接近350 km2，主要勘探目的層為新近系沙灣組及石炭系，沙灣組多為斷塊圈閉油氣藏，石炭系為裂縫型油氣藏，而研究區(qū)石炭系頂面構(gòu)造也可以作為沙灣組底面構(gòu)造，因此在實(shí)際研究中，選取石炭系頂面作為標(biāo)準(zhǔn)層上下開時(shí)窗進(jìn)行屬性分析，利用所提出的方法開展沙灣組斷層識(shí)別及石炭系裂縫發(fā)育帶預(yù)測，預(yù)測效果良好。

1 三參數(shù)小波變換基本原理

小波變換方法目前應(yīng)用比較廣泛，高靜懷等提出三參數(shù)小波變換算法［5-10］，該方法增加了3個(gè)控制參數(shù)，相對(duì)增加了小波分析的適用性，通過3個(gè)參數(shù)的互相調(diào)節(jié)，提高了應(yīng)用的靈活性，與最佳匹配地震子波小波變換（BMSW）等其他小波變換相比，其具有更好的收斂特點(diǎn)。該小波時(shí)域表達(dá)式為：

對(duì)（1）式作傅氏變換，得到頻率域公式：

假設(shè)任意給定的信號(hào)s(t)∈L2(R)，則三參數(shù)小波變換式為：

為驗(yàn)證建立的三參數(shù)小波變換的時(shí)-頻優(yōu)勢，將三參數(shù)小波變換與BMSW 的時(shí)-頻特征進(jìn)行對(duì)比。取Λ=(1,0.5,0)，計(jì)算得到這2 種方法的時(shí)-頻與振幅的變化關(guān)系（圖1），圖1a 和1b 為BMSW 計(jì)算得到的時(shí)間與振幅和頻率與振幅的關(guān)系，其時(shí)-頻特征出現(xiàn)了多個(gè)峰值，如果用這種出現(xiàn)多個(gè)峰值的小波作為基本小波對(duì)信號(hào)進(jìn)行分析時(shí)就會(huì)產(chǎn)生一些假象，而圖1c 和1d 為三參數(shù)小波變換得到的時(shí)-頻特征圖，其時(shí)-頻特征僅出現(xiàn)一個(gè)峰值，表明三參數(shù)小波變換獲得的結(jié)果在時(shí)-頻域具有相對(duì)唯一性，即表明三參數(shù)小波適于分析不同頻帶分量的信號(hào)。斷層及裂縫發(fā)育帶的地震資料含有頻率及振幅快速變化的分量，BMSW 小波具有模糊斷層及裂縫的時(shí)-頻響應(yīng)特性，難以對(duì)斷層及裂縫發(fā)育邊界做出正確的識(shí)別，而三參數(shù)小波變換卻能較好地解決這一問題。

2 技術(shù)原理

常規(guī)邊緣檢測采用梯度極大值技術(shù)，具有算法簡單、程序易于實(shí)現(xiàn)、計(jì)算量較小等特點(diǎn)，被廣泛應(yīng)用于地震資料信噪比較高情況下的河道、砂體、斷層等地質(zhì)體邊界的識(shí)別中。但對(duì)于地震資料品質(zhì)較差地區(qū)的斷層識(shí)別及裂縫發(fā)育帶預(yù)測難以取得良好的效果，究其原因是采用絕對(duì)梯度值，導(dǎo)致裂縫邊界的弱反射信息被強(qiáng)地層背景邊界信息所淹沒，在地層邊界信息的干擾下，難以有效檢測小斷層及微裂縫邊界。而利用三參數(shù)小波變換計(jì)算不同頻帶及不同尺度的地震數(shù)據(jù)體，在此基礎(chǔ)上求取地層傾角和方位角，然后作為約束，利用索伯算子應(yīng)用于小斷層及微裂縫邊界檢測。三參數(shù)小波變換與結(jié)構(gòu)導(dǎo)向梯度聯(lián)合約束的多尺度邊緣檢測方法實(shí)現(xiàn)步驟包括：①利用（6）式計(jì)算出不同尺度的分頻地震數(shù)據(jù)體。②對(duì)這些不同尺度的分頻數(shù)據(jù)體，計(jì)算地層傾角和方位角。③根據(jù)索伯算子，針對(duì)不同尺度的分頻數(shù)據(jù)體，計(jì)算結(jié)構(gòu)導(dǎo)向梯度屬性來進(jìn)行地質(zhì)體的邊緣檢測。

圖1 三參數(shù)小波變換與BMSW的時(shí)-頻特征對(duì)比Fig.1 Comparison between time-frequency characteristics of three-parameter wavelet transform and BMSW

通常情況下地層傾角測井可以測量出地層傾角和方位角，而傾角和方位角是地震幾何屬性中必備的起重要作用的地質(zhì)分析參數(shù)，如可以利用地層傾角和方位角確定構(gòu)造特征、分析斷層、不整合、層理、沙壩、礁灘、鹽丘等的構(gòu)造變形及識(shí)別地層裂縫和破碎帶等。因此本文利用傾角和方位角來計(jì)算結(jié)構(gòu)導(dǎo)向梯度屬性體，通過（7）—（10）式推導(dǎo)出計(jì)算傾角和方位角體的（11）和（12）式，將計(jì)算結(jié)果帶入（13）式進(jìn)行計(jì)算，最終獲得結(jié)構(gòu)導(dǎo)向梯度屬性體。

圖2為地震數(shù)據(jù)體中的一個(gè)小塊體的傾角和方位角幾何模型，根據(jù)圖中的幾何關(guān)系，傾角和方位角表達(dá)式分別為：

圖2 傾角和方位角模型與差分計(jì)算的網(wǎng)格單元Fig.2 Model of dip and azimuth and grid elements of differential calculation

但實(shí)際計(jì)算中很難從三維地震數(shù)據(jù)體中獲得視傾角，為此采用3×3 網(wǎng)格單元，由一階導(dǎo)數(shù)的定義，用差分計(jì)算的方法分別寫為：

利用（11）和（12）式計(jì)算出三維地震數(shù)據(jù)體的傾角體和方位角體。

在此基礎(chǔ)上，建立結(jié)構(gòu)導(dǎo)向的索伯算子，其表達(dá)式為：

將地震數(shù)據(jù)體分成多個(gè)小塊體，將這些小塊體當(dāng)成一個(gè)目標(biāo)點(diǎn)來進(jìn)行計(jì)算，如通常選取小塊體為3 道×3 道或者5 道×5 道等，根據(jù)實(shí)際情況可以任意選取。通過（7）式分別求取每個(gè)小塊體x方向和y方向的索伯梯度值，并對(duì)2 個(gè)方向的梯度值進(jìn)行平方求和，然后沿著法線方向?qū)λ行K體求得的索伯梯度相加求和，最后做歸一化能量均衡處理。其中：

利用（6）式，結(jié)合（7）—（18）式即可實(shí)現(xiàn)三參數(shù)小波變換與結(jié)構(gòu)導(dǎo)向梯度聯(lián)合約束的多尺度邊緣檢測方法，進(jìn)而實(shí)現(xiàn)多尺度三維地震資料的邊緣檢測。

3 應(yīng)用實(shí)例

準(zhǔn)噶爾盆地西緣車排子地區(qū)排691井區(qū)發(fā)育古近系沙灣組斷塊圈閉油氣藏及石炭系頂部風(fēng)化殼裂縫型油氣藏，為了較好地描述該區(qū)沙灣組斷塊圈閉油藏及石炭系裂縫油藏的有利分布區(qū)［11-17］，利用邊緣檢測技術(shù)對(duì)研究區(qū)2套目的層的斷裂進(jìn)行識(shí)別及對(duì)裂縫發(fā)育帶進(jìn)行預(yù)測。

圖3 排691井區(qū)沿石炭系頂向下50 ms時(shí)窗常規(guī)邊緣檢測相干圖Fig.3 Coherence map of conventional edge detection below time window of 50 ms along top of the Carboniferous in Block Well P691

以25 m×25 m 網(wǎng)格精度精細(xì)解釋了研究區(qū)三維地震資料的石炭系頂面反射層，并沿該反射層上下各開50 ms時(shí)窗，基本可以涵蓋2套目的層段。利用該數(shù)據(jù)體通過常規(guī)邊緣檢測技術(shù)計(jì)算獲得的石炭系頂面的相干圖（圖3）顯示，相干整體效果較差，平面上斷層邊界模糊導(dǎo)致特征難以識(shí)別，僅有幾條較大的斷層有一定的顯示，無法描述研究區(qū)斷裂體系的完整平面展布特征，不利于對(duì)沙灣組斷塊圈閉的描述。同時(shí)石炭系裂縫發(fā)育特征的識(shí)別也不夠理想，反映裂縫發(fā)育的黑色條帶也模糊不清，無法反映石炭系風(fēng)化殼裂縫型油氣藏的平面發(fā)育特征。為此，該地區(qū)嘗試應(yīng)用多尺度邊緣檢測技術(shù)進(jìn)行斷層識(shí)別及裂縫發(fā)育帶預(yù)測，首先利用三參數(shù)小波變換計(jì)算得到大中小3 個(gè)尺度的地震數(shù)據(jù)體，分別計(jì)算傾角體和方位角體，利用（13）式計(jì)算出3 個(gè)尺度的結(jié)構(gòu)導(dǎo)向梯度屬性。由計(jì)算得到的3個(gè)尺度的結(jié)構(gòu)導(dǎo)向梯度屬性邊緣檢測效果（圖4）可見，整體上看，利用多尺度邊緣檢測技術(shù)得到的不同尺度的邊緣檢測效果圖斷層及裂縫的邊界都比較清晰，噪音干擾較少，尤其是對(duì)于斷層的刻畫較有利，解釋人員可以較好地對(duì)平面斷層的展布進(jìn)行組合，勾繪出研究區(qū)不同級(jí)別斷層的分布特征。從細(xì)節(jié)上來看，大尺度低頻率的預(yù)測圖反映的大斷裂特征更加清晰，研究區(qū)主要斷層清晰可辨，尤其是斷層的展布方向、斷層的組合、斷層的延伸及斷層的平面排列特征等較好的刻畫出來（圖4a），有利于沙灣組斷塊圈閉描述。中尺度主頻的預(yù)測圖反映較小斷層展布更加明顯（圖4b），可以將中等發(fā)育斷層及裂縫發(fā)育區(qū)識(shí)別出來。作為大尺度和小尺度的一種補(bǔ)充，可以完整描述出研究區(qū)的斷層分布形態(tài)，尤其是斷層組合過程中的某兩條斷層的接觸關(guān)系，中尺度主頻的預(yù)測結(jié)果可以作為參考和輔證。而小尺度則反映裂縫的集中發(fā)育帶。小尺度高主頻的預(yù)測圖將研究區(qū)裂縫發(fā)育特征展示得比較清晰，裂縫帶主要發(fā)育在研究區(qū)中部區(qū)帶（圖4c），因此可以利用該圖進(jìn)行石炭系裂縫發(fā)育有利區(qū)帶預(yù)測。這樣將3個(gè)尺度的邊緣檢測結(jié)果進(jìn)行組合解釋，一方面可以較好地解釋出研究區(qū)不同尺度的斷層以及斷裂組合，同時(shí)也可以明確研究區(qū)的裂縫分布特征，從而預(yù)測裂縫發(fā)育有利區(qū)，為研究區(qū)斷塊圈閉油氣藏及裂縫性油氣藏的分布區(qū)預(yù)測提供較好的基礎(chǔ)研究資料。

圖4 排691井區(qū)三參數(shù)小波變換與結(jié)構(gòu)導(dǎo)向梯度聯(lián)合約束的多尺度邊緣檢測結(jié)果Fig.4 Multi-scale edge detection results under joint constraint of three-parameter wavelet transform and structure-oriented gradient in Block Well P691

4 結(jié)論

三參數(shù)小波變化可以將地震資料進(jìn)行分頻得到不同頻帶的地震數(shù)據(jù)體，進(jìn)而計(jì)算得到不同頻率數(shù)據(jù)體的傾角體和方位角體，再利用改進(jìn)的索伯算子計(jì)算得到結(jié)構(gòu)導(dǎo)向梯度屬性，最后建立了三參數(shù)小波變換與結(jié)構(gòu)導(dǎo)向梯度聯(lián)合約束的多尺度邊緣檢測技術(shù)。結(jié)果表明，該方法理論技術(shù)是可行的，在斷層及裂縫發(fā)育帶的檢測方面具有較好的適用性，尤其是不同尺度下的邊緣檢測方法突破了傳統(tǒng)方法的局限性，具多尺度優(yōu)勢，可以充分挖掘地震資料包含的豐富的地質(zhì)信息，為斷塊圈閉油氣藏及裂縫型油藏勘探提供了新技術(shù)手段，以期為勝利油田西部新區(qū)及東部老區(qū)新層系勘探提供指導(dǎo)。

符號(hào)解釋

a——尺度因子；

b——平移因子；

D——小塊體分析窗口所有樣點(diǎn)在三參數(shù)小波變換后不同頻率的地震數(shù)據(jù)中的取值；

i——橫向測線偏移量；

j——縱向測線偏移量；

k——時(shí)間域約束系數(shù)；

K——上塊體個(gè)數(shù)，2K+1為塊體總數(shù)；

L2——平方可積函數(shù)空間；

Mx，My——x和y方向的索伯算子模板；

n——上下滑動(dòng)時(shí)窗數(shù)量；

p(Λ),k(Λ),q(Λ)——待定函數(shù)；

R——實(shí)數(shù)域集合；

s(t)——任意的信號(hào)；

Sφ——三參數(shù)小波變換式；

Sobel（x，y，t）——結(jié)構(gòu)導(dǎo)向的索伯算子；

t——時(shí)間，ms；

Wx——振幅在x方向的差分值；

Wy——振幅在y方向的差分值；

x——橫向測線值；

y——縱向測線值；

z——三維地震數(shù)據(jù)道振幅；

z1—z9——差分格式網(wǎng)格點(diǎn)的振幅；

β——能量延遲因子；

Δt——采樣率；

Δx，Δy——網(wǎng)格點(diǎn)間的在x和y方向上的距離，即為道間距；

θ——陰影地層反射界面的傾角，（°）；

θx，θy——x和y方向的視傾角，（°）；

Λ——σ,τ,β的一組集合；

σ——小波的調(diào)制頻率；

τ——能量衰減因子；

τx——用傾角體和方位角體計(jì)算的沿x方向的延遲時(shí)間，ms；

τy——用傾角體和方位角體計(jì)算的沿y方向的延遲時(shí)間，ms；

(x,y,t)——橫向測線方向相同層位的延遲時(shí)間，ms；

(x,y,t)——縱向測線方向相同層位的延遲時(shí)間，ms；

φ——陰影地層反射界面的方位角，（°）；

φ(t;Λ)——時(shí)域小波函數(shù)；

φ?(ω;Λ)——頻率域小波函數(shù)；

φ*——復(fù)共軛；

ω——頻率。