IIR数字滤波器设计I

1、IIR数字滤波器设计（）一、 实验课程任务与要求数字信号处理课程是电子、通信、计算机、自动化、信息处理等专业的重要基础课。本课程以信号与系统、工程数学为基础，要求学生掌握时域离散信号和系统的基本理论、基本分析方法以及FFT、数字滤波器、谱分析等数字信号处理技术。数字信号处理是一门理论与实践联系紧密的课程，所以本课程安排一周课程设计，以帮助学生掌握数字信号处理技术，提高学生分析问题和解决问题的能力，并通过设计培养学生的创新意识。本课程设计的基本要求如下：1 学会用MATLAB语言编写数字信号处理的程序，通过上机实习加深对课堂所学知识的理解；2 上机前应按照要求把设计内容准备好，即编好程序及需要改

2、变的参数，能预计出可能出现的结果；3 观察实验结果，得出结论；4 设计结束时提交设计报告。二、设计题目IIR数字滤波器设计（）、三、设计要求（1）用单位冲激响应不变变换法进行设计。（2）给出详细的滤波器设计说明书。（3）给出经过运行是正确的程序清单并加上详细的注释。（4）画出所设计滤波器的幅度特性和相位特性。四 设计原理1、冲激响应不变法冲激响应不变法5的变换原理是使数字滤波器的单位冲激响应序列模仿模拟滤波器的单位冲激响应序列。将模拟滤波器的单位冲激响应加以等简隔抽样，使正好等于的抽样值，即满足 式(3.1)其中T是抽样周期。如果令是的拉普拉斯变换，为的Z变换，根据抽样序列的Z变换和模拟信号的

3、拉普拉斯变换的关系，可得 式(3.2)可以看出，利用冲激响应不变法将模拟滤波器变换数字滤波器，实际上是先对模拟滤波器的系统函数做周期延拓，再经过的映射变换，从而得到数字滤波器的系统函数。假设在S平面上，S在轴上取值，Z在Z平面内的单位圆周上取值，则可得到数字滤波器的频率响应和模拟滤波器的频率响应间的对应关系： 式(3.3)如图3.1所示，在S平面上每一条宽为的横条都将重叠地映射到整个Z平面上，而每一条的左半边映射到Z平面单位圆内，右半边映射到Z平面单位圆外，S平面虚轴映射Z平面单位圆上，虚轴上每一段长为的线段都映射到Z平面单位圆上一周。由于S平面每一横条都要重叠地映射到Z平面上，这正好反映和的

4、周期延拓函数的变换关系，故有冲激响应不变法并不是简单的从S平面映射到Z平面。图3.1 冲激响应不变法的映射关系因为不是简单的一一映射关系，且对于任何一个实际的模拟滤波器，它的频率响应是不可能真正带限的。因而将不可避免的出现频率的重叠，即混叠失真。数字滤波器的频率响应不能重现模拟滤波器的频率响应，只有当模拟滤波器的频率响应在超过重叠频率后的衰减很大时，混叠失真才会很小，此时才能够满足设计要求。综上所述，冲激响应不变法具有以下特点：(1) 模拟频率和数字频率的转换是线性的，并保持了模拟滤波器的时域瞬态特性。(2) 当模拟滤波器的频率响应不是严格限带时，则用冲激不变法设计出的数字滤波器在频率出现混叠

5、失真。(3) 由于(2)而使这种设计方法受到限制，即当不严格限带或g(t)变化不稳定。五、 设计内容 设计一个数字巴特沃斯低通滤波器，设计指标如下： 采样时间间隔。六 设计步骤（1） 把数字频率转换为模拟频率：，。（2） 计算巴特沃斯模拟滤波器的截止频率和阶数。（3） 设计巴特沃斯模拟低通滤波器，给出参数和（此处使用了MATLAB中的buttap(N)函数。（4） 把模拟滤波器用单位冲激响应不变变换法转换成数字滤波器（此处使用了MATLAB中的residuez函数）。（5） 变直接形式为并联形式，并给出结构图。（6） 画出幅度特性和相位特性。七、设计步骤框图 数字滤波器的设计步骤如图1所示 滤

