1. 给出五种常用小波基的时域和频域波形图。
与标准的傅里叶变换相比,小波分析中使用到的小波函数具有不唯一性,即小波函数ψ(t) 具有多样性。小波分析在工程应用中,一个十分重要的问题就是最优小波基的选择问题,因为用不同的小波基分析同一个问题会产生不同的结果。目前我们主要是通过用小波分析方法处理信号的结果与理论结果的误差来判定小波基的好坏,由此决定小波基。常用小波基有Haar 小波、Daubechies(dbN)小波、Mexican Hat(mexh)小波、Morlet 小波、Meyer 小波等5种。 (1)Haar 小波
Haar 函数是小波分析中最早用到的一个具有紧支撑的正交小波函数,也是最简答的一个小波函数,它是支撑域在t
∈[0,1]范围内的单个矩形波。 Haar
⎧10≤t ≤1
⎪12
函数的定义如下:ψ(t)=⎨-1≤t ≤1
2⎪⎩0其他
Haar 小波在时域上是不连续的,所以作为基本小波性能不是特别好。但它也有自己的优点,如: ①计算简单;
②ψ(t)不但与ψ(2j t) [j∈z]正交,而且与自己的整数位移正交。 因此,在a 族。
ψ(t ) 的傅里叶变换是:
=2j 的多分辨率系统中Haar 小波构成一组最简单的正交归一的小波
ψ(Ω)=j
4Ω
sin 2() e -j Ω/2 Ωa
Haar 小波的时域和频域波形图
haar 时域
1.5
1
0.5
-0.5
-1
-1.5
00.5
11.5
t
i=20; wav = 'haar';
[phi,g1,xval] = wavefun(wav,i); subplot(1,2,1);
plot(xval,g1,'-r','LineWidth',1.5); xlabel('t') title('haar 时域'); g2=fft(g1); g3=abs(g2);
subplot(1,2,2);plot(g3); xlabel('f') title('haar 频域')
5
haar 频域
7
6
5
4
3
2
1
05
1015f
x 10
5
(2)Daubechies(dbN)小波
Daubechies 小波是世界著名的小波分析学者Inrid Daubechies 构造的小波函数,简写为dbN ,N 是小波的阶数。小波ψ(t)和尺度函数
φ(t)中的支撑区为
除N =1外,dbN 不具有对称性(即非线性相位)。2N -1,ψ(t)的消失矩为N 。dbN 没有明确的表达式(除N
=1外),但转换函数h 的平方模是明确的。
Daubechies 小波系是由法国学者Daubechies 提出的一系列二进制小波的总称,在Matlab 中记为dbN ,N 为小波的序号,N 值取2,3,…,10。该小波没有明确的解析表达式,小波函数φ与尺度函数Φ的有效支撑长度为2N-1. 当N 取1时便成为Haar 小波。 令
p (y)=∑C
k =0
N -1
N -1+k
k
y
k
,其中C k
2
N -1+k
为二项式的系数,则有
m 0(ω)
=(c o s
2
ω
) p (s i ) 22
2
ω
1
式中,m (ω) =∑h e
2N -1k =0
k
-jk ω
。
Daubechies 小波具有以下特点:
(1)在时域是有限支撑的,即ψ(t)长度有限。 (2)在频域
ψ(ω) 在ω=0处有N 阶零点。
k
ψ(t)ψ(t-k)dt =δ。
(4)小波函数ψ(t)可以由所谓“尺度函数”φ(t)求出来。尺度函数φ(t)为
(3)ψ(t)和它的整数位移正交归一,即⎰低通函数,长度有限,支撑域在t=0~(2N-1)范围内。
Daubechies 小波的时域和频域波形图
db4 时域
1.5
1
0.5
-0.5
-10
2
46
8
t
i=10; wname = 'db4';
[phi,g1,xval] = wavefun(wname,i); subplot(1,2,1);
plot(xval,g1,'-r','LineWidth',1.5); xlabel('t') title('db4 时域'); g2=fft(g1); g3=abs(g2);
subplot(1,2,2);plot(g3,'-r','LineWidth',1.5); xlabel('f') title('db4 频域')
db4 频域
1000900
800700
600500
400300
200100
00
2000
40006000
8000
f
注意
Daubechies 小波常用来分解和重构信号,作为滤波器使用。波形如下:
分解低通滤波
器
0.80.60.4
0.20-0.2
-0.5-110.5
分解高通滤波
器
02468
重构低通滤波
器
10.5
10.50
-0.5
2
4
6
8
-10
重构高通滤波器
-0.5
2468
wname = 'db4'; % 计算该小波的4个滤波器
[Lo_D,Hi_D,Lo_R,Hi_R] = wfilters(wname); subplot(2,2,1); stem(Lo_D); title('分解低通滤波器'); subplot(2,2,2); stem(Hi_D); title('分解高通滤波器'); subplot(2,2,3); stem(Lo_R); title('重构低通滤波器'); subplot(2,2,4); stem(Hi_R); title('重构高通滤波器');
(3)Mexican Hat(mexh)小波
Mexican Hat函数为Gauss 函数的二阶导数:
ψ(t)=(1-t ) e
2
-
2
2
ψ(ω) =2ωe
2
-2
2
因为它的形状像墨西哥帽的截面,所以也称为墨西哥帽函数。
Mexican Hat(mexh)小波的时域和频域波形图
Mexihat 时域
1
15
mexihat 频域
0.5
10
05
-0.5-10
-5
0t
510
50f
100
d=-6; h=6; n=100;
[g1,x]=mexihat(d,h,n);
subplot(2,2,1);plot(x,g1,'-r','LineWidth',1.5); xlabel('t')
title('Mexihat 时域'); g2=fft(g1); g3=(abs(g2)); subplot(2,2,2);plot(g3); xlabel('f')
title('mexihat 频域');
特点:墨西哥帽函数在时间域与频率域都有很好的局部化,并且满足
⎰ψ(t)dt =0。由于它不存在尺度函数,所以小波函数不具有正交性。
R
它是高斯包络下的单频率副正弦函数:
ψ(t ) =C e c o s (5x)
2
-2
C 是重构时的归一化常数。Morlet 小波没有尺度函数
φ(t),而且是非正交分解。
Morlet 频域
Morlet 小波的时域波形图和频域波形图
Morlet 时域
10.5
0-0.5-1
-10
d=-6; h=6; n=100;
[g1,x]=morlet(d,h,n);
subplot(2,2,1);plot(x,g1,'-r','LineWidth',1.5); xlabel('t')
title('Morlet 时域'); g2=fft(g1); g3=abs(g2);
subplot(2,2,2);plot(g3); xlabel('f')
title('Morlet 频域')
15
10
5
-5
0t
510
50f
100
Meyer 小波的小波函数和尺度函数都是在频率域中进行定义的,其定义为:
π32π4π⎧
⎪(2π) e sin (2v(2πω-1)) 3≤ω≤3⎪π34π8π⎪
ψ(ω) =⎨(2π) e cos (v(-1)) ≤ω≤
22π33⎪
⎡2π8π⎤⎪0ω∉⎢, ⎥
⎪⎣33⎦⎩
-1
2
jw 2
-12
jw 2
其中,v(a)为构造Meyer 小波的辅助函数,具有
v (a)=a (35-84a +70a -20a ) a ∈[0,1]
4
2
3
⎧(2π) ⎪π32π4π⎪
φ(ω) =⎨(2π) cos (v(-1)) ≤≤
22π33⎪
04π⎪>⎩
3
--12
12
ω≤
2π3
Meyer 小波不是紧支撑的,但它收敛的速度很快:
(t)≤C (1+t )
n
2
-n
ψ(t)无限可微。
Meyer 小波的时域和频域波形图
Meyer 时域
2
1
-1t
meyer 频域
1510
5
f
d = -6; h = 6; n = 128;
[psi,x] = meyer(d,h,n,'psi');
subplot(2,1,1), plot(x,psi,'-r','LineWidth',1.5) xlabel('t')
title('Meyer 时域'); PSI=fft(psi); PSII=abs(PSI);
subplot(2,1,2), plot(PSII); xlabel('f') title('meyer 频域')
2、在信号 x(t )=sin(2π*30t)+cos(2π*50t)加上噪音后分别进行FFT 和CWT 变换。
解:引入随机噪声randn(1,N)
原信号x(t)波形图
2
50
x(t)的fft 变换图
x (t )
-2
50t
x(t)加噪声后fft 变换图
1000
50f morlet
100
50
64
21
20
406080time (or space) b
尺度为2
100
50f
尺度为1
100
10
-1
s c a l e s a
50
100
20
-2
050100
N=100;fs=1000; n=0:N-1; t=n/fs;
x=sin(60*pi*t)+cos(100*pi*t); %原信号 subplot(3,2,1);
plot(x,'-r','LineWidth',1.5); xlabel('t') ylabel('x(t)')
title('原信号x(t)波形图') F1=fft(x); m1=abs(F1);
subplot(3,2,2);plot(m1);
xlabel('f')
title('x(t)的fft 变换图')
x1=randn(1,N); %加入噪声
x2=x+x1;
F2=fft(x2);
m2=abs(F2);
subplot(3,2,3);plot(m2);
xlabel('f')
title('x(t)加噪声后fft 变换图')
scale=[1 2 4 6]; %设置尺度
subplot(3,2,4);x3=cwt(x2,scale,'morl','plot');
title('morlet'); %加噪声后CWT 变换结果图 subplot(3,2,5);plot(x3(1,:));title('尺度为1');
subplot(3,2,6);plot(x3(2,:));title('尺度为2');
1. 给出五种常用小波基的时域和频域波形图。
与标准的傅里叶变换相比,小波分析中使用到的小波函数具有不唯一性,即小波函数ψ(t) 具有多样性。小波分析在工程应用中,一个十分重要的问题就是最优小波基的选择问题,因为用不同的小波基分析同一个问题会产生不同的结果。目前我们主要是通过用小波分析方法处理信号的结果与理论结果的误差来判定小波基的好坏,由此决定小波基。常用小波基有Haar 小波、Daubechies(dbN)小波、Mexican Hat(mexh)小波、Morlet 小波、Meyer 小波等5种。 (1)Haar 小波
Haar 函数是小波分析中最早用到的一个具有紧支撑的正交小波函数,也是最简答的一个小波函数,它是支撑域在t
∈[0,1]范围内的单个矩形波。 Haar
⎧10≤t ≤1
⎪12
函数的定义如下:ψ(t)=⎨-1≤t ≤1
2⎪⎩0其他
Haar 小波在时域上是不连续的,所以作为基本小波性能不是特别好。但它也有自己的优点,如: ①计算简单;
②ψ(t)不但与ψ(2j t) [j∈z]正交,而且与自己的整数位移正交。 因此,在a 族。
ψ(t ) 的傅里叶变换是:
=2j 的多分辨率系统中Haar 小波构成一组最简单的正交归一的小波
ψ(Ω)=j
4Ω
sin 2() e -j Ω/2 Ωa
Haar 小波的时域和频域波形图
haar 时域
1.5
1
0.5
-0.5
-1
-1.5
00.5
11.5
t
i=20; wav = 'haar';
[phi,g1,xval] = wavefun(wav,i); subplot(1,2,1);
plot(xval,g1,'-r','LineWidth',1.5); xlabel('t') title('haar 时域'); g2=fft(g1); g3=abs(g2);
subplot(1,2,2);plot(g3); xlabel('f') title('haar 频域')
5
haar 频域
7
6
5
4
3
2
1
05
1015f
x 10
5
(2)Daubechies(dbN)小波
Daubechies 小波是世界著名的小波分析学者Inrid Daubechies 构造的小波函数,简写为dbN ,N 是小波的阶数。小波ψ(t)和尺度函数
φ(t)中的支撑区为
除N =1外,dbN 不具有对称性(即非线性相位)。2N -1,ψ(t)的消失矩为N 。dbN 没有明确的表达式(除N
=1外),但转换函数h 的平方模是明确的。
Daubechies 小波系是由法国学者Daubechies 提出的一系列二进制小波的总称,在Matlab 中记为dbN ,N 为小波的序号,N 值取2,3,…,10。该小波没有明确的解析表达式,小波函数φ与尺度函数Φ的有效支撑长度为2N-1. 当N 取1时便成为Haar 小波。 令
p (y)=∑C
k =0
N -1
N -1+k
k
y
k
,其中C k
2
N -1+k
为二项式的系数,则有
m 0(ω)
=(c o s
2
ω
) p (s i ) 22
2
ω
1
式中,m (ω) =∑h e
2N -1k =0
k
-jk ω
。
Daubechies 小波具有以下特点:
(1)在时域是有限支撑的,即ψ(t)长度有限。 (2)在频域
ψ(ω) 在ω=0处有N 阶零点。
k
ψ(t)ψ(t-k)dt =δ。
(4)小波函数ψ(t)可以由所谓“尺度函数”φ(t)求出来。尺度函数φ(t)为
(3)ψ(t)和它的整数位移正交归一,即⎰低通函数,长度有限,支撑域在t=0~(2N-1)范围内。
Daubechies 小波的时域和频域波形图
db4 时域
1.5
1
0.5
-0.5
-10
2
46
8
t
i=10; wname = 'db4';
[phi,g1,xval] = wavefun(wname,i); subplot(1,2,1);
plot(xval,g1,'-r','LineWidth',1.5); xlabel('t') title('db4 时域'); g2=fft(g1); g3=abs(g2);
subplot(1,2,2);plot(g3,'-r','LineWidth',1.5); xlabel('f') title('db4 频域')
db4 频域
1000900
800700
600500
400300
200100
00
2000
40006000
8000
f
注意
Daubechies 小波常用来分解和重构信号,作为滤波器使用。波形如下:
分解低通滤波
器
0.80.60.4
0.20-0.2
-0.5-110.5
分解高通滤波
器
02468
重构低通滤波
器
10.5
10.50
-0.5
2
4
6
8
-10
重构高通滤波器
-0.5
2468
wname = 'db4'; % 计算该小波的4个滤波器
[Lo_D,Hi_D,Lo_R,Hi_R] = wfilters(wname); subplot(2,2,1); stem(Lo_D); title('分解低通滤波器'); subplot(2,2,2); stem(Hi_D); title('分解高通滤波器'); subplot(2,2,3); stem(Lo_R); title('重构低通滤波器'); subplot(2,2,4); stem(Hi_R); title('重构高通滤波器');
(3)Mexican Hat(mexh)小波
Mexican Hat函数为Gauss 函数的二阶导数:
ψ(t)=(1-t ) e
2
-
2
2
ψ(ω) =2ωe
2
-2
2
因为它的形状像墨西哥帽的截面,所以也称为墨西哥帽函数。
Mexican Hat(mexh)小波的时域和频域波形图
Mexihat 时域
1
15
mexihat 频域
0.5
10
05
-0.5-10
-5
0t
510
50f
100
d=-6; h=6; n=100;
[g1,x]=mexihat(d,h,n);
subplot(2,2,1);plot(x,g1,'-r','LineWidth',1.5); xlabel('t')
title('Mexihat 时域'); g2=fft(g1); g3=(abs(g2)); subplot(2,2,2);plot(g3); xlabel('f')
title('mexihat 频域');
特点:墨西哥帽函数在时间域与频率域都有很好的局部化,并且满足
⎰ψ(t)dt =0。由于它不存在尺度函数,所以小波函数不具有正交性。
R
它是高斯包络下的单频率副正弦函数:
ψ(t ) =C e c o s (5x)
2
-2
C 是重构时的归一化常数。Morlet 小波没有尺度函数
φ(t),而且是非正交分解。
Morlet 频域
Morlet 小波的时域波形图和频域波形图
Morlet 时域
10.5
0-0.5-1
-10
d=-6; h=6; n=100;
[g1,x]=morlet(d,h,n);
subplot(2,2,1);plot(x,g1,'-r','LineWidth',1.5); xlabel('t')
title('Morlet 时域'); g2=fft(g1); g3=abs(g2);
subplot(2,2,2);plot(g3); xlabel('f')
title('Morlet 频域')
15
10
5
-5
0t
510
50f
100
Meyer 小波的小波函数和尺度函数都是在频率域中进行定义的,其定义为:
π32π4π⎧
⎪(2π) e sin (2v(2πω-1)) 3≤ω≤3⎪π34π8π⎪
ψ(ω) =⎨(2π) e cos (v(-1)) ≤ω≤
22π33⎪
⎡2π8π⎤⎪0ω∉⎢, ⎥
⎪⎣33⎦⎩
-1
2
jw 2
-12
jw 2
其中,v(a)为构造Meyer 小波的辅助函数,具有
v (a)=a (35-84a +70a -20a ) a ∈[0,1]
4
2
3
⎧(2π) ⎪π32π4π⎪
φ(ω) =⎨(2π) cos (v(-1)) ≤≤
22π33⎪
04π⎪>⎩
3
--12
12
ω≤
2π3
Meyer 小波不是紧支撑的,但它收敛的速度很快:
(t)≤C (1+t )
n
2
-n
ψ(t)无限可微。
Meyer 小波的时域和频域波形图
Meyer 时域
2
1
-1t
meyer 频域
1510
5
f
d = -6; h = 6; n = 128;
[psi,x] = meyer(d,h,n,'psi');
subplot(2,1,1), plot(x,psi,'-r','LineWidth',1.5) xlabel('t')
title('Meyer 时域'); PSI=fft(psi); PSII=abs(PSI);
subplot(2,1,2), plot(PSII); xlabel('f') title('meyer 频域')
2、在信号 x(t )=sin(2π*30t)+cos(2π*50t)加上噪音后分别进行FFT 和CWT 变换。
解:引入随机噪声randn(1,N)
原信号x(t)波形图
2
50
x(t)的fft 变换图
x (t )
-2
50t
x(t)加噪声后fft 变换图
1000
50f morlet
100
50
64
21
20
406080time (or space) b
尺度为2
100
50f
尺度为1
100
10
-1
s c a l e s a
50
100
20
-2
050100
N=100;fs=1000; n=0:N-1; t=n/fs;
x=sin(60*pi*t)+cos(100*pi*t); %原信号 subplot(3,2,1);
plot(x,'-r','LineWidth',1.5); xlabel('t') ylabel('x(t)')
title('原信号x(t)波形图') F1=fft(x); m1=abs(F1);
subplot(3,2,2);plot(m1);
xlabel('f')
title('x(t)的fft 变换图')
x1=randn(1,N); %加入噪声
x2=x+x1;
F2=fft(x2);
m2=abs(F2);
subplot(3,2,3);plot(m2);
xlabel('f')
title('x(t)加噪声后fft 变换图')
scale=[1 2 4 6]; %设置尺度
subplot(3,2,4);x3=cwt(x2,scale,'morl','plot');
title('morlet'); %加噪声后CWT 变换结果图 subplot(3,2,5);plot(x3(1,:));title('尺度为1');
subplot(3,2,6);plot(x3(2,:));title('尺度为2');