IIR数字滤波器设计和应用

...IIR数字滤波器设计和应用一、实验目的1、熟悉IIR数字滤波器的设计方法。2、掌握模拟滤波器的matlab实现。3、熟悉用脉冲响应不变法设计IIR数字滤波器的原理与方法。4、熟悉用双线性变换法设计IIR数字滤波器的原理与方法。5、掌握数字滤波器的计算机仿真方法。6、通过观察对实际心电图信号的滤波作用，获得数字滤波的感性知识。二、实验内容和要求1、用Matlab语言分别设计巴特奥斯低通滤波器和切比雪夫低通滤波器，其技术指标为：通带截止频率)5(2KHzp，通带最大衰减dBP3；阻带起始频率)10(2KHzs，阻带最小衰减dBs30；要求：求出他们的零点、极点、阶数、增益等，并画出图形作比较。...（1）巴特奥斯低通滤波器：%DesignaButterworthAnalogbandpassfilter%Desiredperformentsofthefilterwp=5000*2*pi;ws=10000*2*pi;Rp=3;Rs=30;%computeoderandcutofffrequency[N,Wn]=buttord(wp,ws,Rp,Rs,'s');NFc=Wn/(2*pi)%Computethefilter[b,a]=butter(N,Wn,'s');%Outputw=linspace(1,3000,1000)*2*pi;H=freqs(b,a,w);magH=abs(H);phaH=unwrap(angle(H));plot(w/(2*pi),20*log10(magH));xlabel('Frequency(Hz)');ylabel('Magnidute(dB)');title('ButterworthAnalogbandpassfilter')gridon...运算结果：N=5Fc=5.0124e+003（2）切比雪夫低通滤波器：%Chebyshev¢ñlowpassfilter%Desiredperformentswp=5000*2*pi;ws=10000*2*pi;Rp=3;Rs=30;%ComputeOrderNebs=sqrt(10^(Rp/10)-1);A=10^(Rs/20);Wc=wpWr=ws/wp;g=sqrt(A*A-1)/ebs;N1=log10(g+sqrt(g*g-1))/log10(Wr+sqrt(Wr*Wr-1));N=ceil(N1)...运算结果：Wc=3.1416e+004N=42、切比雪夫低通滤波器，其技术指标为：通带截止频率MHzfp3，通带最大衰减dBP1.0；阻带起始频率MHzfs12，阻带最小衰减dBs60；要求：求出他们的零点、极点、阶数、增益等，并画出图形作比较。3、用脉冲响应不变法把模拟滤波器651)(2ssssHa转换为数字滤波器H(z)，其中采样周期为T=0.1s。%Impulseinvariancemethod%Fromanalogtodigitalfilterconversion%Modelofanalogfilterb=[11];a=[156];Ts=0.1;[R,Ps,K]=residue(b,a);Pz=exp(Ps*Ts);disp('Usedirectprinciple£º')[bz,az]=residue(R,Pz,K)disp('Usefunction¡°IMPINVAR¡±£º')[bz1,az1]=impinvar(b,a,1/Ts)运算结果Usedirectprinciple：bz=1.0000-0.8966az=1.0000-1.55950.6065Usefunction“IMPINVAR”：bz1=0.1000-0.0897az1=1.0000-1.55950.60654、利用双线性变换法设计一个巴特奥斯低通IIR数字滤波器。(1)设计的指标参数为：在通带内频率低于0.2π时，最大衰减小于1dB，在阻带内[0.3π，π]频率区间上，最小衰减大于15dB....(2)以0.02π为采样间隔，打印出数字滤波器在频率区间［0,π/2］上的幅频响应特性曲线。(3)用所设计的滤波器对实际心电图信号采样序列(在本实验后面给出)进行仿真滤波处理，并分别打印出滤波前后的心电图信号波形图，观察总结滤波作用与效果。例如假设最后求出的数字滤波器系统函数H(z)为：)215.0904.01()3583.00106.11)(705.0268.11()1(0007378.0)(21212161zzzzzzzzH)1...................................................................................................).........(31式（zHkk式）（2........................................................................3,2,1,z1)21()(2121kzCBzzAzHkkk式中A=0.09036B1=1.2686,C1=-0.7051B2=1.0106,C2=-0.3583B3=0.9044,C3=-0.2155由(1式)和(2式)可见，滤波器H(z)由三个二阶滤波器H1(z),H2(z)和H3(z)级联组成，如图1所示。图1滤波器H(z)的组成x(n)H1(z)H2(z)H3(z)y1(n)y2(n)y3(n)＝y(n)H(z)...(2)编写滤波器仿真程序,计算H(z)对心电图信号采样序列x(n)的响应序列y(n)。设yk(n)为第k级二阶滤波器Hk(z)的输出序列，yk-1(n)为输入序列，如图1所示。由(2式)可得到差分方程：yk(n)=Ayk-1(n)+2Ayk-1(n-1)+Ayk-1(n-2)+Bkyk(n-1)+Ckyk(n-2)(3式)(3)在通用计算机上运行仿真滤波程序，并调用通用绘图子程序，完成实验内容(2)和(3)。人体心电图信号在测量过程中往往受到工业高频干扰，所以必须经过低通滤波处理后，才能作为判断心脏功能的有用信息。下面给出一实际心电图信号采样序列样本x(n)，其中存在高频干扰。在实验中，以x(n)作为输入序列，滤除其中的干扰成分。程序：x=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,-4,-2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];k=1;closeall;figure(1)subplot(2,2,1);n=0:55;stem(n,x,'.');axis([056-10050]);holdon;n=0:60;m=zeros(61);plot(n,m);xlabel('n');ylabel('x(n)');title('心电图信号采样序列X（n）');B=[0.090362*0.090360.09036];A=[1.2686-0.7051];A1=[1.0106-0.3583];A2=[0.9044-0.2155];...while(k=3)y=filter(B,A,x);%Thefunctionistofilte(Â˲¨)thesingalxx=y;ifk==2;A=A1;endifk==3A=A2;endk=k+1;endsubplot(2,2,3)n=0:55;stem(n,y,'.');axis([056-155]);holdon;n=0:60;m=zeros(61);plot(n,m);xlabel('n');ylabel('y(n)');title('三级滤波后的心电图信号');%求数字滤波器的幅频特性A=[0.09036,0.18072,0.09036];B1=[1,-1.2686,0.7051];B2=[1,-1.0106,0.3583];B3=[1,-0.9044,0.2155];[H1,w]=freqz(A,B1,100);[H2,w]=freqz(A,B2,100);[H3,w]=freqz(A,B3,100);H4=H1.*(H2);H=H4.*(H3);mag=abs(H);db=20*log10((mag+eps)/max(mag));subplot(2,2,2)plot(w/pi,db);axis([0,0.5,-50,10]);title('滤波器的幅频响应曲线');...三、思考题1.答：二者一致。2.答：（1）优点：一定是递归结构、对频率分量的选择性好（零极点可同时起作用）、相同性能下阶次较低、运算误差大，有可能出现极限环振荡、设计有封闭形式的公式，一次完成、对计算手段的要求较低;(2)缺点：一般是非线性的、系统不一定稳定、不能用FFT作快速卷积、有噪声反馈，噪声大。...3.答：脉冲响应不变法具有较好的线性特性，但会产生频率混叠。双线性变换法靠频率的严重非线性关系得到s平面与z平面的单值，固不会产生频率混叠，但其线性关系较差。4.答：T的值对结果有影响，因为整个s平面频率轴上的频率范围压缩的±π/T之间。