基于多尺度加窗希爾伯特變換的地震資料體邊緣檢測

李 斌,陳學(xué)華,賀振華,許 迪,李依佳,裴小剛,唐湘蓉

(1.成都理工大學(xué)油氣藏地質(zhì)及開發(fā)工程國家重點實驗室,四川成都610059;2.成都理工大學(xué)地球探測與信息技術(shù)教育部重點實驗室,四川成都610059)

基于多尺度加窗希爾伯特變換的地震資料體邊緣檢測

李 斌1,2,陳學(xué)華1,2,賀振華1,2,許 迪2,李依佳2,裴小剛2,唐湘蓉2

(1.成都理工大學(xué)油氣藏地質(zhì)及開發(fā)工程國家重點實驗室,四川成都610059;2.成都理工大學(xué)地球探測與信息技術(shù)教育部重點實驗室,四川成都610059)

利用地震資料識別斷層或裂縫發(fā)育帶等非連續(xù)性異常信息是裂縫性油氣藏儲層預(yù)測的關(guān)鍵所在。為此,提出了一種基于多尺度加窗希爾伯特變換的體邊緣檢測新方法。在二維(水平)加窗希爾伯特變換的基礎(chǔ)上,考慮了非連續(xù)性地質(zhì)異常的多尺度三維特征,通過同時調(diào)節(jié)水平和深度方向上的計算孔徑來提取實際三維地質(zhì)異常體的邊緣信息,彌補了二維加窗希爾伯特變換僅提取水平方向邊緣信息的明顯不足。實際地震資料體邊緣檢測試處理結(jié)果表明,該方法能夠完整地刻畫多尺度非連續(xù)性地質(zhì)異常體的特征,有效突出裂縫發(fā)育帶的邊緣位置及斷層的走向,具有較為顯著的實用價值。

體邊緣檢測;廣義希爾伯特變換;多尺度特征;水平孔徑;深度孔徑

勘探實踐表明,含油氣儲層裂縫發(fā)育非常普遍,而大量空隙或裂縫組成的邊緣特征信息是可檢測的。傳統(tǒng)的邊緣檢測方法,如利用Robert算子、Sobel算子、Canny算子等的微分法和擬合法等,由于其計算簡單、定位精度高等特點,已被應(yīng)用于同相軸追蹤[1]。隨著可視化和圖像處理技術(shù)的發(fā)展,源于圖像和信號分析領(lǐng)域的邊緣檢測技術(shù)也在地震裂縫檢測中獲得了廣泛的應(yīng)用[2]。如Bahorich等[3]1995年提出的地震相干數(shù)據(jù)體分析技術(shù)；賀振華等[4]1999年提出的基于地下介質(zhì)橫向變化的地震多尺度邊緣檢測技術(shù)。2001年,Robert[5]首次將曲率作為一種地震屬性應(yīng)用到地震資料解釋中,給出了層面屬性的定義和計算方法,此為第一代曲率屬性;2003年,Luo等[6]提出了一種廣義希爾伯特變換(Generalized Hilbert Transform,GHT),并將其應(yīng)用于地震資料中河道信息的成像;黨志敏等[7]分析了GHT在含噪信號邊緣檢測中的應(yīng)用效果;陳學(xué)華等[8]于2011年提出基于廣義S變換的分頻裂縫邊緣檢測方法。

由于傅里葉變換是一種全局性變換,不能有效突出非平穩(wěn)信號的局部特征[9-10]；而圖像是二維非平穩(wěn)信號,建立在傅里葉變換基礎(chǔ)上的傳統(tǒng)希爾伯特變換雖然具有唯一性,可以檢測到圖像邊緣信息,但在高精度和多尺度方面不能得到滿意的結(jié)果。三維地震資料中斷層或裂縫發(fā)育帶等不連續(xù)信息在內(nèi)部結(jié)構(gòu)上存在復(fù)雜多樣性,在空間上存在多尺度特征,如果僅考慮二維(水平)加窗希爾伯特變換,則難以同時刻畫地質(zhì)異常在水平和深度方向上的多尺度特征,無法給出其完整的地質(zhì)信息[11]。

