实验八离散系统的Z域分析

实验八离散系统的Z 域分析

一、目的

（1）掌握利用MATLAB 绘制系统零极点图的方法（2）掌握离散时间系统的零极点分析方法

（3）掌握用MATALB 实现离散系统频率特性分析的方法（4）掌握逆Z 变换概念及MATLAB 实现方法

二、离散系统零极点

线性时不变离散系统可用线性常系数差分方程描述，即

()()N M

i j a y n i b x n j ==-=-∑∑ （8-1）

其中()y k 为系统的输出序列，()x k 为输入序列。

将式（8-1）两边进行Z 变换的

()

()()

()

j N

i i b z

Y z B z H z X z A z a z

-=-==

∑∑ （8-2）将式（8-2）因式分解后有：

()

()()

j N

i z q H z C

z p ==-=-∏∏ （8-3）

其中C 为常数，(1,2,

,)j q j M =为()H z 的M 个零点，(1,2,

,)i p i N =为()H z 的N 个

极点。

系统函数()H z 的零极点分布完全决定了系统的特性，若某系统函数的零极点已知，则系统函数便可确定下来。

因此，系统函数的零极点分布对离散系统特性的分析具有非常重要意义。通过对系统函数零极点的分析，可以分析离散系统以下几个方面的特性：

● 系统单位样值响应()h n 的时域特性； ● 离散系统的稳定性； ● 离散系统的频率特性；

三、离散系统零极点图及零极点分析 1．零极点图的绘制

设离散系统的系统函数为

()

()()

B z H z A z =

则系统的零极点可用MATLAB 的多项式求根函数roots()来实现，调用格式为：

p=roots(A)

其中A 为待根求多项式的系数构成的行矩阵，返回向量p 则是包含多项式所有根的列向

量。如多项式为231

()48

B z z z =++，则求该多项式根的MATLAB 命令为为：

A=[1 3/4 1/8];

P=roots(A) 运行结果为： P =

-0.5000 -0.2500

需注意的是，在求系统函数零极点时，系统函数可能有两种形式：一种是分子、分母多项式均按z 的降幂次序排列；另一种是分子、分母多项式均按1z -的升幂次序排列。这两种方式在构造多项式系数向量时稍有不同。

（1）()H z 按z 的降幂次序排列：系数向量一定要由多项式最高次幂开始，一直到常数项，缺项要用0补齐；如

34322()3221

z z

H z z z z z +=++++

其分子、分母多项式系数向量分别为A=[1 0 2 0]、B=[1 3 2 2 1]。

（2）()H z 按1z -的升幂次序排列：分子和分母多项式系数向量的维数一定要相同，不足的要用0补齐，否则0z =的零点或极点就可能被漏掉。如

1212()11124

z H z z z ---+=++

其分子、分母多项式系数向量分别为A=[1 2 0]、B=[1 1/2 1/4]。

用roots()求得()H z 的零极点后，就可以用plot()函数绘制出系统的零极点图。下面是求系统零极点，并绘制其零极点图的MATLAB 实用函数ljdt()，同时还绘制出了单位圆。

function ljdt(A,B)

% The function to draw the pole-zero diagram for discrete system 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

例1：绘制如下系统函数的零极点

（1）

3510 (

)

z z z

H z

z z z

-+-

（2）

10.5

()

H z

z z

解：MATLAB命令如下

（1）A=[1 -3 7 -5];

B=[3 -5 10 0];

ljdt(A,B)

绘制的零极点图如图8-1（a）所示。

（2）A=[1 3/4 1/8];

B=[1 -0.5 0];

ljdt(A,B)

绘制的零极点图如图8-1（b）所示。

2．离散系统零极点分析

（1）离散系统零极点分布与系统稳定性

《信号与系统》课程已讲到离散系统稳定的条件为：

●时域条件：离散系统稳定的充要条件为()

h n

∞

=-∞

<∞

∑，即系统单位样值响应绝对可和；

●Z域条件：离散系统稳定的充要条件为系统函数()

H z的所有极点均位于Z平面

的单位圆内。

对于三阶以下的低阶系统，可以利用求根公式求出系统函数的极点，从而判断系统的稳定性，但对于高阶系统，手工求解则显得十分困难，这时可以利用MATLAB来实现。实现方法是调用前述的函数ljdt()绘出系统的零极点图，然后根据极点的位置判断系统的稳定性。

例2：系统函数如例1所示，判断两个系统的稳定性。

（a）（b）

图8-1 离散系统的零极点图

