离散时间LTI系统的时域分析(1)

第3章 离散时间LTI系统的时域分析

3.1 实验目的

学会运用MATLAB求解离散时间系统的零状态响应； 学会运用MATLAB求解离散时间系统的单位取样响应； 学会运用MATLAB求解离散时间系统的卷积和。

3.2 实验原理及实例分析

3.2.1 离散时间系统的响应

离散时间LTI系统可用线性常系数差分方程来描述，即

ay(n i) bx(n j) （3-1）

i

j

i 0

j 0

NM

其中，ai（i 0，1， ，N）和bj（j 0，1， ，M）为实常数。

MATLAB中函数filter可对式（13-1）的差分方程在指定时间范围内的输入序列所产生的响应进行求解。函数filter的语句格式为

y=filter(b,a,x)

其中，x为输入的离散序列；y为输出的离散序列；y的长度与x的长度一样；b与a分别为差分方程右端与左端的系数向量。

【实例3-1】 已知某LTI系统的差分方程为

3y(n) 4y(n 1) 2y(n 2) x(n) 2x(n 1)

试用MATLAB命令绘出当激励信号为x(n) (1/2)u(n)时，该系统的零状态响应。

解：MATLAB源程序为

>>a=[3 -4 2]; >>b=[1 2]; >>n=0:30; >>x=(1/2).^n; >>y=filter(b,a,x); >>stem(n,y,'fill'),grid on

>>xlabel('n'),title('系统响应y(n)')

n

程序运行结果如图3-1所示。

离散时间LTI系统的时域分析(1)

图3-1 实例3-1系统的零状态响应

3.2.2 离散时间系统的单位取样响应

系统的单位取样响应定义为系统在 (n)激励下系统的零状态响应，用h(n)表示。MATLAB求解单位取样响应可利用函数filter，并将激励设为前面所定义的impDT函数。例如，求解实例13-1中系统的单位取样响应时，MATLAB源程序为

>>a=[3 -4 2]; >>b=[1 2]; >>n=0:30; >>x=impDT(n); >>h=filter(b,a,x); >>stem(n,h,'fill'),grid on

>>xlabel('n'),title('系统单位取样响应h(n)')

程序运行结果如图3-2所示。

图3-2 实例13-1的系统单位取样响应

离散时间LTI系统的时域分析(1)

MATLAB另一种求单位取样响应的方法是利用控制系统工具箱提供的函数impz来实现。impz函数的常用语句格式为

impz(b,a,N)

其中，参数N通常为正整数，代表计算单位取样响应的样值个数。

【实例3-2】 已知某LTI系统的差分方程为

3y(n) 4y(n 1) 2y(n 2) x(n) 2x(n 1)

利用MATLAB的impz函数绘出该系统的单位取样响应。

解：MATLAB源程序为

>>a=[3 -4 2]; >>b=[1 2]; >>n=0:30;

>>impz(b,a,30),grid on

>>title('系统单位取样响应h(n)')

程序运行结果如图3-3所示，比较图3-2和图3-3，不难发现结果相同。

图3-3 系统单位取样响应

3.2.3 离散时间信号的卷积和运算

由于系统的零状态响应是激励与系统的单位取样响应的卷积，因此卷积运算在离散时间信号处理领域被广泛应用。离散时间信号的卷积定义为

y(n) x(n)*h(n)

m

x(m)h(n m) （3-2）

可见，离散时间信号的卷积运算是求和运算，因而常称为“卷积和”。

MATLAB求离散时间信号卷积和的命令为conv，其语句格式为

y=conv(x,h)

离散时间LTI系统的时域分析(1)

其中，x与h表示离散时间信号值的向量；y为卷积结果，它默认序列从n＝0开始。但是如果序列是从一负值开始，即

x(n):nx1 n nx2 h(n):nh1 n nh2

如果nx1<0或nh1<0就不能直接采用conv函数。其卷积结果序列为{y(n):nx1 nh1 n nx2 nh2}，这样就可以构成一个新的卷积函数conv_m。如下所示：

