本發(fā)明涉及天然氣地球物理勘探領(lǐng)域，是一種利用地震數(shù)據(jù)來提取地震波瞬時質(zhì)心頻率，進而利用地震瞬時質(zhì)心頻率可靠地指示天然氣的存在。

背景技術(shù)：

在天然氣勘探中，如何從地震信號的瞬時屬性中提取更多的有用信息，尤其是利用頻率異常信息，并結(jié)合地質(zhì)、測井等資料，來尋找有意義的天然氣儲集帶是石油物探研究人員一直以來的追求目標，同時也是難點問題。含流體的巖石地層會造成地震波在傳播過程中發(fā)生能量損失，在含氣層的內(nèi)部及其下部，地震波的能量會發(fā)生明顯的高頻衰減。譜分解技術(shù)是目前利用地震信號高頻衰減異常從地震反射數(shù)據(jù)進行地質(zhì)解釋及油氣指示的一種常用的有效的含氣性預(yù)測技術(shù)。地震數(shù)據(jù)體經(jīng)過譜分解可以產(chǎn)生一系列分頻剖面，每個分頻剖面都是所有地震數(shù)據(jù)體在某個特定頻率或頻率段處瞬時振幅的一種反映。不同頻率處的振幅剖面能夠體現(xiàn)不同尺度處的地質(zhì)體的不同的響應(yīng)特征，而含氣性信息可能在某些頻率處的分頻剖面中得到強化反映，更容易被識別。通常含氣區(qū)域會體現(xiàn)出“低頻強能量，高頻弱能量”的衰減特征。但是譜分解技術(shù)需要利用一系列地震分頻剖面進行分析，然后選取合理的分頻剖面進行解釋，工作量較大。

此外，譜分解技術(shù)的核心是利用時頻分析方法如短時傅里葉變換、s變換、小波變換等。相比于傳統(tǒng)的地震信號時頻分析方法，適用于非線性非平穩(wěn)信號處理的基于經(jīng)驗?zāi)B(tài)分解(empiricalmodedecomposition,emd)的時頻分析方法具有一些獨特的優(yōu)勢，如不受測不準原理的限制，時頻分辨率更高等。目前基于emd的時頻分析技術(shù)已經(jīng)被廣泛應(yīng)用于天然氣檢測中，并取得了很好的效果。emd可以把地震信號從高頻到低頻分解為有限個單頻或窄帶本征模態(tài)函數(shù)(intrinsicmodefunction,imf)的和。huang等人(1998)認為只有本征模態(tài)函數(shù)的瞬時頻率才有物理意義。雖然基于emd時頻分析方法的含氣性檢測技術(shù)較基于傳統(tǒng)時頻分析方法的含氣性檢測技術(shù)具有更好的時空分辨率和精確度，但是由于emd方法本身存在模態(tài)混疊效應(yīng)等問題，會導(dǎo)致emd分解出的imf缺乏物理意義和地質(zhì)意義，從而降低含氣性檢測的準確性。

完備聚合經(jīng)驗?zāi)B(tài)分解(completeensembleempiricalmodedecomposition，ceemd)是torres等人(2011)提出的一種克服emd方法中模態(tài)混疊效應(yīng)的改進方法，它可以把地震信號分解的更徹底，誤差在數(shù)值上可以忽略，得到的imf分量更具有物理意義和地質(zhì)意義。

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

本發(fā)明的目的在于解決上述現(xiàn)有技術(shù)中存在的難題，提供一種新的地層含氣性信息的頻率異常檢測方法，采用的是基于ceemd的地震瞬時質(zhì)心頻率技術(shù)。通過定義一種新的地震屬性：地震瞬時質(zhì)心頻率，進行天然氣儲層含氣性檢測，提高現(xiàn)有技術(shù)的精確度。

為實現(xiàn)上述目的，本發(fā)明采用的技術(shù)方案是：一種利用地震瞬時質(zhì)心頻率進行含氣性檢測的方法，包括如下步驟：

a、綜合利用地質(zhì)、測井及合成地震記錄等資料準確標定目標層；

