时频分析在地震数据处理中的应用

时频分析在地震数据处理中的应用

摘要：经典的傅立叶分析只适用于分析平稳信号，而不适用于非平稳信号。为了分析非平稳信号，我们采用时频分析方法。时频分析能够清楚的揭示信号的时变频谱特征，是对时变、非平稳信号进行分析与处理的有力工具。本文介绍了时频分析中常用的STFT和Wigner－Ville 分布和S变换的特点，并利用Matlab对一个地震信号进行时频分析的实现。

关键词：时频分析；短时傅立叶变换（STFT）；S变换；Wigner－Ville分布

引言

在传统的信号处理领域，基于Fourier 变换的信号频域表示及其能量的频域分布揭示了信号在频域的特征，它们在传统的信号分析与处理的发展史上发挥了极其重要的作用。但是，Fourier 变换是一种整体变换，即对信号的表征要么完全在时域，要么完全在频域，作为频域表示的功率谱并不能告诉我们其中某种频率分量出现在什么时候及其变化情况。然而，在许多实际应用场合，信号是非平稳的，其统计量（如相关函数、功率谱等）是时变函数。这时，只了解信号在时域或频域的全局特性是远远不够的，最希望得到的乃是信号频谱随时间变化的情况。为此，需要使用时间和频率的联合函数来表示信号，这种表示简称为信号的时频表示。

时频分析方法旨在通过构造一种时间和频率的密度函数，将一个一维的时间信号以二维的时间一频率函数形式表示出来，以揭示信号中所包含的频率分量及其随时间的变化特性。这使我们不但能够同时掌握非平稳信号的时域及频域信息，而且可以清楚地了解非平稳信号的频率是如何随时间变化的。

通过时频分析方法技术对地震信号进行分析处理可在获得地震信号的瞬时频率、瞬时相位、瞬时振幅等瞬时参数的同时获得时频谱图等重要时频域信息，实现对地震信号的边缘检测、属性提取等。

时频分析是非平稳信号分析处理领域的重要方法，时频分布的基本任务是建立一个函数，要求这个函数能够同时用时间和频率来描述信号的能量密度。如果

有了这样的一个分布，就可以计算某一确定的频率和时间范围内能量的百分率、

计算某一特定时刻的频率密度、计算该分布的整体和局部的各阶矩。即寻找一个联合密

∆∆在时间t和频率f，在度函数P(t，f)，使p(t，f)=在时间t和频率f的强度，或者P(t，f)t f

∆∆内的部分能量。

时一频单元t f

地震信号的时频分析

一般将时频分析方法分为线性和非线性两种。典型的线性时频表示有短时傅叶变换(简记为STFT)、Gabor展开和小波变换(Wavelet Transformation，简记为WT)等。非线性时频方法是一种二次时频表示方法(也称为双线性)，最典型的是WVD(Wigner-Ville Distribution)和Cohen类。

采用的地震信号如下图所示：

图1 时间域的信号

图2 用FFT方法求取的能量谱

1、短时傅里叶变换（STFT）

传统的傅立叶变换只在频率域具备了局部分析的能力，而在时间域不具有这种能力。要得到满足同时获得时间和频率的局部分析能力的要求，一种最基本的方法就是:取出信号在所关心时刻附近的一小段，而忽略信号的其它部分，对其作傅立叶变换，即可得到这一特定时刻的频率分量。因为所取的时间长度与整个信号相比很短，所以将这种方法称为短时傅立叶变换(STFT)，它是时频分析中最简单的形式。短时傅立叶变换(STFT)的基本思想:用窗函数来截取信号，假定信号在窗内是平稳的，采用傅立叶变换来分析窗内信号，以便确定在那个时间存在的频率，然后沿着信号移动窗函数，得到信号频率随时间的变化关系，这就得到了我们所需要的时频分布。

STFT t f可以视为信号s(t)在该时刻STFT的物理意义在于，对于一定的分析时刻t, (,)

s