function[y,ny]=conv_m(x,nx,h,nh)

ny1=nx(1)+nh(1);ny2=nx(length(x))+nh(length(h)); ny=[ny1:ny2]; y=conv(x,h)

值得注意的是用MATLAB进行卷积和运算时，无法实现无限的累加，只能计算时限信号的卷积。

【实例3-3】 已知某系统的单位取样响应为h n 0.8n u n u n 8 ，试用MATLAB求当激励信号为x(n) u(n) u(n 4)时，系统的零状态响应。

解：MATLAB中可通过卷积求解零状态响应，即x(n)*h(n)。由题意可知，描述h(n)向量的长度至少为8，描述x(n)向量的长度至少为4，因此为了图形完整美观，我们将h(n)向量和x(n)向量加上一些附加的零值。MATLAB源程序为

nx=-1:5; %x(n)向量显示范围(添加了附加的零值) nh=-2:10; %h(n)向量显示范围(添加了附加的零值) x=uDT(nx)-uDT(nx-4);

h=0.8.^nh.*(uDT(nh)-uDT(nh-8)); [y,ny]=conv_m(x,nx,h,nh); subplot(311)

stem(nx,x,'fill'),grid on xlabel('n'),title('x(n)') axis([-4 16 0 3]) subplot(312)

stem(nh,h','fill'),grid on xlabel('n'),title('h(n)') axis([-4 16 0 3]) subplot(313)

stem(ny,y,'fill'),grid on xlabel('n'),title('y(n)=x(n)*h(n)') axis([-4 16 0 3])

离散时间LTI系统的时域分析(1)

程序运行结果如图3-5所示。

图3-5 利用卷积和法求解系统的零状态响应

3.3 编程练习

1.

试用MATLAB命令求解以下离散时间系统的单位取样响应。 （1）3y(n) 4y(n 1) y(n 2) x(n) x(n 1) （2） 2.

已知某系统的单位取样响应为h n () u n u n 10 ，试用MATLAB求当激励

n

5

y(n) 6y(n 1) 10y(n 2) x(n) 2

78

信号为x(n) u(n) u(n 5)时，系统的零状态响应。

离散时间LTI系统的时域分析(1)

第4章 z变换及离散时间LTI系统的z域

分析

4.1 实验目的

学会运用MATLAB求离散时间信号的z变换和z反变换； 学会运用MATLAB分析离散时间系统的系统函数的零极点；

学会运用MATLAB分析系统函数的零极点分布与其时域特性的关系； 学会运用MATLAB进行离散时间系统的频率特性分析。

4.2 实验原理及实例分析

4.2.1 z正反变换

序列x n 的z变换定义为

X z Z x n

n

x n z

n

（4-1）

其中，符号Z表示取z变换，z是复变量。相应地，单边z变换定义为

X z Z x n x n z n （4-2）

n 0

MATLAB符号数学工具箱提供了计算离散时间信号单边z变换的函数ztrans和z反变换函数iztrans，其语句格式分别为

Z=ztrans(x) x=iztrans(z)

上式中的x和Z分别为时域表达式和z域表达式的符号表示，可通过sym函数来定义。

【实例4-1】 试用ztrans函数求下列函数的z变换。

（1）x(n) acos( n)u(n)； （2）x(n) [2解：（1）z变换MATLAB源程序为

>>x=sym('a^n*cos(pi*n)'); >>Z=ztrans(x);

>>simplify(Z) %对Z进行简化运算 ans=

z/(z+a)

n

n 1

( 2)n 1]u(n)。

（2）z变换MATLAB源程序为

>>x=sym('2^(n-1)-(-2)^(n-1)');

离散时间LTI系统的时域分析(1)

>>Z=ztrans(x); >>simplify(Z) ans=

z^2/(z-2)/(z+2)

【实例4-2】 试用iztrans函数求下列函数的z反变换。

8z 19z(2z2 11z 12)

（1）X(z) 2 （2）X(z)

z 5z 6(z 1)(z 2)3

解：（1）z反变换MATLAB源程序为

>>Z=sym('(8*z-19)/(z^2-5*z+6)'); >>x=iztrans(Z); >>simplify(x) ans=

-19/6*charfcn[0](n)+5*3^(n-1)+3*2^(n-1)

其中，charfcn[0](n)是 (n)函数在MATLAB符号工具箱中的表示，反变换后的函数形式为

x(n)

19

(n) (5 3n 1 3 2n 1)u(n)。 6

（2）z反变换MATLAB源程序为

>>Z=sym('z*(2*z^2-11*z+12)/(z-1)/(z-2)^3'); >>x=iztrans(Z); >>simplify(x) ans=

-3+3*2^n-1/4*2^n*n-1/4*2^n*n^2

其函数形式为x(n) ( 3 3 2

n

1n12n

n2 n2)u(n)。 44

如果信号的z域表示式X(z)是有理函数，进行z反变换的另一个方法是对X(z)进行部分分式展开，然后求各简单分式的z反变换。设X(z)的有理分式表示为

b0 b1z 1 b2z 2 bmz mB(z)

（4-3） X(z) 1 2 n

A(z)1 a1z a2z anz

MATLAB信号处理工具箱提供了一个对X(z)进行部分分式展开的函数residuez，其语句格式为

[R,P,K]=residuez(B,A)

其中，B，A分别表示X(z)的分子与分母多项式的系数向量；R为部分分式的系数向量；P为极点向量；K为多项式的系数。若X(z)为有理真分式，则K为零。

Y(z)

从

B(z)r(1)r(2)r(N) 1

... k(1) k(2)z ... 1 1 1

A(z)1 p(1)z1 p(2)z1 p(N)z

得

出

其

时

域

信

号

为

：

而

离散时间LTI系统的时域分析(1)

y(n) r(1)[p(1)]nu(n) r(2)[p(2)]nu(n) ... r(N)[p(N)]nu(n) k(1) (n) k(2) (n 1) ...

【实例4-3】 试用MATLAB命令对函数X(z) 展开，并求出其z反变换。

解：MATLAB源程序为

>>B=[18]; >>A=[18,3,-4,-1]; >>[R,P,K]=residuez(B,A) R= 0.3600 0.2400 0.4000 P= 0.5000 -0.3333 -0.3333 K= []

18进行部分分式

18 3z 1 4z 2 z 3

从运行结果可知，p2 p3，表示系统有一个二重极点。所以，X(z)的部分分式展开为

X(z)

因此，其z反变换为

0.36 0.240.4

1 0.5z 11 0.3333z 1（1 0.3333z 1)2