b、對各個地震道進行ceemd分解產(chǎn)生不同的imf分量，對每個imf分量，結(jié)合希爾伯特變換計算瞬時振幅和瞬時頻率，再利用地震瞬時質(zhì)心頻率計算公式計算各個imf分量的瞬時質(zhì)心頻率，通過相關(guān)加權(quán)系數(shù)對各個imf分量計算的瞬時質(zhì)心頻率進行加權(quán)求和，獲得該地震道的瞬時質(zhì)心頻率值，逐道計算目標層范圍內(nèi)的地震數(shù)據(jù)的瞬時質(zhì)心頻率，獲得該工區(qū)目標層范圍內(nèi)的地震瞬時質(zhì)心頻率特征體；

c、結(jié)合地質(zhì)、測井和試油氣信息相關(guān)資料，利用井旁地震瞬時質(zhì)心頻率特征剖面確定不同巖性、流體等引起的地震瞬時質(zhì)心頻率屬性體區(qū)別，再將其特征外推到無井區(qū)域，從而進行巖性、含氣性檢測。

進一步的，步驟b具體為：

b1、將各個地震道進行ceemd分解，得到imf分量

1)確定ceemd分解過程中添加的高斯白噪聲、集合成員的個數(shù)和最佳最大篩選迭代次數(shù)，添加的白噪聲的幅值a和集合成員的個數(shù)n應(yīng)滿足：

其中，εn是最終的誤差標準偏差，白噪聲的幅值a取值范圍為0.1～0.4，集合成員的個數(shù)n＞100；

ceemd方法使用特定的噪聲添加在每次分解的過程中，計算唯一的剩余量以獲得每個imf分量，最終分解獲得的imf分量是完全的，具有數(shù)值可忽略的誤差；通常，添加的白噪聲幅度不宜過大，否則會在分解過程中引入虛假的imf分量；添加的白噪聲分量的幅度也不宜太小，否則它會導(dǎo)致模態(tài)混疊現(xiàn)象不易消除。

2)ceemd分解：在ceemd分解過程中，使用emd經(jīng)過i次分解信號x(t)+ε0wi(t)(i＝1,2,…,i)以獲得第一個模態(tài)，并計算：

其中，wi(t)(i＝1,2,…,i)是不同的白噪聲序列，x(t)是原始地震信號；

在第一個階段(k＝1)中，計算第一個余量r1(t)：

分解r1(t)+ε1e1(wi(t))(i＝1,2,…,i),直到產(chǎn)生它們的第一個imf分量，然后定義第二個imf分量為：

ej(·)表示產(chǎn)生第j個模態(tài)，

對于k＝2,…,k,計算第k次的余量：

分解rk(t)+εkek(wi(t))(i＝1,2,…,i)直到獲得它們的第一個emd模態(tài)，然后定義第(k+1)個模態(tài)為：

重復(fù)前述步驟計算第k個模態(tài)；

b2、計算每個地震道分解產(chǎn)生的imf分量的瞬時屬性

對每一條地震道提取出來的各個模態(tài)分量c(t)，分別利用希爾伯特變換用下式提取瞬時幅度a(t)和瞬時頻率ω(t)：

其中，h[·]表示hilbert變換，p為柯西主值；為了避免(7)式中計算瞬時頻率時相位解卷繞導(dǎo)致的模糊度，采用下式計算瞬時頻率ω(t)：

b3、計算各條地震道的瞬時質(zhì)心頻率

瞬時質(zhì)心頻率f定義為：

其中，ωk為每個imf分量在每個時間采樣點處的瞬時頻率；ak為每個imf分量在每個時間采樣點處的瞬時振幅；

本發(fā)明中，瞬時質(zhì)心頻率f可以更有效的反映頻率異常信息，這里，首先計算每條地震道各個imf分量計算獲得的瞬時質(zhì)心頻率，然后，采用相關(guān)加權(quán)系數(shù)對該條地震道計算的各個imf分量對應(yīng)的瞬時質(zhì)心頻率進行加權(quán)求和，結(jié)果作為該地震道的瞬時質(zhì)心頻率，相關(guān)加權(quán)方案如下：