為此,本文提出一種基于時間域多尺度加窗希爾伯特變換的體邊緣檢測新方法,將二維(水平)加窗邊緣檢測拓展到三維(體)加窗邊緣檢測,同時在水平和深度方向上采用不同的孔徑計算地質(zhì)異常體邊緣信息,突出不同尺度下不連續(xù)性信息的完整異常特征。實際地震資料的體邊緣檢測試處理結(jié)果證明了新方法的有效性。

1 方法原理及實現(xiàn)流程

1.1 加窗希爾伯特變換的基本原理

希爾伯特變換是一種全通濾波器,原始信號經(jīng)過希爾伯特變換后可以得到相應(yīng)的解析信號。對于離散信號x(t),希爾伯特變換的濾波過程可以表述為:

(1)

式中:x(t)為輸入信號;y(t)為希爾伯特變換后的輸出信號;h(t)為濾波因子。

因此,時間域的希爾伯特變換可表示為一個褶積關(guān)系：

(2)

式中:x(n)為原始信號;y(n)為時間域希爾伯特變換后的信號;h(n)為希爾伯特因子,是一個無窮序列,并且隨n的增加而逐漸衰減。h(n)表達式為:

(3)

其中,m為整數(shù)。

理想情況下濾波因子無限長,濾波后的結(jié)果準(zhǔn)確,但在實際應(yīng)用中,數(shù)字濾波所能處理的濾波因子只能是有限長的,所以將濾波因子截斷為有限長。濾波因子截斷后會出現(xiàn)吉普斯現(xiàn)象,使濾波后的結(jié)果發(fā)生畸變,對濾波結(jié)果產(chǎn)生干擾。

為了解決常規(guī)希爾伯特變換對噪聲敏感、抗噪能力差的問題,Luo等[6]提出了廣義希爾伯特變換,引入了窗函數(shù)和階數(shù),從兩個方面對傳統(tǒng)的希爾伯特變換進行了擴展。廣義希爾伯特變換可以看作是加窗希爾伯特變換,通過使用兩端逐漸衰減的窗函數(shù),可以使濾波因子逐漸變?yōu)榱?從而減小截斷效應(yīng),達到理想的效果。

加窗處理通過局部增益加權(quán),使得有效邊緣信息更加突出,和非邊緣信息產(chǎn)生明顯差異。選取不同窗函數(shù)會有不同的效果[12],窗函數(shù)長度因子的選擇是控制水平方向多尺度分辨率的關(guān)鍵,大尺度和小尺度因子能分別刻畫出不同大小的地質(zhì)體特征,突顯出不同尺度的不連續(xù)信息。

本文的仿真模型和實際數(shù)據(jù)處理選用漢寧窗(即升余弦窗)函數(shù),其表達式為:

(4)

式中:N為自然數(shù),即控制水平方向分辨率的多尺度因子(時窗因子)。N越大,尺度越大,分辨率越低;N越小,尺度越小,分辨率越高。當(dāng)N=31時,升余弦窗函數(shù)如圖1所示。

圖1 漢寧窗函數(shù)

將升余弦窗與希爾伯特因子相乘：

(5)

其中,h′(n)為加窗希爾伯特算子,隨著時窗因子的改變,h′(n)也發(fā)生改變。加窗希爾伯特變換可表示為：

(6)

其中,y(n)為加窗希爾伯特變換后的信號,隨著h′(n)的變化而變化,所以也具有多尺度特征。

圖2為階梯形邊緣和斜坡形邊緣模型的一維模擬信號;圖3是利用不同長度窗函數(shù)計算的希爾伯特變換結(jié)果,可見邊緣處的絕對幅值相對極大,位置對應(yīng)準(zhǔn)確,而無邊緣信息的部分絕對幅值相對極小。通過比較可以看出,加窗希爾伯特變換準(zhǔn)確地提取出了模擬信號的邊緣信息,且隨著窗長度N的增加,邊緣變粗,分辨率降低,反之亦然。因此,窗長度N可用于調(diào)節(jié)邊緣檢測結(jié)果的分辨率(與不同地質(zhì)體的尺度大小相對應(yīng)),檢測具有多尺度信息的不連續(xù)性特征。

圖2 階梯形邊緣和斜坡形邊緣模型的一維模擬信號

圖3 邊緣模型一維信號的加窗希爾伯特變換邊緣提取

1.2 基于加窗希爾伯特變換的多尺度體邊緣檢測

