1、数字信号处理实验五FIR滤波器设计汇总实验五 FIR滤波器设计一实验内容(1)认真复习 FIR 数字滤波器的基本概念,线性相位 FIR 滤波器的条件和特点、幅度函数特点、零点位置的基本特点与性质;窗函数设计法的基本概念与方法,各种窗函数的性能和设计步骤,线性相位 FIR 低通、高通、带通和带阻滤波器的设计方法,频率采样设计法的基本概念和线性相位的实现方法。 (2)掌握几种线性相位的特点,熟悉和掌握矩形窗、三角形窗、汉宁窗、海明窗、布莱克曼窗、 凯塞窗设计 IIR 数字滤波器的方法, 熟悉和掌握频率抽样设计法的线性相位的设计方法,并对各种线性相位的频率抽样法的设计给出调整和改进。(3)熟悉利用
2、MATLAB 进行各类 FIR 数字滤波器的设计方法。 二实验内容a. 设线性相位 FIR 滤波器单位抽样响应分别为 h(n )= -4,1, -1,- 2,5,6,5, -2, -1,1, -4h(n)= -4,1, -1,- 2,5,6,6,5,- 2, -1,1,- 4h(n)= -4,1,- 1, -2,5,0,- 5,2,1,- 1,4h(n)= -4,1, -1, -2,5,6,- 6, -5,2,1, -1,4分别求出滤波器的幅度频率响应 H(),系统函数 H(z)以及零极点分布,并绘制相应的波形和分布图。在matlab中新建函数amplres,代码如下:functionA,w,
3、type,tao=amplres(h)N=length(h);tao=(N-1)/2;L=floor(N-1)/2);n=1:L+1;w=0:500*2*pi/500;if all(abs(h(n)-h(N-n+1)1e-10)A=2*h(n)*cos(N+1)/2-n)*w)-mod(N,2)*h(L+1);type=2-mod(N,2);elseif all(abs(h(n)+h(N-n+1)a1type = 1rz = -0.9807 + 0.1956i -0.9807 - 0.1956i -0.5578 + 0.8300i -0.5578 - 0.8300i 0.4052 + 1.23
4、74i 0.4052 - 1.2374i 1.2169 + 0.0000i 0.8218 + 0.0000i 0.2390 + 0.7299i 0.2390 - 0.7299ians = -0.9807 + 0.1956i -0.9807 - 0.1956i -0.5578 + 0.8300i -0.5578 - 0.8300i 0.2390 + 0.7299i 0.2390 - 0.7299i 0.8218 + 0.0000i 1.2169 + 0.0000i对第二个单位抽样响应,在matlab中新建函数a2,代码如下:h2=-4,1,-1,-2,5,6,6,5,-2,-1,1,-4;M=l
5、ength(h2);n=0:M-1;A,w,type,tao=amplres(h2);typesubplot(2,1,1),stem(n,h2);title(冲激响应h2);ylabel(h(n);xlabel(n);subplot(2,1,2),plot(w/pi,A);ylabel(A);xlabel(pi);title(幅频响应);figurerz=roots(h2)for i=1:8r(i)=1/rz(i);endrzplane(h2,1);title(h2零极点图);生成结果如下:a2type = 2rz = 1.3120 + 0.0000i -1.0000 + 0.0000i -0
6、.8868 + 0.4622i -0.8868 - 0.4622i 0.4778 + 1.1851i 0.4778 - 1.1851i -0.2957 + 0.9553i -0.2957 - 0.9553i 0.7622 + 0.0000i 0.2926 + 0.7258i 0.2926 - 0.7258ians = 0.7622 + 0.0000i -1.0000 + 0.0000i -0.8868 + 0.4622i -0.8868 - 0.4622i 0.2926 + 0.7258i 0.2926 - 0.7258i -0.2957 + 0.9553i -0.2957 - 0.9553i
7、对第三个单位抽样响应,在matlab中新建函数a3,代码如下:h3=-4,1,-1,-2,5,0,-5,2,1,-1,4;M=length(h3);n=0:M-1;A,w,type,tao=amplres(h3);typesubplot(2,1,1),stem(n,h3);title(冲激响应h3);ylabel(h(n);xlabel(n);subplot(2,1,2),plot(w/pi,A);ylabel(A);xlabel(pi);title(幅频响应);figurerz=roots(h3)for i=1:8r(i)=1/rz(i);endrzplane(h3,1);title(h3零
8、极点图);生成结果如下: a3type = 3rz = -1.0000 + 0.0000i -0.8732 + 0.4874i -0.8732 - 0.4874i 0.1010 + 1.2041i 0.1010 - 1.2041i 1.0000 + 0.0000i 0.8280 + 0.5607i 0.8280 - 0.5607i 0.0692 + 0.8247i 0.0692 - 0.8247ians = -1.0000 + 0.0000i -0.8732 + 0.4874i -0.8732 - 0.4874i 0.0692 + 0.8247i 0.0692 - 0.8247i 1.0000
9、 + 0.0000i 0.8280 + 0.5607i 0.8280 - 0.5607i对第四个单位抽样响应,在matlab中新建函数a4,代码如下:h4=-4,1,-1,-2,5,6,-6,-5,2,1,-1,4;M=length(h4);n=0:M-1;A,w,type,tao=amplres(h4);typesubplot(2,1,1),stem(n,h4);title(冲激响应h4);ylabel(h(n);xlabel(n);subplot(2,1,2),plot(w/pi,A);ylabel(A);xlabel(pi);title(幅频响应);figurerz=roots(h4)f
10、or i=1:8r(i)=1/rz(i);endrzplane(h4,1);title(h4零极点图);生成结果如下: a4type = 4rz = 0.2631 + 1.3394i 0.2631 - 1.3394i -0.9505 + 0.3108i -0.9505 - 0.3108i -0.7309 + 0.6825i -0.7309 - 0.6825i 0.9021 + 0.4315i 0.9021 - 0.4315i 1.0000 + 0.0000i 0.1412 + 0.7189i 0.1412 - 0.7189ians = 0.1412 + 0.7189i 0.1412 - 0.7
11、189i -0.9505 + 0.3108i -0.9505 - 0.3108i -0.7309 + 0.6825i -0.7309 - 0.6825i 0.9021 + 0.4315i 0.9021 - 0.4315ib. 设计 FIR 数字低通滤波器,技术指标为:p=0.2,st=0.3,1=0.25dB,2=50dB。 (1) 通过技术指标,选择一种窗函数进行设计; (2) 求滤波器的单位抽样响应、频率响应,并绘制波形。 (3) 选择凯塞窗函数设计该滤波器,并绘制相应的波形图。在matlab中新建函数ideal_lp,代码如下:function hd=ideal_lp(wc,M);% I
12、deal 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;n = 0:1:(M-1);m = n - alpha +eps; % add smallest number to avoi divided by zerohd = sin(wc*m)./(pi*m);end在matlab中新建
13、函数freqz_m,代码如下:function db,mag,pha,grd,w = freqz_m(b,a);H,w = freqz(b,a,1000,whole); H = (H(1:1:501); w = (w(1:1:501); mag = abs(H); db = 20*log10(mag+eps)/max(mag); pha = angle(H); grd = grpdelay(b,a,w);在matlab中新建函数b1,代码如下:%用Hamming窗函数设计FIR数字滤波器wp=0.2*pi;ws=0.3*piN=61n=0:1:N-1wc=(ws+wp)/2;%理想低通滤波器h
14、d=ideal_lp(wc,N);%理想低通的冲激响应w_ham=(hamming(N);h=hd.*w_ham;%FIR滤波器冲激响应db,mag,pha,grd,w=freqz_m(h,1);delta_w=2*pi/1000;Rp=-(min(db(1:1:wp/delta_w+1)%实际的通带衰减As=-round(max(db(ws/delta_w+1:1:501)%实际的最小阻带衰减subplot(221);stem(n,hd);title(理想冲激响应)axis(0 N-1 -0.1 0.3);xlabel(n);ylabel(hd(n)subplot(222);stem(n,w
15、_ham);title(hamming窗)axis(0 N-1 0 1.1);xlabel(n);ylabel(w(n);subplot(223);stem(n,h);title(实际冲激响应)axis(0 N-1 -0.1 0.3);xlabel(n);ylabel(h(n)subplot(224);plot(w/pi,db);axis(0 0.8 -100 0);xlabel(以PI为单位的频率);ylabel(对数幅度/db);实验结果如下:在matlab中新建函数b2,代码如下:%用Kaiser窗函数设计FIR数字滤波器wp=0.2*pi;ws=0.3*pi;As=50tr_width
16、=ws-wpN=ceil(As-7.95)/(14.36*tr_width/(2*pi)+1)+1n=0:1:N-1beta=0.1102*(As-8.7)wc=(wp+ws)/2%理想低通的截止频率hd=ideal_lp(wc,N)w_kai=(kaiser(N,beta)h=hd.*w_kaidb,mag,pha,grd,w=freqz_m(h,1)delta_w=2*pi/1000Rp=-(min(db(1:1:wp/delta_w+1)%实际的通带衰减As=-round(max(db(ws/delta_w+1:1:501)%实际的最小阻带衰减subplot(211);plot(w/pi
17、,db);title(凯森窗幅度响应(dB);gridaxis(0 0.5 -100 0)ylabel(对数幅度/db);xlabel(以/pi为单位的频率)subplot(212);plot(w/pi,pha);title(相位响应);gridaxis(0 0.5 -4 4)ylabel(相位);xlabel(以/pi为单位的频率)实验结果如下:c. 设计 FIR 数字带通滤波器,技术指标为: 下阻带边缘:st1=0.2,s1=60dB,下通带边缘:p1=0.35,p1=1dB; 上通带边缘:p2=0.65,p1=1dB,上阻带边缘:st2=0.8,s2=60dB; (1) 通过技术指标,选
18、择一种窗函数进行设计; (2) 求滤波器的单位抽样响应、频率响应,并绘制波形。在matlab中新建函数c,代码如下:ws1=0.2*pi; wp1=0.35*pi;ws2=0.8*pi;wp2=0.65*pi;Ap=60;Rp=1;tr_width=min(wp1-ws1),(ws2-wp2);M=ceil(11*pi/tr_width);n=0:1:M-1;wc1=(ws1+wp1)/2;wc2=(wp2+ws2)/2;hd=ideal_lp(wc2,M)-ideal_lp(wc1,M);w_bla=(blackman(M);h=hd.*w_bla;H,W=freqz(h,1);subplo
19、t(2,2,1);stem(n,hd);title(理想脉冲抽样);subplot(2,2,2);stem(n,w_bla);title(布莱克曼窗);subplot(2,2,3);stem(n,h);title(实际脉冲抽样);subplot(2,2,4);plot(W/pi,20*log10(abs(H);title(幅度响应(db));生成结果如下:d. 设计 FIR 数字带通滤波器,技术指标为: 下阻带边缘:st1=0.2,s1=60dB,下通带边缘:p1=0.4,p1=1dB; 上通带边缘:p2=0.6,p1=1dB,上阻带边缘:st2=0.8,s2=60dB; (1) 通过技术指标
20、,选择一种窗函数进行设计; (2) 求滤波器的单位抽样响应、频率响应,并绘制波形。在matlab中新建函数d,代码如下:ws1=0.2*pi; wp1=0.4*pi;ws2=0.8*pi;wp2=0.6*pi;Ap=60;Rp=1;tr_width=min(wp1-ws1),(ws2-wp2);M=ceil(11*pi/tr_width);n=0:1:M-1;wc1=(ws1+wp1)/2;wc2=(wp2+ws2)/2;hd=ideal_lp(wc2,M)-ideal_lp(wc1,M);w_bla=( hamming(M);h=hd.*w_bla;H,W=freqz(h,1);subplo
21、t(2,2,1);stem(n,hd);title(理想脉冲抽样);subplot(2,2,2);stem(n,w_bla);title(海明窗);subplot(2,2,3);stem(n,h);title(实际脉冲抽样);subplot(2,2,4);plot(W/pi,20*log10(abs(H);title(幅度响应(db));生成结果如下:e. 设计 FIR 数字带通滤波器,技术指标为: 下阻带边缘:st1=0.2,s1=20dB,下通带边缘:p1=0.4,p1=1dB; 上通带边缘:p2=0.6,p1=1dB,上阻带边缘:st2=0.8,s2=20dB; (1) 通过技术指标,选
22、择一种窗函数进行设计; (2) 求滤波器的单位抽样响应、频率响应,并绘制波形。在matlab中新建函数e,代码如下: ws1=0.2*pi;wp1=0.4*pi;ws2=0.8*pi;wp2=0.6*pi;Ap=20;Rp=1;tr_width=min(wp1-ws1),(ws2-wp2);M=ceil(11*pi/tr_width);n=0:1:M-1;wc1=(ws1+wp1)/2;wc2=(wp2+ws2)/2;hd=ideal_lp(wc2,M)-ideal_lp(wc1,M);w_bla=( boxcar(M);h=hd.*w_bla;H,W=freqz(h,1);subplot(2
23、,2,1);stem(n,hd);title(理想脉冲抽样);subplot(2,2,2);stem(n,w_bla);title(矩形窗);subplot(2,2,3);stem(n,h);title(实际脉冲抽样);subplot(2,2,4);plot(W/pi,20*log10(abs(H);title(幅度响应(db));生成结果如下:f. 设计 FIR 数字高通滤波器,技术指标为:通带截止频率为p=15/27,阻带截止频率为st=11/27,通带最大衰减为1=2.5dB,阻带最小衰减为2=55dB。 (1) 通过技术指标,选择一种窗函数进行设计; (2) 求滤波器的单位抽样响应、频
24、率响应,并绘制波形。系统一:在matlab中新建函数f,代码如下: As=55;ws=11*pi/27;wp=15*pi/27; tr_width=wp-ws;%计算过渡带 M=ceil(As-7.95)*2*pi/(14.36*tr_width)+1)+1; %按凯泽窗计算滤波器的长度disp(滤波器的长度,num2str(M);beta=0.1102*(As-8.7); %计算凯泽窗的beta值n=0:1:M-1;disp(线性相位滤波器,num2str(beta);w_kai=(kaiser(M,beta);%求凯泽窗函数wc=(ws+wp)/2;hd=ideal_lp(pi,M)-id
25、eal_lp(wc,M); %求理想脉冲响应h=hd.*w_kai; db,mag,pha,grd,w=freqz_m(h,1);delta_w=2*pi/1000;Rp=-(min(db(wp/delta_w+1:1:501);disp(实际通带波动为,num2str(Rp);As=-round(max(db(1:1:ws/delta_w+1);disp(最小阻带衰减为,num2str(As);subplot(2,2,1);stem(n,hd);title(理想脉冲响应);axis(0 M-1 -0.4 0.8);ylabel(hd(n);subplot(2,2,2);stem(n,w_ka
26、i);title(凯泽窗); axis(0 M-1 0 1.1);ylabel(wd(n);subplot(2,2,3);stem(n,h);title(实际脉冲响应);axis(0 M-1 -0.4 0.8);xlabel(n);ylabel(h(n);subplot(2,2,4);plot(w/pi,db);title(幅度响应/dB);axis(0 1 -100 10);grid;xlabel(以pi为单位的频率);ylabel(分贝数/dB); 生成结果如下:滤波器的长度47线性相位滤波器5.1023实际通带波动为0.025017最小阻带衰减为?55滤波器的长度47线性相位滤波器5.1
27、023实际通带波动为0.025017最小阻带衰减为55g. 设计 FIR 数字高通滤波器,技术指标为:通带截止频率为p=0.6,阻带截止频率为st=0.4,通带最大衰减为1=0.25dB,阻带最小衰减为2=40dB。 (1) 通过技术指标,选择一种窗函数进行设计; (2) 求滤波器的单位抽样响应、频率响应,并绘制波形。 在matlab中新建函数c,代码如下:function db, mag, pha, w=freqz_m2(b,a);%Modified version of freqz subroutineH,w=freqz(b,a,1000,whole);H=(H(1:501); w=(w(1:501);mag=abs(H);db=20*log10(mag+eps)/max(mag);pha=angle(H);在matlab中新建函数g,代码如下:Wp=0.6*pi;Ws=0.4*pi;tr_width=Wp-Ws;M=ceil(6.2*pi/tr_width);n=0:1:M-1;Wc=(Ws+Wp)/
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1