利用各條地震道分解后生成的imf分量與原始地震道的相關(guān)系數(shù)r的大小對各個imf分量獲得的瞬時質(zhì)心頻率進行加權(quán)，相關(guān)加權(quán)系數(shù)rc定義為：

式(10)對具有強相關(guān)的imf分量計算得到的瞬時質(zhì)心頻率保持不變，對具有較強相關(guān)的imf分量計算得到的瞬時質(zhì)心頻率進行10^-1的衰減，對具有弱相關(guān)的imf分量計算得到的瞬時質(zhì)心頻率進行10^-2的衰減；該操作可以加強該地震道的主要貢獻成分，同時減弱次要貢獻分量的信息成分。

對計算出的最終該條地震道的瞬時質(zhì)心頻率進行歸一化，采用如下歸一化計算公式：

x2＝(x-min(x))/(max(x)-min(x))……(11)

其中，x為該條地震道的瞬時質(zhì)心頻率，x2為該條地震道的歸一化瞬時質(zhì)心頻率，將結(jié)果歸一化到區(qū)間[0,1]范圍內(nèi)，min(·)表示取數(shù)據(jù)的最小值，max(·)表示取數(shù)據(jù)的最大值。

本發(fā)明的有益技術(shù)效果是：本發(fā)明的基于ceemd的地震瞬時質(zhì)心頻率提取方法，具有如下特點，主要表現(xiàn)為：

(1)使用了克服模態(tài)混疊效應(yīng)的ceemd方法進行地震道的分解。通過ceemd方法獲得的各個imf分量是具有不同頻帶寬度的窄帶信號，且較常規(guī)emd方法更具物理意義和地質(zhì)意義。準確的imf分量保證了通過它計算獲得的瞬時頻率具有物理意義和更為準確、明確的地質(zhì)意義。

(2)對地震道進行ceemd分解產(chǎn)生的各個不同imf分量分別結(jié)合希爾伯特變換提取瞬時振幅和瞬時頻率，再計算瞬時質(zhì)心頻率，通過相關(guān)加權(quán)系數(shù)對各個imf分量計算的瞬時質(zhì)心頻率進行加權(quán)求和，弱化了地震信號中的次要貢獻成分，突出了主要貢獻成分，能夠更有效地反映地震信號的微弱變化，加強了天然氣地震響應(yīng)特征，有利于含氣性檢測。

(3)地震瞬時質(zhì)心頻率屬性體的計算可適用于二維或三維數(shù)據(jù)的剖面分析、沿層切片或者時間切片等分析，計算方式靈活多樣。

(4)基于ceemd的地震瞬時質(zhì)心頻率屬性體估計方法較常規(guī)譜分解技術(shù)和衰減分析技術(shù)更為簡單方便，同時具有較高的分辨率。

附圖說明

為了更清楚地說明本發(fā)明實施例或現(xiàn)有技術(shù)中的技術(shù)方案，下面將對實施例或現(xiàn)有技術(shù)描述中所需要使用的附圖作簡單地介紹，顯而易見地，下面描述中的附圖僅僅是本發(fā)明的一些實施例，對于本領(lǐng)域普通技術(shù)人員來講，在不付出創(chuàng)造性勞動性的前提下，還可以根據(jù)這些附圖獲得其他的附圖。

圖1為含氣地質(zhì)衰減模型；

圖2為模型的地震響應(yīng)；

圖3為模型的基于ceemd的地震瞬時質(zhì)心頻率剖面；

圖4為四川盆地某地的須家河組一條二維疊后偏移過含氣井剖面；

圖5為利用本技術(shù)計算的該地震剖面基于ceemd的地震瞬時質(zhì)心頻率屬性剖面；

圖6為利用基于短時傅里葉變換的譜分解方法計算的該地震剖面的分頻剖面的對比結(jié)果；

圖7為利用基于小波變換的譜分解方法計算的該地震剖面的分頻剖面的對比結(jié)果；

