小波分析作业(含自编程序)

地球科学学院小波分析课程作业

课程名称：小波分析

指导老师：

学生姓名：

学号：

几种时频分析方法

1 短时傅里叶变换

为了研究信号在局部时间范围内的瞬时频率特性，1946年，D.GABOR 引进了短时傅氏变换或窗口傅氏变换的概念，其基本原理是取一个称为(t)g 窗口的函数，使它在有限的区间范围外恒等于零或趋于零。设任意信号(t)f ，并假设该信号在一个以时间τ为中心，且范围有限的窗口函数)-(t τg 内是稳定的，这样，窗口内函数)-(t)g(t τf 的傅氏变换就定义为短时傅氏变换，表示为

dt e T t i STFT ωτωτ--)-f(t)g(t ),(?

∞

∞=

STFT 是通过滑动时窗来计算其频谱，因而它的时间分辨率和频率分辨率受Heisenberg 测不准原理约束。因此利用短窗口有较高的时间分辨率，但是频率分辨率差。 2 小波变换

常见的小波变换有连续、二进制以及离散小波变换等。在连续小波变换中，仅要求小波函数满足容许条件即可，这使得在选择小波函数时具有很大的自由度。对任意地震信号函数)((t)2R L f ∈，其连续小波变换定义为 )f(t)dt a b -t (a 1(t)f(t)dt b)(a,-*-b a,??∞∞∞∞

==?ψw T 式中，ａ为尺度因子，ｂ为平移参数，函数ψ（ｔ）称为母小波。小波分析具有可调的时频窗口，被广泛地应用于地震信号处理中，但是也存在着一定的局限性，主要表现在难以选择小波基、固定的基函数、恒定的多分辨率，信号的能量—时间—频率分布也很难定量给出。 3 S 变换

为了解决短时傅氏变换只能以一种分辨率进行时频分析及小波变换不能直接与频率对应的缺陷，1996年美国地球物理学家Stockwell 在前人的基础上提出了S 变换。S 变换中，基本小波是由简谐波与高斯函数的乘积构成的，基本小波中的简谐波在时间域仅作伸缩变换，而高斯函数则进行伸缩和平移。这一点与连续小波变换不同，在连续小波变换中，简谐波与高斯函数进行同样的伸缩和平移。信号(t)f 一维连续正变换表达式如下 dt e e f S ift f πτπτ2-2)-(t --22|f |(t)f),(?∞

∞=

式中，ｆ为频率，ｔ是时间，τ控制时间轴上高斯窗的位置。由于Ｓ变换采用宽度可变的高斯窗函数，傅氏变换仍然是高斯函数，所以可以达到很好的时频聚集性能。它综合了短时傅氏变换和小波变换的优点，避免了它们的不足：频率的倒数决定了Ｓ变换中的高斯窗的尺度 大小，使信号的Ｓ变换的时频谱的分辨率与频率（即尺度）有关，从而克服了短时傅氏变换不能调节分析窗口频率的问题。同时具有了小波变换的多分辨优点，而且含有相位因子，这是小波变换所不具备的特性。由于Ｓ变换分辨率可自适应调节，能够兼顾高低频分量，保持低频部分较高的分辨率，且不存在交叉项。

小波分析在地震资料去噪中的应用

姓名：学号：

本文研究了利用小波分析在地震去噪中的应用，讨论了如何选择小波基及去噪中的阈值问题。

（1）对于如何选取小波基，我们通过对地震数据进行小波分析后的重构信号与原始信号的误差大小，来选取最优小波基。

（2）对于阈值问题的讨论，将信号进行小波分解, 即求函数与各小波基函数之间相关系数, 也就是求小波变换值。矢量分解去噪实质是基于正交变换思想的相关分析方法。有效信号的小波变换值一般要比干扰波的小波变换值大得多。设dj是第j 个尺度下的信号与小波基函数的相关系数给定阈值，当小波系数dj小于阈值, 则置零, 即将记录中的噪声去掉，按式重建得到去噪后的记录。为了说明这种方法在地震数据中的应用，本文用分解重构、强制消噪方法和默认阈值消噪方法分别对合成地震记录进行了消噪处理。

以下为本文实现消噪过程的matlab程序。

clc;clear;

m=1024;

t0=0.25/1000;

fid=fopen('MB20.DAT','rb');

status=fseek(fid,512,'bof');

B=fread(fid,[m 1],'integer*4');

bb=B'+rand(1,m)*10000;

% 画出原始图像

figure(1);

subplot(411);plot(B);

title('原始信号及其傅里叶变换结果');

% 将原始信号进行傅里叶变换

y1=fft(B,1024);

p1=abs(y1);

ff1=(0:511)/(1024*t0);

subplot(412);plot(ff1(1:512),p1(1:512));ylabel('振幅');title('振幅谱');xlabel('频率');

subplot(413);plot(bb);

yy1=fft(bb,1024);

pp1=abs(yy1);

ff1=(0:511)/(1024*t0);

title('加干扰信号及其傅里叶变换结果');

subplot(414);plot(ff1(1:512),pp1(1:512));ylabel('振幅');title('振幅谱');xlabel('频率');

% ================================================================ % 下面将原始声信号用db6小波进行5层分解

[c,l]=wavedec(bb,5,'db6');

% 提取各层的低频系数并进行重构

ca5=appcoef(c,l,'db6',5);a5=wrcoef('a',c,l,'db6',5);

ca4=appcoef(c,l,'db6',4);a4=wrcoef('a',c,l,'db6',4);