趙妮, 李瑞, 陳立軍
(1.油氣藏地質(zhì)開發(fā)工程國家重點實驗室(成都理工大學), 四川 成都 610059; 2.西安工業(yè)大學計算機科學與工程學院, 陜西 西安 710032; 3.延長油田股份有限公司勘探開發(fā)技術(shù)研究中心, 陜西 延安 716000)
層序地層學是研究旋回式、成因上有聯(lián)系的、以侵蝕面或者與其可以對比的整合面為界的年代地層格架以及沉積層序內(nèi)部地層、巖相分布模式的地層學分支學科。一般采用露頭、巖心、測井和地震資料進行層序地層劃分和層序地層界面識別,在露頭和巖心資料不足、地震資料分辨率有限的情況下,測井資料在層序界面劃分方面就顯得尤為重要。隨著高分辨率層序地層學理論的提出,小波分析、Hilbert-Huang變換等新方法也應(yīng)用于層序地層劃分。Hilbert-Huang變換是基于Hilbert譜的算法,包括經(jīng)驗模態(tài)分解(Empirical Mode Decomposition,EMD)和希爾伯特譜分析2部分。經(jīng)驗模態(tài)分解即EMD方法針對非平穩(wěn)信號提取局部均值曲線,將復雜信號分解為一系列分量函數(shù),以得到具有明顯物理意義的瞬時頻率,可以很好地處理非平穩(wěn)、非線性信號[1-2]。本文將EMD方法應(yīng)用于測井數(shù)據(jù)處理當中,為測井層序界面識別提供了一種有效的信號處理手段。
EMD分解是基于數(shù)據(jù)本身的一種自適應(yīng)變換方法,其本質(zhì)在于利用信號的特征時間尺度獲得信號的固有波動模式。與小波變換和傅里葉變換相比,它可以很好地處理非平穩(wěn)、非線性信號。該方法建立在以下假設(shè)之上:①信號至少有1個極大值和1個極小值點;②特征時間尺度由極值點之間的時間間隔定義;③對于沒有極值點但有變形點的情況,可以通過微分獲取極值點再通過積分獲得分解結(jié)果。
EMD分解流程見圖1。其中,固有模態(tài)函數(shù)應(yīng)滿足2個條件(即流程中的IMF條件)[3]:①極值點數(shù)目和過0點數(shù)目相等或相差1個;②在任意點處,最大值包絡(luò)線與最小值包絡(luò)線的均值為0。第1個條件保證了數(shù)據(jù)的局部最大值為正,極小值為負;第2個條件剔除了瞬時頻率中由于波形不對稱引起的震蕩。
圖1 EMD分解流程
從圖1可見,EMD分解是將原信號x(t)分解成一系列從高頻到低頻的固有模態(tài)函數(shù)IMFi(t)和一個趨勢殘余項R(t),即
(1)
傅里葉變換和小波分析是測井信號分析中常用的方法。傅里葉變換是經(jīng)典信號頻譜分析的基礎(chǔ),將原信號分解為若干正弦信號的加權(quán)和,其中每1個正弦信號對應(yīng)1個固定的頻率和固定的幅值,適合于分析頻率不隨時間變化的平穩(wěn)信號。測井數(shù)據(jù)是非平穩(wěn)的空間信號,經(jīng)傅里葉變換后得到的頻率信號在頻域上有最大分辨率,但其本身不包含深度定位信息,也就是說可以找出最大頻率,但不能確定其深度位置。小波變換是傅里葉變換的發(fā)展,可以根據(jù)頻率高低自動調(diào)節(jié)窗口大小,具有多分辨分析功能。離散小波分解是通過2個互補的濾波器產(chǎn)生低頻的近似信號和高頻的細節(jié)信號,在應(yīng)用中通常認為低頻部分是最重要的。小波變換的局限性在于,①小波基函數(shù)的選取會影響分解信號的精度;②每次提取的分量不是完整的信息,確切的說是信號經(jīng)過濾波器后的部分信息。而EMD分解是對整個信號進行的分解,每一個分量中都包含了信號完整的信息,反映信號的瞬時頻率特征,這也是用EMD分解劃分旋回優(yōu)于小波變換的主要原因。
在測井數(shù)據(jù)采集過程中,各種干擾因素不可避免地影響信號的分辨率。對相關(guān)性高的敏感指標進行線性組合是提高測井解釋準確性的有效方法。
(1) 針對某一區(qū)塊選擇取心井測井數(shù)據(jù),通過多指標的雷達圖比較法,并與巖心進行對比,選取能反映巖性變化特征的敏感指標,再對選出的敏感指標求取相關(guān)系數(shù),選擇相關(guān)性較好的指標進行組合。比如,自然電位和自然伽馬的組合、聲波時差和密度的組合等。對多個線性相關(guān)的指標進行組合是對地層特征信息的強化。由于燥聲是隨機的,疊加信號突出了有用信息,降低了干擾的影響。
(2) 不同測井指標的量綱往往不相同。為了將不同量綱的指標組合在一起,需要進行指標的歸一化處理,再將無量綱的數(shù)據(jù)進行組合,構(gòu)成新的分析指標。
2.2.1 數(shù)據(jù)處理原則
EMD分解可以將原始信號分解成不同頻率的內(nèi)蘊函數(shù)IMF(i)和一個趨勢殘余項R(t),每個IMF(i)為信號在某一特征尺度參數(shù)上的模態(tài),其中相對頻率較高的分量反映了主要受天文因素控制的短期旋回特征,低頻分量則反映了受構(gòu)造作用控制的較為穩(wěn)定的中長期旋回特征;R(t)為信號本身的趨勢或儀器漂零[4]。利用這一規(guī)律,對EMD分解后得到的固有模態(tài)函數(shù)IMF(i)展開相關(guān)性研究。①計算各階IMF(i)與原信號的相關(guān)系數(shù)。往往那些與原信號相關(guān)性好的IMF(i)也是揭示穩(wěn)定沉積環(huán)境的最佳信號分量。②再分析各階IMF(i)之間的相關(guān)性,將相關(guān)性好的IMF(i)進行疊加組合。相關(guān)性好說明有一致性的變化規(guī)律,將其疊加可以加強反映地質(zhì)特征的有用信息,削弱干擾因素,從而提高識別地層界面的精確度。
2.2.2 EMD分解響應(yīng)特征
測井信號經(jīng)過EMD分解后,可以從不同尺度上提取出原測井信號中的時頻特征,反映出不同的沉積旋回級次。若能與各級層序界面建立一定的對應(yīng)關(guān)系,則可以作為測井層序地層劃分的依據(jù)。
當河流在動力學的作用下達到平衡時,將此時的平衡剖面抽象為一個勢能面,稱之為基準面[5-6]。地表或沉積界面不斷地通過沉積或侵蝕作用向靠近基準面的方向運動以便達到新的平衡,這樣的一個旋回升降過程稱為成因?qū)有?。層序地層的這種旋回性在測井資料上表現(xiàn)為粒度和巖性的韻律性,反映到測井曲線上表現(xiàn)為不同的幅度和頻率特征,而曲線突變點往往是層序界面[7-9]。據(jù)此對測井指標的EMD分解曲線響應(yīng)特征在各種沉積環(huán)境中的形態(tài)表現(xiàn)作鏡像圖,概括出5種基本形態(tài)類型。不同沉積環(huán)境的測井曲線形態(tài)大都是這5種基本形態(tài)的組合或演化(見表1)。
表1 測井指標EMD分解曲線響應(yīng)特征的沉積環(huán)境基本類型
應(yīng)用EMD分解識別地層基準面旋回的方法對四川盆地川西坳陷中段須家河組×井層序地層進行了劃分。
自然伽馬反映沉積物的粒度、分選及泥質(zhì)含量,可以用來判斷巖性、進行地層對比及估算泥質(zhì)含量[10]。通常情況下,自然伽馬在砂巖中讀數(shù)為低值,在泥巖中讀數(shù)為高值,但若砂巖中存在諸如云母、鋯石等放射性礦物,會直接影響自然伽馬的讀值。自然電位常用于砂泥巖剖面中劃分巖性、進行地層對比、確定滲透層界面位置、估計地層中的泥質(zhì)含量、展開沉積相研究等方面[11]。通過分析,這2個測井參數(shù)均是適合該地區(qū)劃分層序的靈敏指標,并且曲線指示具有很好的對應(yīng)性,可以進行組合。
(1) 分析自然伽馬和自然電位指標的相關(guān)性,得到它們之間的相關(guān)系數(shù)值是0.438 3,說明它們的相關(guān)性較好。
(2) 分別對自然伽馬和自然電位進行歸一化處理,選用線性歸一化公式,即
y=(x-xmin)/(xmax―xmin)
(2)
式中,x為指標測值;xmax、xmin分別為x的最大值和最小值;y為對x歸一化后的值。
(3) 對歸一化后的自然伽馬和自然電位進行疊加,取均值,組合成新的指標記為GR-SP。
(4) 分析GR-SP指標與自然伽馬和自然電位的相關(guān)性,得到GR-SP與自然伽馬的相關(guān)系數(shù)為0.936 3;GR-SP與自然電位的相關(guān)系數(shù)為0.725 9,說明新指標與自然伽馬和自然電位都具有較高的相關(guān)性,可以替代原始指標,具有應(yīng)用價值。
對四川盆地川西坳陷中段須家河組×井的GR-SP數(shù)據(jù)進行EMD分解得到本征模態(tài)函數(shù)IMF(i)(見圖2)。
圖2 GR-SP指標EMD分解圖
從圖2可以看到,IMF(1)至IMF(4),顯示出強烈的震蕩,但高幅值出現(xiàn)的位置又不盡相同,這說明高頻分量中存在噪音信號的影響。隨著各階IMF信號頻率的降低,信號中穩(wěn)定的變化趨勢越來越強,尤其是IMF(5)、IMF(6)、IMF(7)、IMF(9)明顯具有一致的變化規(guī)律,這種一致性反映出地層相對穩(wěn)定的沉積環(huán)境特征。
旋回級次是對地層特征在不同尺度上的劃分。EMD各階分量中的高頻分量顯示出包含細節(jié)的多少,而這些細節(jié)分量中包含有用的信息也包含了噪音的干擾。沉積旋回是一個相對穩(wěn)定的變化過程,更應(yīng)當注重信號中的中低頻成分,那些高頻信號的加入只是起到添加劑的作用,這一點與小波分析的觀點是一致的。
設(shè)Ri為GR-SP與其EMD分解函數(shù)IMF(i)的相關(guān)系數(shù),利用MATLAB提供的相關(guān)函數(shù)corrcoef()可求出原信號與固有模態(tài)函數(shù)IMF(i)的相關(guān)系數(shù)值(見表2)。
表2 GR-SP與其EMD分解函數(shù)IMF(i)的相關(guān)系數(shù)值
從表2可知,EMD分解函數(shù)IMF(9)、IMF(5)、IMF(6)、IMF(7)與GR-SP均有較好的相關(guān)性。這一結(jié)果同前面EMD分解圖得出的結(jié)論是一致的。同理,再計算出IMF(i)函數(shù)之間的相關(guān)系數(shù)值(見表3)。
從表3中可以得出相關(guān)系數(shù)大于0.1的有IMF(6)與IMF(7),為0.118 184;IMF(5)與IMF(9)的相關(guān)系數(shù)為0.106 921,說明IMF(6)與IMF(7)、IMF(5)與IMF(9)具有較好的相關(guān)性。相關(guān)系數(shù)接近0.1的還有IMF(1)與IMF(2)、IMF(3)與IMF(4)、IMF(5)與IMF(6),可以作為參考組合。將相關(guān)性好的EMD分解信號進行組合,目的是加強反映地層特征的固有信息,降低干擾的影響,比較結(jié)果說明在所有組合信號中只有那些同時與原信號相關(guān)性高的組合才能更好的揭示出地層的旋回特征。
表3 EMD分解本征模態(tài)函數(shù)IMF(i)之間的相關(guān)系數(shù)值
已知×井信號取樣間隔為0.1 m,IMF(6)+IMF(7)的瞬時頻率范圍為0.003 4~0.516 1,進一步算出對應(yīng)可分辨的深度為0.2~29.4 m;同理,可求出IMF(5)+IMF(9)的瞬時頻率范圍為(9.931 1×10-4)~0.503 4,其對應(yīng)可分辨的深度為0.2~100.69 m。不同級次的旋回通常以分布時限來劃分,體現(xiàn)到層次實體規(guī)模上,一般短期旋回厚度在數(shù)米級至數(shù)十米級,中期旋回厚度在數(shù)十米至近百米,長期旋回厚度在近百米至數(shù)百米。據(jù)此,可將IMF(6)+IMF(7)作為短旋回的參考依據(jù);IMF(5)+IMF(9)作為中長期的參考依據(jù)。具體旋回劃分結(jié)果見圖3。
圖3 EMD分解法識別的×井短期和中長期旋回層序地層柱狀圖
×井位于四川盆地川西坳陷中段孝泉-豐谷構(gòu)造帶。圖3分析的是須家河組二段地層的基準面旋回特征,巖性以細、中粒巖屑砂巖、頁巖為主,屬于三角洲前緣亞相。依據(jù)表1指示的曲線形態(tài)特征對EMD分解后的不同頻率組合信號進行分析,可以明顯識別短期和中長期旋回特征?,F(xiàn)對沉積微相分析如下。①水下分支河道:巖性為灰色細-中粗粒巖屑砂巖,表現(xiàn)為正粒序旋回特征;發(fā)育交錯層理、斜層理、塊狀層理,可見沖刷-充填構(gòu)造。②河口壩:巖性為分選好的細砂或粉砂組成,發(fā)育交錯層理,具有反粒序旋回特征。③支流間灣:巖性以黏土為主,含少量的粉砂及細砂,具有水平層理。④席狀砂:巖性為分選好的細砂或粉砂組成,發(fā)育交錯層理,具有反粒序旋回特征。⑤遠砂壩:位于河口壩前緣較遠部位,巖性主要為粉砂,砂體內(nèi)部發(fā)育交錯層理,具有反粒序旋回特征。
(1) 測井指標從不同角度反映了地層的巖性、構(gòu)造、層序及儲層的物性特征,針對地區(qū)優(yōu)選靈敏指標加以組合可以提高測井解釋的準確性。
(2) 自然電位和自然伽馬是測井資料中能反映地質(zhì)構(gòu)造和巖性的靈敏指標,采用數(shù)學方法對其進行組合,新指標是對地層固有特征的強化,可以提高分層效果。
(3) EMD分解是一種全新的時頻分析方法,不僅保證了信號分解后的非平穩(wěn)特征,還避免了加窗傅里葉變換受時窗限制以及小波變換受基函數(shù)影響的不足,可以更準確地描述信號時變特征。
(4) 對多個測井組合指標進行EMD分解,探討各階本征模態(tài)函數(shù)及其與原指標的相關(guān)性,尋找有價值的組合規(guī)律,創(chuàng)建IMF(i)組合信號,有利于更好地提取異常,改善識別地層基準面旋回的效果。
參考文獻:
[1] 李彥軍, 胡祥云. 何展翔HHT變換在地球物理中的應(yīng)用現(xiàn)狀及前景 [J]. 工程地球物理學報, 2010, 7(5): 537-543
[2] Huang N E, Zheng Shen, Long S R, et al. The Empirical Mode Decomposition and the Hilbert Spectrum for Nonlinear Non-stationary Time Series Analysis [C]∥Proc R Soc London. 1998, 45(4): 903-995.
[3] 鄭天翔, 楊力華. 經(jīng)驗模式分解算法的探討和改進 [J]. 中山大學學報: 自然科學版, 2007, 46(1): 1-6.
[4] 張婷, 聶春燕. 基于EMD時頻分析的聲波測井油層的信息提取 [J]. 長春大學學報, 2011, 21(8): 15-18.
[5] 鄭榮才, 尹世民, 彭軍. 基準面旋回結(jié)構(gòu)與疊加樣式的沉積動力學分 [J]. 沉積學報, 2001, 19(2): 259-255.
[6] 鄭榮才, 彭軍, 吳朝容. 陸相盆地基準面旋回的級次劃分和研究意義 [J]. 沉積學報, 2000, 18(3): 249-255.
[7] 郭小燕, 王成, 于雷. 測井資料識別層序的方法及問題討論 [J]. 西部探礦工程, 2011, 2: 67-70.
[8] 周祺, 鄭榮才, 李鳳杰, 等. 測井曲線在陸相層序地層界面識別中的應(yīng)用 [J]. 大慶石油地質(zhì)與開發(fā), 2008, 27(4): 135-138.
[9] 李新虎. 測井曲線拐點在測井層序地層分析中的應(yīng)用研究 [J]. 天然氣地球科學, 2006, 17(6): 815-819.
[10] 曾婷, 桂志先, 何加成, 等. 儲層預測中有關(guān)測井參數(shù)的分析及應(yīng)用 [J]. 工程地球物理學報, 2010, 7(3): 296-299.
[11] 李君文, 陳洪德, 田景春. 鄂爾多斯盆地石千峰組高分辨率層序地層學特征 [J]. 大慶石油地質(zhì)與開發(fā), 2007, 26(1): 45-49.