本發(fā)明涉及地震勘探中斷裂檢測的技術領域，具體涉及一種基于地震信號奇異性的隱蔽性斷裂檢測方法。

背景技術：

地震裂縫檢測是一種既經(jīng)濟又方便的一種裂縫檢測方法，Odling等人(Odling N E,Gillespie P,Bourgine B,et al.Variations in fracture system geometry and their implications for fluid flow in fractured hydrocarbon reservoirs[J].Petroleum Geoscience,1999,5(4):373-384.)經(jīng)研究發(fā)現(xiàn)裂縫是重要的儲集空間和流體滲濾通道。在長期的科研中逐漸形成了多種利用地震資料來進行斷裂檢測方法，Marfurt等人(Marfurt K J.,Farmer S l.,Bahorich M S.,et al.3-D seismic attributes using a semblance-based coherency algorithm.Geophysics,1998,63(4):1150-1165.)提出的多道分析窗內(nèi)的多道互相關、Gersztenkorn和Marfurt(Gersztenkorn A,Marfurt K J.Eigenstructure-based coherence computations as an aid to 3-D structural and stratigraphic mapping.Geophysics,1999,64(5):1468-1479.)基于多道本征分解、Marfurt(Blumentritt C H,Marfurt K J,Sullivan E C.Volume-based curvature computations illuminate fracture orientations-Early to mid-Paleozoic,Central Basin Platform,west Texas.Geophysics,2006,71(5):B159-B166.)傾角和方位角估計、Blumentritt等(Blumentritt C H,Marfurt K J,Sullivan E C.Volume-based curvature computations illuminate fracture orientations-Early to mid-Paleozoic,Central Basin Platform,west Texas.Geophysics,2006,71(5):B159-B166.)、Al-Dossary和Marfurt(Al-Dossary S,Marfurt K J.3D volumetric multispectral estimates of reflector curvature and rotation.Geophysics,2006,71(5):P41-P51.)、Chopra和Marfurt(Chopra S,Marfurt K J.Integration of coherence and volumetric curvature images.Leading Edge,2010,29(1):1092-1107.)曲率屬性、(Chopra S,Marfurt K J.Integration of coherence and volumetric curvature images.Leading Edge,2010,29(1):1092-1107.)相干體和曲率屬性為主的斷裂檢測技術，上述疊后斷裂檢測方法在淺層分辨和信噪比較高，受噪聲影響大。此外，Chopra和Marfurt(Chopra S,Marfurt K J.Integration of coherence and volumetric curvature images.Leading Edge,2010,29(1):1092-1107.)結合地震屬性，聯(lián)合相干、曲率開展多屬性綜合斷裂檢測。單一的疊后斷裂檢測很難取得很好的應用效果，地震多屬性可以有效避免因深層巖石相對較致密，速度高，由此導致地震反射信號較弱，同相軸雜亂，地震剖面分辨低的斷裂檢測。

而針對隱蔽性斷裂檢測檢測，在構造層面上，先后有Steen等人(SteenSverdrup E,Hanssen T H.Predicting the distribution of small faults in a hydrocarbon reservoir by combining outcrop,seismic and well data[J].Geological Society London Special Publications,1998,147(1):27-50.)就根據(jù)地層傾角、方位角以及傾角變化率，該方法得到的結果可結合測井解釋和模型數(shù)據(jù)實現(xiàn)油氣儲層中隱蔽性斷裂檢測，但是針對地震反射弱的地層難以檢測；在地震屬性技術上，Lawrence(Lawrence P.Seismic attributes in the characterization of small‐scale reservoir faults in Abqaiq Field.Leading Edge,1998,17(4):521-525.)根據(jù)地震屬性，以及Neves等人(Neves F A,Zahrani M S,Bremkamp S W.Detection of potential fractures and small faults using seismic attributes.Leading Edge,2004,23(9):903-906.)結合譜分解技術和地震屬性技術檢測隱蔽性斷裂。雖然該方法可以提供平面的隱蔽性斷裂檢測屬性圖，仍不能準確地解決剖面隱蔽性斷裂的斷點位置；在現(xiàn)代信息學科上，Duchesne等人(Duchesne M J,Halliday E J,Barrie J V.Analyzing seismic imagery in the time–amplitude and time–frequency domains to determine fluid nature and migration pathways:A case study from the Queen Charlotte Basin,offshore British Columbia.Journal of Applied Geophysics,2011,73(2):111–120.)通過時頻分析技術及Basir等人(Basir H M,Javaherian A,Yaraki M T.Multi-attribute ant-tracking and neural network for fault detection:a case study of an Iranian oilfield.Journal of Geophysics&Engineering,2013,10(1):15009-15018.)提出多屬性螞蟻體和神經(jīng)網(wǎng)絡技術來識別小斷層。雖然該方法可以提供平面的隱蔽性斷裂檢測屬性圖，仍不能準確地解決剖面隱蔽性斷裂的斷點位置。

