IIR数字滤波器课程设计.doc

《IIR数字滤波器课程设计.doc》由会员分享，可在线阅读，更多相关《IIR数字滤波器课程设计.doc（75页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。

1、Four short words sum up what has lifted most successful individuals above the crowd: a little bit more.-author-dateIIR数字滤波器课程设计课 程 设 计 报 告数字信号处理课 程 设 计 报 告基于MATLAB的IIR数字滤波器设计专业班级: 电信工程1302班 学 号: 311308000626 学生姓名: 王 海 龙 指导教师: 王 科 平 2016年 7 月目 录摘要.3一、 课程设计任务及要求.4 1.本次设计的目的.4 2.本次设计的要求.4二、 课程设计原理.4 1.

2、脉冲响应不变法原理.4 2.双向性变换法原理.5三、IIR数字滤波器设计内容 .5 1.总体方法分析.5 2.脉冲相应不变法.6 3.双线性变换法.7四、IIR数字滤波器设计过程.9 1.设计步骤.9 2.程序流程框图.11 3.MATLAB程序.11 4.调试分析过程描述.19 5结果分析.19五、结论.22六、参考文献.23 摘 要在当今社会，数字信号处理技术飞速发展，它不但自成一门学科，更是以不同的方式影响和渗透到其他学科的研究中，它变得与我们的生活联系越来越紧密，不断改变着我们的生产生活方式，因此受到人们越来越多的关注。数字滤波器是对数字信号实现滤波的线性时不变系统。数字滤波实质上是一

3、种运算过程，实现对信号的运算处理。输入数字信号（数字序列）通过特定的运算转变为输出的数字序列。描述离散系统输出与输入关系的卷积和差分方程只是给数字信号滤波器提供运算规则，使其按照这个规则完成对输入数据的处理。时域离散系统的频域特性:Y(ejw)=X(ejw)H(ejw),其中、分别是数字滤波器的输出序列和输入序列的频域特性（或称为频谱特性）,H(ejw)是数字滤波器的单位取样响应的频谱，又称为数字滤波器的频域响应。数字滤波器是具有一定传输选择特性的数字信号处理装置,其输入、输出均为数字信号,实质上是一个由有限精度算法实现的线性时不变离散系统。IIR 数字滤波器的特征是，具有无限持续时间冲激响应

4、，需要用递归模型来实现,其差分方程为：y(n)=i x(n-i)+i y(n-i)系统函数为：H（z）=( Z-r )/( 1+k Z-k) 设计IIR滤波器的任务就是寻求一个物理上可实现的系统函数H(z)，使其频率响应H(z)满足所希望得到的频域指标。本次课程设计分别用脉冲响应不变法、双向性变换法设计IIR低通、高通、带通、带阻滤波器滤波器。并在MATLAB环境下实现了IIR 数字滤波器的设计和仿真。其主要内容概括为：首先对滤波器的原理和设计进行了介绍；接着描述了IIR数字滤波器的基本概念，其中包括系统的描述、系统的传递函数、系统的模型；接着简单介绍MATLAB，并对数字滤波器在MATLAB

5、环境下如何实现进行了介绍；重点描述了IIR数字滤波器的设计过程，最后对IIR滤波器进行仿真。关键词：数字滤波器 频域特性 脉冲响应 双向性变换法 MATLAB一、课程设计任务及要求1.本次设计的目的 1）学会MATLAB的使用，掌握MATLAB的程序设计方法； 2）掌握数字信号处理的基本概念、基本理论和基本方法； 3）掌握MATLAB设计IIR滤波器； 4）学会用MATLAB对信号进行分析和处理。2.本次设计的要求 1）分别用脉冲响应不变法、双向性变换法设计IIR低通、高通、带通、带阻滤波器滤波器；2）分别画出其幅频特性、相频特性图；3）IIR滤波器的各项指标:低通：通带截止频率wc=2s r

