陳亞東，戴華林

(天津城建大學(xué) a. 計算機(jī)與信息工程學(xué)院；b. 計算中心，天津 300384)

信息科學(xué)與技術(shù)

基于短時Fourier變換的瑞雷面波相速度提取

陳亞東a，戴華林b

(天津城建大學(xué) a. 計算機(jī)與信息工程學(xué)院；b. 計算中心，天津 300384)

瑞雷面波勘探的主要思想是利用瑞雷面波的頻散特性來反映淺層地質(zhì)問題.在分析常規(guī)Fourier變換提取瑞雷面波相速度的局限性的基礎(chǔ)上，提出采用短時Fourier變換方法對實(shí)際資料進(jìn)行處理，結(jié)果表明：短時Fourier變換提取瑞雷面波相速度的方法具有較高的準(zhǔn)確性和可靠性，有很好的應(yīng)用效果和發(fā)展前景.

瑞雷面波；頻散特性；時-頻分析；短時Fourier變換

瑞雷面波勘探主要是根據(jù)瑞雷面波的頻散特性，通過對原始數(shù)據(jù)的分析，并進(jìn)一步反演，以獲得淺層、超淺層的地質(zhì)結(jié)構(gòu)及介質(zhì)物理參數(shù)．近些年來，該方法發(fā)展迅速，方興未艾，廣泛應(yīng)用于交通、水電、橋梁和地基勘察等方面，已逐步發(fā)展成為一套獨(dú)立的、完整的勘探方法[1]．瑞雷面波勘探的基本工具是常規(guī)Fourier變換，但這種方法只能反映信號的整體特性，不能反映信號的局部特征，尤其是不能同時觀測信號在時間域和頻率域的變化情況．這在信號分析中就面臨一對最基本的矛盾：時域和頻域的局部化矛盾．筆者提出采用短時Fourier變換方法，以實(shí)現(xiàn)對瑞雷面波相速度的有效提?。?/p>

1 瑞雷面波勘探的理論基礎(chǔ)

瑞雷面波勘探最根本的依據(jù)就是瑞雷面波在層狀介質(zhì)中傳播時具有頻散特性．彈性波理論證明，在自由空間和彈性半無限介質(zhì)中，瑞雷面波傳播時的法向應(yīng)力和切向應(yīng)力都為零，此時，瑞雷面波無頻散特性，但這只是理想狀態(tài)．在實(shí)際中，地表總有一層薄的疏松覆蓋層，這就導(dǎo)致了法向應(yīng)力不再為零，而切向應(yīng)力仍為零．正是由于法向應(yīng)力的改變，使得瑞雷面波在層狀介質(zhì)中傳播時的速度發(fā)生了變化[2]，表達(dá)式為

從式(1)可以看出，瑞雷面波的相速度VR是頻率ω的函數(shù)，因此疏松層的存在導(dǎo)致了瑞雷面波的頻散．瑞雷面波的頻散特性為人們提供了提取瑞雷面波相速度的可能性，進(jìn)而，不同相速度的單色波，看到的是不同深度的介質(zhì)性質(zhì)，組合起來，即可以清晰地得到近地表附近的介質(zhì)分布情況．因此，工程中利用瑞雷面波進(jìn)行淺層、超淺層勘探有著積極的意義．

2 瑞雷面波相速度的提取方法

2.1 常規(guī)Fourier變換方法

將接收到的信號進(jìn)行一系列的譜分析，得到兩信號的自功率譜、互功率譜及相干函數(shù)；由互功率譜得到兩信號在波傳播過程中由于時間滯后所產(chǎn)生的相位差，即求得實(shí)測的頻散曲線，再由成層地基動力學(xué)理論及反分析法，求得剪切波速度和深度的關(guān)系．其計算過程簡述如下[3]．

(1)時間域向頻率域的轉(zhuǎn)換，即得到N道信號X1(t),X2(t),X3(t),…,Xn(t)的頻譜[4]

(2)自功率譜、互功率譜和相位差的求?。孟噜彽赖念l譜可求出Xn-1(t)、Xn(t)的自功率譜

因此，Xn(t)、X(n-1)(t)的互功率譜Sn(n-1)(f)可計算求得

由互功率譜式(4)可求出相位差Δφ．

(3)相速度的確定．利用公式

