Matlab小波变换与傅立叶变换用于检测阶跃信号

Matlab 小波变换与傅立叶变换用于测阶跃信号
1 从傅立叶变换到小波变换
小波分析属于时频分析的一种，传统的信号分析是建立在傅立叶变换的基础上的，由于傅立叶分析使用的是一种全局的变换，要么完全在时域，要么完全在时域，要么完全在频域，因此无法表述信号的时频局域性质，而这种性质恰恰是非平稳信号最根本和最关键的性质。

为了分析和处理非平稳信号，人们对傅立叶分析进行了推广乃至根本性的革命，提出并发展了一系列新的信号分析理论：短时傅立叶变换、Gabor 变换、时频分析、小波变换、分数阶傅立叶变换、线调频小波变换、循环统计量理论和调幅-调频信号分析等。

其中，短时傅立叶变换和小波变换也是应传统的傅立叶变换不能够满足信号处理的要求而产生的。

短时傅立叶变换分析的基本思想是：假定非平稳信号在分析窗函数g(t)的一个短时间间隔内是平稳（伪平稳）的，并移动分析窗函数，使
)
()(τ-t g t f 在不同的有限时
间宽度内是平稳信号，从而计算出各个不同时刻的功率谱。

但从本质上讲，短时傅立叶变换是一种单一分辨率的信号分析方法，因为它使用一个固定的短时窗函数。

因而短时傅立叶变换在信号分析上还是存在着不可逾越的缺陷。

小波变换是一种信号的时间—尺度分析方法，它具有多分辨率分析的特点，而且在时频两域都具有表征信号局部特征的能力，是一种窗口大小固定不变但其形状可改变，时间窗和频率窗都可以改变的时频局部化分析方法。

即在低频部分具有较高的频率分辨率，在高频部分具有较高的时间分辨率和较低的频率分辨率，很适合于探测正常信号中夹带的瞬态反常现象并展示其成分，所以被誉为分析信号的显微镜，利用连续小波变换进行动态系统故障检测与诊断具有良好的效果。

2 Fourier 变换
Fourier 变换由下列公式定义： 正变换
⎰
∞
∞
--∧
=
dt
e
t f f t
j ωω)()(
(2-1)
逆变换
⎰
∞
∞
-∧
=
dt
e
f t f t
j ωωπ
)(21)(
对于确定信号和平稳随机信号，傅里叶变换时信号分析和信号处理技术的理论基础，有着非凡的意义，起着重大作用。

傅里叶变换把时间域与频率域联系起来，)(ω∧
f 具有明确的物理含义，通过
研究
)
(ω∧
f 来研究)(t f ，许多在时间域内难以看清的问题，在频率中往往表现得
非常清楚。

但正是由于傅里叶变换的域变换特性，)(t f 与)(ω∧
f 彼此之间是整体刻画，
不能够反映各自在局部区域上的特征，因此不能用于局部分析。

作为变换积分核
的t
j e
ω±的幅值在任何情况下均为1，即t
j e
ω±=1，因此，频谱
)
(ω∧
f 在任一频率处
的值是由时间过程)(t f 在整个时间域
（∞∞-，）上的贡献决定的；反之，过程
)
(t f 在某一时刻的状态也是由)
(ω∧
f 在整个频率域（∞∞-，）上的贡献决定的。

如果
要知道所分析的信号在突变时刻的频率成分，那么傅里叶变换是无能为力的，因为傅里叶变换的积分作用平滑了非平稳信号的突变成分。

简言之，傅里叶变换能提取出函数在整个频率轴上的频率信息，却不能反映信号在局部时间范围内的特征。

然而，对于变频信号如音乐、地震信号、雷达回波等。

此时所关心的恰恰是信号在局部时间范围（特别是突变时刻）内的信号特征（一般是频率成分）。

例如，在音乐和语音信号中，人们关心的是在什么位置出现什么样的反射波。

对非平稳信号用傅里叶变换进行分析，不能提供完全的信息，也即通过傅里叶变换可以知道信号所含有的频率信息，但无法知道这些频率信息究竟出现在哪些时间段上。

可见，若要提取局部时间段（或瞬间）的频率特征信号，傅里叶变换已经不再适用了。

3 短时Fourier 变换
为了克服傅里叶分析的局限性，使其对非平稳信号也能作较好的分析，通过对信号在时域上加一个窗函数g (t-τ),使其对信号
)(t f 进行乘积运算以实现
在τ附近的开窗盒平移，再对加窗的信号进行傅里叶分析，这就是短时傅里叶变换，或者称为加窗傅里叶变换。

短时傅里叶变换定义如下：
dt
e
t g t f t
j f s ⎰
∞
∞
---=
ωττω)()()
,(
其中，窗口函数f (t )一般取为光滑的低通函数，保证)(τ-t g 只在τ点附近局部地测量了频率分量ω的幅度值，得到信号在t =τ时刻附近的频率信息。

这是时间频率局部化的一种标准技术。

若采用高斯函数作为窗口函数，其相应的傅里叶变换仍旧是高斯函数，从而保证短时傅里叶变换在时域和频域内均有局部化功能。

如果选取的窗口函数在时域和频域都具有良好的局部性质（如呈指数衰减的高斯函数），此时短时傅里叶变换能够同时在频域和时域内提取关于信号的精确信息。

但短时傅里叶变换存在其固有的局限，其时间频率窗口是固定不变的，一旦
窗口函数g (t )选定，其时频分辨率也就确定了，并且不随频率ω和时间τ而变化。

也就是说，它对所有的频率都使用同样的窗口。

我们若想提高时间分辨率，就要把窗口缩的很窄，但这样势必会降低频率分辨率。

由测不准原理可知，不可能在时间和频率上均有任意高的分辨率，因为时间和频率的最高分辨率受下式的制约：
π
ω41≥
∆⋅∆t
式中t ∆和ω∆分别代表时间域和频率域的窗口宽度。

这表明，任一方分辨率的提高都意味着另一方分辨率的降低。

4 小波变换
上节的的分析表明，短时傅里叶变换的问题的症结在于使用了固定的窗口，而对实际时变信号的分析需要时频窗口具有自适应性：对于高频谱的信息，时间间隔要相对地小以给出较高的精度；对于低频谱的信息，时间间隔要相对地宽以给出完全的信息。

换句话说，重要的是要有一个灵活可变的时间-频率窗，能够在“高中心频率”时自动变窄，而在“低中心频率”时自动变宽。

而小波函数就是为此而设计的。

设f (t )是平方可积函数，记作f (t )）（R L 2
∈，)
(t ψ为基本小波或母小波，它一
般是时域上以t =0为中心的带通函数，在时域和频域都具有局部化（紧支撑），
且均值为零，即
0)(=⎰
∞
∞
-dt t ψ
如果)
(t ψ
满足容许性条件
∞
<=
⎰
∞
∧
2
)
(ωω
ωψψd C
则
⎰
∞
∞
--=
=dt
a
b t a
t f f b a Wf b
a )(
1)
(,),(,ψψ
称为f (t )的小波变换。

其中
b a ,ψ= )(
1a
b t a
-ψ, a >0, b ∈R
a 为尺度因子，
b 为位移因子。

小波变换的特点有：
(1) 时频局域性、多分辨分析、数学显微镜；
(2) 自适应窗口滤波：低频宽、高频窄； (3) 适用于去噪、滤波、边缘检测等。

5 连续小波变换
对任意连续函数或信号f (t )进行小波变换，如果基函数b a ,ψ(t )的两个参数a 和b 均为连续变量，则被称为连续变量，则被称为连续小波变换(continuous wavelet transform,CWT)。

连续小波变换的形式可以写为
=
),(b a Wf ⎰
∞
∞
--dt a
b t a
t f )(
1)
(ψ
6 二进小波变换
为了将小波变换应用于信号分析的实践，需要将小波变换进行离散化，就像将傅里叶变换离散化一样。

将小波变换离散化就是对变换参数进行离散化。

本节先讨论尺度因子a 的离散化问题。

在小波变换中，常令a 取离散值 j
j
a
2
= (
Z
j ∈)
称之为二进离散化。

这时，小波基函数的形式为
)2
(
2)(2
/,2j
j b
b t t j
-=-ψψ
对应的小波变换为
=
),2(b Wf j
⎰
∞
∞
--dt b t t f j
j
)2
(
2
1)
(ψ
j 称为分辨率级别。

7 离散小波变换
下面再将二进小波变换中的平移因子也离散化。

令b =j n 2，则可得离散小波变换：
=
)2,2(j
j
n Wf
⎰
∞
∞
--dt n t t f j
j
j
)2
2(
2
1)
(ψ
对小波变换进行离散化的一般情形是小波基函数b a ,ψ的尺度因子a 和位移因子b 都只限定在某些离散点上取值。

一种最通常的离散方法如下：尺度因子按幂级数进行离散化，位移因子在尺度内均匀离散化，在尺度间具有幂次关系，即有：
j
a a
=， j
a n
b b
0=；,0,100>>b a Z
n j ∈,。

此时的小波基函数表示为如下形式：
)
()(
1)(002
/0
00
,nb t a a a a nb t a t j
j j
j
j
n
j -=-=
--ψψψ
任意函数f (t)离散小波变换为
)
(),(),(,t t f n j Wf
n
j ψ
=
⎰∞
∞-=dt
t t f n
j )()(,ψ
⎰
∞
∞
----=dt
nb t a t f a j
j )()(002
/0
ψ
8程序
原始时域信号是一个阶跃信号，为了确定阶跃信号的突变点，分别采用傅立叶变换和db5小波进行处理。

程序如下： clear
load wstep ; whos;
figure(1); plot(wstep)
xlabel('时间');ylabel('幅值'); title('频率突变信号'); figure(2); f=fft(wstep); plot(abs(f));
title('傅立叶变换后的示意图'); figure(3);
[d,a]=wavedec(wstep,3,'db5'); a3=wrcoef('a',d,a,'db5',3); d3=wrcoef('d',d,a,'db5',3); d2=wrcoef('d',d,a,'db5',2); d1=wrcoef('d',d,a,'db5',1);
subplot(411);plot(a3);ylabel('近似信号a3'); title('小波分解后示意图');
subplot(412);plot(d3);ylabel('细节信号'); subplot(413);plot(d2);ylabel('细节信号'); subplot(414);plot(d1);ylabel('细节信号'); xlabel('信号'); 程序的运行结果：
图1原始信号示意图
原始信号载入后由矩阵表示，其中矩阵大小为1*1000，矩阵名为wstep。

矩阵是以双精度表示，相应的图像显示如图1所示。

图2 傅立叶分解后的示意图
利用傅立叶变换对原始信号进行处理，可以得到如图2所示的图像。

从图2可以看出：信号经过傅立叶变换后能够清楚确定出原始信号包含的频率值的大小，但是对于确定频率突变点的位置，傅立叶变换却没有这种能力。

图3 小波分解后的示意图
利用小波变换对原始信号进行处理，可以得到如图3所示的小波分解示意图。

从图中可见，db5小波分解后的3层高频系数重构图形可清楚的确定突变点的位置。

从图3中同样可能看出，第1层分解的d1高频系数重构的图像比d2、d3
高频系数重构的图像更清楚的确定了信号突变点的位置。

命令行输出结果如下：
Name Size Bytes Class Attributes
wstep 1x1000 8000 double。