x(n) [0.36 (0.5)n 0.24 ( 0.3333)n 0.4(n 1)( 0.3333)n]u(n)

4.2.2 系统函数的零极点分析

离散时间系统的系统函数定义为系统零状态响应的z变换与激励的z变换之比，即

H(z)

如果系统函数H(z)的有理函数表示式为

Y(z)

（4-4） X(z)

b1zm b2zm 1 bmz bm 1

（4-5） H(z)

a1zn a2zn 1 anz an 1

那么，在MATLAB中系统函数的零极点就可通过函数roots得到，也可借助函数tf2zp得到。 1， roots的格式语句为：p=roots(A),其中A为待求根的多项式的系数构成的行向量，返回向量p则包含该多项式所有的根位置列向量。 2，tf2zp的语句格式为

离散时间LTI系统的时域分析(1)

[Z,P,K]=tf2zp(B,A)

其中，B与A分别表示H(z)的分子与分母多项式的系数向量。它的作用是将H(z)的有理分式表示式转换为零极点增益形式，即

H(z) k

(z z1)(z z2) (z zm)

（4-6）

(z p1)(z p2) (z pn)

z 0.32

2

z z 0.16

【实例4-4】 已知一离散因果LTI系统的系统函数为

H(z)

试用MATLAB命令求该系统的零极点。

解：1，用roots函数求系统的零极点： a=[1 0.32]; b=[1 1 0.16]; r=roots(a) r =

-0.3200 p=roots(b) p =