可計算出離散點(diǎn)f的瑞雷面波平均速度，其中，VR為瑞雷面波平均速度；f為頻率；Δx為相鄰道之間的距離，一般情況下，Δx應(yīng)小于一個波長．

(4)深度-速度函數(shù)的確定．在國內(nèi)外的瑞雷面波勘探中，確定深度和速度一般采用所謂的半波長解釋法[5]，即勘探深度是波長的一半，即

(5)繪制頻散曲線．根據(jù)以上分析計算，即可得到深度-速度函數(shù)(H-VR)，進(jìn)一步可以得到測試點(diǎn)的頻散曲線，如圖1所示．實(shí)際資料見圖1a．該記錄為12道，1,024個采樣點(diǎn)，道間距為2,m，采樣周期為0.000,5,s，偏移距為零；首先對原始記錄進(jìn)行濾波，去除反射波的干擾，圖1b為濾波后結(jié)果；再將信號從時間域轉(zhuǎn)換到頻率域，得到信號的f-k圖，在有效頻率范圍內(nèi)拾取每一f極大值及相應(yīng)的k值，即頻散曲線的屏幕拾取，見圖1c；據(jù)此計算相速度，最后繪制頻散曲線圖1d．

圖1 常規(guī)Fourier變換處理實(shí)際資料

由此可見該方法存在以下缺陷：①利用互相關(guān)譜求相位差，受到Fourier 變換本身的局限，由于Fourier變換的相位變化范圍為［-π,π］之間，這要求兩個檢波點(diǎn)距離不得大于接收到的波的最小波長，否則計算所得相位差可能與實(shí)際相差2π的整數(shù)倍，這就造成Fourier 變換提取的相位差可能是不真實(shí)的；②對于頻率很低的單色波，在滿足上述檢波器距離的要求下，其相位差可能很小，甚至趨于零，而相位差Δφ在式(5)中作為分母，這樣，一個極小的誤差Δφ都會帶來極大的速度誤差．

2.2 短時Fourier變換方法

針對常規(guī)Fourier變換本身的局限性，人們相繼提出了短時Fourier變換、小波分析以及Gabor展開等[6]．與常規(guī)Fourier變換不同，短時Fourier變換可

以同時獲得時間(空間)域和頻率域的分辨率，它既可以確定信號的頻率成分，又可以確定某一頻率在時間域的分布規(guī)律．其基本思想是將非平穩(wěn)信號看成分段平穩(wěn)信號，那么通過在原始信號上加固定窗，然后由窗口在原始信號上用平移的方法，就可以得到原始信號的時頻譜[7]．

對于一個給定時間寬度的窗函數(shù)r(t)，沿時間軸滑動，則信號s(t＇)的短時Fourier變換定義為[8]