由于二維(水平)邊緣檢測是針對目標(biāo)層段的單一樣點計算的,反映的只是斷層或裂縫發(fā)育帶在這一層位橫向上的變化,為水平方向不連續(xù)信息的多尺度顯示,能夠提取的邊緣特征是極其有限的。而斷層或裂縫發(fā)育帶等不連續(xù)信息的邊緣特征在深度方向也是有延伸的,存在深度方向的多尺度特征，深度方向選取的孔徑大小不同,突顯出的地質(zhì)體信息各不相同。所以,我們提出基于多尺度加窗希爾伯特變換的體邊緣檢測方法,同時在水平方向和深度方向上選取不同孔徑進行三維多尺度體邊緣檢測[13],以保證所提取儲層內(nèi)地質(zhì)異常體信息的完整性。

同時沿水平方向和深度方向采用不同孔徑進行多尺度體邊緣檢測的計算公式為:

(7)

式中:z(n,k,t)表示多尺度體邊緣檢測結(jié)果;i表示深度方向以目標(biāo)層位為基準(zhǔn)向上、下延拓的樣點數(shù)(正值向下延拓,負(fù)值向上延拓);-k和k為深度方向延拓樣點數(shù)的上、下限(深度方向孔徑樣點數(shù)為2k+1);dt+i·Δt,x(x,y)和dt+i·Δt,y(x,y)分別表示沿x和y方向?qū)Φ卣饠?shù)據(jù)進行處理(t表示目標(biāo)層位,Δt表示相鄰樣點之間的采樣間隔);h′(n)為(5)式所示的加窗希爾伯特算子。

1.3 多尺度體邊緣檢測方法實現(xiàn)步驟

基于多尺度加窗希爾伯特變換體邊緣檢測方法的計算流程如下:

1) 根據(jù)給定目標(biāo)層位,從三維地震數(shù)據(jù)中提取目標(biāo)層段的地震數(shù)據(jù)d(x,y),其中x為工區(qū)的總測線數(shù),y為工區(qū)單條測線的道數(shù);

2) 對d(x,y)進行保邊平滑濾波(edge preserving smoothing,EPS)[14]預(yù)處理,以保證斷層、裂縫等不連續(xù)性信息的地震反射去除噪聲以及其它干擾的影響,得到預(yù)處理結(jié)果dEPS(x,y);

3) 選取不同的計算孔徑(水平方向的窗長N及深度方向的樣點數(shù)i),利用公式(7)對dEPS(x,y)進行多尺度加窗希爾伯特變換體邊緣檢測計算,得到邊緣檢測切片集z(n,k,t);

4) 對體邊緣檢測結(jié)果z(n,k,t)進行顯示。

2 實際地震資料處理效果分析

以南海某海區(qū)三維地震資料為例進行方法測試。該區(qū)域目標(biāo)層段裂縫發(fā)育,不連續(xù)信息明顯,圖4是提取的研究工區(qū)目標(biāo)層位沿層切片。

采用本文提出的多尺度加窗體邊緣檢測方法對目標(biāo)層位沿層切片進行邊緣檢測試處理,圖5是水平孔徑固定為11,深度孔徑分別為1，5和11的計算結(jié)果。其中深度孔徑為1的計算結(jié)果(圖5a)實際上即為傳統(tǒng)的二維(水平)加窗希爾伯特變換檢測結(jié)果。通過對比可以看出,無論是在分辨率還是信噪比方面,本文提出的三維(體)加窗多尺度邊緣檢測結(jié)果(圖5b,圖5c)均明顯優(yōu)于二維(水平)加窗邊緣檢測結(jié)果(圖5a)。由于考慮了深度方向的異常信息,三維(體)加窗希爾伯特變換結(jié)果突出了不同尺度下不連續(xù)性信息的完整異常特征,能夠更有效地描述裂縫發(fā)育帶的分布及斷層的位置和走向,由此證明了基于多尺度加窗希爾伯特變換體邊緣檢測方法的有效性。

圖4 南海某海區(qū)三維地震資料目標(biāo)層位沿層切片

通過圖5也可以看出,采用本文方法選取不同的深度孔徑可以得到不同的檢測效果。深度孔徑較小時(圖5b)“背景”較為干凈,大的裂縫發(fā)育帶等不連續(xù)性信息比較清晰,分辨率相對較高;隨著深度孔徑變大(圖5c),裂縫、斷層的連續(xù)性和延展性更加完整,能更好地反映斷層和裂縫的走向與分布情況。