圖8為四川盆地某地的須家河組另一條二維疊后偏移過氣水井剖面；

圖9為利用本技術(shù)計算的該疊后偏移過氣水井剖面的地震瞬時質(zhì)心頻率屬性剖面。

具體實施方式

下面將結(jié)合本發(fā)明實施例中的附圖，對本發(fā)明實施例中的技術(shù)方案進行清楚、完整地描述，顯然，所描述的實施例僅僅是本發(fā)明一部分實施例，而不是全部的實施例。基于本發(fā)明中的實施例，本領(lǐng)域普通技術(shù)人員在沒有作出創(chuàng)造性勞動前提下所獲得的所有其他實施例，都屬于本發(fā)明保護的范圍。

實施例1

一種利用地震瞬時質(zhì)心頻率進行含氣性檢測的方法，包括如下步驟：

a、綜合利用地質(zhì)、測井及合成地震記錄等資料準確標定目標層；

b、對各個地震道進行ceemd分解產(chǎn)生不同的imf分量，對每個imf分量，結(jié)合希爾伯特變換計算瞬時振幅和瞬時頻率，再利用地震瞬時質(zhì)心頻率計算公式計算各個imf分量的瞬時質(zhì)心頻率，通過相關(guān)加權(quán)系數(shù)對各個imf分量計算的瞬時質(zhì)心頻率進行加權(quán)求和，獲得該地震道的瞬時質(zhì)心頻率值，逐道計算目標層范圍內(nèi)的地震數(shù)據(jù)的瞬時質(zhì)心頻率，獲得該工區(qū)目標層范圍內(nèi)的地震瞬時質(zhì)心頻率特征體；具體為：

b1、將各個地震道進行ceemd分解，得到imf分量

1)確定ceemd分解過程中添加的高斯白噪聲、集合成員的個數(shù)和最佳最大篩選迭代次數(shù)，添加的白噪聲的幅值a和集合成員的個數(shù)n應(yīng)滿足：

其中，εn是最終的誤差標準偏差，白噪聲的幅值a取值范圍為0.1～0.4，集合成員的個數(shù)n＞100；

ceemd方法使用特定的噪聲添加在每次分解的過程中，計算唯一的剩余量以獲得每個imf分量，最終分解獲得的imf分量是完全的，具有數(shù)值可忽略的誤差；通常，添加的白噪聲幅度不宜過大，否則會在分解過程中引入虛假的imf分量；添加的白噪聲分量的幅度也不宜太小，否則它會導(dǎo)致模態(tài)混疊現(xiàn)象不易消除。

2)ceemd分解：在ceemd分解過程中，使用emd經(jīng)過i次分解信號x(t)+ε0wi(t)(i＝1,2,…,i)以獲得第一個模態(tài)，并計算：

其中，wi(t)(i＝1,2,…,i)是不同的白噪聲序列，x(t)是原始地震信號；

在第一個階段(k＝1)中，計算第一個余量r1(t)：

分解r1(t)+ε1e1(wi(t))(i＝1,2,…,i),直到產(chǎn)生它們的第一個imf分量，然后定義第二個imf分量為：

ej(·)表示產(chǎn)生第j個模態(tài)，

對于k＝2,…,k,計算第k次的余量：

分解rk(t)+εkek(wi(t))(i＝1,2,…,i)直到獲得它們的第一個emd模態(tài)，然后定義第(k+1)個模態(tài)為：

重復(fù)前述步驟計算第k個模態(tài)；

b2、計算每個地震道分解產(chǎn)生的imf分量的瞬時屬性

對每一條地震道提取出來的各個模態(tài)分量c(t)，分別利用希爾伯特變換用下式提取瞬時幅度a(t)和瞬時頻率ω(t)：

其中，h[·]表示hilbert變換，p為柯西主值；為了避免(7)式中計算瞬時頻率時相位解卷繞導(dǎo)致的模糊度，采用下式計算瞬時頻率ω(t)：

b3、計算各條地震道的瞬時質(zhì)心頻率

瞬時質(zhì)心頻率f定義為：

