数字信号处理|| 用窗函数法设计FIR数字滤波器

一、实验目的

（1）加深对窗函数法设计FIR数字滤波器的基本原理的理解。
（2）学习用MATLAB语言的窗函数法编写设计FIR数字滤波器的程序。
（3）了解MATLAB有关窗函数法设计的常用子函数。

二、实验涉及的MATLAB子函数

（1）boxcar
功能：矩形窗。
格式：w＝boxcar（n）
说明：boxcar（n）函数可产生一长度为n的矩形窗函数。

（2）triang
功能：三角窗。
格式：w＝triang（n）
说明：triang（n）函数可得到n点的三角窗函数。三角窗系数为：当n为奇数时，

当n为偶数时， 

（3）bartlett
功能：bartlett（巴特利特）窗。
格式：w＝Bartlett（n）
说明：bartlett（n）可得到n点的bartlett窗函数。bartlett窗函数系数为

（4）hamming
功能：hamming（哈明）窗。
格式：w＝hamming（n）
说明：hamming（n）可产生n点的hamming窗。hamming窗函数系数为

（5）hanning
功能：hanning（汉宁）窗。
格式：w＝hanning（n）
说明：hanning（n）可产生n点的hanning窗。hanning窗函数系数为

（6）blackman
功能：blackman（布莱克曼）窗。
格式：w＝blackman（n）
说明：blackman（n）可产生n点的blackman窗。blackman窗函数系数为

 与等长度的hamming和hanning窗相比，blackman窗的主瓣稍宽，旁瓣稍低。

（7）chebwin
功能：chebyshev（切比雪夫）窗。
格式：w＝chebwin（n，r）
说明：chebwin（n，r）可产生n点的chebyshev窗函数，其傅里叶变换后的旁瓣波纹低于主瓣r dB。注意：当n为偶数时，窗函数的长度为n＋1。
（8）kaiser
功能：kaiser（凯塞）窗。
格式：w＝kaiser（n，beta）
说明：kaiser（n，beta）可产生n点的kaiser窗函数，其中，beta为影响窗函数旁瓣的b参数，其最小的旁瓣抑制a与b之间的关系为

增加b可使主瓣变宽、旁瓣的幅度降低。

（9）fir1
功能：基于窗函数的FIR数字滤波器设计——标准频率响应，以经典方法实现加窗线性相位FIR滤波器设计，可设计出标准的低通、带通、高通和带阻滤波器。
格式：
※b＝fir1（n，Wn）；设计截止频率为Wn的hamming（哈明）加窗线性相位滤波器，滤波器系数包含在b中。当0≤Wn≤1（Wn＝1相应于0.5fs）时，可得到n阶低通FIR滤波器。
当Wn＝［W1W2］时，fir1函数可得到带通滤波器，其通带为w1<w<w2。

※b＝fir1（n，Wn，¢ftype¢）；可设计高通和带阻滤波器，由ftype决定：
当ftype＝high时，设计高通FIR滤波器；
当ftype＝stop时，设计带阻FIR滤波器。
在设计高通和带阻滤波器时，fir1函数总是使用偶对称N为奇数（即第一类线性相位FIR滤波器）的结构，因此当输入的阶次为偶数时，fir1函数会自动加1。※b＝fir1（n，Wn，Window）；利用列矢量Window中指定的窗函数进行滤波器设计，Window长度为n＋1。如果不指定Window参数，则fir1函数采用hamming窗。
※b＝fir1（n，Wn，¢ftype¢，Window）；可利用ftype和Window参数设计各种加窗的滤波器。由fir1函数设计的FIR滤波器的群延迟为n/2。

三、实验原理

（1）运用窗函数法设计FIR数字滤波器
FIR数字滤波器的系统函数为

 这个公式也可以看成是离散LSI系统的系统函数：

分母a0为1，其余ak全都为0时的一个特例。由于极点全部集中在0点，稳定和线性相位特性是FIR滤波器突出的优点，因此在实际中具有更广泛的使用价值。
FIR滤波器的设计任务是选择有限长度的h（n），使传输函数H（ejw）满足技术要求。主要设计方法有窗函数法、频率采样法和切比雪夫等波纹逼近法等。本实验主要介绍用窗函数法设计FIR数字滤波器。用窗函数法设计FIR数字滤波器的基本步骤是：
（1）根据过渡带和阻带衰减设计指标选择窗函数的类型，估算滤波器的阶数N。
（2）由数字滤波器的理想频率响应H（ejw）求出其单位冲激响应hd（n）。
对于理想的数字低通滤波器频率响应，有下列的子程序可以实现（该程序名为ideal－lp.m）：

functionhd＝ideal－lp（wc，N）
%hd＝点0到N－1之间的理想脉冲响应
%wc＝截止频率（弧度）
%N＝理想滤波器的长度
tao＝（N－1）/2；
n＝［0：（N－1）］；
m＝n－tao＋eps；  %加一个小数以避免0作除数
hd＝sin（wc*m）./（pi*m）；

其它选频滤波器则可以由低通频响特性合成。如一个通带在wc1～wc2之间的带通滤波器，在给定N值的条件下，可以用下列程序实现：