6、波器技术指标直接型转换成并联型指标参数变换巴特沃斯模拟滤波器离散化数字巴特沃斯低通滤波器巴特沃斯模拟低通滤波器设计画出幅度特性和相位特性 图1八、Matlab程序代码%主程序wp=0.2*pi;ws=0.3*pi;Rp=1;%通带波动;As=15;%阻带衰减T=1;%采样周期Omigrp=wp*T;Omigrs=ws*T;t1=1,2*0.2588,0.25882+0.96592; b=1;a=conv(conv(t1,t1),t1);N,wc=buttord(wp,ws,Rp,As,s)%计算巴特沃斯数字滤波器的阶数N和截止频率wcB,A=butter(N,wc,s)z,p,k=buttap

7、(N)%模拟低通原型零、极点系数和增益因子bz,az=impinvar(b,a)%用冲激响应不变法将模拟滤波器转化为数字滤波器,采样频率默认1Hzwz=0:pi/512:pi;hz1=freqz(bz,az,wz);%巴特沃斯模拟低通滤波器频率响应C,B,A=dir2par(b,a)%直接型转换成并联型%绘图subplot(1,2,1); plot(wz/pi,abs(hz1)/hz1(1); title(幅度响应)xlabel(); ylabel(|H|);axis(0,1,0,1);grid;subplot(1,2,2); plot(wz/pi,hz1/pi); title(相位响应)xl

8、abel(); ylabel(单位:pi); axis(0,1,-1,2);grid;%直接型转换成并联型子程序function C,B,A=dir2par(b,a);M=length(b);N=length(a);r1,p1,C=residuez(b,a);p=cplxpair(p1,*eps);I=cplxcomp(p1,p);r=r1(I);K=floor(N/2);B=zeros(K,2);A=zeros(K,3);if K*2=N;for i=1:2:(N-2) Brow=r(i:1:(i+1),:); Arow=p(i:1:(i+1),:); Brow,Arow=residuez(

9、Brow,Arow,); B(fix(i+1)/2),:)=real(Brow); A(fix(i+1)/2),:)=real(Arow);EndBrow,Arow=residuez(r(N-1),p(N-1),);B(K,:)=real(Brow) 0;A(K,:)=real(Arow) 0;Elsefor i=1:2:(N-1) Brow=r(i:1:(i+1),:); Arow=p(i:1:(i+1),:); Brow,Arow=residuez(Brow,Arow,); B(fix(i+1)/2),:)=real(Brow); A(fix(i+1)/2),:)=real(Arow);E

10、ndEnd%比较两个含同样标量元素但（可能）有不同下标的复数对及其相应留数向量子程序function I=cplxcomp(p1,p2);I=;for j=1:length(p2)for i=1:length(p1) if(abs(p1(i)-p2(j)0.0001)I=I,i; EndEndEnd运行结果Rp = 1N = 6wc = 0.7087B = 0 0 0 0 0 0 0.1266A = 1.0000 2.7380 3.7484 3.2533 1.8824 0.6905 0.1266z = p = -0.2588 + 0.9659i -0.2588 - 0.9659i -0.707

11、1 + 0.7071i -0.7071 - 0.7071i -0.9659 + 0.2588i -0.9659 - 0.2588ik = 1bz = -0.0000 0.0060 0.1049 0.1946 0.0625 0.0021az = 1.0000 -2.6340 4.1006 -3.8164 2.4437 -0.9355 0.2117C = B = 0.4019 0.1523 0.3135 -0.0693 0.2846 0.0829A = 1.0000 0.5176 0.9999 1.0000 0.5176 0.9999 1.0000 0.5176 0.9999结构图九 仿真波形 总结：通过本次