数字低通滤波器的窗函数设计

齐齐哈尔大学开放实验室项目论文题目窗函数在数字滤波器中的应用研究院（系）通信与电子工程学院专业班级电子信息工程061班学生姓名张文霞指导教师题原成绩2007年5月13日随着信息时代的到来，数字信号处理已经成为当今一门及其重要的学科和技术，并且在通信、语音、图像、自动控制等众多领域得到了广泛的应用，在数字信号处理中，数字滤波器占有及其重要的地位，它具有精度高，可靠性好、灵活性大等特点，现代数字滤波器可以用软件和硬件两种方式来实现，软件方式实现的优点是可以通过滤波器参数的改变去调整滤波器的性能。MATLAB是一种面向科学和工程计算的语言，它集数值分析，矩阵运算、信号处理和图形显示于一体，具有编程效率高、调试手段丰富、扩充能力强等特点，MATLAB的信号处理工具箱具有强大的函数功能，它不仅可以用来设计数字滤波器，还可以使设计达到最优化，是数字滤波器设计的强有力工具。1.FIR滤波器简介FIR数字滤波器设计最简单的方法是窗函数法，通常也称为傅立叶级数法。它是在时域进行的，因而必须由理想滤波器的频率响应)(jwdeH推导出其单位冲激响应)(nhd，在设计一个FIR数字滤波器的单位冲激响应)(nh去逼近)(nhd。根据冲激响应的时域特性，数字滤波器可分为无限长冲激响应（IIR）和有限长冲激响应滤波器（FIR），FIR的突出优点是：系统总是稳定的、易于实现线性相位、允许设计多通带（或多阻带）滤波器，但与IIR相比，在满足同样阻带衰减的情况下需要的阶数较高，滤波器的阶数越高，占用的运算时间越多，因此在满足指标要求的情况下应尽量减少滤波器的阶数。2．窗函数设计法的步骤窗函数设计法是一种通过截短和计权的方法使无限长非因果序列成为有限长脉冲应响应序列的设计方法，通常在设计滤波器之前，应该先根据具体的工程应用确定滤波器的技术指标，在大多数实际应用中，数字滤波器常常被用来实现选频操作，所以指标的形式一般为在频域中以分贝值给出的相对幅度响应和相位响应。用窗函数法设计FIR滤波器的步骤如下：（1）根据过渡带宽及阻带衰减要求，选择窗函数的类型并估计窗口长度N（或阶数M=N-1），窗函数类型可根据最小阻带衰减As独立选择，因为窗口长度N对最小阻带衰减As没有影响，在确定窗函数类型以后，可根据过渡带宽小于给定指标确定所拟用的窗函数的窗口长度N，设待求滤波器的过渡带宽为Δw，它与窗口长度N近似成反比，窗函数类型确定后，其计算公式也确定了，不过这些公式是近似的，得出的窗口长度还要在计算中逐步修正，原则是在保证阻带衰减满足要求的情况下，尽量选择较小的N，在N和窗函数类型确定后，即可调用MATLAB中的窗函数求出窗函数wd（n）。（2）根据待求滤波器的理想频率响应求出理想单位脉冲响应hd（n），如果给出待求滤波器频率应为Hd（ejw），则理想的单位脉冲响应可以用下面的傅里叶反变换式求出：dweeHnhjwjwdd)(21)(在一般情况下，hd（n）是不能用封闭公式表示的，需要采用数值方法表示；从w=0到w=2π采样N点，采用离散傅里叶反变换（IDFT）即可求出。（3）计算滤波器的单位脉冲响应h（n），它是理想单位脉冲响应和窗函数的乘积，即h（n）=hd（n）×wd（n），在MATLAB中用点乘命令表示为h=hd×wd。（4）验算技术指标是否满足要求，为了计算数字滤波器在频域中的特性，可调用freqz子程序，如果不满足要求，可根据具体情况，调整窗函数类型或长度，直到满足要求为止。使用窗函数法设计时要满足以下两个条件：（1）窗谱主瓣尽可能地窄，以获得较陡的过渡带；（2）尽量减少窗谱的最大旁瓣的相对幅度，也就是使能量尽量集中于主瓣，减小峰肩和纹波，进行增加阻带的衰减。3．几种常用的窗3.1矩形窗矩形窗函数的时域形式可以表现为：其他,010,1)()()(NnnRnwN3.2.汉宁窗（升余弦窗）=0.5利用付氏变换的移位特性，汉宁窗频谱的幅度函数W（ω）可用矩形窗的幅度函数表示为：三部分矩形窗频谱相加，使旁瓣互相抵消，能量集中在主瓣，旁瓣大大减小，主瓣宽度增加1倍。3.3.海明窗（改进的升余弦窗）对汉宁窗的改进，在主瓣宽度（对应第一零点的宽度）相同的情况下，旁瓣进一步减小，可使99.96%的能量集中在主瓣内。3.4.布莱克曼窗（三阶升余弦窗）增加一二次谐波余弦分量，可进一步降低旁瓣，但主瓣宽度进一步增加，增加N可减少过渡带。频谱的幅度函数为：4．几种数字滤波器的窗函数设计4.1数字低通滤波器的窗函数设计假设理想低通数字滤波器的频率响应)(jwdeH，其幅频特性为1)(jwdeH。那么该滤波器的：ccjwjwd，，01)(式中cw表示截止频率（rad），表示采样延迟。由此可得理想数字低通滤波器的单位冲激响应)(nhd为：ccwwcjwnjwjwnjwddnnwdweedweeHnh)()(sin121)(21)(（1）根据式（1），可知)(nhd是一个以α为中心的偶对称的无限长非因果序列，如果截取一段n=0～N-1的hd(n)作为h(n)，则为要保证所得到是线性相位FIR滤波器，延时α应为h(n)长度N的一半,α=(N-1)/2（2）其中WR(n)=RN(n)（为了表示一致）设W(ejω)为该窗口函数的频谱:用幅度函数和相函数来表示，则有。其线性相位部分e-jωα则是表示延时一半长度α=(N-1)/2，对频响起作用的是它的幅度函数图矩形窗函数及其幅度函数理想频响也可以写成幅度函数和相位函数的表示形式Hd(ejω)=Hd(ω)e-jωα其中幅度函数为二个信号时域乘积对应于频域卷积，所以有如果也以幅度函数H(ω)和相位函数来表示H（ejω）,H（ejω）=H(ω)e-jw则实际FIR滤波器的幅度函数H（ω）为正好是理想滤波器幅度函数与窗函数幅度函数的卷积。为了满足工程上的需要，可以通过改变窗函数的形状来改善滤波器的幅度频率特性，而窗函数的选择原则是：1)具有较低的旁瓣幅度，尤其是第一旁瓣的幅度；2)旁瓣的幅度下降的速率要快，以利于增加阻带的衰减；3)主瓣的宽度要窄，这样可以得到比较窄的过渡带。通常上述的几点难以同时满足。实际中设计FIR数字滤波器往往要求是线性相位的，因此要求)(nw满足线性相位的条件，即要求w(n)满足：)1()(nNwnw所以，窗函数不仅有截短的作用，而且能够起到平滑的作用，在很多领域得到了应用。4.2利用海明窗设计Ⅱ型数字低通滤波器：技术指标bBAbBAwwspsp50,25.0,4.0,2.0设计程序如下：clearall;Wp=0.2*pi;Ws=0.4*pi;tr_width=Ws-Wp;N=ceil(6.6*pi/tr_width)+1n=0:1:N-1;Wc=(Ws+Wp)/2;hd=ideal_lp(Wc,N);w_ham=(hamming(N))';h=hd.*w_ham;[db,mag,pha,w]=freqz_m2(h,[1]);delta_w=2*pi/1000;Ap=-(min(db(1:1:Wp/delta_w+1)))As=-round(max(db(Ws/delta_w+1:1:501)))subplot(221)stem(n,hd)title('理想单位脉冲响应hd（n）')subplot(222)stem(n,w_ham)title('海明窗w(n)')subplot(223)stem(n,h)title('实际单位脉冲响应hd(n)')subplot(224)plot(w/pi,db)title('幅度响应（dB）')axis([0,1,-100,10])运行结果：N=34Ap=0.0477As=52图（1）Ⅱ型数字低通滤波器冲激响应与幅度响应4.3数字高通滤波器的窗函数设计4.3.1无相移时的高通FIR滤波器的设计幅频特性为1)(jwdeH，而且相移的理想高通滤波器的频率响应)(jwdeH表达式为：ccjwjwd,0,1)(式中cw表示截止频率（rad），表示采样延迟。此时，该理想数字高通滤波器的单位冲激响应)(nHd为：nnwndweedweedweeHnhcwjwnjwwjwnjwjwnjwddcc)()(sin1121)(21)(（3）根据式（3）可知，)(nHd也为无限非因果序列，关于对称。若利用长度为N的窗函数)(nwN对)(nhd进行截取，则截取后)(nhd的序列的表达式如式（2）所示，它关于21N偶对称。当N为奇数时，所设计的FIR数字高通滤波器为I型滤波器；当N为偶数时，为II型滤波器。4.3.2利用汉宁窗设计I型数字高通滤波器其中技术指标为：bBAbBAwwsPcp50,25.0,4.0,6.0设计程序如下:clearall;Wp=0.6*pi;Ws=0.4*pi;tr_width=Wp-Ws;N=ceil(6.2*pi/tr_width)n=0:1:N-1;Wc=(Ws+Wp)/2;hd=ideal_hpl(Wc,N);w_han=(hanning(N))';h=hd.*w_han;[db,mag,pha,w]=freqz_m2(h,[1]);delta_w=2*pi/1000;Ap=-(min(db(Wp/delta_w+1:1:501)))As=-round(max(db(1:1:Ws/delta_w+1)))subplot(221)stem(n,hd)title('理想单位脉冲响应hd(n)')subplot(222)stem(n,w_han)title('汉宁窗w(n)')subplot(223)stem(n,h)title('实际单位脉冲响应hd（n）')subplot(224)plot(w/pi,db)title('幅度响应（dB）')axis([0,1,-100,10])运行结果：N=31Ap=0.0887As=44I型FIR数字高通滤波器冲激响应与幅度响应（三）数字带通滤波器的窗函数设计无相移时的带通FIR滤波器的设计幅频特性为1)(jwdeH，而且相移的理想带通滤波器的频率响应)(jwdeH表达式为：chclchcljwjwd,,0,1)(式中clw表示下截止频率（rad），chw表示上截止频率（rad），表示采样延迟。此时，该理想数字带通滤波器的单位冲激响应)(nHd为：nnwnwdweedweedweeHnhclchwwjwnj)()(sin1121)(21)(（4）根据式（4）可知，)(nHd也为无限非因果序列，关于对称。若利用长度为N的窗函数)(nwN对)(nhd进行截取，则截取后)(nhd的序列的表达式如式（2）所示，它关于21N偶对称。当N为奇数时，所设计的FIR数字带通滤波器为I型滤波器；当N为偶数时，为II型滤波器。利用布拉克曼窗设计II型数字高通滤波器其中技术指标为：bBAbBAbBAbBA设计程序如下：clearall;Wpl=0.4*pi;Wph=0.6*pi;Wsl=0.2*pi;Wsh=0.8*pi;tr_width=min((Wpl-Wsl),(Wsh-Wph));N=ceil(11*pi/tr_width)+1n=0:1:N-1;Wcl=(Wsl+Wpl)/2;Wch=(Wsh+Wph)/2;hd=ideal_bpl(Wcl,Wch,N);w_bman=(blackman(N))';h=hd.*w_bman;[db,mag,pha,w]=freqz_m2(h,[1]);delt