hd＝ideal－lp（wc2，N）－ideal－lp（wc1，N）；

由于第一类线性相位滤波器（类型Ⅰ）能进行低通、高通、带通、带阻滤波器的设计，因此，本实验所有滤波器均采用第一类线性相位滤波器。

（3）计算数字滤波器的单位冲激响应h（n）＝w（n）hd（n）。
（4）检查设计出的滤波器是否满足技术指标。
如果不满足技术指标，则需要重新选择或调整窗函数的类型，估算滤波器的阶数N。再重复前面的四个步骤，直到满足指标为止。
常用的窗函数有矩形窗、三角形窗、汉宁窗、哈明窗、切比雪夫窗、布莱克曼窗、凯塞窗等，MATLAB均有相应的子函数可以调用。
另外，MATLAB信号处理工具箱还提供了fir1子函数，可以用于窗函数法设计FIR滤波器。

四、实验任务

（1）阅读并输入实验原理中介绍的例题程序，观察输出的数据和图形，结合基本原理理解每一条语句的含义。
（2）认真阅读实验原理，明确本次实验任务，读懂例题程序，了解实验方法。
（3）列写调试通过的实验程序，打印或描绘实验程序产生的曲线图形。

五、实验程序及运行结果

（1）各种窗函数特性的比较
1）在同一图形坐标上显示矩形窗、三角形窗、汉宁窗、哈明窗、布莱克曼窗、凯塞窗的特性曲线。
MATLAB 程序如下：

N = 64;         % 输入N
beta = 7.865;   % 凯塞窗的β值
n = 0:N-1;     % 创建一个从0到N-1的序列
wbo = boxcar（N）;    % 矩形窗
wtr = triang（N）;   % 三角形窗
whn = hanning（N）;  % 汉宁窗
whm = hamming（N）;  % 哈明窗
wbl = blackman（N）; % 布莱克曼窗
wka = kaiser（N, beta）; % 凯塞窗
% 在同一界面上作图
plot（n, [wbo, wtr, whn, whm, wbl, wka]）;
% 限定y轴显示范围为0到1.1
axis（[0, N, 0, 1.1]）;
% 添加图例
legend（'矩形', '三角形', '汉宁', '哈明', '布莱克曼', '凯塞'）;

 运行结果：

为了便于滤波器设计，表24-1给出了六种窗函数的特性参数。
表24-1 六种窗函数的特性参数表

（2）用窗函数法设计FIR数字低通滤波器
1）用矩形窗设计一个FIR数字低通滤波器，要求：N＝64，截止频率为wc＝0.4p，描绘理想和实际滤波器的脉冲响应、窗函数及滤波器的幅频响应曲线。
MATLAB程序如下：

%矩形窗设计FIR低通滤波器
wc=0.4*pi;            %输入设计指标
N=64;n=0:N-1;
hd=ideal_lp（wc,N）;    %建立理想低通滤波器
windows=（boxcar（N））'; %使用矩形窗，并将列向量变为行向量
b=hd.*windows;        %求FIR系统函数系数
[H,w]=freqz（b,1）;     %求解频率特性
dbH=20*log10（（abs（H）+eps）/max（abs（H）））; %化为分贝值
subplot（2,2,1）,stem（n,hd,'filled','k'）;
axis（[0,N,1.1*min（hd）,1.1*max（hd）]）;title（'理想脉冲响应'）;
xlabel（'n'）;ylabel（'hd（n）'）;
subplot（2,2,2）,stem（n,windows,'filled','k'）;
axis（[0,N,0,1.1]）;title（'窗函数特性'）;
xlabel（'n'）;ylabel（'wd（n）'）;
subplot（2,2,3）,stem（n,b,'filled','k'）;
axis（[0,N,1.1*min（b）,1.1*max（b）]）;title（'实际脉冲响应'）;
xlabel（'n'）;ylabel（'h（n）'）;
subplot（2,2,4）,plot（w/pi,dbH,'k'）;
axis（[0,1,-80,10]）;title（'幅度频率响应'）;
xlabel（'频率（单位:\pi）'）;ylabel（'H（e^{j\omega}）'）;
set（gca,'XTickMode','manual','XTick',[0,wc/pi,1]）;
set（gca,'YTickMode','manual','YTick',[-50,-20,-3,0]）;grid

function hd=ideal_lp（wc,M）;
% Ideal LowPass filter computation
% --------------------------------
% [hd] = ideal_lp（wc,M）;
% hd = ideal impulse response between 0 to M-1
% wc = cutoff frequency in radians
% M = length of the ideal filter
alpha = （M-1）/2;

运行结果：

六、实验心得

 通过用窗函数法设计FIR数字滤波器实验，我学习到了窗函数法设计FIR数字滤波器的基本原理，并掌握了使用MATLAB进行FIR滤波器设计的编程方法。了解了多种窗函数如矩形窗、汉宁窗、哈明窗等对滤波器性能的影响，并通过fir1函数设计出了满足特定要求的FIR滤波器。实验中对不同窗函数特性的比较和实际滤波器设计的实践，使我收获满满。