鄭露露 林建民 龐立臣 祝捍皓
摘要:
通過Hilbert-Huang變換對淺水環(huán)境下氣槍震源激發(fā)子波特征進行了研究,結(jié)果表明:氣槍震源首先激發(fā)高頻壓力脈沖,其頻帶相對較寬,優(yōu)勢頻率為40~80 Hz,主頻隨時間降低;隨后高壓氣體在水下形成的氣泡振蕩激發(fā)低頻氣泡脈動,其頻帶相對較窄,優(yōu)勢頻率為4~12 Hz,且水深較淺處的氣泡振蕩激發(fā)的氣泡脈動主頻相對較低。研究結(jié)果有助于提高對淺水環(huán)境下氣槍震源激發(fā)的復(fù)雜近場信號的認(rèn)識。
關(guān)鍵詞:Hilbert-Huang變換;氣槍震源;高頻壓力脈沖;低頻氣泡脈動
中圖分類號:P3153文獻標(biāo)識碼:A文章編號:1000-0666(2016)03-0466-07[HJ]
0引言
地震波是唯一能夠穿透地球內(nèi)部的振動,迄今為止,關(guān)于地球內(nèi)部的結(jié)構(gòu)、組成、過程和狀態(tài)的研究等多來自地震波的信息(陳颙,朱日祥,2005)。震源是產(chǎn)生地震波的源頭,其所產(chǎn)生的地震波信號的品質(zhì)將直接影響地震探測的效果。氣槍震源自1964年由美國Bolt公司的卡爾米思克發(fā)明以來,逐漸被發(fā)展成為具有綠色環(huán)保、可重復(fù)性好、可控性強、子波可測等卓越性能的水中震源,而被廣泛應(yīng)用于海洋地震勘探(其中以地下油氣、礦藏等資源的局部尺度探測為主)(Calvert,2004)。近年來,氣槍震源也逐漸被應(yīng)用于區(qū)域尺度的地震探測,出現(xiàn)了利用海洋激發(fā)氣槍、陸地接收信號的科學(xué)嘗試,如美國的LARSE(Los Angeles Region Seismic Experiment)計劃(Lutter et al, 1999;Fuis et al, 2003)、新西蘭的SIGHT(South Island Geophysical Transect)計劃(Okaya et al, 2002;Stern et al, 2002;Avendonk et al, 2004)以及國內(nèi)的南海北部海陸聯(lián)測試驗(趙明輝等,2004;丘學(xué)林等,2007)。由于水中炸藥爆破存在可控性差、激發(fā)可重復(fù)性低以及作業(yè)費用高等缺點,氣槍震源亦開始逐漸成為陸地水域激發(fā)地震波的選擇。國內(nèi)在陸上水庫進行的氣槍震源激發(fā)實驗表明:大容量氣槍震源是一種優(yōu)良的水中震源,具有激發(fā)能量大、地震波轉(zhuǎn)換率高、傳播距離遠等特點,是陸上激發(fā)地震波的一種新方法(陳颙等,2007a,b;Chen et al, 2008;林建民等,2008a,b;楊微等,2013;劉自鳳等,2015;蘇金波等,2015;王彬等,2015)。
通過大量的理論研究與實驗觀測,人們已經(jīng)很好地掌握了深海氣槍源的震源特性和以此組建氣槍陣列的技術(shù)(Mayn,Quay, 1971;Caldewll, Dragoset, 2000;Dragoset, 2000;羅桂純等,2007)。深海中氣槍作業(yè)水深多在30 m以上,如LARSE 計劃中氣槍的沉放水深大于100 m(Fuis et al,2003),周圍水體較大,震源相對比較獨立,激發(fā)信號受水底、水面反射等影響較小,激發(fā)過程相對比較簡單,這使得震源特性幾乎完全由氣槍本身決定。但水深較淺的陸地水域水體相對較小,使用氣槍震源必須綜合考慮氣槍震源與水域的相互作用、耦合等諸多問題,其激發(fā)地震波受激發(fā)條件和激發(fā)環(huán)境多種因素(水深、沉放深度、水底地形、邊界及沉積等)的影響,因此,其激發(fā)特性與規(guī)律和深海氣槍源存在很大差異。盡管水庫氣槍激發(fā)作為陸地震源的有效性已經(jīng)得到驗證(陳颙等,2007b),但至今針對淺水環(huán)境下氣槍震源激發(fā)特性及激發(fā)地震波規(guī)律的研究甚少。林建民等(2010)曾重點分析了沉放深度、工作壓力等參數(shù)對水庫氣槍震源激發(fā)波形的影響,但尚未從氣槍震源信號本身特征出發(fā)加以分析。
氣槍震源信號是典型的非線性、非平穩(wěn)信號,即統(tǒng)計量(相關(guān)函數(shù)、功率譜等)隨時間而變化。Huang等(1998)提出了一種分析非平穩(wěn)非線性時間序列的新方法,即Hilbert-Huang變換(Hilbert-Huang Transform, 簡稱HHT)。該方法通過經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition, 簡稱EMD),根據(jù)信號本身的時間尺度特征,將信號分解為有限個可以視為窄帶信號的本征模態(tài)函數(shù)(Intrinsic Mode Function, 簡稱IMF),并通過各IMF的Hilbert變換獲得具有明確定義和物理意義的瞬時頻率,建立了以本征模態(tài)函數(shù)為基本時域信號、以瞬時頻率作為基本量表征信號變化的新時頻分析方法。本文將利用HHT方法研究水庫氣槍震源激發(fā)信號時頻特征,并獲得瞬時頻率等瞬時地震屬性,提高對淺水環(huán)境下氣槍震源激發(fā)的復(fù)雜近場信號的認(rèn)識。
11經(jīng)驗?zāi)B(tài)分解
經(jīng)驗?zāi)B(tài)分解的關(guān)鍵是將一個非線性非平穩(wěn)時間序列分解為有限個IMF分量和一個趨勢項。由于每個時間序列都是由多個時間尺度的振蕩波構(gòu)成,因此,該方法從本質(zhì)上是將一個信號序列進行平穩(wěn)化的過程,其處理過程是將信號中的這些固有的、內(nèi)在的波動或趨勢項按照不同的時間尺度逐級分解開來,產(chǎn)生一系列具有不同時間尺度的IMF。
該方法可通過包絡(luò)線方法進行逐級篩選,實現(xiàn)時間序列的經(jīng)驗?zāi)B(tài)分解,求出全部有限個IMF,主要步驟包括:
(1)對于原始序列X (t)的局部極大值和極小值,用三次樣條函數(shù)插值形成極大值包絡(luò)Xmax(t)和極小值包絡(luò)Xmin(t);
(2)求取極大值包絡(luò)與極小值包絡(luò)的瞬時均值線
m1(t)=12[Xmax(t)+Xmin(t)]. (1)
(3)考察原始序列X(t)與瞬時均值m1(t)差值,即去掉低頻的新序列h1(t)=X(t)-m1(t),判定其是否滿足以下IMF要求,如果滿足,則將h1(t)作為第一個IMF分量C1;如果不滿足,則將h1(t)作為原始序列重復(fù)以上篩選過程,依次篩選,直至第k次篩選結(jié)果符合要求,即得到第一個IMF分量C1。
這樣得到的C1應(yīng)是時間序列中特征時間尺度最小的高頻IMF分量。然后將C1從原序列 X(t)中分離出去,得到一個剩余序列r1=X(t) -C1,由于r1中可能含有比 C1的特征時間尺度大的IMF,因此又將r1作為新的序列實施上述(1)~(3)篩選過程,以求得第二個IMF分量C2。如此重復(fù),求出一系列IMF,直至剩余序列rn(趨勢項)成為一個單調(diào)函數(shù),從中再也不能分解出IMF分量時,整個分解過程結(jié)束。此時,原始序列可表示成為n個IMF和1個趨勢項之和,即
X(t)=∑nj=1Cj(t)+rn(t).(2)
其中,IMF分量Cj需滿足以下條件:
(1)在整個資料集里,極值點的數(shù)目和穿零點的數(shù)目必須相等或者最多相差1個;
(2)由局部極大值所構(gòu)成的包絡(luò)線和由局部極小值所構(gòu)成的包絡(luò)線的平均值為零。
從上述EMD過程可以看出,該方法的優(yōu)點在于:直觀、直接、基于數(shù)據(jù)本身特征以及自適應(yīng)性強。
其中,P表示Cauchy主值。數(shù)學(xué)上,式(3)表示X(t)與1t的卷積,因此,Hilbert變換強調(diào)X(t)的局部性。在1t的基礎(chǔ)上,構(gòu)造解析信號Z(t)為
Z(t)=X(t)+iY(t)=a(t)eiθ(t). (4)
式中,a(t)為瞬時振幅,θ(t)為瞬時相位,對其求導(dǎo)可得瞬時頻率ω(t)為
ω(t)=dθ(t)dt.(5)
為使瞬時頻率ω(t)具有實際意義,信號X(t)必須為窄帶信號,其函數(shù)圖像關(guān)于零均值線是局部對稱的,并且穿零點和極值的數(shù)目是一樣的。若以上條件不滿足,所求出的瞬時頻率就可能出現(xiàn)異常值或負(fù)數(shù)。
通過EMD所分解得到的IMF為窄帶信號,符合作Hilbert變換求取瞬時頻率的要求,因此,原始信號可以用其各個IMF分量的瞬時頻率ωj(t)和瞬時振幅aj(t)近似地表示成
H(ω,t)=RE∑nj=1aj(t)eiωj(t)dt
.(6)
式(6)即為Hilbert譜,其中,RE表示求取實部運算。ωj(t)和 aj(t)都隨時間而變化,體現(xiàn)了非平穩(wěn)序列的本質(zhì)特征,突破了傳統(tǒng)傅立葉分析中定常振幅和固定頻率的限制。Hilbert譜將信號的振幅表示成為瞬時頻率和時間的函數(shù),描述了其在整個頻率域上隨時間和頻率的變化規(guī)律。
2氣槍震源信號時頻特征及瞬時屬性
本文所研究信號為容量2 000 in3氣槍震源(Bolt 1500LL型)在水深18 m的陸上水庫環(huán)境下的激發(fā)信號,氣槍沉放深度10 m,工作壓力(氣槍中壓縮氣體的壓力)10 MPa,距離氣槍激發(fā)點水平距離50 m、深8 m位置處布設(shè)水聽器(采樣率200 Hz,頻帶50 Hz~40 kHz)所記錄的近場子波信號(林建民等,2010)。海洋石油勘探通常采用氣槍陣列作為震源,其子波測量點與氣槍源之間的距離通常在幾百米以上,主要是為了使陣列中各單槍激發(fā)信號能同時(在一個采樣時間間隔內(nèi))到達水聽器,而本文工作中采用單槍激發(fā),且重點研究近場信號特征。
圖1為氣槍激發(fā)信號的原始波形及其EMD結(jié)果,圖1a是水聽器記錄的實際波形,圖1b是分解得到的7個本征模態(tài)函數(shù)(IMF1~IMF7),其中IMF7是剩余項。觀察各IMF分量時域波形,不同IMF分量包含了不同的時間特征尺度,分別以不同的分辨率顯示了氣槍信號特征,并且這種分辨率是自適應(yīng)的。IMF1~IMF7分別是從高頻至低頻分量逐級分解獲得,其中,IMF1和IMF2主要對應(yīng)氣槍震源信號的高頻壓力脈沖,平穩(wěn)性相對較差;而IMF3和IMF4主要對應(yīng)后續(xù)的低頻氣泡脈沖,平穩(wěn)性比較好。
若直接利用原始波形進行Hilbert變換求瞬時頻率,則因其不完全滿足Hilbert變換的前提假設(shè)條件而產(chǎn)生異常值,不利于真正利用瞬時頻率掌握源信號頻率的時間變化,如圖3所示。由圖2可知,盡管IMF3的主要頻率范圍對應(yīng)低頻氣泡脈沖,但是在高頻壓力脈沖位置仍有部分微小能量分布,導(dǎo)致在瞬時頻率曲線上出現(xiàn)高頻成分分布,這也說明不同IMF之間有頻率交叉,并不是上述的簡單頻域帶通濾波。
氣槍震源激發(fā)信號分高頻壓力脈沖和低頻氣泡脈沖,前者為高壓氣體瞬間釋放到水中激發(fā)壓力波,因具有較高頻率通常用于淺層石油勘探,后者為高壓氣體釋放到水中之后形成后續(xù)氣泡振蕩所產(chǎn)生,主頻相對較低,適合長炮檢距的深部結(jié)構(gòu)探測。從圖2可以看出,IMF3在02 s附近瞬時頻率迅速由高頻轉(zhuǎn)為低頻,對應(yīng)壓力脈沖向氣泡脈沖過渡,即高壓氣體瞬間釋放到水中激發(fā)壓力波之后形成的氣泡在02 s附近開始振蕩。但為掌握氣槍震源信號整體頻率隨時間的變化規(guī)律以及信號振幅或能量隨頻率和時間的變化分布情況,則需通過Hilbert變換的譜分析獲得Hilbert時頻譜,將振幅表示成為瞬時頻率和時間的函數(shù)。圖4給出了本文所分析氣槍信號的波形、功率譜密度分析及相應(yīng)的Hilbert時頻譜,能清晰地刻畫出氣槍信號的時頻特征。由于Hilbert時頻譜包含了詳細的IMF分量頻率變化的信息,能準(zhǔn)確反應(yīng)信號的時頻特性,可以清晰地分辨出信號主頻的突變及所對應(yīng)的時刻,因此,具有更高的時間—頻率域分辨率。
3不同時頻分析方法討論
目前,用于分析地震信號的最常用的、且較為成熟的方法是傅立葉變換。但由于傅立葉變換嚴(yán)格要求系統(tǒng)的線性、數(shù)據(jù)的周期性或平穩(wěn)性(即傅立葉頻率代表著信號的周期性),在處理非平穩(wěn)、非線性數(shù)據(jù)時,其實際應(yīng)用受到很大的限制。如本文中氣槍信號為非平穩(wěn)信號,其特點之一就是無周期性,若按傅立葉變換方法定義其頻率進行頻譜分析將缺乏物理基礎(chǔ)。因此,引入具有真實物理意義的瞬時頻率ω(t)對于研究氣槍信號瞬態(tài)特性具有重要的意義。
通??衫肏ilbert變換求取信號具有實際物理意義的瞬時屬性,但要求被分析信號是窄帶或平穩(wěn)的,而且該方法對噪聲很敏感。而氣槍震源信號等實際地震信號既非平穩(wěn)又含有噪聲,直接對其進行Hilbert變換求取瞬時頻率等屬性參數(shù),將缺乏物理意義甚至失真(圖3)。HHT的根本創(chuàng)新就在于先通過EMD將信號進行平穩(wěn)化,將信號中固有的、內(nèi)在的波動成分按照不同的時間尺度逐級分解成窄帶信號分量IMF,使其滿足利用Hilbert變換求取具有實際物理意義瞬時屬性的要求,再分別計算每個分量對應(yīng)的瞬時頻率等瞬時屬性。瞬時頻率作為時間的函數(shù),能敏銳地識別出待分析信號中的多尺度嵌套結(jié)構(gòu)(圖2)。
傅立葉變換是一種純頻域分析方法,是將信號從時域轉(zhuǎn)換到頻域的整體變換,無法提供時域信息。若要獲得某一時刻的頻率特性,則通常采用時頻分析方法,這也是針對非平穩(wěn)信號較為常用的處理方法,比較典型的有短時傅立葉變換、小波變換等(李春峰,Linea, 2005;陳雨紅等,2006)。但這些方法基本上都是將信號序列在某種基底函數(shù)系上展開,然后分析展開系數(shù)及各分量的特征。如短時傅立葉變換實際上是加窗傅立葉變換,通常假設(shè)信號在窗函數(shù)的有效持續(xù)時間內(nèi)是平穩(wěn)的,但此條件通常無法滿足或只是近似滿足,因而它的時間和頻率分辨率都很低。圖5為所分析的氣槍震源信號通過短時傅立葉變換所獲得的時頻特征。小波分析雖然通過一系列可伸縮平移的小波函數(shù)實現(xiàn)了信號時頻局部化分析目的,但在本質(zhì)上仍然是具有柔性時頻窗的加窗傅立葉變換(崔錦泰,2006),它也要求小波窗內(nèi)的信號是短時平穩(wěn)的,還沒有根本擺脫傅立葉變換的局限性。而且,小波基選擇的多樣性和使用過程中的不可變更性也影響了它在非平穩(wěn)信號分析中的適應(yīng)能力。同時,小波基的有限長會造成信號的能量泄漏,因而難以對信號作精確的時頻分析。
與傳統(tǒng)時頻分析方法最大的不同就是不采用固定的基底分解,而是基于信號本身的局部特征時間尺度進行分解獲得IMF分量,大大提高了其自適應(yīng)性,有效克服了傅立葉變換中測不準(zhǔn)原理的限制,在客觀性和分辨率方面具有明顯的優(yōu)勢,更適用于非平穩(wěn)、非線性信號處理。該方法能用分解得到的IMF分量及1個剩余項來揭示信號序列的振蕩結(jié)構(gòu)特征和非平穩(wěn)性,并最終由每個IMF分量經(jīng)過Hilbert變換得到時域和頻域均具有較高分辨率的Hilbert譜(圖4c),具有明確的物理意義,且可以準(zhǔn)確給出原信號序列及其IMF分量的主要振幅變化所對應(yīng)的頻率和時間,即反映物理過程能量—頻率—時間的分布。
由于Hilbert時頻譜是通過計算IMF分量的瞬時頻率和瞬時幅度而直接得到,因此,其時間和頻率分辨率不同于其他時頻分析方法受Heisenberg不確定原理的限制,而是相互獨立的。相比于圖5所示基于短時傅里葉變換所得的時頻特征,圖4c所示的Hilbert時頻譜具有更好的局部時頻特性表現(xiàn)能力,從圖中可以觀察到氣槍信號能量隨時間和頻率的精細變化趨勢。從012 s左右開始,高頻壓力脈沖激發(fā),由于信號脈沖較陡,其頻帶覆蓋范圍較寬,優(yōu)勢頻段則主要集中在40~80 Hz。同時,得益于Hilbert譜的高時頻分辨率,還能觀察到壓力脈沖(圖4c中012~018 s之間原始波形中的4個脈沖峰)中后面的脈沖比前面的脈沖主頻低,即壓力脈沖主頻隨時間遞減,從80 Hz降低到40~60 Hz,然后直接過渡到后續(xù)的低頻氣泡脈沖。其原因可能是后面的脈沖是由壓力脈沖在水面及水底的反射或壓力脈沖本身的振蕩產(chǎn)生,其主頻在水介質(zhì)中傳播及兩相界面上反射過程中下降或自身振蕩過程中下降。從02 s附近開始,釋放的高壓氣體在水下形成的氣泡開始振蕩激發(fā)低頻氣泡脈沖,持續(xù)時間達05 s以上,其優(yōu)勢頻率集中在4~12 Hz的低頻范圍。從圖4b的分析結(jié)果可得,低頻氣泡脈沖在6 Hz和10 Hz附近有2個尖峰,但缺少了時域的分辨信息。而從圖4c中可以推測,10 Hz的尖峰能量可能主要來自于第一次氣泡振蕩,而6 Hz的尖峰能量主要來自于第二次氣泡振蕩。該實際計算結(jié)果與林建民等(2010)的分析結(jié)果中公式(1)相符合,由于氣泡脈沖激發(fā)過程對應(yīng)著逐漸上浮的水下氣泡振蕩過程,第二次氣泡振蕩水深位置相對于第一次氣泡振蕩更淺,因此,公式(1)中P00更小,且R0更大,導(dǎo)致氣泡振蕩周期T更大,即對應(yīng)氣泡脈沖主頻更低,所以,后續(xù)的氣泡脈沖主頻要低于前面的氣泡脈沖(圖4)。此外,在012 s壓力脈沖到達之前,存在微小振幅的振蕩波形,這可能由采集過程中濾波因素引起,它對后續(xù)的Hilbert譜計算可能引入一定的干擾。
本文對氣槍震源信號的處理分析結(jié)果顯示了HHT方法的有效性,但同時也注意到該方法在應(yīng)用過程中對噪聲的敏感性以及邊界問題的影響,如圖4中t=0附近出現(xiàn)近微小能量的近100 Hz高頻能量,這個是由邊界效應(yīng)引起的干擾。如何進一步消除信號中夾雜噪聲及邊界效應(yīng)的影響,對于進一步提高HHT的應(yīng)用效果至關(guān)重要,這是該方法需要進一步完善的地方,也是在應(yīng)用過程中需要進一步研究的問題。
4結(jié)論
HHT是近年發(fā)展起來的處理非線性、非平穩(wěn)信號的高效時頻分析方法,該方法通過EMD將信號分解為IMF窄帶分量,并計算其具有實際物理意義的瞬時頻率,掌握信號頻率的時間變化特征。在此基礎(chǔ)上,進一步通過Hilbert譜分析獲得在時域和頻域均具較高分辨率的時頻譜,更為清晰地表征信號的時頻分布特征。氣槍震源信號為典型的非平穩(wěn)信號,若直接利用基于傅立葉變換的各種處理方法,所獲得的頻率信息物理意義不明確,時頻分析結(jié)果分辨率不高。本文利用HHT對水庫氣槍震源近場信號實際水聽器記錄數(shù)據(jù)進行了分析,相對于以往工作,在氣槍震源信號在不同時間尺度上振蕩特征、具真實物理意義瞬時頻率時間變化特征、時頻域高分辨率能量分布特征及原因等方面獲得了更多認(rèn)識,有助于對淺水環(huán)境下氣槍震源激發(fā)的復(fù)雜近場信號的理解。
感謝中國地震局地球物理研究所水庫氣槍試驗及其他參與單位對本文工作的支持。
參考文獻:
陳颙, 朱日祥. 2005. 設(shè)立“地下明燈研究計劃”的建議[J]. 地球科學(xué)進展, 20(5): 485-489.
陳颙, 李宜靖. 2007a. 地震波雷達研究展望: 用人工震源探測大陸地殼結(jié)構(gòu)[J]. 中國科學(xué)技術(shù)大學(xué)學(xué)報, 37(8): 813-819.
陳颙, 張先康, 丘學(xué)林等. 2007b. 陸地人工激發(fā)地震波的一種新方法[J]. 科學(xué)通報, 52(11): 1317-1321.
陳雨紅, 楊長春, 曹齊放等. 2006. 幾種時頻分析方法比較[J]. 地球物理學(xué)進展, 21(4): 1180-1185.
崔錦泰. 2006. 小波分析導(dǎo)論[M]. 西安: 西安交通大學(xué)出版社.
李春峰, Liner C. 2005. 基于小波多尺度分析的奇性指數(shù): 一種新地震屬性[J]. 地球物理學(xué)報, 48(4): 882-888.
林建民, 王寶善, 葛洪魁等. 2008a. 大容量氣槍震源特征及深部介質(zhì)中傳播的震相分析[J]. 地球物理學(xué)報, 51(1): 206-212.
林建民. 2008b. 基于人工震源的長偏移距地震信號檢測和探測研究[D]. 合肥: 中國科學(xué)技術(shù)大學(xué).
林建民, 王寶善, 葛洪魁等. 2010. 大容量氣槍震源子波激發(fā)特性分析[J]. 地球物理學(xué)報, 53(2): 342-349.
劉自鳳, 蘇有錦, 王寶善等. 2015. 賓川主動源地震波走時變化分析方法研究[J]. 地震研究,38(4):591-597.
羅桂純, 葛洪魁, 王寶善等. 2007. 氣槍震源激發(fā)模式及應(yīng)用[J]. 中國地震, 23(3): 225-232.
丘學(xué)林, 陳颙, 朱日祥等. 2007. 大容量氣槍震源在海陸聯(lián)測中的應(yīng)用: 南海北部試驗結(jié)果分析[J]. 科學(xué)通報, 52(1): 1-7.
蘇金波, 王寶善, 王海濤等. 2015. 利用大容量氣槍震源資料研究北天山地區(qū)介質(zhì)衰減特征[J]. 地震研究,38(4): 598-605.
王彬, 吳國華, 蘇有錦等. 2015. 賓川地震信號發(fā)射臺的選址、建設(shè)及初步觀測結(jié)果[J]. 地震研究,38(1): 1-6.
楊微, 王寶善, 葛洪魁等. 2013. 大容量氣槍震源主動探測技術(shù)系統(tǒng)及試驗研究[J]. 中國地震,29(4): 339-410.
趙明輝, 丘學(xué)林, 夏戡原等. 2004. 南海東北部海陸聯(lián)測地震數(shù)據(jù)處理及初步結(jié)果[J]. 熱帶海洋學(xué)報, 23(1): 58-63.