资源描述
第二章_离散信号频谱的窗谱校正方法.txt41滴水能穿石,只因为它永远打击同一点。42火柴如果躲避燃烧的痛苦,它的一生都将黯淡无光。华中理工大学博士学位论文
第二章离散信号频谱的窗谱校正方法
----基本理论
§2.1 引言
利用DFT可以对离散信号进行频谱分析,但是计算工作量相当大,因此,在快速算
法没用发明之前,DFT并没用多大的实际意义。直到1965年,Cooley-Tukey在《计算数
学》杂志上首先提出FFT算法之后,
DFT才得到广泛的应用。这一快速算法的出现对数字
信号分析领域的发展起到了极大的推动作用。从此以后,它作为频谱分析的基础得到了
广泛的应用
[75,76,77,78]
。
由于计算机只能对信号的有限多个样本进行计算,信号的
FFT谱分析也只能在时域
信号的有限区间内进行,这就不可避免地存在由于时域截断(加矩形窗)而产生泄漏[61],
使谱峰值减小,精度降低,求得的信号相位更是面目全非。在数字信号处理中,由
DFT
或
FFT得到的幅值谱是离散谱,是信号与窗函数频谱卷积后,按频率分辨率
Δfs
=
fs
/ N
( fs为信号采样频率,N为分析信号样本长度)等间隔频域抽样的结果(如图
21
所示)[78]。
A幅值
f
图2-1 频谱抽样的离散谱线
如果周期信号的频率正好
表2-1 离散频谱幅值、相位和频率误差表
落在某一谱线上,经FFT后得
矩形窗 Hanning窗
Hamming
窗
幅值误差(%)
0--36.4 0--15.3 0--18.3
相位误差(0)
±90 ±90 ±90
频率误差(Hz)
±05Δfs. ±05Δfs. ±05Δfs.
到的频率、幅值和相位是准确
的。在一般情况下,信号频率
落于两条相邻谱线之间,由于
谱线不在主瓣中心,由峰值谱
线反映的频率和幅值都不准
10
华中理工大学博士学位论文
确,相位误差更大。从理论上分析,加矩形窗时,最大误差可达
36 4
.%,即使加其他窗
时,也不能完全消除这一影响,在加Hanning窗时,只进行幅值恢复时的最大幅值误差仍
高达15 3
.%,相位误差将更大,表2-1是离散频谱只进行幅值恢复,不进行其他处理时幅
值、相位和频率误差
[141]
。
§2.2 单频谐波的频谱分析误差产生原因
无限长信号
x()t(如振动信号、噪声信号等)的频谱分析所采用的方法为对信号进
行截取,然后再对截取得到的有限长度信号进行频谱分析。窗函数
wt()的作用就相当于
对无限长的信号开一窗口,从窗口中取出一段数据,从而完成信号的截取。窗函数都是
选择实偶函数,并在时域上将窗函数的中心放于被分析的那段信号的中心。
加窗信号的傅氏变换为:
Fx.
()] ()() 2π
dt.............................................................(2.1)
[ () twtxtwteTTjft=....∫(∞)
.∞
其中,wt
T()由对称窗wt()在时间上平移T/2得到,即
T() t.T/)
wt=w( 2 ................................................................................................(2.2)
设wt()的傅氏变换(如图2-2a)
Fwt
=W( f
[( )] ) ......................................................................................................(2.3)
根据傅氏变换的奇偶性质,当wt()是实偶函数时,Wf
()此时也为实偶函数。又由
傅氏变换的时移特性可知(如图2-2b),
Fw
t
=Wfe
() jfT
π
[ T( )] ...........................................................................................(2.4)
设有一周期信号
x()t=ACos(2π.f0 .t+.),则其傅氏变换结果为(如图2-2c):
A .j.
Aj.
() 0 .f0 ) ...............................................................(2.5)
Xf =.e .δ(f +f) +.e .δ(f
22
根据卷积定理,加窗后的谐波信号
x()twt.T()的傅氏变换可表示为(图2-2d):
Xf
() =Fxt
w
t
T( )] =Fxt
[()] Fw
t
( )]
T
[() ..[ T
(这里“*”表示卷积)
AA
.jj
.
..
π
..
j
fT
={ .e
.δ( +
) +.e
.δ( ..Wfe
ff
ff
)}{ () }
00
22
A....
+
)+.] Aj[πTf
f
0).
j[πTf
f
....
.
(0 (]
=.Wf
fe
+
) +.
( .
)
( Wf
fe
.
...............(2.6)
00
22
由此可知,在加窗信号的傅氏分析中,当
f≠f0时,将存在泄漏情况。此时的幅值
及相位分别为:
A
Y=.Wf
f
( .
0 ) ....................................................................................................(2.7)
2
Φ=..
.
Tf( .f
+.......................................................................................(2.8)
π
0)
11
华中理工大学博士学位论文
对窗长度T
=
N
/fs
作归一化处理,则T
=
1,且令Δf
=
f
.
f0代入上面两式可得:
A
(.
)
Φ=.π.Δf
+.
......................................................................................................(2.10)
Y
=
2
.Wf
...........................................................................................................(2.9)
W
( f
)
w
() T
t
t
f
-T/2
T/2
1.1/Τ0.2/Τ2/Τ1/Τ
时域波形窗谱模函数
(a)对称矩形窗函数的时域波形和频谱模函数
wt
T
()
1
t
0T
T
.2/Τ.1/Τ0 2/Τ1/ΤW()fTf
时域波形窗谱模函数
(b)实际窗函数的时域波形和频谱模函数
T0
Ax(t)
t
X
( f
)
A/2
A/2
f
-f
0 f
00
时域波形傅里叶变换模函数
(c)单频率谐波的时域波形和频谱模函数
A
xtT( )
0 T t
Y
n
yKyK.1
k k+2k+1k-1k-2 f0
0
时域波形离散频谱模函数
(d)单频率谐波离散频谱模函数
图2-2 单频率谐波离散频谱的误差产生原因
12
华中理工大学博士学位论文
显然,当
f=f0时,Y=
A
,Φ=.不存在泄漏情况,得到的幅值、相位和频率都是
2
准确无误差的。在大多数情况下,当
f≠f0时。由加窗信号的傅氏分析得到的频率
f、
幅值Y和相位Φ并不是真实值,且有旁瓣产生,这就是所谓的离散频谱的栅栏效应、梳
状效应、能量泄漏和假频等(如图
2-2d所示)。当信号真实频率位于两个相邻离散谱线中
间时,即
fK.1 =f0 .Δfs/2 ,
fK
=f0 +Δfs/2 (这里
Δfs为频率分辩率)时,求得的信号幅值、相位和频率的误差最大。
§2.3 离散频谱信号的窗谱校正方法
假设加窗信号的频谱主瓣中心为
f(即为信号的真实频率),信号幅值为
A0;加窗
信号FFT结果的频谱中,最高的频谱频率(0) 为 f,高度为Y;次高谱线频率为 f2,高度为
Y2。显然,
f和
f相差仅为一个频率分辨率(1) ,对此归(1) 一化后,即有 f1 =f2 ±1。当
f1 >f2时取“+(1) ”号(2) ;当
f1 <f2时取“—”号。加窗信号的频谱由窗频谱对信号真实频
谱调制而得到,因而可利用窗函数的频谱图形对傅氏分析结果进行校正,以求出真实的
频率、幅值和相位
[65][67]
。
(1) 频率校正
0
WWxyf+1
f
f+1fΔΔΔΔ0
yyxyf+1
f
f+1ff0
A
图2-3 窗函数的频谱函数图2-4 频率和幅值的校正
频率校正即求出信号幅值谱图上窗谱主瓣中心所处的横坐标
f0。设窗函数的频谱函
数为
Wfff,其窗谱函数为
(),
W()对称于
Y轴(如图
2-3所示)。对于任一
ΔWf(Δ),相应的信号离散频谱幅值为
y;对于处于Y轴右侧的
Δf+1谱线,其窗谱函
数为Wf1,相应的信号离散频谱为
yf+1(1、
yf和
yf+1
(Δ+)(Δ+)
幅值(f) ;由WfΔ)、Wf
可求出
Δf的值,即求出了频率修正量
Δffff
=.0。由于W()的函数表达式为已知,
可构成一函数:
WfΔ) yf
(
mFf
=
(Δ) ==
..........................................................................(2.11)
Wf+1) y
(.
f+1
13
华中理工大学博士学位论文
m为窗谱主瓣中心两旁相隔为1的两谱线幅值比值,为
f的函数,求上式的反函数:
ΔfG=
m
() ...............................................................................................................(2.12)
=
这样便可得到频率修正量
Δf
ff.
0。
在实际计算中,主瓣中心位于信号真实频率
f0处。在图2-4中
yf和
yf+1是幅值谱主
瓣内谱峰左侧和右侧的谱线,将
m=
yf
代入式(2.12)可求出谱线修正量
Δf。频率校正
y
f+1
结果为∶
f0 =( f..
) s/
ff
N.............................................................................................(2.13)
N为分析长度,
fs为采样频率。可将
ΔfG()m称为频率校正函数,对于不同的窗
=
函数,其Gv
()是不同的。如果在同样的采样频率和同样的信号分析样本长度的情况下,
对加窗信号的频域抽样位置也就确定下来了,由于信号频率
f0的位置固定,即信号频域
峰值位置固定,则可得不管信号加了何种窗,所产生的频率误差都是一样的。
(2) 幅值校正
设窗函数的频谱模函数为Wf
(),则图2-4中主瓣函数为∶
yAW( .
0)
=.
xf
...................................................................................................(2.14)
这就是信号频谱与窗函数卷积的结果,式中,
A为真实幅值,对应主瓣中心
f0,现
将
yyf, xf代入式得∶
==
y
AW( ff
.
f=.
0 ) ................................................................................................(2.15)
.=
f
y
式中,
ff
0 Δ,故可解出幅值A的值。
A=
f
.............................................................................................................(2.16)
Wf
(.
)
(3)相位校正
频谱分析所用窗函数都不是对称于Y轴,要向右平移
N/2点,其频谱函数相对于Y
轴来说有一相移因子
e
.jNω/2,其相位角∶
φ=.Nω/ 2 ............................................................................................................(2.17)
ω与谱线
f的关系为∶
ω=2π
f/ N...........................................................................................................(2.18)
将上式代入式(2.17)可得∶
φ=.π
f.................................................................................................................(2.19)
这表明窗函数的相位是线性相位(图2-5所示)。
信号频谱函数与窗函数的频谱作复卷积时为复数相乘,此时相位相加。由图2-5可看
出,频率误差为半个谱线时,相位误差为最大可达900,显然经
FFT得到的信号实部与虚
14
华中理工大学博士学位论文
部所得到的相位如不加校正是不能用的。由式
(2.12)得到的频率修正量后,相位修正量
为∶
Δφ
=.Δf
π
.............................................................................................................(2.20)
当加窗信号经FFT后,得到的相位.'时,则其真实相位角为∶
φ=+.' Δφ
................................................................................................................(2.21)
从这可看出,如果所加窗为对称窗,不论窗函
φ
数是简单的还是复杂的,对加窗信号进行傅氏分析
900 时,泄漏所产生的相位误差与频率误差都成线性关
系。因此,对不同的对称窗,加窗信号的幅值误差
.0.50.5
f
可能不同,但由频率校正可知,所有加窗信号的相
位误差是相同的。
-900
从上面的分析可看出上述的校正方法的实质是
图2-5 窗函数的相位利用归一化后窗函数幅值谱的主瓣中心两旁差值为
一个频率分辨率
Δfs的两条谱线的高度建立一个以
修正频率为变量的方程,解出修正频率,进而对其相应的信号幅值和相位进行修正。这
一修正方程建立在所加窗函数的频谱公式(即窗谱图形)之上,因此称这一修正方法为
窗谱校正方法。
下面以矩形窗和汉宁(Hanning)窗的较正为例进行窗谱校正方法的具体分析:
1、矩形窗
在频谱分析中,矩形窗的定义
[78]
为:
Wn
=
1,n
=
012, 1 .................................................................................(2.22)
() ,, ...
, N
.
其频谱函数为:
Nω
Sin
N
.1
.
j
ω
2
ω
.......................................................................................(2.23)
W()=
e
2
ω
Sin
2
则频谱函数的模为式(2.24),其波形图如图2-6(a)所示
Nω
Sin
2
ω
...................................................................................................(2.24)
W0()=
ω
Sin
2
相应的相位为:
N
.
1
.=.
ω
.........................................................................................................(2.25)
2
15
华中理工大学博士学位论文
从矩形窗频谱模函数波形图可看出主瓣宽度为
4π
,为归一化后的频率分辨率
Δfs的
N
两倍。由此可知在主瓣内谱线的位置及条数有如下两种情况:
(1)只有一条谱线。此时谱线正好落在主瓣中心上,显然没有泄漏发生。
(2)主瓣内有两条谱线,且落在主瓣中心的两边,这时有泄漏发生。对于这一情况作
如下推导:
N2πN.2πN4πN.4πN0
ωW0 ()ω012-1-21yyxy1
2
ΔxΔx-1
(a) 矩形窗频谱模图形 (b)归一化矩形窗频谱模图形
图
2-6 矩形窗频谱模图形与其校正方法
由于信号加窗截取后的DFT结果为信号频谱和窗函数频谱在频域中的卷积。相当于
将窗函数频谱主瓣中心平移至信号的各频率点后,在频率分辨率整数倍的横坐标点上抽
样而得到DFT结果
[78]
。利用这可进一步简化以下的推导:
由于谱线分辨率为2Nπ
,令ω=
2
N
π
x,代入式(2.24)即可得:
Sinπ.
x
Wx
.....................................................................................................(2.26)
() =
π.
x
Sin
N
由式(2.24)可得图2-6(b),因此可设谱线处于
Δx和
Δx
.1二点,且已知相应的谱线高
度为
y1, y2则有:
Sin(π
.Δx)
y
=
A
..................................................................................................(2.27)
1 π
.Δx
Sin()
N
Sin[(π.
Δx
.
1)]
y
=
A
.........................................................................................(2.28)
2 π.
(Δx
.
1)
Sin[]
N
由式(2.27)(2.28)可求得
π
y2 .
Sin
N
.1 N
Δx
=
tg
.....................................................................................(2.29)
ππ
y
+
y
Cos
12 N
当N>>1时,简化式(2.29)可得下面的近似式子:
16
华中理工大学博士学位论文
y
Δx
=
2 ............................................................................................................(2.30)
y1 +y2
如果信号中有一频率分量,则可看为将窗谱平移至如图
2-7所示的位置,此时的Δx和
Δx.1点变为K,K-1点,则信号的实际频率为:
y
x
=K
.Δx
=K
.
2 .....................................................................................(2.31) 0 y1 +y2
0 k-1yyxy1
2
ΔxΔx-1)(
K( )
f0
图2-7 信号将窗谱平移图形
将式(2.30)代入式(2.27)可得这一频率分量的幅值为:
y
..Δx
π
A0 =
1 ...............................................................................................(2.32)
NSin
(π.Δx
.
)
相应相位的校正结果为:
N
.12π
N
.12π
N
.1
.=.
..x
=.
..(K
.Δx) =.
+π.
.Δx
..........(2.33)
00 K
2 N
2 NN
因此对于矩形窗可用式(2.30)(2.31)(2.32)(2.33)分别对频率、相应的幅值及相应的相位
进行校正。
2、Hanning窗
Hanning窗的定义
[78]
为:
12π
()1 201 ...
....................................(2.34)
Wn
=
[ +Cos( .nn)], =.
N
/,.
,, , N
/2
2 N
其频谱函数为:
1 12π
12π
W() =
D() +
D(ω.
) +
D
ωω
(ω+
) ...............................................(2.35)
24 N
4 N
Nω
Sin
ω
式中
D()=
ω
2 ej
2,为了简化推算,在
N>>1时,我们忽略三项之间的相位差
ω
Sin
2
2Nπ,则W()
ω的模可近似为三项模函数之和。使用类似于上面矩形窗的讨论方法可得:
17
华中理工大学博士学位论文
Sinπ.x
1 Sin[(x
.)] 1 Sin[(x
+)] Sinπ.x
1 π
1 π
11
Wx
=.
+.=.
2
() +.
.x
4 π(x
.1) 4 π(x
+1) π.x
( .x
2 π
21 )
...................................................................................................................................(2.36)
其函数图形如图2-8所示:
yxyy12ΔxΔx-1
图2-8 汉宁窗频谱模图形与其校正方法
由图2-8可知通常在频谱的主瓣内有四条谱线,且主瓣中心处于两个连续的最高谱线
之间。
设
Δx和
Δx
.1处谱线所对应的幅值高度为
y1和
y2,则有如下方程组:
.
Sinπ.Δx
1
y1 =.
2
.
..
π
21
.Δx
( .Δx
)
..................................................................(2.37)
Sinπ(Δx
)1
..1
.y2 =.
2
.
π..(Δx
1) 21[ .(Δx
.1) ]
.
由此方程组可求得:
2 y2 .y1
Δx
=
............................................................................................................(2.38)
y1 +y2
与矩形窗类似,按照式(2.30)(2.31)(2.32)(2.33)的推导过程可得汉宁窗校正公式为:
2 y2 .y1
=.
x0 K
....................................................................................................(2.39)
y
+y2
Δxy(1)
π..12
A0 =+21.Δx
( ) .....................................................................................(2.40)
Sinπ.Δx
N
.12π
N
.1
.=.
.
(K
.Δx) =.
π
0 K
+.
.Δx
..................................................(2.41)
2 NN
以上是矩形窗和汉宁窗的校正公式推导,对于其它窗函数因其公式过于复杂,难以
找出类似的校正公式。从这一角度来看,基于公式推导进行频谱校正有其局限性。
§2.4 仿真计算
18
华中理工大学博士学位论文
在上一节中已推导出矩形窗和汉宁窗的校正公式,下面的仿真计算就利用这些校正
公式来进行信号频谱的频率、幅值和相位修正。
设一模拟信号在无噪声情况下表达为:
10π
15π
18π
y
=
05 .
Cos
(2 .
25 +
) +
30 .
Cos
(2 .
50 +
) +
0 75 .
Cos
(
. π
t
. π
t
.2π.100t
+
)
180 180 180
30π
60π
+05 .Cos(2π.150t
+
) +
0 15 .Cos(2π.200t
+
) ...................................(2.42)
..
180 180
信号的采样频率固定为1000HZ,分别对信号加矩形窗和Hanning窗进行2048点FFT,
采用各窗的相应的校正公式进行计算,所得结果分别列于表2-2和表2-3中。从结果可看
出,信号频谱经修正后可达较高的精度。
表2-2 信号加矩形窗的频域参数校正前后对照表
谱峰频率谱峰幅值谱峰相位
矩形窗
不加校正
矩形窗
公式校正
真实
值
矩形窗
不加校正
矩形窗
公式校正
真实
值
矩形窗
不加校正
矩形窗
公式校正
真实
值
24.902344 24.990749 25.0 0.478713 0.505535 0.50 46.604000 14.014318 10.0
49.804688 50.000698 50.0 2.265033 2.999154 3.0 86.963432 14.706041 15.0
100.097656 100.191962 100.0 0.707923 0.753301 0.75 -18.818087 15.946660 18.0
149.902344 150.004197 150.0 0.462865 0.497733 0.50 65.201889 27.654678 30.0
200.195312 200.011447 200.0 0.118115 0.150937 0.30 -14.291993 53.487831 60.0
表2-3 信号加Hanning窗的频域参数校正前后对照表
谱峰频率谱峰幅值谱峰相位
Hanning
不加校正
Hanning
公式校正
真实
值
Hanning
不加校正
Hanning
公式校正
真实
值
Hanning
不加校正
Hanning
公式校正
真实
值
24.902344 25.000009 25.0 0.487229 0.499997 0.50 45.999359 9.996161 10.0
49.804688 49.999999 50.0 2.702956 3.000000 3.0 86.999680 15.000054 15.0
100.097656 99.999998 100.0 0.730850 0.750000 0.75 -18.000568 18.000162 18.0
149.902344 149.999997 150.0 0.487235 0.500000 0.50 65.999077 30.000103 30.0
200.195312 199.999997 200.0 0.135147 0.150000 0.30 -12.001210 60.000042 60.0
§2.6 本章小结
FFT作为频谱分析的基础已得到了广泛的应用。如何提高谱分析精度现已成为进一步
推广应用FFT的关键问题。本章在分析了离散频谱型信号分析时误差产生的机理的基础
19
华中理工大学博士学位论文
上,研究了离散频谱信号的实用校正方法
-----窗谱校正方法,这一方法的实质是利用归一
化后窗谱主瓣峰值两旁频率差值为一个频谱分辨率的两条谱线的幅值建立一个以修正频
率为变量的方程,解出修正频率,进而对其相应的幅值和相位进行修正。由加窗信号的
频谱分析的校正公式推导得到如下结论∶
(1) 对于加窗离散频谱信号的频谱分析,可以利用所加窗函数的频谱模图形对其频谱
分析结果进行校正,以克服信号频谱分析时的栅栏效应和能量泄漏误差。
(2) 在同样的采样频率和同样的信号分析样本长度的情况下,不管信号加了何种对称
窗,所产生的频率误差都是一样的。
(3) 不论窗函数是简单的还是复杂的,对加窗信号进行傅氏分析时,泄漏所产生的相
位误差与频率误差都成线性关系。
20
展开阅读全文