-0.8000 -0.2000

2，用tf2zp函数求系统的零极点，MATLAB源程序为

>>B=[1,0.32]; >>A=[1,1,0.16]; >>[R,P,K]=tf2zp(B,A) R= -0.3200 P= -0.8000 -0.2000 K= 1

因此，零点为z 0.32，极点为p1 0.8与p2 0.2。

若要获得系统函数H(z)的零极点分布图，可直接应用zplane函数，其语句格式为

zplane(B,A)

其中，B与A分别表示H(z)的分子和分母多项式的系数向量。值注意的是：求系统函数零极点时，离散系统的系统函数可能有两种形式，一种是分子分母多项式按z的降幂次序排列，另一种是按z的升幂次序排列。若是以z的降幂次序排列，则系数向量一定要由多项

1

离散时间LTI系统的时域分析(1)

式的最高幂次开始，一直到常数项，缺项要用0补齐；若以z的升幂次序排列，则分子分母多项式系数向量的维数一定要相同，不足的要用0补齐，否则零点或极点有可能被漏掉。

【实例4-5】 已知一离散因果LTI系统的系统函数为

1

z2 0.361 0.36z 2

H(z) 2

z 1.52z 0.681 1.52z 1 0.68z 2

试用MATLAB命令绘出该系统的零极点分布图。

解：用zplane函数求系统的零极点，MATLAB源程序为

>>B=[1,0,-0.36]; >>A=[1,-1.52,0.68]; >>zplane(B,A),grid on >>legend('零点','极点') >>title('零极点分布图')

程序运行结果如图14-1所示。可见，该因果系统的极点全部在单位圆内，故系统是稳定的。

图4-1 零极点分布图

4.2.3 系统函数的零极点分布与其时域特性的关系

与拉氏变换在连续系统中的作用类似，在离散系统中，z变换建立了时域函数h(n)与z域函数H(z)之间的对应关系。因此，z变换的函数H(z)从形式可以反映h(n)的部分内在性质。我们仍旧通过讨论H(z)的一阶极点情况，来说明系统函数的零极点分布与系统时域特性的关系。

【实例14-6】 试用MATLAB命令画出现下列系统函数的零极点分布图、以及对应的时域单位取样响应h(n)的波形，并分析系统函数的极点对时域波形的影响。

（1）H1(z)

zzz

（2）H2(z) （3）H

3(z) 2

z 0.8z 0.8z 1.2z 0.72

离散时间LTI系统的时域分析(1)

zzz （5）H5(z) 2 （6）H6(s) z 1z 1.2z 1.6z 1

z

（7）H7(z) 2

z 2z 1.36

（4）H4(z) 解：MATLAB源程序为

>>b1=[1,0]; >>a1=[1,-0.8]; >>subplot(121) >>zplane(b1,a1)

>>title('极点在单位圆内的正实数') >>subplot(122)

>>impz(b1,a1,30);grid on; >>figure >>b2=[1,0]; >>a2=[1,0.8]; >>subplot(121) >>zplane(b2,a2)

>>title('极点在单位圆内的负实数') >>subplot(122)

>>impz(b2,a2,30);grid on; >>figure >>b3=[0,1,0]; >>a3=[1,-1.2,0.72]; >>subplot(121) >>zplane(b3,a3)

>>title('极点在单位圆内的共轭复数') >>subplot(122)

>>impz(b3,a3,30);grid on; >>figure >>b4=[1,0]; >>a4=[1,-1]; >>subplot(121) >>zplane(b4,a4)

>>title('极点在单位圆上为实数1') >>subplot(122) >>impz(b4,a4);grid on; >>figure >>b5=[0,1,0]; >>a5=[1,-1.6,1]; >>subplot(121) >>zplane(b5,a5)

>>title('极点在单位圆上的共轭复数') >>subplot(122)

>>impz(b5,a5,30);grid on;

离散时间LTI系统的时域分析(1)

>>figure >>b6=[1,0]; >>a6=[1,-1.2]; >>subplot(121) >>zplane(b6,a6)