6、adk /2,阻带截止频率为8KHZ=Wp，通带衰减pR小于3dB，阻带衰减大于 15dB,采样频率20000Hz;高通：通带截止频率为2.5KHZ,通带衰减不大于2dB,阻带上限截止频率为1.5KHZ，阻带衰减不小于15dB;带通：中心频率为p0=0.5,通带截止频率p1=0.4,p2=0.6;通带最大衰减p=3dB;阻带最小衰减s=15dB;阻带截止频率s2=0.7;带阻：抽样频率为10KHZ,在-2dB衰减处边带频率是1.5KHZ,4KHZ,在-13dB处边带频率为2KHZ和3KHZ。二、课程设计原理1.脉冲响应不变法原理 脉冲响应不变法是实现模拟滤波器数字化的一种直观而常用的方法，它特

7、别适合于对滤波器的时域特性有一定要求的场合。 脉冲响应不变法是从滤波器的脉冲响应出发，使数字滤波器的单位脉冲响应序列h(n) 模仿模拟滤波器的冲击响应ha(t), 使h(n)正好等于ha(t)的采样值，即h(n)=ha(nT)T为采样周期。 如以Ha(s)及H(z)分别表示ha(t)的拉氏变换及h(n)的z变换，即 Ha(s)=Lha(t) H(z)=Zh(n)则根据采样序列z变换与模拟信号拉氏变换的关系，可知：采用脉冲响应不变法将模拟滤波器变换为数字滤波器时，它所完成的S平面到Z平面的变换，正是以前讨论的拉氏变换到Z变换的标准变换关系，即首先对Ha(s)作周期延拓，然后再经过z=e的映射关系

8、映射到Z平面上。脉冲响应不变法映射关系见图2。2.双向性变换法： 脉冲响应不变法的主要缺点是频谱交叠产生的混淆，这是从S平面到Z平面的标准变换z=e的多值对应关系导致的,为了克服这一缺点，设想变换分为两步。1) 将整个S平面压缩到S1平面的一条横带里。2) 通过标准变换关系将此横带变换到整个Z平面上去。由此建立S平面与Z平面一一对应的单值关系，消除多值性，也就消除了混淆现象。 图1 双线性换法映射关系图 双线性换法的主要优点是S平面与Z平面一单值对应，S平面的虚轴(整个j)对应于Z平面单位圆的一周，S平面的=0处对应于Z平面的=0处，对应即数字滤波器的频率响应终 止于折迭频率处，所以双线性变换

9、不存在混迭效应。三、 IIR数字滤波器设计内容1.总体方法分析IIR数字滤波器是一种离散时间系统，其系统函数为：假设MN，当MN时,系统函数可以看作一个IIR的子系统和一个(M-N)的FIR子系统的级联。IIR数字滤波器的设计实际上是求解滤波器的系数和 ，它是数学上的一种逼近问题，即在规定意义上（通常采用最小均方误差准则）去逼近系统的特性。如果在S平面上去逼近，就得到模拟滤波器；如果在z平面上去逼近，就得到数字滤波器。2.脉冲相应不变法脉冲响应不变法是从滤波器的脉冲响应出发，使数字滤波器的单位脉冲响应序列h(n)模仿模拟滤波器的冲激响应ha(t)，即将ha(t)进行等间隔采样，使h(n)正好等

10、于ha(t)的采样值，满足:h(n)=ha(nT) 式中,T是采样周期。 如果令Ha(s)是ha(t)的拉普拉斯变换，H(z)为h(n)的Z变换，利用采样序列的Z变换与模拟信号的拉普拉斯变换的关系得 (1-1) 则可看出，脉冲响应不变法将模拟滤波器的S平面变换成数字滤波器的Z平面，这个从s到z的变换z=esT是从S平面变换到Z平面的标准变换关系式。 图2 脉冲响应不变法的映射关系 由（1-1）式，数字滤波器的频率响应和模拟滤波器的频率响应间的关系为 (1-2)这就是说，数字滤波器的频率响应是模拟滤波器频率响应的周期延拓。正如采样定理所讨论的，只有当模拟滤波器的频率响应是限带的，且带限于折叠频率

11、以内时，即 (1-3)才能使数字滤波器的频率响应在折叠频率以内重现模拟滤波器的频率响应，而不产生混叠失真，即 |w| （1-4） 但是，任何一个实际的模拟滤波器频率响应都不是严格限带的，变换后就会产生周期延拓分量的频谱交叠，即产生频率响应的混叠失真，如图7-4所示。这时数字滤波器的频响就不同于原模拟滤波器的频响，而带有一定的失真。当模拟滤波器的频率响应在折叠频率以上处衰减越大、越快时，变换后频率响应混叠失真就越小。这时，采用脉冲响应不变法设计的数字滤波器才能得到良好的效果。 图3 脉冲响应不变法中的频响混叠现象 对某一模拟滤波器的单位冲激响应ha(t)进行采样，采样频率为fs，若使fs增加，即

