基于相空間濾波策略的地震瞬時屬性提取方法
【技術領域】
[0001] 本發(fā)明屬于地球物理勘探中的信號處理領域，涉及地震資料的瞬時屬性提取，尤 其涉及一種基于相空間濾波策略的地震瞬時屬性提取方法。
【背景技術】
[0002] 地震屬性是指由地震數(shù)據(jù)經(jīng)過數(shù)學變換導出的有關地震波的幾何形態(tài)、運動學特 征、動力學特征和統(tǒng)計學特征，這些特征能夠從視覺上反映儲層的形態(tài)及其含油氣性。20世 紀70年代，地震屬性分析開始被引入地震解釋中。20世紀90年代以來，由于儲層描述和三 維數(shù)據(jù)體解釋的需要，地震屬性分析技術迅速發(fā)展，地震屬性在地層構造解釋、儲層巖性和 物性特征描述、油氣藏預測與動態(tài)監(jiān)視方面得到了廣泛應用。目前還沒有一個公認的地震 屬性分類，按照物理意義，地震屬性可以分為時間、振幅、頻率、相干、衰減等幾大類，按照屬 性拾取方法，每類又可以分為瞬時屬性、單道時窗屬性、多道時窗屬性。其中瞬時屬性是在 地震波到達的位置上拾取的屬性，包括瞬時振幅、瞬時相位、瞬時頻率、瞬時帶寬等等。近年 來，地震屬性的數(shù)量突增，然而，瞬時屬性仍然是地震數(shù)據(jù)地質解釋的支柱，也是目前商業(yè) 處理軟件必備的技術模塊。本發(fā)明主要研宄地震資料瞬時屬性的提取。
[0003] 在信號處理領域，以瞬時振幅、瞬時相位、瞬時頻率為代表的信號瞬時屬性的概念 由來已久，在Gabor和Ville的工作之后，大量學者在相關領域做了研宄，Boashash對這 些工作做了綜述。Taner等人于1979年引入復地震道分析，提出通過復地震道求取瞬時屬 性的方法，并給出了這些屬性的物理意義以及在地震解釋中的應用。瞬時振幅與相鄰層的 巖性變化及油氣聚集有關。瞬時相位反映界面的不連續(xù)性、斷層、不整合面和層序邊界等。 瞬時頻率的變化可以有效刻畫地層的厚度和巖性變化，指示油氣的分布等。Robertson和 Nogami將瞬時屬性用于孔隙砂巖薄層厚度估計，Chopra和Marfurt利用瞬時屬性進行不連 續(xù)性、斷層和橫向不連續(xù)性檢測。Liu和Marfurt利用瞬時屬性檢測和刻畫曲流河的分布并 確定其厚度。Zeng利用瞬時頻率異常來指示薄層。Gao等人利用瞬時頻率進行地震資料Q 值估計。
[0004] 常用的提取地震瞬時屬性的方法可以分為以下三類：
[0005] (I)Hilbert變換方法。對地震道做Hilbert變換，轉化為復地震道，其中實部為 原地震道數(shù)據(jù)，虛部為其Hilbert變換。得到復地震道之后，就可以在每個采樣點計算振 幅、相位和頻率等屬性，即瞬時屬性。在Taner的工作之后，Hilbert變換法被廣泛用于計 算地震瞬時屬性，至今大部分商業(yè)軟件仍采用該方法。很多學者對該方法進行了發(fā)展和改 進，Barnes于1996年提出二維復地震道分析的概念。Luo等人于2003年提出廣義Hilbert 變換并給出其在地球物理方面的應用。Barnes于2007年提出加權瞬時頻率的概念。Lu和 Zhang將Hilbert變換方法進行推廣，提出基于自適應濾波器計算瞬時頻率的方法。
[0006] (2)基于時頻分析的方法。時頻分析方法利用信號的時頻分布來求取瞬時屬性。 Boashash等人提出基于時頻分析的自適應瞬時頻率估計方法。Stankovic等人利用自適應 窗時頻分布計算瞬時頻率。高靜懷等人提出了在相空間計算地震瞬時屬性的方法。Marfurt 等人提出基于窄帶譜分析的地震資料瞬時屬性提取方法。Steeghs和Drijkoningen提出 基于二次時頻分布的地震層序分析以及屬性提取方法。Huang等人提出基于經(jīng)驗模式分解 (EMD)計算瞬時頻率的方法。Han等人將完全總體經(jīng)驗模式分解（CEEMD)方法用于提取地 震資料的瞬時頻率。
[0007] (3)基于反演的方法。Fomel等人提出利用反演的方法來獲得局部屬性，相比于瞬 時頻率，局部頻率物理意義更加明確，且應用效果明顯。Liu等人提出基于反演的方法，計算 地震資料的瞬時頻率。
[0008] 常用的基于Hilbert變換估計瞬時參數(shù)的方法對噪聲很敏感，且由于濾波器的截 斷效應，使得計算出的瞬時屬性精度低。