>>title('极点在单位圆外的正实数') >>subplot(122)

>>impz(b6,a6,30);grid on; >>figure >>b7=[0,1,0]; >>a7=[1,-2,1.36]; >>subplot(121) >>zplane(b7,a7)

>>title('极点在单位圆外的共轭复数') >>subplot(122)

>>impz(b7,a7,30);grid on;

程序运行结果分别如图14-2的（a）、（b）、（c）、（d）、（e）、（f）、（g）所示。

(a)

离散时间LTI系统的时域分析(1)

(b)

(c)

(d)

(e)

离散时间LTI系统的时域分析(1)

(f)

(g)

图4-2 系统函数的零极点分布与其时域特性的关系

从图14-2可知，当极点位于单位圆内时，h(n)为衰减序列；当极点位于单位圆上时，

h(n)为等幅序列；当极点位于单位圆外时，h(n)为增幅序列。若h(n)有一阶实数极点，

则h(n)为指数序列；若h(n)有一阶共轭极点，则h(n)为指数振荡序列；若h(n)的极点位于虚轴左边，则h(n)序列按一正一负的规律交替变化。

4.2.4 离散时间LTI系统的频率特性分析

对于因果稳定的离散时间系统，如果激励序列为正弦序列x(n) Asin(n )u(n)，则系

离散时间LTI系统的时域分析(1)

统的稳态响应为yss(n) A|H(ej )|sin[n ( )]u(n)。其中，H(ej )通常是复数。离散时间系统的频率响应定义为

H(ej ) |H(ej )|ej ( ) （14-7）

其中，|H(ej )|称为离散时间系统的幅频特性； ( )称为离散时间系统的相频特性；

H(ej )是以 s（ s

2

，若零T 1， s 2 ）为周期的周期函数。因此，只要分T

析H(ej )在| | 范围内的情况，便可分析出系统的整个频率特性。

MATLAB提供了求离散时间系统频响特性的函数freqz，调用freqz的格式主要有两种。一种形式为

[H,w]=freqz(B,A,N) 其中，B与A分别表示H(z)的分子和分母多项式的系数向量，值得注意的是，B，A采用z

1

的升次幂形式；N为正整数，默认值为512；返回值w包含[0, ]范围内的N个频率等分点；返回值H则是离散时间系统频率响应H(ej )在0~ 范围内N个频率处的值。另一种形式为

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

与第一种方式不同之处在于角频率的范围由[0, ]扩展到[0,2 ]。

z2 0.96z 0.9028

【实例4-6】 用MATLAB命令绘制系统H(z) 2的频率响应曲线。

z 1.56z 0.8109

解：利用函数freqz计算出H(e

j

)，然后利用函数abs和angle分别求出幅频特性与相

频特性，最后利用plot命令绘出曲线。MATLAB源程序为

>>b=[1 -0.96 0.9028]; >>a=[1 -1.56 0.8109];

>>[H,w]=freqz(b,a,400,'whole'); >>Hm=abs(H); >>Hp=angle(H); >>subplot(211) >>plot(w,Hm),grid on

>>xlabel('\omega(rad/s)'),ylabel('Magnitude') >>title('离散系统幅频特性曲线') >>subplot(212) >>plot(w,Hp),grid on

>>xlabel('\omega(rad/s)'),ylabel('Phase') >>title('离散系统相频特性曲线')

程序运行结果如图4-2所示。

离散时间LTI系统的时域分析(1)

图4-3 离散系统频响特性曲线

4.3 编程练习

1.

2z4 16z3 44z2 56z 32

试用MATLAB的residuez函数，求出X(z) 的部分分

3z4 3z3 15z2 18z 12

式展开和。

试用MATLAB画出下列因果系统的系统函数零极点分布图，并判断系统的稳定性。

2.

2z2 1.6z 0.9

（1）H(z) 3 2

z 2.5z 1.96z 0.48

（2）H(z)

z 1

432

z 0.9z 0.65z 0.873z

3.

z2

试用MATLAB绘制系统H(z) 的频率响应曲线。

31z2 z

48