其中，ωk為每個imf分量在每個時間采樣點處的瞬時頻率；ak為每個imf分量在每個時間采樣點處的瞬時振幅；

本發(fā)明中，瞬時質(zhì)心頻率f可以更有效的反映頻率異常信息，這里，首先計算每條地震道各個imf分量計算獲得的瞬時質(zhì)心頻率，然后，采用相關(guān)加權(quán)系數(shù)對該條地震道計算的各個imf分量對應(yīng)的瞬時質(zhì)心頻率進行加權(quán)求和，結(jié)果作為該地震道的瞬時質(zhì)心頻率，相關(guān)加權(quán)方案如下：

利用各條地震道分解后生成的imf分量與原始地震道的相關(guān)系數(shù)r的大小對各個imf分量獲得的瞬時質(zhì)心頻率進行加權(quán)，相關(guān)加權(quán)系數(shù)rc定義為：

式(10)對具有強相關(guān)的imf分量計算得到的瞬時質(zhì)心頻率保持不變，對具有較強相關(guān)的imf分量計算得到的瞬時質(zhì)心頻率進行10^-1的衰減，對具有弱相關(guān)的imf分量計算得到的瞬時質(zhì)心頻率進行10^-2的衰減；該操作可以加強該地震道的主要貢獻成分，同時減弱次要貢獻分量的信息成分。

對計算出的最終該條地震道的瞬時質(zhì)心頻率進行歸一化，采用如下歸一化計算公式：

x2＝(x-min(x))/(max(x)-min(x))……(11)

其中，x為該條地震道的瞬時質(zhì)心頻率，x2為該條地震道的歸一化瞬時質(zhì)心頻率，將結(jié)果歸一化到區(qū)間[0,1]范圍內(nèi)，min(·)表示取數(shù)據(jù)的最小值，max(·)表示取數(shù)據(jù)的最大值。

c、結(jié)合地質(zhì)、測井和試油氣信息相關(guān)資料，利用井旁地震瞬時質(zhì)心頻率特征剖面確定不同巖性、流體等引起的地震瞬時質(zhì)心頻率屬性體區(qū)別，再將其特征外推到無井區(qū)域，從而進行巖性、含氣性檢測。

應(yīng)用實施例1

(1)利用地質(zhì)、測井、速度等資料結(jié)合彌散黏滯方程一個建立地質(zhì)衰減模型。表1是依據(jù)鄂爾多斯盆地某氣田的資料建立的含氣衰減地質(zhì)模型的參數(shù)。其中，層④為含氣層，層③為干層(地質(zhì)模型如圖1所示)。

表1模型的地質(zhì)參數(shù)

其中，vp是縱波速度，ρ是密度，ζ是彌散系數(shù)，η是黏滯系數(shù)。q是衰減系數(shù)。

模型含氣層厚度為75m，子波頻率為40hz，采樣頻率為512hz。

(2)圖2為模型的地質(zhì)響應(yīng)。

(3)圖3為模型對應(yīng)的基于ceemd的瞬時質(zhì)心頻率剖面。如圖所示，本技術(shù)很好的檢測到了含氣區(qū)域。

(4)圖4為四川盆地某地的須家河組一條二維疊后偏移過井剖面。該區(qū)域為致密砂巖儲層。圖中上部的橢圓所示區(qū)域為含氣區(qū)域，下部的橢圓所示區(qū)域為含水區(qū)域。含氣井命名為well。

(5)利用本技術(shù)計算的基于ceemd的地震瞬時質(zhì)心頻率剖面如圖5所示。從圖中可以看到，在上部分的橢圓所示的含氣區(qū)域中存在強烈的瞬時質(zhì)心頻率異常值。在下部分的橢圓所示的含水區(qū)域有較強的瞬時質(zhì)心頻率異常值。排除巖性等其他地層影響因素，所提方法給出了一個含烴類的統(tǒng)計性解釋結(jié)果。同時，所提方法對于氣、水的衰減程度在檢測上有所差異，含氣區(qū)域的異常值更大，水層異常值較小。氣層、水層區(qū)分較大。圖像時間分辨率和空間分辨率較強。

