第8章 小波分析及其MATLAB 实现
本文档节选自:
Matlab 在数学建模中的应用,卓金武等编著, 北航出版社,2011年4月出版
8.1小波分析基本理论
8.1.1 Fourier 变换的局限性
8.1.1.1变换的含义
我们把那些定义域和因变域都不是数值或常量的函数称为变换或算子,它们是定义域和值域本身为函数集的函数,如傅里叶变换(Fourier Transform )和拉普拉斯变换(Laplace Transform ),其定义域是时间的函数,而因变域是频率的函数。简单地说,变换的基本思想仍然是映射,变换是函数的函数。
8.1.1.2 Fourier 变换的局限性
信号分析的主要目的是寻一种简单有效的信号变换方法,使信号包含的重要特征能显示出来。在小波变换兴起之前,Fourier 级数展开和Fourier 分析是刻画函数空间、求解微分方程、进行数学计算的有效方法
和有效的数学工具。从物理直观上看,一个周期性振动的量可以看成是具有简单频率的简谐振动的叠加。Fourier 级数展开则是这一物理过程的数学描述,Fourier 变化和Fourier 逆变换公式如下:
函数)()(1R L t f ∈的连续Fourier 变换定义为t t f f
t d e )()(ˆi -⎰
+∞
-=ωω
)(ˆωf 的Fourier 逆变换定义为  ωωωd e )(ˆπ
21)(i ⎰+∞∞-=t f t f  从公式上看,Fourier 变换是域变换,它把时间域和频率域联系起来,在时间域上难以
观察到的现象和规律,在频域往往能十分清楚地显示出来。频谱分析的本质就是对)(ˆωf
的加工、分析和滤波等处理。Fourier 变换是平稳信号分析的最重要的工具。然而在实际运用中,所遇到
的信号大多数并不平稳,而是时变频率信号,这时人们需要知道信号在突变时刻所对应的频率成分,显然Fourier 变换的积分作用平滑了非平稳过程的突变部分,作为积分核t
ω-i e
的幅值在任何情况下均为1,因此,在频谱)(ˆωf
的任一频点值是由时间过程)(t f 在整个时间域上)(∞+-∞,
上的贡献决定的;反之,时间过程)(t f 在某一时刻的状态也是由)(ˆωf
在整个频率域)(∞+-∞,上贡献决定的。也就是说,Fourier 变换是整体变换,只能反映信号的整体特征,而对信号的局部(奇异性)反映不敏感,对描述信号激烈震荡的细节无能为力[1]。
1946年诺贝尔奖获得者D.Gabor 引入了短时傅里叶变换(Short-time Fourier Transform),虽然短时傅里叶变换随着窗口在时间轴t 上滑动可以抠取频谱上的所有信息,但从本质上讲,短时傅里叶变换是一种单一分辨率的信号分析方法,因为它使用一个固定的短时窗函数(常用Gauss 函数),窗口大小缺乏自适应功能,窗口位置可以移动,但是形状不能改变。
下面以声波为例来说明傅里叶变换的局限性。声音的特征与振幅、时间、频率以及波形有关联。音强对
应于振幅大小;音长对应于声波持续时间;音高对应于频率高低;音质不同时波形有别。Fourier 变换无法反映信号在哪一时刻有高音,在哪一时刻有低音,因此所有的音符都挤在了一起。小波变换有效地克服了Fourier 变换的这一缺点,信号变换到小波域后,小波不仅能检测到高音与低音,而且还能将高音与低音发生的位置与原始信号相对应。可以非常直观地知道什么时候发出了什么频率的音符,如图8-1所示。
图8-1 声音信号经过Fourier 变换和Wavelet 变换
8.1.2 伸缩平移和小波变换
从事石油信号处理的法国工程师J.Morlet 在研究地下岩石油层分布时,于1974年首先提出了小波变换的概念,并且J.Morlet 构建了光滑柔性的、时频域局部性能都较好的Morlet 小波函数。数学家Meyer ,Mallat ,Daubechies ,K.Chui 等人的工作为小波分析的诞生和发展奠定了基础[2]
。同时计算机的发展,也为小波分析的成长提供了有利土壤,有些不能用解析式表示的小波,借助于计算机得到了广泛的应用。小波分析是纯数学、应用数学和工程技术的完美结合。小波分析的核心是小波函数的构建和多尺度分析。小波函数的主要特质是快速衰减性和震荡性,其子函数都是相互正交。这里正交的概念不是狭义上的垂直,而是指不能用任意1-N 个子小波来表征第N 个子小波(假如母小波函数经过伸缩平移得到N 个子小波)。就像平面向量永远无法表征立体向量一样。
8.1.2.1 小波的定义
小波(Wavelet)是由英文单词“Wave ”和“小”的后缀“-let ”构成,表示的是一种长度有限,平均值为0的波形。小波函数的确切定义为:设)(t ψ为一平方可积函数,即
)()(2R L t ∈ψ,若其Fourier 变换)(ˆωψ
满足条件: ⎰∞<=R
d C ωωωψ
ψ||)(ˆ2
(8-1)
则称)(t ψ为一个基本小波、母小波或者容许小波,我们称式(8-1)为小波函数的可容许条件。)(2
R L 表示满足
+∞<R
dt t f 2
)(的函数空间。更一般地,)(R L p 表示满足
+∞<R
p
dt t f )(的函数空间。小波)(t ψ有两个基本特点:
一是“小”,在)(2R L 空间内,我们常选取紧支集或近似紧支集(具有时域的局限性)且具有正则性(具有频域的局限性)的实数或复数做为小波母函数,这样的小波母函数在时频域都会具有较好的局部特性,即快速衰减性。
二是正负交替的“波动性”,也即直流分量为零,即震荡性。
如下图8-2中图(a )正弦函数不是小波函数因为它不具有快速衰减性,是无限延伸的;而图(b )也不是小波函数因为它既不具有正负交替的震荡性也不具有衰减性;图(c )和图(d )为小波函数,分别是4阶复Gaussian 小波的实部和虚部部分的函数图像。
图8-2 波形示意图
8.1.2.2母小波伸缩平移
将小波母函数)(t ψ进行伸缩和平移,就可以得到子函数)(,t b a ψ:                  )(
1)(,a
b
t a
t b a -=
ψψ,  R b a ∈,,0>a            (8-2) 我们称)(,t b a ψ为依赖于参数a ,b 的小波家族子函数。上式中,a 为伸缩因子,用来控制家族子小波)(,t b a ψ图像的“体形”(“胖瘦”和“高矮”),b 为平移因子,用来控制小波子函数的中心位置。由于尺度因子a 和平移因子b 是连续变化的值,因此称)(,t b a ψ为连续小波函数基。它们是由同一母函数)(t ψ经伸缩和平移得到的一组函数序列。下图8-3为Morlet 小波经伸缩平移后得到的几何图像。
-10
010
x
y
图(a )
-10
010
x
y
图(b )
-5
5
x
y短时傅里叶变换matlab程序
图(c )
-5
5
x
y
图(d )
图8-3 Morlet 小波母函数经过伸缩平移在时频域的图像
关于式(8-2)两点说明:
(1)由同一母小波经过伸缩和平移得到的小波家族子函数序列均是相互正交的。 (2))
(
1)(,a
b t a
t b a -=ψψ而不是
)()(,a b t t b a -=ψψ是为了使所有子小波以及其母函数具有相等的能量。
8.1.2.3连续小波变换
定义符号>∙∙<,为内积运算符。假设在区间),(21t t 定义两个连续时间实能量信号)(t x 和)(t y ,则此两个信号的内积表达式为:⎰
>=
<2
1
)()()(),(t t dt t y t x t x t y 。连续信号)(t f 在
小波基下展开称为连续小波变换(Continue Wavelet Transform ,CWT ),其表达式为:
->=
=<R
b a f dt a
b
t t f a
t t f b a WT )(
)
(1)(),(),(*,ψψ 其中,
)(*a
b
t -ψ是)(a b t -ψ的共轭函数。所谓“共轭”就是实部相同虚部相反。显然,当)(,t b a ψ为实小波函数时,有:)()(*
a
b
t a b t -=-ψψ。 8.1.3小波变换入门和多尺度分析
8.1.3.1小波变换的入门----平均和细节
小波知识艰深难懂,一般初学者难以掌握小波分析的核心思想, 孙延奎教授编著的文献
[3]
将小波本质解说的简明透彻,通俗易懂,读者可以系统参阅一下。
设{}21,x x 是一个由两个元素组成的信号,定义这两个元素的平均和细节为:
2
/)(2/)(2121x x d x x a -=+=
则可将{}d a ,作为原信号另一种表示,且原信号{}21,x x 可由{}d a ,恢复如下:
d
a x d a x -=+=21
信号{}d a ,是{}21,x x 的线性变换结果,可以看出,当1x 与2x 非常接近时,d 会很小,这时{}21,x x 可近似地用{}a 来表示,由此实现信号的压缩。重构信号为{}a a ,,误差信号为
{}{}d d a x a x
,,21
=--。
特别地,当21x x =时,0=d ,这时可以实现信号的无损重构。这表明,平均值a 可以看成是原信号的整体信息,而d 可以看成是用a 表示原始信号所丢失的细节信息。用{}a 近似地表示{}21,x x 可以实现信号压缩,而且丢失的很小细节信息对最终信号的重构不会造成影响。
对于较多元素的信号{}4321,,,x x x x ,通过如下平均和细节运算, 平均:
2
/)(2/)(431,1210,1x x a x x a +=+=
细节:
2
/)(2/)(431,1210,1x x d x x d -=-=
可得到原信号的另一种表示{}1,10,1,1,10,1,,d d a a 。若细节信号1,10,1,d d 都很小,则丢失细节信号,可得压缩后的信号为{}
1,10,1,a a 。用类似的方法对压缩后的信号{}
1,10,1,a a 进行压缩:
平均:      2/)(1,10,10,0a a a += 细节:      2/)(1,10,10,0d d d +=
如果0,01,10,1,,d d d 都非常小,则我们可以用{}
0,0a 代替原信号{}4321,,,x x x x 。事实上,
4/)(2/)2/)(2/)((2/)(432143211,10,10,0x x x x x x x x a a a +++=+++=+=
即0,0a 为整个信号所有元素的平均值,它保留了信号最基本的信息。
易知,由{}
0,00,0,d a 可以重构{}1,10,1,a a ,由{}1,10,1,a a 和{}
1,10,1,d d 可以重构
{}4321,,,x x x x ,于是,由{}1,10,10,00,0,,,d d d a 可恢复{}4321,,,x x x x ,故我们到了原信号

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。