式中：“*”表示復(fù)共軛；r(t)是有限支集的函數(shù)．在這個變換中，e-j2πft＇起著頻限的作用，r(t)起著時限的作用．隨著時間t的變化，r(t)所確定的“時間窗”在t軸上移動，使s(t')“逐漸”進(jìn)行分析．由式(7)可以看出，正是窗函數(shù)r(t)的時(間)移(位)和頻(率)移(位)使短時Fourier變換具有了局域特性，它既是時間的函數(shù)，又是頻率的函數(shù)．式(7)即為信號s(t＇)的時頻分布函數(shù)，利用這個時頻分布函數(shù)，可以輕松地求得瑞雷面波的相速度．

(1)對相距為Δx的N個檢波器接收到的信號作時-頻分析，可以得到N個時頻分布圖．

(2)對于每一個固定頻率f，分別在每個時頻分布圖上尋找該頻率能量極值點(diǎn)所對應(yīng)的時間t1,t2,…，tN，進(jìn)而可以計算出該頻率在相鄰兩個時頻分布圖上的時間差Δt=tn-tn-1(0＜n≤N)，這個時差就是該頻率以其相速度在相鄰兩個檢波器之間傳播形成的時間差．

(3)相速度的計算．利用公式[9]

可計算出每一固定頻率f在相鄰兩道傳播的相速度．對于實(shí)測的信號，可能有N道記錄，可以對其作統(tǒng)計分析，得到相速度的統(tǒng)計量．

(4)頻散曲線的繪制．根據(jù)第(3)步的計算結(jié)果，可以繪出瑞雷面波相速度隨頻率的變化曲線，如圖2所示．對上述實(shí)際記錄選取其中第7道和第9道短時Fourier變換，得到的時頻分布見圖2a、圖2b．在兩道信號的時頻圖上尋找相同頻率的極值點(diǎn)位置，如圖2a和圖2b中14,Hz的極點(diǎn)(即圖中白色十字游標(biāo)中心位置)分別位于96，102,ms位置上，則該頻率的單色波傳播在兩個道間形成的時差為6,ms，相速度為2/0.006＝333.4,m/s．對每一道數(shù)據(jù)作短時Fourier變換，得出相速度隨頻率變化的統(tǒng)計量，進(jìn)而繪制頻散曲線圖2c．

顯然，短時Fourier變換處理結(jié)果比常規(guī)Fourier變換的精度要高，漏掉的巖層被顯示出來．由于短時Fourier變換縱橫坐標(biāo)控制精確，其穩(wěn)定性強(qiáng)，該方法在瑞雷面波相速度提取中必將有廣闊的應(yīng)用前景．

圖2 短時Fourier變換處理實(shí)際資料

3 結(jié) 論

短時Fourier變換克服了常規(guī)Fourier變換利用互相關(guān)譜法求相位差，進(jìn)而再求相速度的缺陷．短時Fourier變換求取相速度的精度取決于時間分辨率，而對于合格的面波勘探資料，完全可以保證短時Fourier變換的時間分辨率，從而得到的是確切的相速度．因此，該方法可以極大地提高瑞雷面波相速度

提取的穩(wěn)定性和可靠性，從而提高工程勘探地基土分層的準(zhǔn)確性．

［1］ 楊成林. 瑞雷波勘探[M]. 北京：地質(zhì)出版社，1993.

［2］ 張大洲，熊章強(qiáng)，秦 臻. 基于 Fourier 變換的瑞雷面波分離提取及實(shí)例分析[J]. 中南大學(xué)學(xué)報：自然科學(xué)版，2010，41(2)：643-648.

［3］ 程乾生. 信號數(shù)字處理的數(shù)學(xué)原理[M]. 北京：石油工業(yè)出版社，1979.

［4］ CHEN Xiao. Leaky Rayleigh wave at the interface between inviscid fluid and solid[J]. Chinese Journal of Theoretical and Applied，2010，42(1)：51-55.

［5］ 胡昌華，周 濤. 基于MATLAB的系統(tǒng)分析與設(shè)計——時頻分析[M]. 西安：西安電子科技大學(xué)出版社，2002.

［6］ LIU Bin，DONG Shiyun，XU Binshi，et al. Study on surface stress measurement of laser cladding Fe-based alloy coating based on Rayleigh wave signal analysis [J]. Advanced Materials Research，2012，399/400/401：2,177-2,180.

［7］ PHAM-CHI Vinh，PETER G M. Improved approximations of the Rayleigh wave velocity[J]. Journal of Thermoplastic Composite Materials，2008，21(4)：337-352.

［8］ LI Aibing. Influence of surface stress and atomic defect generation on Rayleigh wave propagation in laser-excited solids[J]. Physica B：Condensed Matter，2013，42(1)：28-33.

［9］ 豐 赟，周竹生，沙 椿. 瑞雷波數(shù)值模擬中的數(shù)值頻散校正策略及實(shí)例分析[J]. 中南大學(xué)學(xué)報：自然科學(xué)版，2012，43(6)：2,231-2,237.

The Rayleigh Surface Wave Phase Velocity Extraction Based on Short-time Fourier Transform

CHEN Ya-donga，DAI Hua-linb
(a. School of Computer and Information Engineering；b. Computing Center，Tianjin Chengjian University，Tianjin 300384，China)

The main idea of Rayleigh wave exploration is making use of frequency dispersion of Rayleigh wave to reveal geology structure of low layer. Based on the limitation analysis of extraction of the Rayleigh wave speed by conventional Fourier transform,authors present a new use actual data processing method based on the short-time Fourier transform. The results show that the extraction of the Rayleigh wave speed by short-time Fourier transform is high accuracy and reliability,which has a good application effect and development prospects.

Rayleigh wave；frequency dispersion；time-frequency analysis；short-time Fourier transform

P631.4

A

2095-719X(2015)04-0299-04

2014-12-02；

2015-04-20

天津市自然科學(xué)基金(10JCYBJC00600)

陳亞東（1974—），女，內(nèi)蒙赤峰人，天津城建大學(xué)副教授，碩士．