選取不同的水平孔徑也會得到不同的檢測效果。圖6是深度孔徑固定為11,水平孔徑分別為5,31和101得到的三維(體)加窗邊緣檢測結(jié)果。對比可見小水平孔徑的優(yōu)勢在于斷層和裂縫顯示清晰、分辨率高及易于主要異常的定位;而大水平孔徑可以突出高信噪比,更好地反映不連續(xù)性地質(zhì)異常的完整特征。

圖5 南海某海區(qū)三維地震資料采用傳統(tǒng)方法與本文方法的邊緣檢測結(jié)果(水平孔徑均為11)

圖6 南海某海區(qū)三維地震資料采用本文方法的體邊緣檢測結(jié)果(深度孔徑均為11)

3 結(jié)束語

在二維(水平)加窗希爾伯特變換的基礎(chǔ)上,本文提出了基于多尺度加窗希爾伯特變換的體邊緣檢測方法。通過同時在水平和深度方向上采用不同的孔徑提取地質(zhì)異常體的邊緣信息,可以在不同尺度刻畫出不連續(xù)性三維地質(zhì)異常體的完整特征,有效地突出裂縫發(fā)育帶的邊緣位置及斷層的位置和走向。實際地震資料的試處理結(jié)果表明,作為一種新的地震資料邊緣檢測方法,本文方法具有較為顯著的實用價值,可望在裂縫性油氣藏的地震儲層預(yù)測中發(fā)揮效用。在水平和深度方向上采用不同孔徑能夠得到不同的檢測效果,所以在實際應(yīng)用中需要對不同孔徑的體邊緣檢測圖像進行對比和綜合分析,選取最佳的水平孔徑和深度孔徑,以提高對地震資料不連續(xù)性信息檢測的完整性和準(zhǔn)確性。

致謝：感謝成都理工大學(xué)優(yōu)秀科研創(chuàng)新團隊培育計劃(KYTD201410)對本文研究工作的支持。

[1] 李紅星,劉財,陶春輝.圖像邊緣檢測方法在地震剖面同相軸自動檢測中的應(yīng)用研究[J].地球物理學(xué)進展,2007,22(5):1607-1610 Li H X,Liu C,Tao C H.The study of application of edge measuring technique to the detection of phase axis of the seismic setion[J].Progress in Geophysics,2007,22(5):1607-1610

[2] 羅夢貞.傳統(tǒng)的圖像邊緣檢測算法的分析與比較[J].電腦知識與技術(shù),2011,7(16):3907-3909 Luo M Z.The analysis and comparison of traditional image edge detection optimization[J].Computer Knowledge and Technology,2011,7(16):3907-3909

[3] Bahorich M,Farmer S.3D seismic discontinuity for faults and stratigraphic features:the coherence cube[J].The Leading Edge,1995,14(10):1053-1058

[4] 賀振華,黃捍東,胡光岷,等.地下介質(zhì)橫向變化的地震多尺度邊緣檢測技術(shù)[J].物探化探技術(shù),1999,24(4):289-294 He Z H,Huang H D,Hu G M,et al.Lateral identification by seismic multi-scale edge detection[J].Computing Techniques for Geophysical and Geochemical,1999,24(4):289-294

[5] Roberts A.Curvature attributes and their application to 3D interpreted horizons[J].First Break,2001,19(2):85-100

[6] Luo Y,Saleh A,Marhoon M,et al.Generalized Hilbert transform and its applications in geophysics[J].The Leading Edge,2003,22(3):198-202

[7] 黨志敏,賀振華,黃德濟.GHT在含噪信號邊緣檢測中的應(yīng)用效果分析[J].石油工業(yè)計算機應(yīng)用,2008,16(2):19-21 Dang Z M,He Z H,Huang D J.GHT application results analysis in noisy signal edge detection analysis[J].Oil Industry Computer Applications,2008,16(2):19-21

[8] 陳學(xué)華,賀振華,文曉濤,等.基于廣義S變換的裂縫分頻邊緣檢測方法[J].吉林大學(xué)學(xué)報(地球科學(xué)版),2011,41(5):1605-1609 Chen X H,He Z H,Wen X T,et al.Fracture multi-frequency edge detection based on generalized S transform[J].Jilin University Journal(Earth Sciences),2011,41(5):1605-1609

