時頻分析在地震數(shù)據(jù)處理中的應(yīng)用

1、時頻分析在地震數(shù)據(jù)處理中的應(yīng)用摘要：經(jīng)典的傅立葉分析只適用于分析平穩(wěn)信號，而不適用于非平穩(wěn)信號。為了分析非平穩(wěn)信號，我們采用時頻分析方法。時頻分析能夠清楚的揭示信號的時變頻譜特征，是對時變、非平穩(wěn)信號進行分析與處理的有力工具。本文介紹了時頻分析中常用的STFT和WignerVille分布和S變換的特點，并利用Matlab對一個地震信號進行時頻分析的實現(xiàn)。關(guān)鍵詞：時頻分析；短時傅立葉變換（STFT）；S變換；WignerVille分布引言在傳統(tǒng)的信號處理領(lǐng)域，基于 Fourier 變換的信號頻域表示及其能量的頻域分布揭示了信號在頻域的特征，它們在傳統(tǒng)的信號分析與處理的發(fā)展史上發(fā)揮了極其重要的作用

2、。但是，F(xiàn)ourier 變換是一種整體變換，即對信號的表征要么完全在時域，要么完全在頻域，作為頻域表示的功率譜并不能告訴我們其中某種頻率分量出現(xiàn)在什么時候及其變化情況。然而，在許多實際應(yīng)用場合，信號是非平穩(wěn)的，其統(tǒng)計量（如相關(guān)函數(shù)、功率譜等）是時變函數(shù)。這時，只了解信號在時域或頻域的全局特性是遠遠不夠的，最希望得到的乃是信號頻譜隨時間變化的情況。為此，需要使用時間和頻率的聯(lián)合函數(shù)來表示信號，這種表示簡稱為信號的時頻表示。 時頻分析方法旨在通過構(gòu)造一種時間和頻率的密度函數(shù)，將一個一維的時間信號以二維的時間一頻率函數(shù)形式表示出來，以揭示信號中所包含的頻率分量及其隨時間的變化特性。這使我們不但能夠同

3、時掌握非平穩(wěn)信號的時域及頻域信息，而且可以清楚地了解非平穩(wěn)信號的頻率是如何隨時間變化的。通過時頻分析方法技術(shù)對地震信號進行分析處理可在獲得地震信號的瞬時頻率、瞬時相位、瞬時振幅等瞬時參數(shù)的同時獲得時頻譜圖等重要時頻域信息，實現(xiàn)對地震信號的邊緣檢測、屬性提取等。時頻分析是非平穩(wěn)信號分析處理領(lǐng)域的重要方法，時頻分布的基本任務(wù)是建立一個函數(shù)，要求這個函數(shù)能夠同時用時間和頻率來描述信號的能量密度。如果有了這樣的一個分布，就可以計算某一確定的頻率和時間范圍內(nèi)能量的百分率、計算某一特定時刻的頻率密度、計算該分布的整體和局部的各階矩。即尋找一個聯(lián)合密度函數(shù)P(t，f)，使p(t，f)=在時間t和頻率f的強度

4、，或者P(t，f)在時間t和頻率f，在時一頻單元內(nèi)的部分能量。地震信號的時頻分析一般將時頻分析方法分為線性和非線性兩種。典型的線性時頻表示有短時傅葉變換(簡記為STFT)、Gabor展開和小波變換(Wavelet Transformation，簡記為WT)等。非線性時頻方法是一種二次時頻表示方法(也稱為雙線性)，最典型的是WVD(Wigner-Ville Distribution)和Cohen類。 采用的地震信號如下圖所示：圖1 時間域的信號圖2 用FFT方法求取的能量譜1、短時傅里葉變換（STFT）傳統(tǒng)的傅立葉變換只在頻率域具備了局部分析的能力，而在時間域不具有這種能力。要得到滿足同時獲得時

5、間和頻率的局部分析能力的要求，一種最基本的方法就是:取出信號在所關(guān)心時刻附近的一小段，而忽略信號的其它部分，對其作傅立葉變換，即可得到這一特定時刻的頻率分量。因為所取的時間長度與整個信號相比很短，所以將這種方法稱為短時傅立葉變換(STFT)，它是時頻分析中最簡單的形式。短時傅立葉變換(STFT)的基本思想:用窗函數(shù)來截取信號，假定信號在窗內(nèi)是平穩(wěn)的，采用傅立葉變換來分析窗內(nèi)信號，以便確定在那個時間存在的頻率，然后沿著信號移動窗函數(shù)，得到信號頻率隨時間的變化關(guān)系，這就得到了我們所需要的時頻分布。STFT的物理意義在于，對于一定的分析時刻t, 可以視為信號s(t)在該時刻的“局部頻譜”，從而整個變