12、令采样时间间隔（T=1/fs）减小，则系统频率响应各周期延拓分量之间相距更远，因而可减小频率响应的混叠效应。3.双线性变换法脉冲响应不变法的主要缺点是产生频率响应的混叠失真。这是因为从S平面到平面是多值的映射关系所造成的。为了克服这一缺点，可以采用非线性频率压缩方法，将整个频率轴上的频率范围压缩到-/T/T之间，再用z=esT转换到Z平面上。这样就使S平面与Z平面建立了一一对应的单值关系，消除了多值变换性就消除了频谱混叠现象，映射关系如图4所示。图4 双线性变换的映射关系为了将S平面的整个虚轴j压缩到S1平面j1轴上的-/T到/T段上，可以通过以下的正切变换实现（1-5）式中,T仍是采样间隔。

13、当1由-/T经过0变化到/T时，由-经过0变化到+，也即映射了整个j轴。将式（1-5）写成将此关系解析延拓到整个S平面和S1平面，令j=s，j1=s1，则得再将S1平面通过以下标准变换关系映射到Z平面 z=es1T从而得到S平面和Z平面的单值映射关系为： （1-6） （1-7） 式（1-6）与式（1-7）是S平面与Z平面之间的单值映射关系，这种变换都是两个线性函数之比，因此称为双线性变换式（1-5）与式（1-6）的双线性变换符合映射变换应满足的两点要求。首先,把z=ej，可得 （1-8）即S平面的虚轴映射到Z平面的单位圆。其次，将s=+j代入式（1-8），得因此由此看出，当0时，|z|0时，|

14、z|1。也就是说，S平面的左半平面映射到Z平面的单位圆内，S平面的右半平面映射到Z平面的单位圆外，S平面的虚轴映射到Z平面的单位圆上。因此，稳定的模拟滤波器经双线性变换后所得的数字滤波器也一定是稳定的。四、IIR数字滤波器设计过程根据以上IIR数字滤波器设计方法，下面运用双线性变换法基于MATLAB设计一个IIR带通滤波器，其中带通的中心频率为p0=0.5，;通带截止频率 图5 部分滤波器设计指标图示 wp1=0.4,p2=0.6;通带最大衰减p=3dB;阻带最小衰减s=15dB;阻带截止频率s2=0.71. 设计步骤 （以带通为例）(1)根据任务,确定性能指标:在设计带通滤波器之前,首先根据

15、工程实际的需要确定滤波器的技术指标，带通滤波器的阻带边界频率关于中心频率p0几何对称，因此ws1=wp0- (ws2-wp0)=0.3通带截止频率wc1=0.4,wc2=0.6;阻带截止频率wr1=0.3，wr2=0.7；阻带最小衰减s=3dB和通带最大衰减p=15dB;(2)用=2/T*tan(w/2)对带通数字滤波器H(z)的数字边界频率预畸变，得到带通模拟滤波器H(s)的边界频率主要是通带截止频率p1,p2;阻带截止频率s1，s2的转换。为了计算简便，对双线性变换法一般T=2s通带截止频率wc1=(2/T)*tan(wp1/2)=tan(0.4/2)=0.7265wc2=(2/T)*ta

16、n(wp2/2)=tan(0.6/2)=1.3764阻带截止频率wr1=(2/T)*tan(ws1/2)=tan(0.3/2)=0.5095wr2=(2/T)*tan(ws2/2)=tan(0.7/2)=1.9626阻带最小衰减s=3dB和通带最大衰减p=15dB;(3)运用低通到带通频率变换公式=(2)-(02)/(B*)将模拟带通滤波器指标转换为模拟低通滤波器指标。B=wc2-wc1=0.6499normwr1=(wr12)-(w02)/(B*wr1)=2.236normwr2=(wr22)-(w02)/(B*wr2)=2.236normwc1=(wc12)-(w02)/(B*wc1)=1