解：由例1绘出的零极点图可以看出两个系统的稳定性分别为：第（1）个系统不稳定；第（2）个系统稳定。

（2）零极点分布与系统单位样值时域特性的关系从《信号与系统》课程中已经得知，离散系统的系统函数()H z 与单位样值响应()h n 是一对Z 变换对；因而，()H z 必然包含了()h n 的固有特性。

离散系统的系统函数可以写成

()

()()

j N

i z q H z C

z p ==-=-∏∏ （8-4）

若系统的N 个极点均为单极点，可将()H z 进行部分分式展开为：

1()N

i i i k z

H z z p ==-∑ （8-5）

由Z 逆变换得：

1()()()N

n i i i h n k p u n ==∑ （8-6）

从式（8-5）和（8-6）可以看出离散系统单位样值响应()h n 的时域特性完全由系统函数()H z 的极点位置决定。从《信号与系统》的学习中已经得出如下规律： ● ()H z 位于Z 平面单位圆内的极点决定了()h n 随时间衰减的信号分量； ● ()H z 位于Z 平面单位圆上的一阶极点决定了()h n 的稳定信号分量；

● ()H z 位于Z 平面单位圆外的极点或单位圆上高于一阶的极点决定了()h n 的随时

间增长的信号分量；

下面以例子证明上述规律的正确性：

例3：已知如下系统的系统函数()H z ，试用MATLAB 分析系统单位样值响应()h n 的时域特性。

（1）1

()1

H z z =-，单位圆上的一阶实极点；

（2）21

()2cos()1

H z z z π

=-+，单位圆上的一阶共轭极点；

（3）2

()(1)z

H z z =-，单位圆上的二阶实极点；（4）1

()0.8

H z z =-，单位圆内的一阶实极点；

（5）2

()(0.5)H z z =-，单位圆内的二阶实极点；（6）1

() 1.2

H z z =-，单位圆外的一阶实极点；

解：利用MATLAB 提供的函数impz()绘制离散系统单位样值响应波形，impz()基本调用方式为（其他方式，请读者参看MATLAB 帮助）：impz(b,a,N)，其中，b 为系统函数分子多项式的系数向量，a 为系统函数分母多项式的系数向量，N 为产生序列的长度；需

要注意的是，b 和a 的维数应相同，不足用0补齐，例如2211

()(1)21

H z z z z =

=--+的

b=[0 0 1]，a=[1 –2 1]。下面是求解个系统单位样值响应的MATLAB 命令：（1）a=[1 -1];

b=[0 1]; impz(b,a,10)

运行结果如图8-2（a ）所示。（2）a=[1 –2*cos(pi/8) 1];

b=[0 0 1]; impz(b,a,50)

运行结果如图8-2（b ）所示。（3）a=[1 -2 1];

b=[0 1 0]; impz(b,a,10)

运行结果如图8-2（c ）所示。（4）a=[1 -0.8];

b=[0 1]; impz(b,a,10)

运行结果如图8-2（d ）所示。（5）a=[1 -1 0.25];

b=[0 0 1]; impz(b,a,10)

运行结果如图8-2（e ）所示。（6）a=[1 -1.2];

b=[0 1]; impz(b,a,10)

运行结果如图8-2（f ）所示。

（a ）

（b ）

图8-2 系统的单位样值响应

四、离散系统频率特性分析

1．离散系统的频率响应()j H e ω

对于某因果稳定离散系统，如果激励序列为正弦序列：

0()sin()()x n A n u n ω= 则，根据《信号与系统》课程给出的结果有，系统的稳态响应为：

()()sin[()]()j ss y n A H e n u n ωω?ω=+

定义离散系统的频率响应为

()()()()j j j j z e H e H z H e e ωωω?ω=== 其中，()j H e ω——称为离散系统的幅频特性；

()?ω——称为离散系统的相频特性；

()j H e ω是以2π为周期的周期函数，只要分析()j H e ω在ωπ≤范围内的情况，便可分析出系统的整个频率特性。

2．用MATLAB 实现离散系统的频率特性分析方法（1）直接法

（c ）

（d ）

（e ）

（f ）

图8-2 系统的单位样值响应（续）

设某因果稳定系统的系统函数()H z ，则系统的频响特性为：

()()()()j j j j z e H e H z H e e ωωω?ω===

MATLAB 提供了专门用于求离散系统频响特性的函数freqz()，调用freqz()的格式有以下两种：

● [H,w]=freqz(B,A,N)

B 和A 分别为离散系统的系统函数分子、分母多项式的系数向量，N 为正整数，