6、換的結(jié)果也就能揭示信號頻譜的變化特性。短時Fourier變換的時頻分辨率受制于窗函數(shù)的形狀和寬度短時傅立葉變換的時間分辨率與分析窗函數(shù)的時間域?qū)挾瘸烧?，而其頻率分辨率與分析窗的頻寬成正比。從而，一個好的時間分辨率需要一個短的窗函數(shù)，而一個好的頻率分辨率需要一個長的窗函數(shù)。因此，短時傅立葉變換不能同時兼顧時間分辨率和頻率分辨率。當選用的窗函數(shù)為Gaussian函數(shù)時，該變換為Gabor變換。圖3 用STFT做時頻分析2、S變換主要對S變換的定義、推導(dǎo)和特性進行詳細闡述，S變換綜合短時傅立葉變換和小波變換的優(yōu)點，又避免了它們的不足:它與傅立葉變換有著直接的聯(lián)系，具有無損可逆性;與短時傅立葉變換和

7、小波變換一樣，也是一種線性時頻表示，因此不存在交義項的干擾;S變換具有多種分辨率，克服了短時傅立葉變換固定分辯率的不足;S變換中含有相位因子，這是小波變換所不具備的特性?？傮w來說，S變換是近幾年發(fā)展起來的一種新的時頻分析方法.S變換結(jié)合了短時傅里葉變換和小波變換的優(yōu)點,具有相位信息,同時該變換與小波變換一樣,其時頻窗可以調(diào)節(jié)大小以適應(yīng)非平穩(wěn)信號的特點.S變換的這些優(yōu)點,使得它在地球物理方面得到廣泛的應(yīng)用。S變換首先是由Stockwell等人提出的，是以Morlet小波為基本小波的連續(xù)小波變換的延展。在S變換中，簡諧波與高斯函數(shù)的乘機構(gòu)成了基本小波，因為簡諧波在時間域可以作伸縮變換，而高斯函數(shù)則

8、進行伸縮和平移。ST也可以認為是CWT的“相位校正”。函數(shù)的S變換表示為：其中：可以用表示的傅里葉變換，其中和是相同的意義。圖4 S變換時頻分析3、WVDWigner-Ville變換是1932年由Wigner首次提出的，并應(yīng)用于量子力學(xué)領(lǐng)域，后來Ville等人將其引入到信號分析處理領(lǐng)域。1966年，Cohen發(fā)現(xiàn)各種發(fā)現(xiàn)各種時頻分析只是Wigner-Ville變換的不同形式，可以統(tǒng)一起來，成為Cohen類雙線性時頻分析。信號的Wigner-Ville變換用公式表示為式中是的解析信號，即是實信號的Hilbert變換。Wigner-Ville變換也可以用解析信號的頻譜來表示從上面兩種不同形式的Wi

