matlab求解零状态零输入响应

1. 已知离散时间系统的差分方程为：
2y(n) - y(n -1) - 3y(n -2)=2x(n) - x(n -1)
x(n)=0.5n u(n) , y(-1)=1,y(-2)=3 , 试用filter 函数求系统的零输入响应、零状态响应和全响应.
解：将差分方程Z 变换得：
12112()[()(1)]3[()(1)(2)]2()[()(1)]Y z z Y z y z Y z z y y X z z X z x -----+--+-+-=-+- (1)
依题意有：x(-1)=0,x(-2)=0,y(-1)=1,y(-2)=3 ,X(z)=
1110.50.5
z z z -=-- 将上式变形如下： 1211(23)()[(1)3(1)3(2)](2)()z z Y z y z y y z X z --------+-+-=- ………..(2) 1211(23)()(2)()[(1)3(1)3(2)]z z Y z z X z y z y y ------=-+-+-+-
1211(23)()(2)()[103]z z Y z z X z z ------=-++ (3)
易得系统函数为H(z)= 12122222323
z z z z z z z -----=---- ① 零输入时
零输入时，x(n)=0,差分方程右边为0，z 变换后应为
121(23)()103z z Y z z -----=+
1
12
103()23z Y z z z ---+=-- =2210323
z z z z +-- =7183
5152
z z z z ++- 将Y(z)进行Z 反变换，得到其零输入响应为： y(n)= 7
183[(1)()]()552
n n u n -+ ② 零状态时
零状态时，将y(-1)=0,y(-2)=0代入上面的式（2）中，得
Y(z)= 112223z z z ------X(z)= 112223z z z ------1110.5z --=22223
z z z --
=23
3 515
2 z z z z
+
+-
将其Z反变换，得到零状态响应为：
y(n)=
233 [(1)()]() 552
n n u n -+
③全响应
与上面同理，y(-1)=1,y(-2)=3将上面式(3)变形得：
Y(z)=
2
2
123
23
z z
z z
+
--
=
921
3
515
2
z z
z z
+
+-
Z反变换得全响应为
Y(n)=
921
[]()
3
515
2
z z
u n z z
+
+-
程序代码：
%第二章Z变换第2.12题程序
clear all;close all;
num=[2 -1 0]; %系统函数分子的系数
den=[2 -1 -3]; %系统函数分母的系数
n=0:50;
nl=length(n);
%求零输入响应
y01=[1 3]; %y的初始状态
x01=[0 0]; %x 的初始状态
x1=zeros(1,nl);
zi1=filtic(num,den,y01,x01); %为filter函数准备初始值y1=filter(num,den,x1,zi1); %求零输入响应
subplot(311);
stem(n,y1,'r.');
title('零输入响应');
grid on;
%求零状态响应
y02=[0 0];
x02=[0 0];
x2=0.5.^n;
zi2=filtic(num,den,y02,x02);
y2=filter(num,den,x2,zi2);
subplot(312);
stem(n,y2,'r.');
title('零状态响应');
grid on;
%求全响应
y03=[1 3];
x03=[0 0];
x3=0.5.^n;
zi3=filtic(num,den,y03,x03);
y3=filter(num,den,x1,zi3);
subplot(313);
stem(n,y3,'r.');
title('全响应');
grid on;
运行结果如下：
2. 已知离散系统的系统函数分别为
(1)
2
3
21
()
21
z z
H z
z
--
=
-
(2)
3
1
()
1
z
H z
z
+
=
-
(3)
2
32
2
()
2241
z
H z
z z z
+
=
+-+
(4)
3
32
()
0.20.30.4
z
H z
z z z
=
+++
试用MATLAB实现下列分析过程：
①求出系统的零极点位置；
②绘出系统的零极点图，根据零极点图判断系统的稳定性；
③绘出系统单位响应的时域波形，并分析系统稳定性与系统单位响应时域特性的关系。