17、normwc2=(wc22)-(w02)/(B*wc2)=1得出，normwc=1，normwr=2.236模拟低通滤波器指标：normwc=1，normwr=2.236，p=3dB，s=15dB(4)设计模拟低通原型滤波器。用模拟低通滤波器设计方法得到模拟低通滤波器的传输函数Ha(s);借助巴特沃斯(Butterworth)滤波器、切比雪夫(Chebyshev)滤波器、椭圆(Cauer)滤波器、贝塞尔(Bessel)滤波器等。(5)调用lp2bp函数将模拟低通滤波器转化为模拟带通滤波器。(6)利用双线性变换法将模拟带通滤波器Ha(s)转换成数字带通滤波器H(z).图6 四种数字滤波器2.程序

18、流程框图开始读入数字滤波器技术指标将指标转换成归一化模拟低通滤波器的指标设计归一化的模拟低通滤波器阶数N和3db截止频率模拟域频率变换，将G(P)变换成模拟带通滤波器H(s)用双线性变换法将H(s)转换成数字带通滤波器H(z)输入信号后显示相关结果结束图7 程序流程图3.MATLAB程序MATLAB 程序如下:带通滤波器：clearwp0=0.5*pi;wp1=0.4*pi;wp2=0.6*pi;Ap=3;ws2=0.7*pi;As=15;T=2; %数字带通滤波器技术指标ws1=wp0-(ws2-wp0); %计算带通滤波器的阻带下截止频率wc1=(2/T)*tan(wp1/2);wc2=(

19、2/T)*tan(wp2/2);wr1=(2/T)*tan(ws1/2);wr2=(2/T)*tan(ws2/2);w0=(2/T)*tan(wp0/2); %频率预畸变B=wc2-wc1; %带通滤波器的通带宽度 normwr1=(wr12)-(w02)/(B*wr1);normwr2=(wr22)-(w02)/(B*wr2);normwc1=(wc12)-(w02)/(B*wc1);normwc2=(wc22)-(w02)/(B*wc2); %带通到低通的频率变换if abs(normwr1)abs(normwr2) normwr=abs(normwr2)else normwr=abs(n

20、ormwr1)end normwc=1; %将指标转换成归一化模拟低通滤波器的指标N=buttord(normwc,normwr,Ap,As,s); %设计归一化的模拟低通滤波器阶数N和3db截止频率bLP,aLP=butter(N,normwc,s); %计算相应的模拟滤波器系统函数G(p)bBP,aBP=lp2bp(bLP,aLP,w0,B); %模拟域频率变换，将G(P)变换成模拟带通滤波器H(s) b,a=bilinear(bBP,aBP,0.5); %用双线性变换法将H(s)转换成数字带通滤波器H(z)w=linspace (0,2*pi,500);h=freqz(b,a,w);su