(6)圖6為地震剖面的利用基于短時傅里葉變換的譜分解方法計算的分頻剖面的對比結(jié)果。圖6a和6b所示為利用短時傅里葉變換分別提取的低頻分頻剖面和高頻分頻剖面。該數(shù)據(jù)主頻為50hz。這里，低頻分頻剖面為45hz，高頻分頻剖面為65hz。從圖中可以看出，基于短時傅里葉變換的低頻分頻剖面在含氣區(qū)域和含水區(qū)域都檢測到了較強的強振幅異常，在高頻分頻剖面中在含氣區(qū)域和含水區(qū)域強能量減弱，都表現(xiàn)出了“低頻強能量、高頻弱能量”的特征，排除巖性等其他地層影響因素，都給出了一個含烴類的統(tǒng)計解釋結(jié)果。但是分頻剖面時空分辨率較低，對氣、水區(qū)域差異不大。

(7)圖7為利用基于小波變換的譜分解方法計算的地震剖面的分頻剖面的對比結(jié)果。圖7a和7b所示為利用小波變換分別提取的低頻分頻剖面和高頻分頻剖面。這里，低頻分頻剖面為45hz，高頻分頻剖面為65hz。從圖中可以看出，基于小波變換的低頻分頻剖面和高頻分頻剖面在含氣區(qū)域和含水區(qū)域都呈現(xiàn)出了“低頻強能量、高頻弱能量”的特征，都給出了一個含烴類的統(tǒng)計解釋結(jié)果。與圖6基于短時傅里葉變換提取的分頻剖面相比，基于小波變換的分頻剖面的時空分辨率較高，但是從圖6a可以看出，氣、水區(qū)域差異不大。對比圖6、圖7利用常規(guī)譜分解技術(shù)與圖5所提技術(shù)得到的結(jié)果，本技術(shù)給出的結(jié)果時空分辨率較好，對儲層信息體現(xiàn)較好，對氣、水的區(qū)分能力更好。同時，本技術(shù)不需要類似常規(guī)譜分解技術(shù)需要分析大量的分頻剖面。

(8)圖8為四川盆地某地的須家河組另一條二維疊后偏移過氣水井剖面。圖中的橢圓所示區(qū)域為含氣水區(qū)域。

(9)圖9為利用本技術(shù)計算的該疊后偏移過氣水井剖面的地震瞬時質(zhì)心頻率屬性剖面。從圖中可以看到，在含氣水區(qū)域存在較強的瞬時質(zhì)心頻率異常值。排除巖性等其他地層影響因素，所提方法給出了一個含烴類的統(tǒng)計性解釋結(jié)果。

本技術(shù)領(lǐng)域技術(shù)人員可以理解，除非另外定義，這里使用的所有術(shù)語(包括技術(shù)術(shù)語和科學(xué)術(shù)語)具有與本發(fā)明所屬領(lǐng)域中的普通技術(shù)人員的一般理解相同的意義。還應(yīng)該理解的是，諸如通用字典中定義的那些術(shù)語應(yīng)該被理解為具有與現(xiàn)有技術(shù)的上下文中的意義一致的意義，并且除非像這里一樣定義，不會用理想化或過于正式的含義來解釋。

本領(lǐng)域普通技術(shù)人員可以理解：實施例中的裝置中的部件可以按照實施例描述分布于實施例的裝置中，也可以進行相應(yīng)變化位于不同于本實施例的一個或多個裝置中。上述實施例的部件可以合并為一個部件，也可以進一步拆分成多個子部件。

最后所應(yīng)說明的是：以上實施例僅用以說明而非限制本發(fā)明的技術(shù)方案，盡管參照上述實施例對本發(fā)明進行了詳細說明，本領(lǐng)域的普通技術(shù)人員應(yīng)該理解：依然可以對本發(fā)明進行修改或者等同替換，而不脫離本發(fā)明的精神和范圍的任何修改或局部替換，其均應(yīng)涵蓋在本發(fā)明的權(quán)利要求范圍當(dāng)中。