[9] 張賢達,保錚.非平穩(wěn)信號分析與處理[M].北京:國防工業(yè)出版社,1998:20-25 Zhang X D,Bao Z.Non-stationary signal analysis and processing[M].Beijing:Defense Industry Press,1998:20-25

[10] 陳學(xué)華.時頻分布與地震信號譜分析研究[D].成都:成都理工大學(xué)信息工程學(xué)院,2006 Chen X H.Time-frequency distribution and investigation in spectral analysis of seismic signal[D].Chengdu:College of Information Engineering,Chengdu University of Technology,2006

[11] 陳學(xué)華,賀振華,黃德濟.地震資料的高階偽希爾伯特變換邊緣檢測[J].地球物理學(xué)進展,2008,23(4):1106-1110 Chen X H,He Z H,Huang D J.Seismic data edge detection pseudo Hilbert based on higher-order transform[J].Progress in Geophysics,2008,23(4):1106-1110

[12] 劉麗娟,王山山.廣義S變換窗函數(shù)的分析和改進[J].巖性油氣藏,2007,19(2):76-79 Liu L J,Wang S S.Analysis and improvement of window function of generalized S-transform[J].Lithologic Reservoirs,2007,19(2):76-79

[13] 伍鵬,陳學(xué)華,賀振華,等.基于小波分頻的地震資料相干屬性對比分析[J].長江大學(xué)學(xué)報,2011,8(1):39-42 Wu P,Chen X H,He Z H,et al.Comparative analysis of seismic data coherence properties of wavelet-based crossover[J].Yangtze University Journal,2011,8(1):39-42

[14] Luo Y,Maher M.Edge-preserving smoothing and applications[J].The Leading Edge,2002,21(2):136-141

[15] 趙巖,賀振華,黃德濟.基于保邊去噪的高精度相干分析應(yīng)用研究[J].長江大學(xué)學(xué)報,2009,6(4):169-171 Zhao Y,He Z H,Huang D J.Based on analysis of the coherent application of precision edge-preserving denoising[J].Yangtze University Journal,2009,6(4):169-171

(編輯：戴春秋)

Seismic data 3D edge detection based on multi-scale windowed Hilbert transform

Li Bin1,2,Chen Xuehua1,2,He Zhenhua1,2,Xu Di2,Li Yijia2,Pei Xiaogang2,Tang Xiangrong2

(1.StateKeyLaboratoryofOil&GasReservoirGeologyandExploitation,ChengduUniversityofTechnology,Chengdu610059,China; 2.KeyLaboratoryofEarthExplorationandInformationTechnology,MinistryofEducation,ChengduUniversityofTechnology,Chengdu610059,China)

The identification of discontinuity abnormal information in seismic data,such as the seismic responses from faults or fracture zone,is crucial to fractured oil-gas reservoir prediction.Therefore,we proposed a new 3D edge detection method based on multi-scale windowed Hilbert transform.It is based on 2D (horizontal) windowed Hilbert transform，considering the characteristics of 3D multi-scale for discontinuous geological anomalies.The edge information of 3D geological anomalies is calculated by simultaneously adjusting the aperture in the horizontal and depth direction,thus overcame the obvious inadequacy of edge information extracting only in the spatial direction for 2D (horizontal) windowed Hilbert transform.Trial processing results of the actual seismic data for 3D edge detection show that the method can completely describe the characteristics of multi-scale discontinuous geological anomalies,effectively reveal the location and direction of fractures developed zone and faults with more significant practical value.

3D edge detection,generalized Hilbert transform,multi-scale features,horizontal space aperture,depth aperture

2014-06-17;改回日期：2014-10-19。

李斌(1990—)，男，碩士在讀，主要從事地震信號分析與儲層結(jié)構(gòu)地震檢測方面的研究。

國家自然科學(xué)基金項目(41374134，41174114)、四川省杰出青年科學(xué)基金(2013JQ0011)和國家科技重大專項項目(2011ZX05023-005-010)聯(lián)合資助。

P631

A

1000-1441(2015)03-0345-05

10.3969/j.issn.1000-1441.2015.03.014