21、bplot(2,1,2);plot(w,abs(h);grid onxlabel(w(rad)ylabel(|H(jw)|)title(频谱函数)subplot(2,2,1);plot(w,20*log10(abs(h);axis(0,2*pi,-120,20);grid onxlabel(w(rad)ylabel(20*lg|H(jw)|(db)title(20*lg|H(jw)|-w)带阻滤波器：%周期方波时域f1=64;fs=10000; t=0:1/fs:0.1;x1=square(2*pi*f1*t);subplot(321);plot(t,x1);axis(0,0.1,-1.3,1

22、.3);xlabel(t/s);ylabel(x(t);title(时域谱)%周期方波频域L=1250;X1=fft(x1,L);ws=2*pi*fs;w=(0:L-1)*ws/L;subplot(322);stem(0:L-1),abs(X1);xlabel(w/rad);ylabel(X);title(幅度谱);%采用BW和双线性变换法Wp1=0.1*pi,Wp2=0.4*pi;Ws1=0.2*pi,Ws2=0.3*pi;Ap=1;As=50;wp1=2*fs*tan(Wp1/2),wp2=2*fs*tan(Wp2/2);ws1=2*fs*tan(Ws1/2),ws2=2*fs*tan(W

23、s2/2);B=ws2-ws1,w0=sqrt(ws1*ws2);ws=1;wp1=(B*wp1)/(-wp12+w02),wp2=(B*wp2)/(-wp22+w02);wp=max(abs(wp1),abs(wp2);%wp=wp1wp2,ws=ws1ws2;N,wc=buttord(wp,ws,Ap,As,s);num,den=butter(N,wc,s);numt,dent=lp2bs(num,den,w0,B);numd,dend=bilinear(numt,dent,fs);w=linspace(0,pi,1024);H=freqz(numd,dend,w);subplot(324

24、);plot(w/pi,20*log10(abs(H);grid;xlabel(W);ylabel(H/db);title(幅度响应);h,t1=impz(numd,dend);subplot(323);plot(t1,h);axis(0,60,min(h),max(h);xlabel(t/s);ylabel(h);title(单位冲击响应);%输出信号y=filter(numd,dend,x1);subplot(325);plot(t,y);axis(0,0.1,-1.3,1.3);Y=fft(y,L);ws=2*pi*fs;w=(0:L-1)*ws/L;subplot(326);stem(

25、0:L-1),abs(Y);高通滤波器：%周期方波时域f1=64;fs=10000;t=0:1/fs:0.1;x1=square(2*pi*f1*t);subplot(321);plot(t,x1);axis(0,0.1,-1.3,1.3);xlabel(t/s);ylabel(x(t);title(时域谱)%周期方波频域L=1250;X1=fft(x1,L);ws=2*pi*fs;w=(0:L-1)*ws/L;subplot(322);stem(0:L-1),abs(X1);xlabel(w/rad);ylabel(X);title(幅度谱);%采用BW和双线性变换法Wp=0.5*pi,Ws

26、=0.1*pi;Ap=1,As=50;wp=2*fs*tan(Wp/2);ws=2*fs*tan(Ws/2);wp1=1/wp,ws1=1/ws;N,wc=buttord(wp1,ws1,Ap,As,s);num,den=butter(N,wc,s);numt,dent=lp2hp(num,den,1);numd,dend=bilinear(numt,dent,fs);w=linspace(0,pi,1024);H=freqz(numd,dend,w);subplot(324);plot(w/pi,abs(H);xlabel(W);ylabel(H);title(幅度响应);h,t1=impz

27、(numd,dend);subplot(323);plot(t1,h);axis(0,max(t1),min(h),max(h);xlabel(t/s);ylabel(h);title(单位冲击响应);%输出信号y=filter(numd,dend,x1);subplot(325);plot(t,y);axis(0,0.1,-1.3,1.3);Y=fft(y,L);w=(0:L-1)*ws/L;subplot(326);stem(0:L-1),abs(Y);低通滤波器：%周期方波时域f1=64;fs=10000;t=0:1/fs:0.1;x1=square(2*pi*f1*t);subplot

28、(321);plot(t,x1);axis(0,0.1,-1.3,1.3);xlabel(t/s);ylabel(x(t);title(时域谱)%周期方波频域L=1250;X1=fft(x1,L);ws=2*pi*fs;w1=(0:L-1)*ws/pi;subplot(322);stem(0:L-1),abs(X1);xlabel(w/rad);ylabel(X);title(幅度谱);%采用BW和双线性变换法Wp=0.1*pi;Ws=0.5*pi;Ap=1;As=50;wp=2*fs*tan(Wp/2);ws=2*fs*tan(Ws/2);N,wc=buttord(wp,ws,Ap,As,s

29、);num,den=butter(N,wc,s);numd,dend=bilinear(num,den,fs);w=linspace(0,pi,1024);H=freqz(numd,dend,w);subplot(324);plot(w/pi,abs(H);xlabel(W);ylabel(H);title(幅度响应);h,t1=impz(numd,dend);subplot(323);plot(t1,h);axis(0,max(t1),min(h),max(h);xlabel(t/s);ylabel(h);title(单位冲击响应);%输出信号y=filter(numd,dend,x1);s

30、ubplot(325);plot(t,y);axis(0,0.1,-1.3,1.3);Y=fft(y,L);ws=2*pi*fs;w=(0:L-1)*ws/pi;subplot(326);stem(0:L-1),abs(Y);4. 调试分析过程描述 MATLAB是矩阵实验室（Matrix Laboratory）之意。除具备卓越的数值计算能力外，同时它还提供了专业水平的符号计算，文字处理，可视化建模仿真和实时控制等功能。它的指令表达式与数学,工程中常用的形式十分相似,故用MATLAB来解算问题要比用C,FORTRAN等语言完成相同的事情简捷得多。 用MATLAB进行模拟原型的数字滤波器的设计，一

31、般步骤如下: （1）按一定规则将给出的数字滤波器的技术指标转换成模拟低通滤波器的技术指标； （2） 根据转换后的技术指标使用滤波器阶数选择函数，确定最小阶数N和固有频率Wn，根据选用的模拟低通滤波器的类型可分别用: buttord,cheblord，cheb2ord，ellipord等函数； （3） 运用最小阶数N产生模拟滤波器原型，模拟低通滤波器的创建函数有:buttap,cheblap,cheb2ap,ellipap,besselap等； （4） 运用固有频率Wn把模拟低通滤波器原型转换成模拟低通、高通、带通、带阻滤波器，可分别用函数lp2lp,lp2hp,lp2bp,lp2bs； （5）

32、 运用冲激响应不变法或双线性变换法把模拟滤波器转换成数字滤波器，分别用函数impinva和bilinear来实现。 图8 带通滤波器结果图图9 带阻滤波器结果图图10 高通滤波器结图11 低通滤波器结程序运行结果：normwr=2.2361由设计流程计算得normwr=2.236与运行结果相同。低通原型的每一个边界频率都映射为带通滤波器两个相应的边界频率。根据通带截至频率和阻带截至频率与频谱函数曲线比较，满足设计要求。五总结 本次通过这个数字信号处理课程的课程设计，使得自己的能力和知识得到了很大的提升。不但对于上学期学到的数字信号处理内容有了回顾和复习，而且自己也真正的在我们的实验中体会到了理

33、论知识对我们的实验的指导意义，在自己的切身的实践动手的过程中，自己不仅对于上学期学到的内容中设计带通数字滤波器的整个过程流程的思路有了一定的了解，同时我也通过自己的实际动手操作发现原来自己的学习效果是如此的差，自己对于我们学到的知识从来都是一知半解，自己学习不动脑筋，不去思考，尤其是在自己的课程设计遇到差错的时候，自己没有一点耐心，十分烦躁，这一切的原因都是自己学习的不够深入，如果自己学习过程中注重思考，那么遇到的问题也就不会显得那么突然和无可奈何。同时在本次的实验中对于我们学到的双线性变换法，巴特沃斯设计模拟滤波器的运用，也算是比较熟悉了。同时在对数字带通滤波器的设计过程中，自己也重新学习了

34、MATLAB的运行环境，自己也通过复习课本的知识，对于MATLAB了解的更加深入，同时在上网学习过程中学习初步掌握了MATLAB语言在数字信号处理中一些基本库函数的调用和编写基本程序等的应用；了解和熟悉了在滤波器设计中，我们需要具备的设计思路；掌握到滤波器设计的一般原理，对滤波器的工作情况以及最终的设计的结果有了一定的了解；同时自己也掌握到了数字高通滤波器设计的一般步骤；加深了对滤波器设计中产生误差的原因以及双线性变换法优缺点的理解和认识。这些对于自己知识和实践能力的提高有很大的作用。在这次课程设计中，自己主要是负责寻找资料，以及实验报告的排版编写工作，主要是关于在本次课程设计中需要用到的工作

35、原理和设计思路的流程设计，自己和自己的队友分工明确，各取所长。自己也在本次的课程设计中实实在在锻炼了自己，同时自己也好好补习了了自己所欠缺的MATLAB语言知识，对于自己的编程能力有了提高，同时在复习到自己上学期学习的数字信号处理过程中，使理论联系了实际，巩固并深化了对课本基本知识的认识和理解，使理论得以升华。也算是一次不小的收获吧。 六、参考文献 1高西全,丁美玉.数字信号处理.第1版.北京：西安电子科技大学出版社, 2008 2郭仕剑，王宝顺，贺志国，杨可心.MATLAB7.X数字信号处理.人民邮电出版社, 2006 3李正周.Matlab数字信号处理与应用工程.北京：清华大学出版社，2008 4刘舒帆.数字信号处理实验（MATLAB）版.北京：西安电子科技大学出版社. 5张圣勤.MATLAB 7.0使用教程.北京：机械工业出版社，2008 -