技術實現(xiàn)要素：

本發(fā)明的目的在于克服上述背景技術的不足，而提供一種基于地震信號奇異性的隱蔽性斷裂檢測方法，該檢測方法通過小波變換處理后，提取地震信號的奇異值，建立地震信號奇異值與斷點的映射關系，準確地解決剖面隱蔽性斷裂的斷點位置，可以實現(xiàn)一個同相軸內(nèi)、同相軸扭曲、斷距小及弱信號的隱蔽性斷裂。

為實現(xiàn)上述目的，本發(fā)明提供的一種基于地震信號奇異性的隱蔽性斷裂檢測方法，包括以下步驟：

S1：讀入原始地震數(shù)據(jù)；

S2：將步驟S1中讀入的原始地震數(shù)據(jù)進行地震構造分析和區(qū)域沉積特征分析分別確定斷裂大小和地震層位，斷裂大小以垂直斷距(F_i+1-F_i)表示，其中，F(xiàn)_i表示第i層上盤地震反射時間，F(xiàn)_i+1表示第i+1層下盤地震反射時間，F(xiàn)_i單位為ms；地震層位以地震同相軸的連續(xù)性來確定，以T_i表示第i層的地震層位，T_i單位為ms；

S3：將步驟S2中所得的斷裂大小及地震層位進行差值運算，分別求取斷裂大小(F_i+1-F_i)及地震層位厚度(T_i+1-T_i)，并將斷裂大小(F_i+1-F_i)及地震層位厚度(T_i+1-T_i)作為約束參數(shù)；

S4：將步驟S3中的斷裂大小(F_i+1-F_i)及地震層位厚度(T_i+1-T_i)作為約束參數(shù)代入小波變換的公式其中，尺度因子a由地震層位厚度(T_i+1-T_i)確定，平移因子b由斷裂大小(F_i+1-F_i)確定，然后對地震資料的隨機噪聲進行濾波，完成信噪分離，得到濾波后地震數(shù)據(jù)；

S5：將步驟S4得到的濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)，先分別進行傅里葉變換，確定濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)的頻帶寬度及主頻；然后利用小波變換進行時頻分析，確定濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)的能量；若濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)滿足質量控制標準則執(zhí)行步驟S6，若不滿足質量控制標準則返回至步驟S2重新分析T_i和F_i，當最小地震層位厚度(T_i+1-T_i)內(nèi)斷裂大小滿足F_i+1-F_i≠0時作為濾波處理的最終參數(shù)；

S6：根據(jù)小波變換原理對滿足步驟S5質量控制標準的地震數(shù)據(jù)進行小波變換：根據(jù)小波函數(shù)滿足正交性、支撐性、對稱性、消失矩及正則性五個方面來確定小波函數(shù)對尺度因子a和平移因子b按照二進制方式進行離散化，即通過一對共軛濾波器{h(n)}和{g(n)}對信號f進行二進制小波分解：

其中，hj和gj分別由h和g在每對相鄰的樣本間插2j-1個零元得到，從而獲得1～J的小波系數(shù) 當時，第J層的逼近系數(shù)為0；

S7：將經(jīng)過步驟S6小波變換后的地震數(shù)據(jù)進行奇異性處理：確定各層小波系數(shù)奇異性的位置及其相應的值，對于第j層奇異值為 如果且則在點t處的值即為奇異值，記所有的奇異值點時間為t^j₁,t^j₂,…,t^j_Nj，那么相 應的奇異值為：

上述技術方案中，在所述步驟S5中，傅里葉變換具體過程如下：設x(n)為N點的有限長序列，則其傅里葉變換為：

其中利用傅立葉變換可以確定濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)的頻帶寬度及主頻。

上述技術方案中，在所述步驟S5中，質量控制標準為同時滿足以下兩個條件：其一，濾波后地震數(shù)據(jù)與原始地震數(shù)據(jù)經(jīng)傅里葉變換后的頻帶寬度和主頻一致；其二，濾波后地震數(shù)據(jù)與原始地震數(shù)據(jù)進小波變換后的能量一致。

上述技術方案中，在所述步驟S6中，將小波變換定義為：對任意函數(shù)f(x)∈L²(R)，其連續(xù)小波變換基本關系由下式表示：

為滿足一定條件的小波基函數(shù)，其基本關系由下式表示：

式中，尺度因子a由地層厚度參數(shù)(T_i+1-T_i)確定，且a≠0，若a>1，則基函數(shù)相當于將函數(shù)拉伸，使窗口的時寬增大，其頻譜縮窄并向低頻方向移動；若a<1，其頻譜展寬并向高頻方向移動；平移因子b由斷裂大小(F_i+1-F_i)確定，對應于時間軸的平移位置；t為雙程旅行時。

