窗函数的实现与分析.docx
《窗函数的实现与分析.docx》由会员分享,可在线阅读,更多相关《窗函数的实现与分析.docx(22页珍藏版)》请在冰豆网上搜索。
窗函数的实现与分析
摘要I
3.1MATLAB软件简介9
摘要
现代图像、语声、数据通信对线性相位的要普遍的。
正是此原因,使得具有线性相位的FIR数字滤波器得到大力发展和广泛应用。
在实际进行数字信号处理时,往往需要把信号的观察时间限制在一定的时间间隔,只需要选择一段时间信号对其进行分析。
这样,取用有限个数据,即将信号数据截断的过程,就等于将信号进行加窗函数操作。
而这样操作以后,常常会发生频谱分量从其正常频谱扩展开来的现象,即所谓的“频谱泄漏”。
当进行离散傅立叶变换时,时域中的截断是必需的,因此泄漏效应也是离散傅立叶变换所固有的,必须进行抑制。
而要对频谱泄漏进行抑制,可以通过窗函数加权抑制DFT的等效滤波器的振幅特性的副瓣,或用窗函数加权使有限长度的输入信号周期延拓后在边界上尽量减少不连续程度的方法实现。
而在后面的FIR滤波器的设计中,为获得有限长单位取样响应,需要用窗函数截断无限长单位取样响应序列。
另外,在功率谱估计中也要遇到窗函数加权问题。
由此可见,窗函数加权技术在数字信号处理中的重要地位。
1.窗函数
1.1基本概念
在实际进行数字信号处理时,往往需要把信号的观察时间限制在一定的时间间隔,只需要选择一段时间信号对其进行分析。
这样,取用有限个数据,即将信号数据截断的过程,就等于将信号进行加窗函数操作。
而这样操作以后,常常会发生频谱分量从其正常频谱扩展开来的现象,即所谓的“频谱泄漏”。
当进行离散傅立叶变换时,时域中的截断是必需的,因此泄漏效应也是离散傅立叶变换所固有的,必须进行抑制。
而要对频谱泄漏进行抑制,可以通过窗函数加权抑制DFT的等效滤波器的振幅特性的副瓣,或用窗函数加权使有限长度的输入信号周期延拓后在边界上尽量减少不连续程度的方法实现。
而在后面的FIR滤波器的设计中,为获得有限长单位取样响应,需要用窗函数截断无限长单位取样响应序列。
另外,在功率谱估计中也要遇到窗函数加权问题。
窗函数的基本概念。
设x(n)是一个长序列,w(n)是长度为N的窗函数,用w(n)截断x(n),得到N点序列xn(n),即
xn(n)=x(n)w(n)
在频域上则有
由此可见,窗函数w(n)不仅仅会影响原信号x(n)在时域上的波形,而且也会影响到频域的形状。
1.2设计原理
窗函数设计法的基本原理是用有限长单位脉冲响应序列
逼近
。
由于
往往是无限长序列,而且是非因果的,所以用窗函数
将
截断,并进行加权处理,得到:
就作为实际设计的FIR数字滤波器的单位脉冲响应序列,其频率响应函数
为
式中,N为所选窗函数
的长度。
用窗函数法设计的滤波器性能取决于窗函数
的类型及窗口长度N的取值。
设计过程中,要根据对阻带最小衰减和过渡带宽度的要求选择合适的窗函数类型和窗口长度N。
2.窗函数的种类
对时间序列作傅里叶变换时,实际上要作周期延拓,如果取长序列的一段进行处理,还要先作截断。
截断会引起频谱的泄漏问题,这是由于抽样后的离散序列与矩形截断序列相乘,在频域造成两者的频谱卷积形成的。
为了避免频谱泄露对结果的影响,在对非周期信号作时间截断时,除尽量增加截断序列的宽度外,也应选其频谱的旁瓣较小的截断窗函数,以减小泄漏的影响。
在信号处理中窗函数是一种除在给定区间之外取值均为0的实函数。
常用的窗函数很多,例如矩形窗、三角窗、Hanning窗、Hamming窗、Parzen窗、Kaiser窗、Chebyshev窗、Tukey窗、Poisson窗、Caushy窗、Gaussian窗和Blackman窗等等,定义方式不同,性质也不同,可以根据实际需要来选择合适的截断窗函数减轻泄露问题。
窗函数在光谱分析、滤波器设计以及音频数据压缩等方面有广泛的应用。
在这里,应用窗函数对输入数据进行截断(时域加窗),从而得到要处理的数据及其长度。
表2.1MATLAB窗函数
窗
窗函数
矩形窗
Boxcar
巴特利特窗
Barlett
三角窗
Triang
布莱克曼窗
Blackman
海明窗
Hamming
汉宁窗
Hanning
凯塞窗
Kaiser
切比雪夫窗
Chebwin
表2.2各种窗函数的基本参数
窗函数
旁瓣峰值幅度/dB
过渡带宽
阻带最小衰减/dB
矩形窗
-13
4π
-12
三角形窗
-25
8π/N
-25
汉宁窗
-31
8π/N
-44
海明窗
-41
8π/N
-53
不莱克曼窗
-57
12π/N
-74
凯塞窗(α=7.865)
-57
10π/N
-80
这样选定窗函数类型和长度N之后,求出单位脉冲响应
,并按照上式求出
。
是否满足要求,要进行演算。
一般在
尾部加零使长度满足2的整数次幂,以便用FFT计算
。
如果要观察细节,补零点数增多即可。
如果
不满足要求,则要重新选择窗函数类型和长度N,再次验算,直至满足要求。
如果要求线性相位特性,则
还必须满足
根据上式中的正、负号和长度N的奇偶性又将线性相位FIR滤波器分成四类。
要根据所设计的滤波特性正确选择其中一类,例如,要设计线性相位低通特性,可以选择
这一类,而不能选择
这一类。
主程序框图如图2.1所示。
其中幅度特性要求用dB表示。
设
画图时,用
打印幅度特性。
第k点对应的频率
。
为使曲线包络更接近
的幅度特性曲线,DFT变换区间要选大些。
例如窗口长度N=33时,可通过在
末尾补零的方法,使长度变为64,再进行64点DFT,则可以得到更精确的幅度衰减特性曲线。
2.1基本窗函数
数字信号处理领域中所用到的基本窗函数主要有:
矩形窗函数、三角窗函数和巴特利特窗函数。
下面就对这些窗函数展开介绍。
矩形窗(RectangularWindow)函数的时域形式可以表示为:
(公式2-1)
它的频域特性为
公式(2-2)
Boxcar函数:
生成矩形窗
调用方式w=boxcar(n):
输入参数n是窗函数的长度;输出参数w是由窗函数的值组成的n阶向量。
从功能上讲,该函数又等价于w=ones(n,1)。
三角窗是最简单的频谱函数
为非负的一种窗函数。
三角窗函数的时域形式可以表示为:
当n为奇数时
公式(2-3)
当n为偶数时
公式(2-4)
它的频域特性为:
(公式2-5)
三角窗函数的主瓣宽度为
,比矩形窗函数的主瓣宽度增加了一倍,但是它的旁瓣宽度却小得多。
Triang函数:
生成三角窗
调用方式w=triang(n);输入参数n是窗函数的长度;输出参数w是由窗函数的值组成的n阶向量。
三角窗也是两个矩形窗的卷积。
三角窗函数的首尾两个数值通常是不为零的。
当n是偶数时,三角窗的傅立叶变换总是非负数。
2.2广义余弦窗
汉宁窗、海明窗和布莱克曼窗,都可以用一种通用的形式表示,这就是广义余弦窗。
这些窗都是广义余弦窗的特例,汉宁窗又被称为余弦平方窗或升余弦窗,海明窗又被称为改进的升余弦窗,而布莱克曼窗又被称为二阶升余弦窗。
采用这些窗可以有效地降低旁瓣的高度,但是同时会增加主瓣的宽度。
这些窗都是频率为0、2π/(N–1)和4π/(N–1)的余弦曲线的合成,其中
为窗的长度。
通常采用下面的命令来生成这些窗:
(公式2-8)
(公式2-9)
其中,A、B、C适用于自己定义的常数。
根据它们取值的不同,可以形成不同的窗函数,分别是:
汉宁窗A=0.5,B=0.5,C=0;
海明窗A=0.54,B=0.54,C=0;
布莱克曼窗A=0.5,B=0.5,C=0.08;
汉宁窗函数的时域形式可以表示为:
(公式2-10)
它的频域特性为:
(公式2-11)
其中,
为矩形窗函数的幅度频率特性函数。
汉宁窗函数的最大旁瓣值比主瓣值低31dB,但是主瓣宽度比矩形窗函数的主瓣宽度增加了1倍,为8π/N。
hanning函数:
生成汉宁窗
调用方式
(1)w=hanning(n):
输入参数n是窗函数的长度;输出参数w是由窗函数的值组成的n阶向量。
注意:
此函数不返回是零点的窗函数的首尾两个元素。
(2)w=hanning(n,'symmetric'):
与上面相类似。
(3)w=hanning(n,'periodic'):
此函数返回包括为零点的窗函数的首尾两个元素。
海明窗函数的时域形式可以表示为
(公式2-12)
它的频域特性为
(公式2-13)
其中,
为矩形窗函数的幅度频率特性函数。
海明窗函数的最大旁瓣值比主瓣值低41dB,但它和汉宁窗函数的主瓣宽度是一样大的。
Hamming函数:
生成海明窗
调用方式
(1)w=hamming(n):
输入参数n是窗函数的长度;输出参数w是由窗函数的值组成的n阶向量。
(2)w=hamming(n,sflag):
参数sflag用来控制窗函数首尾的两个元素值;其取值为symmetric或periodic;默认值为symmetric。
上面所讨论的几种窗函数,在获得旁瓣抑制的同时却增加了主瓣的宽度。
而凯塞窗定义了一组可调的窗函数,它是由零阶贝塞尔函数构成的,其主瓣能量和旁瓣能量的比例是近乎最大的。
而且,这种窗函数可以在主瓣宽度和旁瓣高度之间自由选择它们的比重,使用户的设计变得非常灵活。
凯塞窗函数的时域形式可表示为
(公式2-16)
其中,
是第1类变形零阶贝塞尔函数,
是窗函数的形状参数,由下式确定:
(公式2-17)
其中,
为凯塞窗函数的主瓣值和旁瓣值之间的差值(dB)。
改变β的取值,可以对主瓣宽度和旁瓣衰减进行自由选择。
β的值越大,窗函数频谱的旁瓣值就越小,而其主瓣宽度就越宽。
Kaiser函数:
生成凯塞窗
调用方式w=kaiser(n,beta):
输入参数n是窗函数的长度;输入参数beta用于控制旁瓣的高度;输出参数w是由窗函数的值组成的n阶向量。
n一定时,beta越大,其频谱的旁瓣就越小,但主瓣宽度相应的增加;当beta一定时,n发生变化,其旁瓣高度不会发生变化。
对于给定的旁瓣高度,切比雪夫窗的主瓣宽度最小。
这是因为它的旁瓣具有相同的高度,也就是具有等波纹性。
切比雪夫窗在边沿的采样点有尖峰。
Chebwin函数:
生成切比雪夫窗
调用方式w=chebwin(n,r):
输入参数n是窗函数的长度;输入参数r用于控制旁瓣的峰值低于主瓣的分贝数。
3.基于matlab的实现
3.1matlab软件简介
MATLAB是美国MathWorks公司出品的商业数学软件,用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境,主要包括MATLAB和Simulink两大部分。
它将数值分析、矩阵计算、科学数据可视化以及非线性动态系统的建模和仿真等诸多强大功能集成在一个易于使用的视窗环境中,代表了当今国际科学计算软件的先进水平。
MATLAB和Mathematica、Maple并称为三大数学软件。
它在数学类科技应用软件中在数值计算方面首屈一指。
MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连matlab开发工作界面接其他编程语言的程序等,主要应用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。
MATLAB具有以下六个特点:
1.编程效率高
用MATLAB编写程序犹如在演算纸上排列出公式与求解问题,MATLAB语言也可通俗地称为演算纸式的科学算法语言。
由于它编写简单,所以编程效率高,易学易懂。
2.用户使用方便
MATLAB语言把编辑、编译、连接和执行融为一体,其调试程序手段丰富,调试速度快,需要学习时间少。
它能在同一画面上进行灵活操作快速排除输入程序中的书写错误、语法错误以至语意错误,从而加快了用户编写、修改和调试程序的速度,可以说在编程和调试过程中它是一种比VB还要简单的语言。
3.扩充能力强
高版本的MATLAB语言有丰富的库函数,在进行复杂的数学运算时可以直接调用,而且MATLAB的库函数同用户文件在形成上一样,所以用户文件也可作为MATLAB的库函数来调用。
因而,用户可以根据自己的需要方便地建立和扩充新的库函数,以便提高MATLAB使用效率和扩充它的功能。
4.语句简单,涵丰富
MATLAB语言中最基本最重要的成分是函数,其一般形式为(a,6,c…)=fun(d,e,f,…),即一个函数由函数名,输入变量d,e,f,…和输出变量a,b,c….组成,同一函数名F,不同数目的输入变量(包括无输入变量)及不同数目的输出变量,代表着不同的含义。
这不仅使MATLAB的库函数功能更丰富,而大大减少了需要的磁盘空间,使得MATLAB编写的M文件简单、短小而高效。
5.高效方便的矩阵和数组运算
MATLAB语言像Basic、Fortran和C语言一样规定了矩阵的一系列运算符,它不需定义数组的维数,并给出矩阵函数、特殊矩阵专门的库函数,使之在求解诸如信号处理、建模、系统识别、控制、优化等领域的问题时,显得大为简捷、高效、方便,这是其它高级语言所不能比拟的。
6.方便的绘图功能
MATLAB的绘图是十分方便的,它有一系列绘图函数(命令),使用时只需调用不同的绘图函数(命令),在图上标出图题、XY轴标注,格绘制也只需调用相应的命令,简单易行。
另外,在调用绘图函数时调整自变量可绘出不变颜色的点、线、复线或多重线。
3.2各窗函数的图形
用MATLAB编程绘制各种窗函数的形状,窗函数的长度为21。
程序如下:
clf;Nwin=21;n=0:
Nwin-1;
figure
(1)
forii=1:
4
switchii
case1
w=boxcar(Nwin);%矩形窗函数
stext='矩形窗函数';
case2
w=hanning(Nwin);%hanning窗函数
stext='hanning窗函数';
case3
w=hamming(Nwin);%hamming窗函数
stext='hamming窗函数';
case4
w=triang(Nwin);%三角窗函数
stext='triang窗函数';
end
posplot=['2,2,'int2str(ii)];
subplot(posplot);
stem(n,w);%绘出窗函数
holdon
plot(n,w,'r');%绘出包络线
xlabel('n');ylabel('w(n)');title(stext);
holdoff;gridon;
end
程序运行结果如下图:
3.1各窗函数的图形
3.3各窗函数的幅频特性
用MATLAB编程,采用512个频率点绘制各窗函数的幅频特性。
程序如下:
clf;Nf=512;%窗函数复数频率特性的数据点数
Nwin=20;%窗函数数据长度
figure
(1)
forii=1:
4
switchii
case1
w=boxcar(Nwin);%矩形窗
stext='矩形窗函数幅频';
case2
w=hanning(Nwin);%汉宁窗
stext='hanning窗函数幅频';
case3
w=hamming(Nwin);%海明窗
stext='hamming窗函数幅频';
case4
w=triang(Nwin);%三角窗
stext='triang窗函数幅频';
end
[y,f]=freqz(w,1,Nf);%求解窗函数的幅频特性
mag=abs(y);%求得窗函数幅频特性
posplot=['2,2,'int2str(ii)];
subplot(posplot);
plot(f/pi,20*log10(mag/max(mag)));%绘制窗函数的幅频特性
xlabel('归一化频率');ylabel('振幅/dB');
title(stext);gridon;
end
程序运行结果如下图:
3.2各窗函数的幅频特性
4.频谱泄露
4.1频谱泄漏原理
对于频率为fs的正弦序列,它的频谱应该只是在fs处有离散谱。
但是,在利用DFT求它的频谱做了截短,结果使信号的频谱不只是在fs处有离散谱,而是在以fs为中心的频带围都有谱线出现,它们可以理解为是从fs频率上“泄露”出去的,这种现象称为频谱“泄露”。
在实际问题中遇到的离散时间序列x(n)通常是无限长序列,因而处理这个序列的时候需要将它截短。
截短相当于将序列乘以窗函数w(n)。
根据频域卷积定理,时域中x(n)和w(n)相乘对应于频域中它们的离散傅立叶变换X(jw)和W(jw)的卷积。
4.2产生机理
窗函数修正了由于信号的非周期性并减小了频谱中由于泄露而带来的测量不准确性。
快速傅里叶变换假定了时间信号是周期无限的。
但在分析时,我们往往只截取其中的一部分,因此需要加窗以减小泄露。
窗函数可以加在时域,也可以加在频域上,但在时域上加窗更为普遍。
截断效应带来了泄漏,窗函数是为了减小这个截断效应,其设计成一组加权系数。
例如,一个窗函数可以定义为:
w(t)=g(t) -T/2w(t)=0 其他
g(t)是窗函数,T是窗函数的时间
待分析的数据x(t)则表示为:
x(t)=w(t)*x(t)'
x(t)'表示原始信号x(t)表示待分析信号。
加窗在时域上表现的是点乘,因此在频域上则表现为卷积。
卷积可以被看成是一个平滑的过程。
这个平滑过程可以被看出是由一组具有特定函数形状的滤波器,因此,原始信号中在某一频率点上的能量会结合滤波器的形状表现出来,从而减小泄漏。
基于这个原理,人们通常在时域上直接加窗。
大多数的信号分析仪一般使用矩形窗(rectangular),汉宁(hann),flattop和其他的一些窗函数。
矩形窗函数:
w(k)=1
汉宁窗:
w(k)=0.5*(1-cos(2*pi*k/(N-1))) 0<=k<=N-1
由于加窗计算中衰减了原始信号的部分能量,因此对于最后的结果还需要加上修正系数。
在线性谱分析中,一般使用幅度系数(amplitude correction),在功率谱中,一般使用能量系数(energy correction)。
4-2频谱泄漏
上图阐述了非整周期信号采样时造成泄露的理论原因。
4.3窗函数的频谱泄漏的抑制方法
如果一个信号由两个不同频率的正弦波叠加而成,则泄露可能使用户无法区分其各自的频率分量。
如果它们的频率相差很多,而其中的一个分量的能量远大于另外一个,则能量小的那个频率可能完全会被从能量大的频率分量泄露过来的能量所覆盖。
若它们的频率相近,则它们双方泄露的能量可能会使它们之间的频率点和它们原有频率点的能量保持一致,就是说,可能会使用户只看到一个频率峰值。
有两种方法可以避免泄露的发生。
第一个方法是采集信号足够长,基本上可以覆盖到整个有效信号的时间跨度。
这种方法经常在瞬态捕捉中被使用到,比如说冲击试验,如果捕捉的时间够长,捕捉到的信号可以一直包括了振动衰减为零的时刻。
在这种情况下,可以不加窗函数。
第2种方法是对于信号的采样时间正好能和信号整周期倍数吻合,即一帧数据的长度正好是原始信号中频率的整数倍。
比如一个正弦波的频率在1000Hz,因此采样频率就需要设到8000Hz..每个正弦周期有8个采样点。
1024点的数据正好涵盖了8个整周期。
在这种情况下,没有窗函数也可以避免泄露的情况。
以下的图片显示了一个1000Hz频率的正弦波的频谱,其中没有泄露。
4-3-1
下图显示了一个1010Hz的正弦波,但是其频率谱中有一个较宽的泄露情况。
该频率谱在1010Hz以外的地方也有着明显的能量值。
可以看到,泄露的能量总是围绕着其原始的峰值。
4-3-2
很多种的窗函数被设计出来以减小这种泄露情况,以下的图就展示了对1010Hz正弦波加floattop窗之后得到的功率谱。
4-3-3
加上floattop窗之后,泄露明显的减小。
正弦频率和噪音可以被很好的区别开。
但是,这个窗函数仍然使得频率点的估计不是那么的精确。
4.4窗函数的选择
如果在测试中可以保证不会有泄露的发生,则不需要用任何的窗函数(在软件中可选择uniform)。
但是如同刚刚讨论的那样,这种情况只是发生在时间足够长的瞬态捕捉和一帧数据中正好包含信号整周期的情况。
如果测试信号有多个频率分量,频谱表现的十分复杂,且测试的目的更多关注频率点而非能量的大小。
在这种情况下,需要选择一个主畔够窄的窗函数,汉宁窗是一个很好的选择。
如果测试的目的更多的关注某周期信号频率点的能量值,比如,更关心其EUpeak,EUpeak-peak,EUrms或者EUrms2,那么其幅度的准确性则更加的重要,可以选择一个主畔稍宽的窗,flattop窗在这样的情况下经常被使用。
对冲击实验的数据进行分析时,因为在数据帧开始段的一些重要信息会被一般的窗函数所衰减,因此可以使用force/exponential窗。
Force窗一移去了数据帧末端的噪声,对激励信号有用。
而exponential窗则确保响应信号在末端的振动衰减为零值。
如果被测信号是随机或者未知的,选择汉宁窗
5.实验结果分析
由上图各窗函数的幅频特性可以看到各种窗函数都具有明显的主瓣和旁瓣。
主瓣频宽和旁瓣的幅值衰减特性决定了窗函数的应用场合。
矩形窗具有最窄的主瓣,但也有最大的旁瓣峰值(第一旁瓣衰减为13dB);不同窗函数在这两方面的特点是不同的,因此应根据具体的问题进行选择。
通常来讲,海明窗和汉宁窗的主瓣,具有较小的旁瓣和较大的衰减速度,是较为常用的窗函数。
本文针对连续信号在频域分析中常常出现的频谱泄漏现象产生的原因进行了分析研究,提出了在频谱分析中减小频谱泄露的方法与途径,给出了仿真结果。
采用此种处理手段可以使信号的频谱分析精度更准确,函数谱更加贴近实际,处理效果更加理想。
6.心得体会
通过本次对窗函数的研究,更深层次的了解了窗函数的设计原理和频谱分析,并学会了matlab软件的使用及分析,熟练掌握了matlab软件的各项功能及用途,对以后的学习及研究运用打下了好的基础。
通过一个星期的钻研及学习,对于窗函数的各项指标有了进一步的掌握,更增强了自己的学习能力和搜集资料整合资料的能力,受益颇多也获益良多,更感在这期间老师及同学给予我的帮助,没有你们,我也无法完成这份报告,如今报告完整了,心里既充实又感慨颇多。
参考文献
[1]丁玉美,高西全.数字信号处理[M].:
电子科技大学,2000.
[2]程佩青.数字信号处理教程第二版[M].:
清华大学,2001.
[3]敏歌.基于窗函数法的FIR数字滤波器的设计[J].师大学学报,2007,35:
72-74.
[4]卫东,绍全.加窗离散傅里叶变换测频分辨率研究[J].电子科技大学学报,2000,27
(2):
157-160.
[5]许珉,鸿博.基于Blankman-harris窗的加窗FFT插值修正算法[J].大学学报,2005,26(4):
99-101.