解：程序代码如下：
%%第二章Z变换第2.13题程序
clear all;close all;
%题(1)
a1=[2 0 0 -1]; %系统函数分母的系数
b1=[0 2 -2 -1]; %系统函数分子的系数
p1=roots(a1), %求极点
pa1=abs(p1), %求极点到坐标原点的距离，看它是否大于1,若有一个大于1, %则系统不稳定；若所有的都小于1，则系统稳定
q1=roots(b1), %求零点
h1=impz(b1,a1); %求单位响应
subplot(421);
zplane(b1,a1);%画零极点图
title('(1)的零极点图');
subplot(425);
stem(h1,'.'); %单位响应的时域波形
grid on;
title('(1)的单位响应的时域波形');
%题(2)
a2=[3 0 0 -1];
b2=[0 0 1 1];
p2=roots(a2),
pa2=abs(p2),
q2=roots(b2),
h2=impz(b2,a2);
subplot(422);
zplane(b1,a1);
title('(2)的零极点图');
subplot(426);
stem(h2,'.');
grid on;
title('(2)的单位响应的时域波形');
%题(3)
a3=[1 2 -4 1];
b3=[0 1 0 2];
p3=roots(a3),
pa3=abs(p3),
q3=roots(b1),
h3=impz(b3,a3);
subplot(423);
zplane(b3,a3);
title('(3)的零极点图');
subplot(427);
stem(h3,'.');
grid on;
title('(3)的单位响应的时域波形');
%题(4)
a4=[1 0 0 0];
b4=[1 0.2 0.3 0.4];
p4=roots(a4),
pa4=abs(p4),
q4=roots(b4),
h4=impz(b4,a4);
subplot(424);
zplane(b1,a1);
title('(1)的零极点图');
subplot(428);
stem(h4,'.');
grid on;
title('(1)的单位响应的时域波形');
运行结果如下：
3. 已知描述离散系统的差分方程为：
y(n) - y(n-1) - y(n-2)=4x(n) - x(n-1) - x(n-2)
试用MATLAB绘出系统的零极点分布图，并绘出系统的幅频和相频特性曲线，分析该系统的作用
解：
程序代码如下：
clear all;close all;
num=[4,-1,-1];
den=[1 -1 -1];
[H,w]=freqz(num,den);
subplot(311);
zplane(num,den);
subplot(312);
plot(w/pi,abs(H));
grid on;
title('幅频响应曲线')
subplot(313);
plot(w/pi,angle(H));
title('相频响应曲线');
grid on;
运行结果如下：
4. 已知因果（单边）离散序列的Z变换分别如下所示，试用MATLAB求出其Z反变换
(1)
2
2
1
()
2
z z
F z
z z
++
=
+-
(2)
2
32
21
()
1
2
z z
F z
z z z
-+
=
++
(3)
2
()
F z=(4)
32
432
21
()
32321
z z z
F z
z z z z
+++
=
++++
解：
程序代码如下：
clear all;close all;
F1=sym('(z^2+z+1)/(z^2+z-2)'); f1=iztrans(F1),
F2=sym('(2*z^2-z+1)/(z^3+z^2+z/2)');
f2=iztrans(F2),
F3=sym('(z^2)/(z^2+sqrtm(2)*z+1)');
f3=iztrans(F3),
F4=sym('(z^3+2*z^2+z+1)/(3*z^4+2*z^3+3*z^2+2*z+1)');
f4=iztrans(F4)
运行结果如下：
f1 =
(-2)^n/2 - kroneckerDelta(n, 0)/2 + 1
δ
注：kroneckerDelta(n, 0)=()n
f2 =
2*kroneckerDelta(n -1, 0) -6*kroneckerDelta(n, 0) + 3*(-1)^n*2^(1 -n)*i*(i + 1)^(n -1) -3*(-1)^n*2^(1 - n)*i*(1 - i)^(n - 1)
f3 =
2*(-1)^n*cos(n*acos(sqrtm(2)/2)) + ((-1)^n*(sqrtm(2)/2 + (sqrtm(2)^2/4 -1)^(1/2))^(n -1))/(2*(sqrtm(2)^2/4 -1)^(1/2)) -((-1)^n*(sqrtm(2)/2 -(1/4*sqrtm(2)^2 -1)^(1/2))^(n -1))/(2*(sqrtm(2)^2/4 - 1)^(1/2))
f4 =
sum(-(r3*r3^n + r3^n + 2*r3^2*r3^n + r3^3*r3^n)/(2*r3^3 + 6*r3^2 + 6*r3 + 4), r3 in RootOf(z1^4 + (2*z1^3)/3 + z1^2 + (2*z1)/3 + 1/3, z1)) + kroneckerDelta(n, 0)
sum( -(r3*r3^n + r3^n + 2*r3^2*r3^n + r3^3*r3^n)/(2*r3^3 + 6*r3^2 + 6*r3 + 4), r3 in RootOf(z1^4 + (2*z1^3)/3 + z1^2 + (2*z1)/3 + 1/3, z1) ) + kroneckerDelta(n, 0)
注：
r3 in RootOf(z1^4 + (2*z1^3)/3 + z1^2 + (2*z1)/3 + 1/3, z1)
就是说r3是关于Z1的方程z1^4 + (2*z1^3)/3 + z1^2 + (2*z1)/3 + 1/3=0的根。

sum( -(r3*r3^n + r3^n + 2*r3^2*r3^n + r3^3*r3^n)/(2*r3^3 + 6*r3^2 + 6*r3 + 4), r3 in RootOf(z1^4 + (2*z1^3)/3 + z1^2 + (2*z1)/3 + 1/3, z1) )就是将上面方程的每个根（即r3的值）代入-(r3*r3^n + r3^n + 2*r3^2*r3^n + r3^3*r3^n)/(2*r3^3 + 6*r3^2 + 6*r3 + 4)，然后相加。