與現(xiàn)有技術相比，本發(fā)明具有如下優(yōu)點：

其一，本發(fā)明通過小波變換處理后，提取地震信號的奇異值，建立地震信號奇異值與斷點的映射關系，準確地解決剖面隱蔽性斷裂的斷點位置，可以實現(xiàn)一個同相軸內(nèi)、同相軸扭曲、斷距小及弱信號的隱蔽性斷裂。

其二，本發(fā)明的檢測方法可以求取任意時刻t處的記所有的奇異值點時間為t^j₁,t^j₂,…,t^j_Nj，那么相應的奇異值為： 而現(xiàn)有技術中不能實現(xiàn)任意時刻的斷點提取，因此斷裂檢測精度高，可以實現(xiàn)一個同相軸內(nèi)、同相軸扭曲、斷距小及弱信號的隱蔽性斷裂；

其三，本發(fā)明推導出經(jīng)過一系列處理后的地震數(shù)據(jù)滿足 且則在點t處的值即為奇異值，因此該檢測方法可以直接提取反映斷裂信息的奇異性信息，建立地震信號奇異值與斷點的映射關系，準確地解決剖面隱蔽性斷裂的斷點位置，消除了其他間接手段因自身方法及其適用性造成的誤差；

附圖說明

圖1是本發(fā)明的實施流程示意圖；

圖2a是本實施例中原始地震資料頻譜分析圖；

圖2b是本實施例中濾波后地震資料頻譜分析圖；

圖3a是本實施例中原始地震資料的局部時頻特征圖；

圖3b是本實施例中濾波后地震資料的局部時頻特征圖

圖4a是檢驗本實施例方法的地質模型；

圖4b是檢驗本實施例方法的地震正演記錄；

圖4c是檢驗本實施例方法的奇異值斷點檢測結果；

圖5a是實際地震剖面；

圖5b是經(jīng)本實施例方法的斷點檢測結果；

圖6a過WP1井的實際地震剖面；

圖6b為過WP1井相干體斷裂檢測剖面；

圖6c為過WP1井經(jīng)本發(fā)明方法的奇異值斷裂檢測剖面。

具體實施方式

下面結合實施例詳細說明本發(fā)明的實施情況，但它們并不構成對本發(fā)明的限定，僅作舉例而已。同時通過說明本發(fā)明的優(yōu)點將變得更加清楚和容易理解。

本實施例的一種基于地震信號奇異性的隱蔽性斷裂檢測方法，包括以下步驟：

S1：讀入原始地震數(shù)據(jù)；

S2：將步驟S1中讀入的原始地震數(shù)據(jù)進行構造分析，首先進行斷裂解釋，確定斷裂的大小，以垂直斷距(F_i+1-F_i)表示(單位為ms，F(xiàn)i表示第i層上盤地震反射時間，F(xiàn)_i+1表示第_i+1層下盤地震反射時間)；然后進行地震層位解釋，根據(jù)地震同相軸的連續(xù)性來確定地震層位Ti(Ti單位為ms，Ti表示第i層)；在斷裂及層位解釋的基礎上，根據(jù)地質背景知識確定具有重要沉積相意義的地震反射特征，確定有沉積相意義的地震反射特征的地震層位厚度(T_i+1-T_i)。

S3：將步驟S2中所得的斷裂大小及地震層位進行差值運算，分別求取斷裂大小(F_i+1-F_i)及地震層位厚度(T_i+1-T_i)，并將斷裂大小(F_i+1-F_i)及地震層位厚度(T_i+1-T_i)作為約束參數(shù)；

S4：將S3中的約束參數(shù)代入小波變換的公式，將小波變換定義為：對任意函數(shù)f(x)∈L²(R)，其連續(xù)小波變換為：

小波變換的基函數(shù)是窗函數(shù)的平移因子b和尺度因子a的結果。

式中：a≠0為尺度因子；b為平移因子；為滿足一定條件的小波基函數(shù)。a的變化可以改變窗口的大小，若尺度因子a>1，則基函數(shù)相當于將函數(shù)拉伸，使窗口的時寬增大，其頻譜縮窄并向低頻方向移動；而a<1，其頻譜展寬并向高頻方向移動。大的尺度因子對應于低頻端，頻率分辨率高、時間分辨率低；小尺度的參數(shù)對應于高頻端，頻率分辨率低、時間分辨率高，這就是小波變換的多分辨率特性。

在公式其中尺度因子a由地震層位厚度(T_i+1-T_i)決定，平移因子b由斷裂大小F_i+1-F_i來確定，對地震資料的隨機噪聲進行濾波，完成信噪分離，得到濾波后地震數(shù)據(jù)；

