湖南大学信号与系统实验九

发布时间：2020-03-02 08:36:57 来源：范文大全收藏本文下载本文手机版

实验九

Z变换

§9.1 作离散时间零极点图目的

学习如何将一个离散时间有理系统函数H(z)的零极点展现一张零极点图。

下面的练习需要用到函数plotdpz，其具体实现如下。函数plotdpz(b,a)画出离散时间系统的零极点。 plotdpz的输入和filter具有相同的格式，并且如果分子和分母多项式不是同阶次的话，plotdpz将自动的对a或b补上适当个数的零。同时，plotdpz在图中还包括单位圆，以及如果在单位圆内有多于一个零点或极点的话，还会指出极点或零点的数目。 function plotdpz(b,a) la=length(a); lb=length(b); if (la>lb)

b=[b zeros(la-lb)]; elseif (lb>la)

a=[a zeros(lb-la)]; end ps=roots(a); zs=roots(b); mx=max(abs([ps\'zs\'.95]))+.05; clg; axis([-mx mx -mx mx]); axis(\'equal\'); hold on; w=[0:.01:2*pi]; plot(cos(w),sin(w),\'.\'); plot([-mx mx],[0 0]); plot([0 0],[-mx mx]); text(0.1,1.1,\'Im\',\'sc\'); text(1.1,.1,\'Re\',\'sc\'); plot(real(ps),imag(ps),\'X\'); plot(real(zs),imag(zs),\'o\'); numz=sum(abs(zs)==0); nump=sum(abs(ps)==0); if numz>1

text(-.1,-.1,num2str(numz)); elseif nump>1

text(-.1,-.1,num2str(nump)); end hold off 基本题

1．用plotdpz画出(9.3)或H(z)的零极点图。解：不用plotdpz,自定义函数： function ljdt(A,B) p=roots(A);%求系统极点 q=roots(B);%求系统零点

p=p\';%将极点列向量转置为行向量 q=q\';%将零点列向量转置为行向量 x=max(abs([p q 1]));%确定纵坐标范围 x=x+0.1; y=x;%确定横坐标范围 clf hold on axis([-x x -y y])%确定坐标轴显示范围 w=0:pi/300:2*pi; t=exp(i*w); plot(t)%画单位圆 axis(\'square\') plot([-x x],[0 0])%画横坐标轴 plot([0 0],[-y y])%画纵坐标轴 text(0.1,x,\'jIm[z]\') text(y,1/10,\'Re[z]\') plot(real(p),imag(p),\'x\')%画极点 plot(real(q),imag(q),\'o\')%画零点

title(\'pole-zero diagram for discrete system\')%标注标题 hold off a=[1 3 2]; b=[1 -1 0]; ljdt(a,b) p=roots(a) q=roots(b) pa=abs(p)

使用plotdpz(): b=[1 -1 0]; a=[1 3 2]; plotdpz(b,a);

2．用plotdpz画出满足差分方程y[n]y[n1]0.5y[n2]x[n]的滤波器的零极点图。解：代码如下：使用自定义函数ljdt():

a=[1 1 0.5]; b=[1 0 0]; ljdt(a,b)

使用plotdpz(): b=[1 0 0]; a=[1 1 0.5]; plotdpz(b,a);

3．用plotdpz画出满足差分方程

y[n]1.25y[n1]0.75y[n2]0.125y[n3]x[n]0.5x[n1] 所描述滤波器的零极点图。解：

使用自定义函数ljdt(): a=[1 -1.25 0.75 -0.125] b=[1 0.5] ljdt(a,b)

使用plotdpz（）： b=[1 0.5 0 0]; a=[1 -1.25 0.75 -0.125]; plotdpz(b,a);

§9.2离散时间频率响应的几何解释基本题

对于这些题将考察一个简单的二阶离散时间系统，其系统函数是

H1(z)110.9z10.81z2 z0.9

1．定义b1和a1包含按filter要求格式的H1(z)分子和分母多项式的系数，按10.1节所介绍的方法求出并画出H1(z)的零极点图。解：代码如下： b1=[1 0 0]; a1=[1 -0.9 0.81]; plotdpz(b1,a1);

2．定义omega=[0:511]*pi/256和unitcirc=exp(j*omega)得到在单位圆上512等分点，在这些点上将要对频率响应H1(ej)求值。定义polevectors1是一个2×512的矩阵，其中每一行包含一些复数，这些复数是由unitcirc的相应列减去一个极点位置得到的。如果ps1是一个列向量，它包含了极点的位置，那么就能用下列命令来完成

>> polevectors1=ones(2,1)*unitcirc-ps1*ones(1,512); 在MATLAB中总是愿意用矩阵和向量运算，因为MATLAB几乎总是计算矩阵运算比用for循环构成的运算要快。用abs和atan2定义polelength1和poleangle1作为polevectors1中每一元素的幅值和相位。解： b1=[1 0 0]; a1=[1 -0.9 0.81]; ps1=roots(a1); omega=[0:511].*pi/256; unitcirc=exp(i.*omega); polevectors1=zeros(2,512); polevectors1=ones(2,1)*unitcirc-ps1*ones(1,512); polelength1=abs(polevectors1) poleangle1=atan2(imag(polevectors1),real(polevectors1)) 3．类似于polevecters1定义zerovectors1，使得它是2×512的矩阵，其中包含从零点位置到unitcirc元素的向量。定义zerolength1和zeroangle1分别是这些向量的幅值和相位。解： b1=[1 0 0]; a1=[1 -0.9 0.81]; zs1=roots(a1); omega=[0:511].*pi/256; unitcirc=exp(i.*omega); zerovectors1=zeros(2,512); zerovectors1=ones(2,1)*unitcirc-ps1*ones(1,512); zerolength1=abs(zerovectors1) zeroangle1=atan2(imag(zerovectors1),real(zerovectors1)) 4．画出polelength1和zerolength1对omega的图。根据这些图，你预计H1(ej)在哪里有最大值和最小值？解： b1=[1 0 0]; a1=[1 -0.9 0.81]; ps1=roots(a1); zs1=roots(b1); omega=[0:511]*pi/256; unitcirc=exp(i*omega); zerovectors1=zeros(2,512); zerovectors1=ones(2,1)*unitcirc-zs1*ones(1,512); zerolength1=abs(zerovectors1); subplot(2,1,1) plot(omega,zerolength1); polevectors1=zeros(2,512); polevectors1=ones(2,1)*unitcirc-ps1*ones(1,512); polelength1=abs(polevectors1); subplot(2,1,2) plot(omega,polelength1);

5．用polelength1和zerolength1计算H1(ej)，并将这个结果存入geonH1mag中。用poleangle1和zeroangle1计算H1(ej)，并将这个结果存入gemoH1phase中。在定义geonH1mag和gemoH1phase时，你或许会发现prod和sum是很有用的。画出从几何上到处的幅值和相位，并将它和经由下面计算所得结果比较： >> H1=freqz(b1,a1,512,\'whole\'); 你在4中估计的峰值频率正确吗？解： b1=[1 0 0]; a1=[1 -0.9 0.81]; ps1=roots(a1); zs1=roots(b1); omega=[0:511]*pi/256; unitcirc=exp(i*omega); zerovectors1=zeros(2,512); zerovectors1=ones(2,1)*unitcirc-zs1*ones(1,512); zerolength1=abs(zerovectors1); zeroangle1=atan2(imag(zerovectors1),real(zerovectors1)); polevectors1=zeros(2,512); polevectors1=ones(2,1)*unitcirc-ps1*ones(1,512); polelength1=abs(polevectors1); poleangle1=atan2(imag(polevectors1),real(polevectors1)); geonH1mag=zerolength1./polelength1; gemoH1phase=zeroangle1-poleangle1; subplot(2,1,1) plot(omega,geonH1mag); subplot(2,1,2) plot(omega,gemoH1phase); freqz(b1,a1,512,\'whole\');

中等题

下面这些题的第一组将考查移动H1(z)一个零点的效果。具体说来，考虑系统函数为

H2(z)10.5z110.9z10.81z2 z0.9

5．求出并画出H2(z)的零极点图。你怎样预计这个系统的polevectors2或者zerovectors2与H1(z)的有什么不同吗？解： b2=[1 -0.5 0]; a2=[1 -0.9 0.81]; plotdpz(b2,a2);

结果分析：这个系统的polevectors2与H1(z)相同，但zerovectors2与H1(z)零点值不同，矩阵第二行的列值发生变化

6．计算H2(z)对polevectors2和zerovectors2以及所有这些向量的幅值和相位。画出幅值和相位对omega的图。你在5中的预计正确吗？解： b2=[1 -0.5 0]; a2=[1 -0.9 0.81]; ps2=roots(a2); zs2=roots(b2); omega=[0:511]*pi/256; unitcirc=exp(i*omega); zerovectors2=zeros(2,512); zerovectors2=ones(2,1)*unitcirc-zs2*ones(1,512); zerolength2=abs(zerovectors2); zeroangle2=atan2(imag(zerovectors2),real(zerovectors2)); polevectors2=zeros(2,512); polevectors2=ones(2,1)*unitcirc-ps2*ones(1,512); polelength2=abs(polevectors2); poleangle2=atan2(imag(polevectors2),real(polevectors2)); geonH1mag=zerolength2./polelength2; gemoH1phase=zeroangle2-poleangle2; subplot(2,1,1) plot(omega,geonH1mag); subplot(2,1,2) plot(omega,gemoH1phase); freqz(b2,a2,512,\'whole\');

结果分析：在5中的预计正确

7．根据这些零点的变化，预计H2(ej)与H1(ej)将会有怎样的不同。用freqz计算H2验证你的答案。解： b2=[1 -0.5 0]; a2=[1 -0.9 0.81]; H2=freqz(b2,a2,512,\'whole\') 结果分析：H2(ej)与H1(ej)相比，相位和幅值都发生了变化以下的题将考虑系统函数为

H332zz(z)132z0.25z0.251122 z0.5

8．求出并画出H3(z)的零极点图。这些零点和极点的位置是如何关联的？解：代码如下： b3=[0.25 -1.732/2 1]; a3=[1 -1.732/2 0.25]; plotdpz(b3,a3);

结果分析：这些零点和极点的位置关于x轴对称

9．类似于你在基本题中所做的方式定义polevectors3和zerovectors3。定义polelength3和zerolength3是这些复数的幅值。在同一组坐标轴上画出所有这些幅值，也就是polelength3和zerolength3每一行的幅值。这些幅值有什么关系？根据这一点，你如何预计频率响应的幅值H3(ej)岁频率变化？利用这些长度计算频率响应幅值，并将它存入gemoH3mag中。画出gemoH3mag对于omega的图。解：

b3=[0.25 -1.732/2 1]; a3=[1 -1.732/2 0.25]; ps3=roots(a3); zs3=roots(b3); omega=[0:511]*pi/256; unitcirc=exp(i*omega); zerovectors3=zeros(2,512); zerovectors3=ones(2,1)*unitcirc-zs3*ones(1,512); zerolength3=abs(zerovectors3); zeroangle3=atan2(imag(zerovectors3),real(zerovectors3)); polevectors3=zeros(2,512); polevectors3=ones(2,1)*unitcirc-ps3*ones(1,512); polelength3=abs(polevectors3); poleangle3=atan2(imag(polevectors3),real(polevectors3)); plot(omega,polelength3); hold on; plot(omega,zerolength3,\'r\'); geonH1mag=zerolength3./polelength3; gemoH1phase=zeroangle3-poleangle3; plot(omega,geonH1mag);

结果分析：polelength3和zerolength3幅值变化趋势一样，但是大小不等，形成闭合曲线；频率响应的幅值H3(ej)的频率也将周期性的单调递减 10．利用freqz计算H3，并验证你在9中的答案。解：

b3=[0.25 -1.732/2 1]; a3=[1 -1.732/2 0.25]; H3=freqz(b2,a2,512,\'whole\') 实验报告：

通过这次实验，学会了如何将一个离散时间有理系统函数的零极点展现一张零极点图。并通过考察一个简单的二阶离散时间系统的系统函数，能够大致预先判断零极点和系统的相关特性如稳定性，因果性等，加深了对Z变换的理解。基本达到本次实验的目的。