9、gner-Ville變換表達式中可以看出，式中不包含任何窗函數(shù)，從而避免了線性時頻變換中實踐分辨率和頻率分辨率的相互牽制，難以兼顧的問題。因為Wigner-Ville變換的時間-帶寬積可達到Heisenberg測不準原理(不確定原理給出的下界，故可以證明，沒有任何一種時頻變換方法的時間-頻率分辨率及聚集性能出其右。在Wigner-Ville變換表達式中，信號出現(xiàn)了兩次，故稱之為雙線性時頻變換。Wigner-Ville變換不是線性的，即兩信號之和的Wigner-Ville變換不等于每一個信號的Wigner-Ville變換之和，其中多出了一個附加項。令則有其中式中前兩項是自由項(auto term

10、s)，第三項是交叉項(cross terms)。交叉項常常導(dǎo)致時頻平面上出現(xiàn)偽影現(xiàn)象：交叉項是實的，混雜于自由項成分當中，并且其幅度是自由項成分的兩倍；交叉項是震蕩型的，每兩個信號分量就會產(chǎn)生一個交叉項。交叉項的存在嚴重地干擾著對Wigner-Ville分布的解釋，當信號變得非常復(fù)雜時，Wigner-Ville分布甚至變得毫無意義。圖5 WVD時頻分析4、小波變換 在短時傅立葉變換和Gabor展開中我們都使用了固定的時間窗函數(shù)，這就引出了時間分辨率和頻率分辨率的概念，時間分辨率和頻率分辨率是一對矛盾。根據(jù)海森堡的測不準原理，即時間窗函數(shù)的長度越長，頻率分辨率就越高，而對于時間分辨率則越差。為了

11、平衡時間分辨率和頻率分辨率這個矛盾，可以采取對存在高頻分量的部分采用高的時間分辨率和低的頻率分辨率，而對于低頻分量則采用高的頻率分辨率和低的時間分辨率的方法，這就是多分辨分析的思想。 小波變換是一種在時間-尺度平面內(nèi)，利用多分辨率分析思想分析非平穩(wěn)信號的方法。所謂小波，就是一個滿足容許條件的一個函數(shù)族 可以看出函數(shù)族是由窗函數(shù)在時間上平移b，在尺度上伸縮a，再乘上歸一化因子后的結(jié)果，所以非平穩(wěn)信號的連續(xù)小波變換定義為 :其中是小波基函數(shù)的共軛。 將小波變換和短時傅立葉變換兩者的基函數(shù)相比較，可以看出，小波變換基函數(shù)的尺度參數(shù)決定了小波變換的多分辨分析特性，即利用時間-尺度聯(lián)合函數(shù)來分析非平穩(wěn)信

12、號的“變焦距”法，以達到分析信號局部特性的目的。 小波變換由于其本身分辨力的優(yōu)良性能，因此一經(jīng)提出，很快就成了非平穩(wěn)信號分析和處理的一大熱點，經(jīng)過近 20 年的發(fā)展，小波變換取得了突破性的發(fā)展，形成了多分辨分析，框架和濾波器組三大完整豐富的小波變換理論體系?，F(xiàn)在小波變換己經(jīng)被廣泛地應(yīng)用在信號的奇異性檢測、計算機視覺、圖像處理、語音分析與合成等等諸多領(lǐng)域、在分形和混沌理論中也有了很多的應(yīng)用。 聲音信號分析：1. 女聲版江南style：圖6 女聲版聲音信號圖7 FFT所得振幅信息圖8 短時傅里葉變換時頻分析圖9 S變換時頻分析 圖10 WVD時頻分析圖11 小波變換時頻分析圖12 Hilbert變

13、換后小波變換時頻分析圖13 Hilbert變換后WVD時頻分析2. 男聲版江南style：圖12 男聲版聲音信號圖13 FFT所得振幅信息圖14 短時傅里葉變換時頻分析圖15 S變換時頻分析圖16 WVD時頻分析圖17 小波變換時頻分析圖18 加入Hilbert變換小波變換時頻分析圖19 加入Hilbert變換WVD時頻分析對比男生與女生版兩段相似的聲音信號發(fā)現(xiàn)女生的頻率段基本在0.35kHz左右，而男生版的聲音信號則有一段明顯的低頻段。再者由于聲音信號內(nèi)存在的其他干擾，可能是的信號的頻率與能量等有差異，但是最基本的頻率信息還是可以看出來的。另外，交叉項是震蕩型的，每兩個信號分量就會產(chǎn)生一個交

14、叉項。交叉項的存在嚴重地干擾著對Wigner-Ville分布的解釋，當信號變得非常復(fù)雜時，Wigner-Ville分布甚至變得毫無意義。就像上面分析的聲音信號，除了看出簡單的對稱性之外，其他的分析意義基本沒有。幾種不同的Wigner-Ville變換原始信號： x1=cos(2*pi*200*t).*(t>=0 & t<=0.8); x2=cos(2*pi*50*t).*(t>=0.2 & t<=1); x=x1+x2;圖20 原始信號圖21 三種WVD時頻分析圖22 WVD時頻分析三維顯示結(jié)論由于地下介質(zhì)的復(fù)雜性，地震信號是一種非平穩(wěn)隨機信號，時頻分析是

15、分析處理地震信號的一種非常有效的方法。前人在此領(lǐng)域已做了大量工作，將很多時頻分析的方法應(yīng)用在地震數(shù)據(jù)處理解釋中。本文對各種時頻分析方法的比較研究和分析其各自特點。并將時頻分析的方法用于地震信號分析處理中，結(jié)果表明:通過時頻分析方法對地震剖面進行處理，可以得到瞬時地震屬性參數(shù)(包括頻率類瞬時屬性參數(shù)、瞬時振幅譜)，這些新的地震剖面對地質(zhì)解釋有一定的輔助作用。另外后面聲音信號的分析說明，時頻分析對于頻率有很好的敏感性，對于區(qū)分不同頻帶的信息有著非常好的分辨率。像不同的樂器，相同物體不同的激發(fā)聲音方式等，都可以通過時頻分析進行解釋。時頻分析以聯(lián)合時頻分布的形式來表示信號的特性，具有很多的優(yōu)點: 1）

16、它克服了傅里葉分析時域和頻域完全分離的缺陷，將時頻兩域聯(lián)合起來對信號進行分析，能同時考慮到兩個方面的性能。 2）彌補了信號的時間能量密度和頻譜能量密度不能充分描述信號的物理特性的缺陷。 3）在時頻相平面上，可以精確地定位在某一時刻出現(xiàn)了哪些頻率分量，以及某一分量出現(xiàn)在哪些時刻。 4）對于不同情況的信號，通過對核函數(shù)施加一些約束條件，就可以設(shè)計符合期望性能的時頻分布，來滿足處理信號的要求。 5）為非平穩(wěn)信號的分析提供了有效的工具，為信號的分析開辟了新的途徑。時頻分析雖然具有很多優(yōu)點，但同時也具有不少缺點: 1）由于雙線性形式的時頻分布是非線性的，使得兩個信號和的時頻分布已不再是兩個信號各自分布的

17、和，即存在交叉項。 2）在時頻域進行去噪時計算復(fù)雜，不易實現(xiàn)。 3）時頻分布雖然反映了信號的能量分布，但不能用信號的“瞬時能量”來解釋某一時刻或某一頻率處的時頻分布。 現(xiàn)實中很多信號，比如語音信號，都是時變非平穩(wěn)的，時變非平穩(wěn)特性是現(xiàn)實信號的普遍規(guī)律，聯(lián)合時頻分析技術(shù)正是應(yīng)現(xiàn)實的科學(xué)和工程應(yīng)用需求而產(chǎn)生和發(fā)展起來的。對于許多信號，僅用時域或頻域里的各種方法去分析往往不能揭示信號內(nèi)部的局部特征和信息，而時頻分析作為一種能將頻譜隨時間的演變關(guān)系明確表現(xiàn)出來的新手段，自然更符合實際應(yīng)用的需要。參考文獻1. 鄒紅星,周小波,李衍達. 時頻分析:回溯與前瞻J. 電子學(xué)報,2000,09:78-84.2.

18、 劉葵,劉招君,朱建偉,張峰. 時頻分析在石油地球物理勘探中的應(yīng)用J. 世界地質(zhì),2000,03:282-285.3. 劉麗娟. 時頻分析技術(shù)及其應(yīng)用D.成都理工大學(xué),2008.4. 陳斌. 時頻分析及在地震信號分析中的應(yīng)用研究D.成都理工大學(xué),2007. 5. 于振江. 一種基于Matlab的語音信號采集與分析系統(tǒng)設(shè)計J. 科技情報開發(fā)與經(jīng)濟,2012,12:112-114.6. 陳家焱,陳冬嬌,張達響. 基于Matlab的聲音信號采集與分析處理J. 計算機與現(xiàn)代化,2005,06:91-92+96.7. 題原,張勁松. 基于MATLAB的語音信號采集和分析系統(tǒng)的可視化設(shè)計J. 齊齊哈爾大學(xué)

19、學(xué)報,2006,06:43-46.8. 蔣濟同,唐世振. 基于MATLAB的振動信號采集與分析系統(tǒng)的研究A. 中國力學(xué)學(xué)會、中國振動工程學(xué)會、中國航空學(xué)會、中國機械工程學(xué)會、中國宇航學(xué)會.第九屆全國振動理論及應(yīng)用學(xué)術(shù)會議論文集C.中國力學(xué)學(xué)會、中國振動工程學(xué)會、中國航空學(xué)會、中國機械工程學(xué)會、中國宇航學(xué)會:,2007:5.9. 宋楊潔. 基于LabVIEW與MATLAB的語言信號的采集與分析D.武漢理工大學(xué),2012.10. 周淵,王炳和,劉斌勝. 基于MATLAB的噪聲信號采集與分析系統(tǒng)的設(shè)計J. 電聲技術(shù),2004,07:52-54.附：1、時頻分析基于Matlab的實現(xiàn)clc;clear

20、;nfft=512;fid=fopen('xinhao1.txt','r');fclose(fid);%x=u(:,200);figure(1);tfrstft(u);figure(2);tfrsp(u);figure(3);tfrwv(u);2、 不同的Wigner-Ville變換clear all;N=2048;fs=2048;n=0:N;t=n*1/2000;x1=cos(2*pi*200*t).*(t>=0 & t<=0.8);x2=cos(2*pi*50*t).*(t>=0.2 & t<=1);x=x1+x2;figure(1);plot(t,x);title('原始信號');xlabel('時間 t');ylabel('幅值 A');grid on;figure(2);subplot(311);tfr,t,f=tfrwv(x',1:N,N);cont