返回量H 则包含了离散系统频响()j H e ω在0~π范围内N 个频率等分点的值，向量w 则包含0~π范围内N 个频率等分点。调用中若N 默认，默认值为512。

● [H,w]=freqz(B,A,N,’whole ’)

该调用格式将计算离散系统在0~2π范围内N 个频率等分点的频率响应()j H e ω的值。

因此，可以先调用freqz()函数计算系统的频率响应，然后利用abs()和angle()函数及plot()函数，即可绘制出系统在0~π或0~2π范围内的频响曲线。例4：绘制如下系统的频响曲线

0.5()z H z z -=

解：MATLAB 命令如下： B=[1 -0.5]; A =[1 0];

[H,w]=freqz(B,A,400,'whole'); Hf=abs(H); Hx=angle(H); clf

figure(1) plot(w,Hf)

title('离散系统幅频特性曲线') figure(2) plot(w,Hx)

title('离散系统相频特性曲线')

运行结果如图8-3所示。

（2）几何矢量法

图8-3 系统的幅频特性曲线和相频特性曲线

利用几何矢量求解示意图如图8-4所示。

j j j j j e q B e ψ

ω-=

i j j i i e p A e θω-=

有：

1212()

()()

()()M N M

j j

j j j j N j i

i B e

H e H e e Ae ψψψωω?ωθθθ+++=+++===∏∏

则系统的幅频特性和相频特性分别为：

()M

j j N i i B

H e A

ω===

∏∏ （8-7）

()M

j i j i ?ωψθ===-∑∑ （8-8）

根据式（8-7）和（8-8），利用MATLAB 来求解频率响应的过程如下： ● 根据系统函数()H z 定义分子、分母多项式系数向量B 和A ； ● 调用前述的ljdt()函数求出()H z 的零极点，并绘出零极点图； ● 定义Z 平面单位圆上的k 个频率分点；

● 求出()H z 所有的零点和极点到这些等分点的距离； ● 求出()H z 所有的零点和极点到这些等分点矢量的相角； ● 根据式（8-7）和（8-8）求出系统的()j H e ω和()?ω；

● 绘制指定范围内系统的幅频曲线和相频曲线；

下面是实现上述过程的实用函数dplxy()。有四个参数：k 为用户定义的频率等分点数目；B 和A 分别为系统函数分子、分母多项式系数向量；r 为程序绘制的频率特性曲线的频率范围（0~r π?）。

function dplxy(k,r,A,B)

%The function to draw the frequency response of discrete system p=roots(A); %求极点 q=roots(B); %求零点 figure(1)

ljdt(A,B) %画零极点图 w=0:r*pi/k:r*pi; y=exp(i*w); %定义单位圆上的k 个频率等分点 N=length(p); %求极点个数 M=length(q); %求零点个数 yp=ones(N,1)*y; %定义行数为极点个数的单位圆向量 yq=ones(M,1)*y; %定义行数为零点个数的单位圆向量 vp=yp-p*ones(1,k+1); %定义极点到单位圆上各点的向量 vq=yq-q*ones(1,k+1); %定义零点到单位圆上各点的向量 Ai=abs(vp); %求出极点到单位圆上各点的向量的模 Bj=abs(vq); %求出零点到单位圆上各点的向量的模

Ci=angle(vp); %求出极点到单位圆上各点的向量的相角Dj=angle(vq); %求出零点到单位圆上各点的向量的相角fai=sum(Dj,1)-sum(Ci,1); %求系统相频响应

H=prod(Bj,1)./prod(Ai,1); %求系统幅频响应

figure(2)

plot(w,H); %绘制幅频特性曲线

title('离散系统幅频特性曲线')

xlabel('角频率')

ylabel('幅度')

figure(3)

plot(w,fai)

title('离散系统的相频特性曲线')

xlabel('角频率')

ylabel('相位')

例5：已知某离散系统的系统函数为：

/4(1)

()

11/4

H z

绘出该系统的零极点图及频响特性。

解：MATLAB命令如下：

A=[1 -1/4];

B=[5/4 -5/4];

dplxy(500,2,A,B)

运行结果如图8-4所示。

五、实验内容

已知离散系统的系统函数分别为：

（1）

()

z z

H z

（2）

()

241

H z

z z z

+-+

试用MATLAB分析：

（1）绘出系统的零极点图，根据零极点图判断系统的稳定性；（2）如果系统稳定，绘出幅频特性和相频特性曲线。

图8-4 离散系统的零极点图、幅频和相频曲线