【發(fā)明內容】

[0009] 本發(fā)明目的在于克服現(xiàn)有技術的不足，提供了一種具有良好的抗噪性能和精度的 基于相空間濾波策略的地震瞬時屬性提取方法。
[0010] 為達到上述目的，本發(fā)明采用以下技術方案：
[0011] 基于相空間濾波策略的地震瞬時屬性提取方法，包括以下步驟：
[0012] 步驟1 :由目標區(qū)時間和頻率范圍確定m和η的區(qū)間；
[0013] 小波原子由基本小波進行伸縮和平移得到，即
[0014]
[0015] 在進行數(shù)值計算時，需要對尺度因子a和平移因子t進行離散化，對尺度因子進行 指數(shù)化離散記為[<}，其中％為離散化步長且a(l> 1，平移因子t的離散化步長為t ^，離散 化后的小波族表示為
[0016]
[0017] 平移因子t所對應的η的范圍由待分析數(shù)據(jù)的時間范圍確定；根據(jù)廣義Morse小 波尺度和頻率的關系，確定尺度因子a的范圍，進而確定m的范圍；
[0018] 步驟2 :由公式（14)求變換域系數(shù)； _9] Cm，n=〈s，Φ m，n〉， (13)
[0020] 其中Cm，n為標架系數(shù)，s為離散化的信號構成的向量，{Ψ π，η}為小波緊標架， 〈s，U表示二者的內積；
[0021] 步驟3 :閾值操作，得到有效信號對應系數(shù)；
[0022] 對系數(shù)運用閾值策略，軟閾值和硬閾值操作分別定義為
[0023]
[0024]
[0025] 其中λ為閾值參數(shù)；
[0026] 步驟3 :由公式⑵計算解析信號；
[0027]
[0028] 其中，s (t)是信號，h(t)是s (t)的Hilbert變換，c (t)是s (t)對應的解析信號；
[0029] 步驟4 :根據(jù)公式（3)-(5)計算瞬時屬性；
[0030] 利用c (t)來計算瞬時振幅、瞬時相位和瞬時頻率：
[0034] 其中Re[c(t)]和Im[c(t)]分別表示c(t)的實部和虛部。[0035] 進一步，步驟1中廣義Morse小波尺度a和頻率f的關系由下式確定：
[0031]
[0032]
[0033]
[0036]
[0037] 進一步，步驟3中閾值參數(shù)λ由下式計算得到：
[0038]
[0039] 其中ρ為百分比。
[0040] 本發(fā)明提出了基于相空間濾波策略的地震瞬時屬性提取方法，在小波變換域自適 應確定有效信號能量分布空間，在此基礎上，利用小波變換與Hilbert變換的關系，提出了 含噪信號瞬時屬性分析的方法。本發(fā)明提出衡量瞬時屬性抗噪性能的定量指標，通過計算 不同信噪比的合成地震記錄的瞬時屬性，對比不同方法的抗噪性能，本文方法的計算結果 具有良好的抗噪性能和精度。將該方法應用于實際資料，能夠更加清晰地刻畫儲層的特征。
[0041] 此外，在數(shù)值計算時，需要對連續(xù)小波變換中的平移因子和尺度因子進行離散化 以得到小波族，本發(fā)明給出一種離散化方式，推導出基于該離散化方式小波標架的上下界 估計及標架判定準則，證明相比常用的Morlet小波，離散化的廣義Morse小波更易構成緊 標架。
【附圖說明】
[0042] 圖I 50Hz的Ricker子波的時間尺度譜圖；
[0043] 圖2合成地震記錄和采用不同方法計算出來的瞬時頻率；
[0044] 圖3不同信噪比的合成地震記錄；
[0045] 圖4采用不同方法計算的含噪合成記錄的瞬時頻率；
[0046] 圖5用不同的小波函數(shù)計算的IOdB含噪合成記錄的瞬時頻率；
[0047] 圖6實際資料算例。
【具體實施方式】
[0048] 以下通過具體實施例和附圖對本發(fā)明方案做具體說明：
[0049] 小波變換法計算瞬時屬性
[0050] 實地震道s (t)的小波變換定義為
[0051]
[0052] 當基本小波Φ (t)為解析小波時，高靜懷等證明了以下定理：
[0053] 如果Φ⑴為解析小波，其實部Φκα)為偶函數(shù)且其傅里葉變換Ψκ(ω)滿足如 下條件
，則對任意的信號外；)e Z2(R)，有下列等式成立
[0054]
[0055] 其中h(t)是s⑴的Hilbert變換，c⑴是s⑴對應的解析信號。該定理建立 了小波變換和Hilbert變換之間的關系。因此，本發(fā)明可以利用c(t)來計算瞬時振幅、瞬 時相位和瞬時頻率：
[0056]
[0057]
[0058]
[0059] 其中Re[c(t)]和Im[c(t)]分別表示c(t)的實部和虛部。
[0060] 小波函數(shù)的選取
[0061] 廣義Morse小波是一類雙參數(shù)解析小波族，最早由Daubechies等人研宄時頻局域 化算子時導出，是時頻聯(lián)合局域化問題的解。其頻域表達式為
[0062]
[0063] 其中υ(ω)為單位階躍函數(shù)，β和γ為小波的參數(shù)，且β > 〇，γ > 〇。α ρ，γ為歸一化常量：
[0064] α 2(θγ/β) β/γ. (27)
[0065] 本發(fā)明選用廣義Morse小波計算信號的瞬時屬性，主要原因有：
[0066] (1)廣義Morse小波是嚴格解析的，而常用的Morlet小波只是在調制頻率足夠大 的時候近似解析，在時間局域化要求比較高的時候，Morlet小波會出現(xiàn)負頻率泄露。更重 要的是，上面定理中要求基本小波為解析小波，小波偏離解析性會造成該定理不再成立，影 響所求瞬時屬性的準確性。而廣義Morse小波嚴格解析，完全滿足定理要求；
[0067] (2)廣義Morse小波具有更高的自由度，可以通過調整參數(shù)，表現(xiàn)出不同的形態(tài)和 特點，以匹配不同的地震信號；
[0068] (3)離散后的廣義Morse小波更容易構成緊標架，使得小波變換的運算量大大降 低，計算小波變換的系數(shù)和從小波系數(shù)重構信號變得很方便，同樣也有利于后面的濾波策 略中對系數(shù)的處理。另外，利用冗余的標架對于噪聲的壓制也非常有利。
[0069] 小波函數(shù)的離散化方法
[0070] 小波原子由基本小波進行伸縮和平移得到，即
[0071]
[0072] 在進行數(shù)值計算時，需要對尺度因子a和平移因子t進行離散化，對尺度因子進行 指數(shù)化離散，記為[<}，其中％為離散化步長且a 1，平移因子t的離散化步長為t ^。這 樣，離散化后的小波族可以表示為
[0073]
[0074] 在本發(fā)明的離散化方法中，h和待分析信號的采樣間隔相同，因此離散化的小波 族滿足平移不變性，可以看做由字典
平移產(chǎn)生。對于具備平移不變性 的字典，有以下定理：
[0075] 如果存在B彡A > 0,使得對于》eR，有
[0076]
[0077] 則該平移不變字典構成一個標架，A和B分別為標架的上界和下界。該定理將字 典構成標架的條件轉化為Φ^χ)的傅里葉變換Φπ(ω)應滿足的條件。
[0078] 對于小波函數(shù)，基于本發(fā)明的離散化方法，有Φ," (似)