的“局部频谱”，从而整个变换的结果也就能揭示信号频谱的变化特性。短时Fourier变换的时频分辨率受制于窗函数的形状和宽度短时傅立叶变换的时间分辨率与分析窗函数的时间域宽度成正比，而其频率分辨率与分析窗的频宽成正比。从而，一个好的时间分辨率需要一个短的窗函数，而一个好的频率分辨率需要一个长的窗函数。因此，短时傅立叶变换不能同时兼顾时间分辨率和频率分辨率。当选用的窗函数为Gaussian函数时，该变换为Gabor变换。

图3 用STFT 做时频分析

2、S 变换

主要对S 变换的定义、推导和特性进行详细阐述，S 变换综合短时傅立叶变换和小波变换的优点，又避免了它们的不足:它与傅立叶变换有着直接的联系，具有无损可逆性;与短时傅立叶变换和小波变换一样，也是一种线性时频表示，因此不存在交义项的干扰;S 变换具有多种分辨率，克服了短时傅立叶变换固定分辩率的不足;S 变换中含有相位因子，这是小波变换所不具备的特性。总体来说，S 变换是近几年发展起来的一种新的时频分析方法.S 变换结合了短时傅里叶变换和小波变换的优点,具有相位信息,同时该变换与小波变换一样,其时频窗可以调节大小以适应非平稳信号的特点.S 变换的这些优点,使得它在地球物理方面得到广泛的应用。

S 变换首先是由Stockwell 等人提出的，是以Morlet 小波为基本小波的连续小波变换的延展。在S 变换中，简谐波与高斯函数的乘机构成了基本小波，因为简谐波在时间域可以作伸缩变换，而高斯函数则进行伸缩和平移。ST 也可以认为是CWT 的“相位校正”。函数()h t 的S 变换表示为：

2(,)()()j ft ST t f h t e d πτωττ+∞

--∞=-⎰ 其中：

2

222())2()(,)()exp(2)2t t t f ST t f h j f d ωσττπττ+∞-∞=--=-⎰ 可以用2

222()exp()W f

πγγ=-表示()t ω的傅里叶变换，其中γ和f 是相同的意义。

图4 S 变换时频分析

3、WVD

Wigner-Ville 变换是1932年由Wigner 首次提出的，并应用于量子力学领域，后来Ville 等人将其引入到信号分析处理领域。1966年，Cohen 发现各种发现各种时频分析只是Wigner-Ville 变换的不同形式，可以统一起来，成为Cohen 类双线性时频分析。

信号()s t 的Wigner-Ville 变换用公式表示为

*(,)().().exp(2)22z W t f z t z t j f d ττ

πττ+∞

-∞=+--⎰ 式中()z t 是()s t 的解析信号，即

()()()()()()[()]s z t s t js t h t s t j d jH s t t μμμ+∞

-∞

=+⊗=++-⎰

[()]H s t 是实信号()s t 的Hilbert 变换。 Wigner-Ville 变换也可以用解析信号的频谱来表示

*(,)().().exp(2)22

z W t f Z f Z f j t d ννπνν+∞-∞=-+-⎰ 从上面两种不同形式的Wigner-Ville 变换表达式中可以看出，式中不包含任何窗函数，从而避免了线性时频变换中实践分辨率和频率分辨率的相互牵制，难以兼顾的问题。因为Wigner-Ville 变换的时间-带宽积可达到Heisenberg 测不准原理(不确定原理给出的下界，故可以证明，没有任何一种时频变换方法的时间-频率分辨率及聚集性能出其右。

在Wigner-Ville 变换表达式中，信号()s t 出现了两次，故称之为双线性时频变换。Wigner-Ville 变换不是线性的，即两信号之和的Wigner-Ville 变换不等于每一个信号的Wigner-Ville 变换之和，其中多出了一个附加项。

令12()()()s t s t s t =+

则有

1212(,)(,)(,)2Re{(,)}s s s s s W t f W t f W t f W t f =++