南昌航空大学电子信息工程学院
数字信号处理 课程实验报告
实验名称:IIR 数字滤波器设计与信号滤波
实验时间: 2018 年 5 月 28 日
指导教师:
班 级 :
150415
学 号 :
姓 名 :
成 绩 :
一、实验目的
1、熟悉用双线性变换法设计IIR数字滤波器的原理与方法;
2、熟悉用脉冲响应不变法设计IIR数字滤波器的原理与方法。
二、实验内容
1、已知低通滤波器的指标为:
通带边缘频率:0.4π,Ap=0.5dB;
阻带边缘频率:0.6π,As=50dB;
a、采用脉冲响应不变法设计巴特沃斯,T=1.画出幅度响应和数字滤波器的脉冲响应h(n);
b、采用脉冲响应不变法设计巴特沃斯,T=1.画出幅度响应和数字滤波器的脉冲响应h(n)。
2、用双线性变换法设计低通滤波器,满足技术指标
wp=0.2π,Ap=0.25dB;
ws=0.4π,As=50dB,
并对方波信号进行滤波,画出滤波前后的波形图并进行简要分析。
3、设计一个数字高通滤波器H(z),它用在结构
xa(t) A/D
H(z)
D/A
ya(t)
中,满足下列要求:
a、采样速率为10kHZ;
b、阻带边缘频率为1.5kHZ,衰减为40dB;
c、通带边缘频率为2kHZ,衰减为3dB;
d、单调的通带和阻带。
4、设计一个带阻滤波器,要求通带上下截止频率为0.8π,0.2π,通带内衰减不大于1dB,
阻带起始频率为0.7π,0.4π,阻带内衰减不小于30dB。
设计巴特沃斯带阻滤波器并画出该数字高通滤波器的幅度响应和脉冲响应。
三、实验原理
利用双线性变换设计 IIR 滤波器(只介绍巴特沃斯数字低通滤波器的设计),首先要设
计出满足指标要求的模拟滤波器的传递函数 Ha(s)然后由 Ha(s)通过双线性变换可得所要设
计的 IIR 滤波器的系统函数 H(z)。
四、实验步骤、相应的实验结果、相应的结果分析
1. clear;close all;
wp=0.4*pi;
ws=0.6*pi;
Ap=0.5;
As=50;
T=1;
Fs=1/T;
OmegaP=wp/T;
OmegaS=ws/T;
[cs,ds]=afd_butt(OmegaP,OmegaS,Ap,As);
[b,a]=impinvar(cs,ds,Fs);
[h,w]=freqz(b,a);
subplot(2,2,1);plot(w/pi,abs(h));
title('幅度响应');
grid;
subplot(2,2,2);plot(w/pi,angle(h));
title('相位响应');
grid;
subplot(2,2,3);
plot(w/pi,20*log(abs(h)));
title('幅度响应dB');grid;
n=[0:1:59];
imp=[1;zeros(59,1)];
y=filter(b,a,imp);
subplot(2,2,4);plot(n,y);
title('脉冲响应');
grid;
实验分析:;用脉冲响应不变法设计低通巴特沃斯滤波器时,因为 w=Ωt,w 与Ω呈线性关系,
所以其相位响应图是呈线性的,如上图的“相位响应”图所示;因为其设计的是低通滤波器,
所以会把其高频的波滤掉,而留下低频的,从而在变化的时候表现的平缓些,而不是特别的
陡峭,如上图的“脉冲响应”图所示,从而实现低通滤波的功能。
clear;close all;
wp=0.4*pi;
ws=0.6*pi;
0.6*pi
Ap=0.5;
As=50;
T=1;
Fs=1/T;
OmegaP=(2/T)*tan(wp/2);
OmegaS=(2/T)*tan(ws/2);
[cs,ds]=afd_butt(OmegaP,OmegaS,Ap,As);
[b,a]=bilinear(cs,ds,Fs);
[h,w]=freqz(b,a);
subplot(2,2,1);plot(w/pi,abs(h));
title('幅度响应');
grid;
subplot(2,2,2);plot(w/pi,angle(h));
title('相位响应');grid;
subplot(2,2,3);plot(w/pi,20*log(abs(h)));
title('幅度响应dB');
axis([0,1,-80,5]);grid;
n=[0:1:59];
imp=[1;zeros(59,1)];
y=filter(b,a,imp);
subplot(2,2,4);plot(n,y);
title('脉冲响应');
grid;
实验分析:用双线性变换法实现低通巴特沃斯滤波器时,Ω=(2/T)tan(w/2),可以知道 w 与Ω
是呈非线性关系的,所以其图形如上图的“相位响应”图所示,呈现的是不是直线,而是曲
线。由于其设计的是低通滤波器,所以把其高频的部分滤掉了,留下低频的部分,从而在变
化的时候表现的平缓些,而不是特别的陡峭,如上图的“脉冲响应”图所示,从而实现低通
滤波的功能。
2. wp=0.2*pi;
ws=0.4*pi;
Ap=0.25;
As=50;
T=1;
Fs=1/T;
OmegaP=(2/T)*tan(wp/2);
OmegaS=(2/T)*tan(ws/2);
N=ceil((log10((10^(Ap/10)-1)/(10^(As/10)-1)))/
(2*log10(OmegaP/OmegaS)));
fprintf('\n***Butterworth Filter Order =%2.0f\n',N);
OmegaC=OmegaP/((10^(Ap/10)-1)^(1/(2*N)));
wn=2*atan((OmegaC*T)/2);
wn=wn/pi;
[b,a]=butter(N,wn);
[h,w]=freqz(b,a);
subplot(2,2,1);
plot(w/pi,abs(h));
title('幅度响应');grid ;
subplot(2,2,3);
plot(w/pi,20*log10(abs(h)));
title('幅度响应dB');
axis([0,1,-80,5]);grid ;
n=[0:1:59];
imp=[1;zeros(59,1)];
y=filter(b,a,imp);
subplot(2,2,4);plot(n,y);
title('脉冲响应');grid ;
n=[0:100];
t=0.001*n;
x=2*sin(2*pi*400*t)+5*sin(2*pi*50*t);
x1=ones(1,101);
y1=filter(b,a,x1);
figure(2);
subplot(2,1,1);
plot(n,x1,'p-');
subplot(2,1,2);
plot(n,y1,'p-');
实验分析:用双线性法设计这个低通滤波器,由实验程序中可知,其 fs=1000HZ,当输入的
是方波函数时,把其高频的部分滤掉,留下低频部分,所以其滤波后的图形在变化时表现的
稍微平缓些,而不是特别的陡峭,如上图的“脉冲响应”图所示,从而实现低通滤波的功能。
3、f=10;fp=2;fs=1.5;
wp=2*pi*fp*(1/f);
ws=2*pi*fs*(1/f);
Ap=3;
As=40;
[N,wn]=buttord(wp/pi,ws/pi,Ap,As);
[b,a]=butter(N,wn,'high');
[h,w]=freqz(b,a);
subplot(2,2,1);plot(w/pi,abs(h));
title('幅度响应');grid ;
subplot(2,2,3);plot(w/pi,20*log10(abs(h)));
title('幅度响应dB');grid ;
n=[0:1:59];imp=[1;zeros(59,1)];
y=filter(b,a,imp);
subplot(2,2,4);plot(n,y);
title('脉冲响应');grid ;
n=[0:1000];t=0.0001*n;
x=2*sin(2*pi*4000*t)+5*sin(2*pi*50*t);
y1=filter(b,a,x);
figure(2);
subplot(2,1,1);
plot(n,x,'p-');
subplot(2,1,2);
plot(n,y1,'p-');
实验分析:这里设计的是一个高通数字滤波器,由实验程序中可知,其 fs=10000HZ,而输
入信号 x=x1+x2,(其中 x1=2*sin(2*pi*4000*t),x2=5*sin(2*pi*50*t))是由两种
不同频率的 sin 函数混合而成的,fs/2=5000HZ,f1=4000HZ,相对来说是高频,可以通
过,f2=50HZ,相对来说是低频,会被滤波器滤掉,所以当 x 经过高通滤波器滤波后,只
会留下 x1 的频谱,x2 会被滤掉;从上图中可以看出,没经滤波前,其波形的幅度变化范
围为 7 到-7,x1、x2 两种频率的波是混合在一起的,经滤波后其波形的幅度范围变,2 到
-2,只剩下了 x1 波形,即 x2 被高通滤波器滤掉了,实现了高通滤波器的功能。
附录:
function [b,a] =afd_butt(Omegp,Omegs,Ap,As);
% b = Ha(s)分子多项式的系数
% a = Ha(s)分母多项式的系数
% Omegp = 通带截止频率 单位为 rad/s; Ωp > 0
% Omegs =阻带起始频率单位为 rad/s; Ωs > Ωp > 0
% Ap = Passband ripple in +dB; (Ap > 0)
% As = Stopband attenuation in +dB; (As > 0)
if Omegp<= 0
error('Passband edge must be larger than 0')
end
if Omegs<=Omegp
end
if (Ap <= 0) | (As < 0)
error('Stopband edge must be larger than Passband edge')
error('PB ripple and/or SB attenuation must be larger than 0')
end
N=ceil((log10((10^(Ap/10)-1)/(10^(As/10)-1)))/(2*log10(Omegp/O
megs)) );
fprintf('\n*** Butterworth Filter Order = %2.0f \n',N)
OmegaC = Omegp/((10^(Ap/10)-1)^(1/(2*N)));
[b,a]=u_buttap(N,OmegaC);
end
function [b,a] = u_buttap(N,Omegac);
%未归一化的Butterworth模拟低通滤波器原型
% [b,a] = u_buttap(N,Omegac);
% b = Ha(s)分子多项式的系数
% a = Ha(s)分母多项式的系数
% N =滤波器的阶数 r
% Omegac = 截止频率 单位 rad/s
[z,p,k] = buttap(N);
p = p*Omegac;
k = k*Omegac^N;
B = real(poly(z));
b = k*B;
b0= k;
a = real(poly(p));
end
五、实验小结