S5：將步驟S4得到的濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)，先分別進行傅里葉變換，傅里葉變換具體過程如下：設x(n)為N點的有限長序列，則其傅里葉變換為：

其中W_N＝e^-j2*π/N，利用傅立葉變換可以確定濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)的頻帶寬度及主頻，利用傅立葉變換(公式1)可以確定地震資料的頻帶寬度及主頻信息，通過對濾波前后地震資料的頻譜分析，可以在全局上判定地震資料處理的好壞。因傅立葉變換只能在整個時間段的頻譜信息，但要控制地震資料濾波處理過程的每個時間段的細節(jié)特征，需要對信號的局部特征進行篩選，然后對局部信號進行頻譜分析。小波變換可以有效地從信號中提取信號局部特征，該方法通過伸縮和平移等運算對信號進行多尺度的細化，達到分析信號的局部頻譜特征，利用小波變換進行時頻分析，根據(jù)小波函數(shù)滿足正交性、支撐性、對稱性、消失矩及正則性五個方面來確定小波函數(shù)對尺度因子a和平移因子b按照二進制方式進行離散化，即通過一對共軛濾波器{h(n)}和{g(n)}對信號f進行二進制小波分解：

其中，hj和gj分別由h和g在每對相鄰的樣本間插2j-1個零元得到，從而獲得1～J的小波系數(shù) 當時，第J層的逼近系數(shù)為0；確定濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)的能量(局部時頻特征)，質量控制標準為同時滿足以下兩個條件：其一，濾波后地震數(shù)據(jù)與原始地震數(shù)據(jù)經(jīng)傅里葉變換后的頻帶寬度和主頻一致；其二，濾波后地震數(shù)據(jù)與 原始地震數(shù)據(jù)進小波變換后的能量一致；地震資料濾波效果的好壞主要從頻帶寬度、主頻和能量三個方面來控制，三者的一致性說明濾波過程是有效的，圖中2a為原始地震資料頻譜特征，圖2b為濾波后地震資料頻譜特征，對比圖2a與圖2b發(fā)現(xiàn)頻帶寬度和主頻特征一致；圖中3a為原始地震資料局部時頻特征，圖3b為濾波后地震資料局部時頻特征，對比圖3a與圖3b發(fā)現(xiàn)兩者能量一致，因此可以根據(jù)頻帶寬度、主頻和能量三個方面來確定濾波效果，則濾波后地震數(shù)據(jù)及原始地震數(shù)據(jù)滿足質量控制標準則執(zhí)行步驟S6；

S6：根據(jù)小波變換原理對對滿足步驟S5質量控制標準的地震數(shù)據(jù)進行小波變換：

根據(jù)小波函數(shù)滿足正交性、支撐性、對稱性、消失矩及正則性5個方面來選擇小波函數(shù)為了減少冗余度，對尺度因子a和平移因子按照二進制方式進行離散化，即一對共軛濾波器{h(n)}和{g(n)}對信號f進行二進制小波分解：

其中h_j和g_j分別由h和g在每對相鄰的樣本間插2j-1個零元得到，從而獲得1～J的小波系數(shù) 當時，第J層的逼近系數(shù)為0。

S7：將經(jīng)步驟S6小波變換后數(shù)據(jù)進行奇異性處理，圖4a為不同斷裂性質的不同斷距的理論地質模型，圖4b為圖4a模型條件下的正演記錄，圖4c為通過本發(fā)明的斷點檢測，對比圖4a地質模型與圖4b 正演記錄發(fā)現(xiàn)，通過人工斷裂解釋難以識別，而通過奇異值可以確定各層小波系數(shù)奇異性的位置及其相應的值(見圖4c)，即在圖4c上可以清楚檢測出斷點，且與圖4a模型斷點位置保持一致。

為了進一步證實本發(fā)明的適用性，在實際地震資料上，圖5a為實際地震資料，從人工地震解釋層位上也很難看出斷點，從圖5b上可以清楚檢測出斷點。

為了進一步說明與其他現(xiàn)有技術的差異，本發(fā)明獨特之處在于可以直接提取地震數(shù)據(jù)任意時刻t的奇異值如果 且則在點t處的值即為奇異值，記所有的奇異值點時間為t^j₁,t^j₂,…,t^j_Nj。那么相應的奇異值為：

圖6a為過WP1井地震剖面，圖6b為過WP1井不同時刻t的相干體斷裂檢測剖面，由此可見雙程旅行時t的影響，圖6c為過WP1井任意時刻t奇異性斷裂檢測剖面，因此從圖6c與圖6b對比可以看出，圖6c斷裂檢測精度比圖6b高。

本說明書未作詳細描述的內(nèi)容屬于本領域專業(yè)技術人員公知的現(xiàn)有技術。