实验七 离散系统分析的MATLAB 实现一、实验目的1、掌握利用MATLAB 绘制系统零极点图的方法；2、掌握离散时间系统的零极点分析方法；3、学习离散系统响应的MATLAB 求解方法；4、掌握用MATALB 实现离散系统频率特性分析的方法；5、深刻理解离散系统的系统函数零极点对系统频响的影响，可以根据 零极点知识设计简单的滤波器。

二、基本原理（一）离散系统零极点线性时不变离散系统可用线性常系数差分方程描述，即()()N Miji j a y n i b x n j ==-=-∑∑ （1）其中()y k 为系统的输出序列，()x k 为输入序列。

将式（1）两边进行Z 变换，00()()()()()Mjjj Nii i b zY z B z H z X z A z a z-=-====∑∑ （2） 将式（2）因式分解后有：11()()()Mjj Nii z q H z Cz p ==-=-∏∏ （3）其中C 为常数，(1,2,,)j q j M =L 为()H z 的M 个零点，(1,2,,)i p i N =L 为()H z 的N 个极点。

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

因此，系统函数的零极点分布对离散系统特性的分析具有非常重要意义。

（二）离散系统零极点图及零极点分析 1、零极点图的绘制设离散系统的系统函数为()()()B z H z A z =则系统的零极点可用MATLAB 的多项式求根函数roots()来实现，调用格式为：p=roots(A) 其中A 为待求根多项式的系数构成的行矩阵，返回向量p 则是包含多项式所有根的列向量。

如多项式为231()48B 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()3221z zH z z z z z +=++++其分子、分母多项式系数向量分别为A=[1 0 2 0]、B=[1 3 2 2 1]。

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

如11212()11124z H z z z ---+=++其分子、分母多项式系数向量分别为A=[1 2 0]、B=[1 1/2 1/4]。

用roots()求得()H z 的零极点后，就可以用plot()函数绘制出系统的零极点图。

下面是求系统零极点，并绘制其零极点图的MATLAB 实用函数ljdt()，同时还绘制出了单位圆。

函数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 onaxis([-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）32323510()375z z zH z z z z -+=-+- （2）11210.5()31148z H z z z ----=++解：MATLAB 命令如下：（1） A=[1 -3 7 -5];B=[3 -5 10 0]; ljdt(A,B)绘制的零极点图如图7-1（a ）所示。

（2） A=[1 3/4 1/8];B=[1 -0.5 0]; ljdt(A,B)绘制的零极点图如图7-1（b ）所示。

图7-1 离散系统的零极点图2、离散系统零极点分析《信号与系统》课程已讲到离散系统稳定的条件为：①时域条件：离散系统稳定的充要条件为()n h n ∞=-∞<∞∑，即系统单位样值响应绝对可和；②Z 域条件：离散系统稳定的充要条件为系统函数()H z 的所有极点均位于Z 平面的单位圆内。

对于三阶以下的低阶系统，可以利用求根公式求出系统函数的极点，从而判断系统的稳定性，但对于高阶系统，手工求解则显得十分困难，这时可以利用MATLAB 来实现。

实现方法是调用前述的函数ljdt()绘出系统的零极点图，然后根据极点的位置判断系统的稳定性。

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

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

（三）离散系统响应的求解除可以使用MATLAB 命令lsim 求解外，还可以使用命令filter 来求解系统响应。

例3：已知系统函数为1111121212()(10.2)(10.6)10.40.12z z H z z z z z------++==-++-，求 （1） 系统的脉冲响应()h n ；（2） 输入()()x n u n =，求系统的零状态响应()y n ；（3） 输入()()x n u n =，初始条件(1)1,(2)2y y -=-=，求系统的完全响应()y n 。

解：（1） 计算前11个时刻的()h n N=11；b=[1,2];a=[1,0.4,-0.12]; x=[1,zeros(1,N-1)]; y=filter(b,a,x)（2） 计算前11个时刻的零状态响应()y n N=11；b=[1,2];a=[1,0.4,-0.12]; x=ones(1,N); y=filter(b,a,x)（3） 计算前11个时刻的完全响应()y n注意filter(b,a,x,zi)中的初始值zi 不是(1)1,(2)2y y -=-=，它可以用命令filtic 来求得。

N=11；b=[1,2];a=[1,0.4,-0.12]; x=ones(1,N);zi=filtic(b,a,[1,2]); y=filter(b,a,x,zi);（四）离散系统频率特性分析 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 实现离散系统的频率特性分析方法设某因果稳定系统的系统函数()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) clffigure(1) plot(w,Hf)title('离散系统幅频特性曲线') figure(2)plot(w,Hx)title('离散系统相频特性曲线')运行结果如图7-2所示。

图7-2 系统的幅频特性曲线和相频特性曲线3、离散系统函数的零极点对系统频域特性的影响例5：已知某系统系统函数为1121()10.5zH zz z---+=-+，求系统的幅频曲线。

解： b=[1,1,0];a=[1,-1,0.5];[H,w]=freqz(b,a);plot(w,abs(H));xlabel('Frequency(rad)'); ylabel('Magnitude');title('Magnitude response');其运行结果如下图7-3所示。

图7-3 系统的幅频特性曲线由图7-3可见，系统的零点1j z e π=-=迫使高频处πΩ=的幅频响应的幅度很小。

三、预习练习1、为了使实验能够顺利地进行，课前对教材中离散系统的Z 域分析的相关内容和实验原理、方法及内容做好充分预习，并预期实验的结果。

2、学习MATLAB 软件，尤其是其中的和本次实验有关的一些函数的使用。