1、1“数字信号处理”实验指导书(三)一、实验课程编码:105003二、实验课程名称:数字信号处理 三、实验项目名称:数字滤波器结构的 MATLAB 实现四、实验目的掌握应用 MATLAB 实现 IIR、FIR 数字滤波器结构的方法,即熟悉应用 MATLAB 实现IIR、FIR 数字滤波器结构的相关函数。五、主要设备安装有 MATLAB 软件的电脑六、实验内容编写 MATLAB 程序,实现以下题目:1已知某四阶 IIR 滤波器的系统函数如下 -12240.283-.4067z.083(z)15zH分别实现系统的级联型结构和并联型结构,并求系统的单位脉冲响应和单位阶跃响应。2. 用直接型结构、级联型
2、结构实现系统函数 )81.09.1)(5.0(364) 2zzz并求系统的单位脉冲响应和单位阶跃响应。3. 用直接型结构、并联型结构实现系统函数 211211.65.43874)( zzzH并求系统的单位脉冲响应和单位阶跃响应。4. 实现下列线性相位 FIR 滤波器的级联型结构,求系统的单位脉冲响应和单位阶跃响应。(1)1234567()0.128.0.0.169.17658Hzzzzz(2)12345678.4.3.4.0.3zzzz 七、实验步骤1、熟悉与离散信号频谱分析相关的 MATLAB 函数(参考附录 1) ;22、通过运行附录 2 中提供的例题,熟悉用 MATLAB 实现数字滤波器
3、结构的基本方法;3、根据“六、实验内容”中各个题目的要求,编写 MATLAB 程序代码,调试程序,分析并保存结果。八、实验结果对实验练习题编写 MATLAB 程序并运行,在计算机上输出仿真结果。附录 1 主要的相关 MATLAB 函数1. yfilter(b,a,x )2. z,p,k tf2zp(num ,den)3. r,p,k residuez (num ,den )附录 2 例题例 1、用直接型实现系统函数为 5432161873)( zzzH的 IIR 数字滤波器,求单位脉冲响应和单位阶跃响应。MATLAB 程序 :%用直接型结构实现 IIR 数字滤波器b=1,-3,11,27,18
4、;a=16,12,2,-4,-1;N=25;h=impz(b,a,N); %直接型单位脉冲响应x=ones(1,5),zeros(1,N-5); %单位阶跃信号y=filter(b,a,x); % 直接型输出信号subplot(1,2,1);stem(h,.);line(0,30,0,0);title(直接型单位脉冲响应);subplot(1,2,2);stem(y,.);line(0,30,0,0);title(直接型单位阶跃响应);运行结果 : 3例 2、用级联型结构实现系统函数为 5432161873)( zzzH的 IIR 数字滤波器,求单位脉冲响应和单位阶跃响应。MATLAB 程序
5、:%用级联型结构实现直接形式的系统函数N=25;%delta=impseq(0,0,N);delta=1,zeros(1,N-1);b=1,-3,11,27,18;a=16,12,2,-4,-1;b0,B,A=dir2cas(b,a);h=casfilter(b0,B,A,delta);x=ones(1,5),zeros(1,N-5); y=casfilter(b0,B,A,x);subplot(1,2,1);stem(h,.);line(0,25,0,0);title(级联型单位脉冲响应);subplot(1,2,2);stem(y,.);line(0,25,0,0);title(级联型单位
6、阶跃响应);运行结果 :4例 3、用并联型结构实现系统函数 5432161873)( zzzHMATLAB 程序 :%用并联型结构实现直接形式的系统函数b=1,-3,11,27,18;a=16,12,2,-4,-1;N=25;%delta=impseq(0,0,N);delta=1,zeros(1,N-1);C,B,A=dir2par(b,a);x=ones(1,5),zeros(1,N-5);h=parfilter(C,B,A,delta); %并联型单位脉冲响应y=parfilter(C,B,A,x); %并联型单位阶跃响应subplot(1,2,1);stem(h,.);line(0,2
7、5,0,0);title(并联型单位脉冲响应);subplot(1,2,2);stem(y,.);line(0,25,0,0);title(并联型单位阶跃响应); 运行结果:例题中的函数:5function b0,B,A=dir2cas(b,a)%直接型结构到级联型结构的转换子程序%Direct-form to Cascade-form conversion in z-plane%(complex pair version)% b0,B,A=dir2cas(b,a);% C=gain coefficient% B=K by 3 matrix of real coefficients conta
8、ining bks% A=K by 3 matrix of real coefficients containing aks% b=numerator polynomial coefficients of direct form% a=denominator polynomial coefficients of direct form% compute gain coefficient b0b0=b(1);b=b/b0;a0=a(1);a=a/a0;b0=b0/a0;M=length(b);N=length(a);if NMb=b zeros(1,N-M);elseif MNa=a zeros
9、(1,M-N);N=M;elseNM=0;endK=floor(N/2);B=zeros(K,3);A=zeros(K,3);if K*2=Nb=b 0;a=a 0;endbroots=cplxpair(roots(b);aroots=cplxpair(roots(a);for i=1:2:2*KBrow=broots(i:1:i+1,:);Brow=real(poly(Brow);B(fix(i+1)/2),:)=Brow;Arow=aroots(i:1:i+1,:);Arow=real(poly(Arow);A(fix(i+1)/2),:)=Arow;endfunction C,B,A=d
10、ir2par(b,a)% 直接型转成并联型子程序M=length(b);N=length(a);r1,p1,C=residuez(b,a);p=cplxpair(p1,10000000*eps);I=cplxcomp(p1,p);r=r1(I);K=floor(N/2);B=zeros(K,2);A=zeros(K,3);if K*2=Nfor i=1:2:(N-2)Brow=r(i:1:(i+1),:);Arow=p(i:1:(i+1),:);Brow,Arow=residuez(Brow,Arow,);6B(fix(i+1)/2),:)=real(Brow);A(fix(i+1)/2),:
11、)=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);endendfunction I=cplxcomp(p1,p2)% I=cplxcomp(p1,p2)% 比较两个包含同样标量元素但(可能)有不同下标的复数对% 本程序必须用 cplxpair 函数之后使用,以便重新排序频率极点向量% 及其相应的留数向量;% p2=cplxpair(p1);%I=;for j=1:1:length(p2)for i=1:1:length(p1)if(abs(p1(i)-p2(j)0.0001)I=I,i;endendendI=I;