河海大学数字信号处理实验五

-1-实验五有限长单位脉冲响应滤波器设计班级：姓名：学号：一、实验目的1.掌握用窗函数法、频率采样法及优化设计法设计FIR滤波器的原理及方法，熟悉相应的MATLAB编程。2.熟悉线性相位FIR滤波器的幅频特性和相频特性。3.了解各种不同窗函数对滤波器性能的影响。二、实验原理window=ones(1,N):产生N点矩形窗，行向量。window=hann(N):产生N点汉宁窗，列向量。window=hanning(N):产生N点非零汉宁窗，列向量。等价于去除hann(N+2)的第一个零元素和最后一个零元素，得到的N点非零窗函数。window=hamming(N):产生N点海明窗，列向量。window=blackman(N):产生N点布莱克曼窗，列向量。window=kaiser(N,beta):产生参数为beta的N点凯塞窗，列向量。[M,Wd,beta,ftype]=kaiserord(f,a,dev,fs):凯塞窗参数估计。f为一组边界频率，最高频率为fs/2。a为f中各个频带的幅度值，通带取1，阻带取0。如果f中有2个元素，则形成3个频带，其中第1个和第3个是通带或阻带，第2个是过渡带，a中也有2个元素，指明第1个和第3个频带是通带还是阻带；如果f中有4个元素，则形成5个频带，其中1，3和5是通带或阻带，2和4是过渡带，a中有3个元素，指明1，3和5是通带还是阻带。dev的维数与a相同，指明每个频带上的波动值。fs为采样频率。M为FIR滤波器的阶数，M=N-1。Wd为归一化边界频率，等于数字边界角频率除以π，或者边界频率除以fs/2。beta就是凯塞窗的参数β。ftype为滤波器的类型。fs=8000;f=[10001500];a=[10];dev=[0.050.01];[M,Wd,beta,ftype]=kaiserord(f,a,dev,fs);h=fir1(M,Wd,ftype,kaiser(M+1,beta));-2-freqz(h);fs=8000;f=[1000130022102410];a=[010];dev=[0.010.050.01];[M,Wd,beta,ftype]=kaiserord(f,a,dev,fs);hh=fir1(M,Wd,ftype,kaiser(M+1,beta));[H,f]=freqz(hh,1,1024,fs);figure;plot(f,abs(H));gridon;00.10.20.30.40.50.60.70.80.91-1500-1000-5000NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-150-100-50050NormalizedFrequency(rad/sample)Magnitude(dB)0500100015002000250030003500400000.20.40.60.811.21.4b=fir1(M,Wd,'ftype',window):用窗函数法求FIR滤波器的系数b(单位脉冲响应)。M为滤波器的阶数，M=N-1。Wd为一组归一化边界频率，通带和阻带间隔分布，无过渡带；只有一个元素，表示低通或高通滤波器；有两个元素表示带通和带阻滤波器；有三个及以上元素，表示多带滤波器。'ftype'表示滤波器类型，'high'表示高通滤波器，'stop'表示带阻滤波器，'DC-0'表示多带滤波器的第一个频带为阻带，'DC-1'表示多带滤波器的第一个频带为通带。window为窗口类型，缺省为海明窗。b=fir2(M,f,m,window):用频率采样法求FIR滤波器的系数b。M为滤波器的阶数，M=N-1。f为一组归一化频率，第一个元素必须为0，最后一个元素必须为1(对应奈奎斯特频率，即采样频率的一半)，中间的元素按升序排列。m的维数与f相同，指明f中每个频率上的理想幅度。window为窗口类型，缺省为海明窗。Fir2可以实现任意幅度特性的滤波器。三、实验内容1、用窗函数法设计一个线性相位FIR高通滤波器，通带边界频率为0.7π，阻带边界频率为0.5π，要求阻带衰减不小于50dB，窗函数从矩形窗、汉宁窗、海明窗和布莱克曼窗中选取，且要求滤波器的阶数最小。给出该滤波器的单位脉冲响应h(n)的解析式，并用MATLAB-3-绘出时域波形和幅频特性。2、用窗函数法设计一个线性相位FIR带通滤波器，通带边界频率为3kHz和7kHz，阻带边界频率为2kHz和9kHz，采样频率为20kHz，要求阻带衰减不小于40dB，窗函数从矩形窗、汉宁窗、海明窗和布莱克曼窗中选取，且要求滤波器的阶数最小。给出该滤波器的单位脉冲响应h(n)的解析式，并用MATLAB绘出时域波形和幅频特性。3、用窗函数法设计一线性相位FIR带阻滤波器，通带边界频率为600Hz和1400Hz，阻带边界频率为700Hz和1200Hz，采样频率为4000Hz，要求阻带衰减不小于50dB，窗函数从矩形窗、汉宁窗、海明窗和布莱克曼窗中选取，且要求滤波器的阶数最小。给出该滤波器的单位脉冲响应h(n)的解析式，并用MATLAB绘出时域波形和幅频特性。4、分别用矩形窗、汉宁窗、海明窗和布莱克曼窗设计FIR低通滤波器，理想滤波器的边界频率为0.8π，N=51，用Fir1函数设计。根据所设计滤波器的阻带最小衰减值确定实际阻带边界频率。5、用频率采样法设计一个线性相位低通滤波器，N=40，通带边界频率为π/4，过渡带设置一个采样点|H(k)|=0.5，给出单位脉冲响应h(n)的序列值，并作出幅频响应。6、用频率采样法设计一个线性相位高通滤波器，通带边界频率为2π/3，过渡带设置一个采样点|H(k)|=0.5，分别求N=31和N=32时的单位脉冲响应h(n)的序列值，并分别作出幅频响应。7、用频率采样法设计一个多带线性相位滤波器，理想幅频特性如下图所示，N=50，过渡带设置一个采样点，过渡点的值设为0.5，窗函数选矩形窗。给出该滤波器的单位脉冲响应h(n)的解析式，并用MATLAB绘出时域波形和幅频特性。00.10.20.30.40.50.60.70.80.9100.20.40.60.811.21.41.61.8/|Hd(ej)|理想滤波器的幅频特性-4-四、实验结果与分析（一）窗函数法：1、解：（1）求数字边界频率:0.7,0.5cr（2）求理想低通滤波器的边界频率：理想低通滤波器的边界频率是实际滤波器幅度函数正负肩峰频率的中心，但这里不知道正负肩峰的频率，因此用通带和阻带边界频率的中心近似:0.62crn（3）选择窗函数和窗口长度：阻带衰减不小于50dB，因此选择海明窗(其阻带最小衰减为-53dB)；滤波器的过渡带宽为0.7π-0.5π=0.2π，因此窗口长度为:6.610.233,162NNN线性相位延迟常数为:N1162（4）根据理想边界频率n和线性相位延迟常数，求理想单位脉冲响应d()hn:dsin()sin[()]()()1nnnnnnhnn（5）窗函数与理想单位脉冲响应相乘，得到线性相位FIR低通滤波器的单位脉冲响应:-5-33dsin(16)sin[0.6(16)]0.540.46cos()16()()()16(16)0.416nnnRnnhnwnhnnn程序：clear;clc;N=33;n=0:N-1;hd=(sin(pi*(n-16))-sin(0.6*pi*(n-16)))./(pi*(n-16));hd(17)=0.4;win=hanning(N);h=win'.*hd;figure;stem(n,h);xlabel('n');ylabel('h(n)');grid;title('FIR高通滤波单位脉冲响应h(n)');[H,w]=freqz(h,1);H=20*log10(abs(H));figure;plot(w/pi,H);axis([01-10010]);xlabel('\omega/\pi');ylabel('幅度/dB');grid;title('FIR高通滤波器，海明窗，N=33');结果：05101520253035-0.4-0.3-0.2-0.100.10.20.30.4nh(n)FIR高通滤波单位脉冲响应h(n)-6-00.10.20.30.40.50.60.70.80.91-100-90-80-70-60-50-40-30-20-10010/幅度/dBFIR高通滤波器，海明窗，N=332、解：（1）求数字边界频率:12120.3,0.7,0.2,0.9ccrr（2）求理想滤波器的边界频率:120.25,0.8nn（3）求理想单位脉冲响应:21d21sin()sin[()]()()nnnnnnnnhnn（4）选择窗函数：阻带最小衰减为-40dB，因此选择汉宁窗(其阻带最小衰减为-44dB)；滤波器的过渡带宽为0.3π-0.2π=0.1π，因此6.210.162,30.52NNN-7-（5）求FIR滤波器的单位脉冲响应h(n):62dsin0.8(30.5)sin[0.25(30.5)]()()()1cos()30.5(30.5)nnnhnwnhnRnn程序：clear;clc;N=62;n=0:N-1;hd=(sin(0.8*pi*(n-30.5))-sin(0.25*pi*(n-30.5)))./(pi*(n-30.5));win=hanning(N);h=win'.*hd;figure;stem(n,h);xlabel('n');ylabel('h(n)');grid;title('FIR带通滤波单位脉冲响应h(n)');[H,w]=freqz(h,1);H=20*log10(abs(H));figure;plot(w/pi,H);axis([01-10010]);xlabel('\omega/\pi');ylabel('幅度/dB');grid;title('FIR带通滤波器，汉宁窗，N=62');结果：010203040506070-0.4-0.3-0.2-0.100.10.20.30.4nh(n)FIR带通滤波单位脉冲响应h(n)-8-00.10.20.30.40.50.60.70.80.91-100-90-80-70-60-50-40-30-20-10010/幅度/dBFIR带通滤波器，汉宁窗，N=623、解：（1）求数字边界频率:12120.3,0.7,0.35,0.6ccrr（2）求理想滤波器的边界频率:120.325,0.65nn（3）求理想单位脉冲响应:12d12sin()sin()sin[()]()()1nnnnnnnnnhnn（4）选择窗函数：阻带最小衰减为-50dB，因此选择汉明窗(其阻带最小衰减为-53dB)；滤波器的过渡带宽为0.35π-0.3π=0.05π，因此-9-6.610.05133,662NNN（5）求FIR滤波器的单位脉冲响应h(n):12ddsin(66)sin(66)sin[(66)]66()()()(66)0.67566nnnnnnhnwnhnnn程序：clear;clc;N=133;n=0:N-1;hd=(sin(pi*(n-66))+sin(0.325*pi*(n-66))-sin(0.65*pi*(n-66)))./(pi*(n-66));hd(67)=0.675;win=hamming(N);