离散时间系统的Matlab实现

离散时间系统的Matlab实现Impz函数功能：求解系统的单位冲击响应调用方式：[h,t]=impz(b,a)：b、a分别为系统传递函数的分子和分母的系数向量。返回系统(b,a)的冲激响应h和相应的时间轴向量t。[h,t]=impz(b,a,n)：返回n点冲激响应。impz(b,a);在当前窗口用stem(t,h)函数出图。【例1】计算线性系统(b,a)的冲激响应：b=[0.20.10.30.10.15];a=[1-1.11.45-0.60.3];impz(b,a,50);dstep功能：求解数字系统的阶跃响应。调用格式：[h,t]=dstep(b,a);求解数字系统的阶跃响应h，取样点数为缺省值。[h,t]=dstep(b,a,n);求解数字系统的阶跃响应h，取样点数为n值。dstep(b,a)；在当前窗口用stairs(t,h)函数出图。【例2】已知一个因果系统的差分方程为6y(n)+2y(n-2)=x(n)+3x(n-1)+3x(n-2)+x(n-3)满足初始条件y(-1)=0,x(-1)=0.求系统的单位冲激响应和阶跃响应。解：将上述方程对y(n)项系数进行归一化，得到其系统函数分子和分母系数a0=1,a1=0,a2=1/3,a3=0b0=1/6,b1=1/2,b2=1/2,b3=1/6用impz()函数的MATLAB程序（取N=32点作图）a=[1,0,1/3,0];b=[1/6,1/2,1/2,1/6];N=32;n=0:N-1;hn=impz(b,a,n);gn=dstep(b,a,n);subplot(1,2,1),stem(n,hn,’filled’);title(‘系统的单位冲击响应’)；ylabel(‘h(n)’);xlabel(‘n’);axis([0,N-1,-1.1*min(hn),1.1*max(hn)]);subplot(1,2,2),stem(n,gn,’k’);title(‘系统的单位阶跃响应’)；ylabel(‘g(n)’);xlabel(‘n’);axis([0,N-1,-1.1*min(gn),1.1*max(gn)]);zplanezplane(z,p)：绘制系统零极点图，“o”表示零点，“x”表示极点。z，p分别为零点和极点向量。zplane(b,a)：b、a分别为系统传递函数的分子和分母系数向量。【例3】计算线性系统(b,a)的零点和极点。b=[0.20.10.30.10.15];a=[1-1.11.45-0.60.3];zplane(b,a);legend('零点','极点');filter功能：对数字系统的输入信号进行滤波处理。因为一个离散系统可以看作是一个滤波器，系统的输出就是输入经过滤波器滤波的结果。调用格式：y=filter(b,a,x)，对于由矢量b,a决定的数字系统（b和a分别表示系统函数H（z）对应的分子项和分母项系数构成的数组，而且分母系数要归一化处理。）当输入信号为x时，对x中的数据进行滤波，结果存于y中，长度取max(na,nb).[y,zf]=filter(b,a,x);除得到结果矢量y外，还得到x的最终状态矢量zf。y=filter(b,a,x,zi)；可在zi中指定x的初始状态。【例4】计算低通滤波器的冲激响应。x=[1zeros(1,100)];[b,a]=cheby1(11,1,.4);y=filter(b,a,x);stem(y);filtic子函数功能：为filter子函数选择初始条件。调用格式：zi=filtic(b,a,y,x)；求给定输入x和y时的初始状态。zi=filtic(b,a,y)；求x=0，给定输入y时的初始状态。其中，x和y分别是表示过去的输入和输出。【例5】已知一个因果系统的差分方程为6y(n)+2y(n-2)=x(n)+3x(n-1)+3x(n-2)+x(n-3)满足初始条件y(-1)=0,x(-1)=0.求系统的单位冲激响应和阶跃响应。用filter()函数的MATLAB程序（取N=32点作图）a=[1,0,1/3,0];b=[1/6,1/2,1/2,1/6];xi=filtic(b,a,0,0);N=32;n=0:N-1;x1=[n==0];%单位冲激信号hn=filter(b,a,x1,xi);subplot(1,2,1),stem(n,hn,’filled’);title(‘系统的单位冲激响应’)；ylabel(‘h(n)’);xlabel(‘n’);axis([0,N-1,-1.1*min(hn),1.1*max(hn)]);x2=[n>=0];%单位阶跃信号gn=filter(b,a,x2,xi);subplot(1,2,2),stem(n,gn,’k’);title(‘系统的单位阶跃响应’)；ylabel(‘g(n)’);xlabel(‘n’);axis([0,N-1,-1.1*min(gn),1.1*max(gn)]);freqz功能：用于求离散时间系统的频率响应函数。

调用格式：1）[h,w]=freqz(b,a,n)。可以得到数字滤波器的n点复频响应值，这n个点均匀地分布之[0,pi]上，并将这n个频点的频率记录在w中，相应的频响值记录在h中，n缺省时取512点。2）[hf]=freqz(b,a,n,Fs)；用于对H(jw)在[0,Fs/2]上等间隔采样n点，采样点频率及相应频响值分别记录在f和h中，由用户指定Fs（以Hz为单位）的值。3）h=freqz(b,a,w)；用于对H(jw)在[0,2*pi]上进行采样，采样频率点由矢量w指定。4）h=freqz(b,a,f,Fs)；用于对H(jw)在[0,Fs]上采样，采样频率点由矢量f指定。【例6】已知离散时间系统的系统函数求系统在0-pi频率范围内，归一化的绝对幅度频率响应，相对幅度频率响应，相位频率响应和零极点分部图。b=[0.2,0.1,0.3,0.1,0.2];a=[1,-1.1,1.5,-0.7,0.3];n=(0:500)*pi/500;%在[0,pi]的范围内取501个采样点[h,w]=freqz(b,a,n);%求系统的频率响应subplot(2,2,1),plot(n/pi,abs(h));grid%作系统的绝对幅度频响图axis([0,1,1.1*min(abs(h)),1.1*max(abs(h))]);ylabel('幅度');title('幅频响应（V）');subplot(2,2,2),plot(n/pi,angle(h));grid%作系统的相位频响图axis([0,1,1.1*min(angle(h)),1.1*max(angle(h))]);ylabel('相位');xlabel('以pi为单位的频率');title('相频响应');db=20*log10(abs(h));subplot(2,2,3),plot(n/pi,db);grid%作系统的相对幅度频响图title('幅频响应（db）');subplot(2,2,4),zplane(b,a);grid%作零极点分布图title('零极点分布');状态空间的转换tf2ss函数：由传递函数到状态空间的转换[A,B,C,D]=TF2SS(NUM,DEN) H(s)=NUM(s)/DEN(s