电子科技大学 信号与系统 软件实验3

电子科技大学
实验报告
学生姓名：xxx 学号：2901305032 指导教师：崔琳莉
一、实验室名称：信号与系统实验室
二、实验项目名称：离散系统的转移函数，零、极点分布和模拟
三、实验原理：
离散系统的时域方程为
∑∑
= =-
= -
M
m m
N
k
k
m
n
x
b
k
n
y
a
]
[
]
[
其变换域分析方法如下：
系统的频率响应为
ω
ω
ω
ω
ω
ω
ω
jN
N
j
jM
M
j
j
j
j
e
a
e
a
a
e
b
e
b
b
e
A
e
B
e
H
-
-
-
-
+
+
+
+
+
+
=
=
...
...
)
(
)
(
)
(
1
1
Z域
)
(
)
(
)
(
]
[
]
[
]
[
]
[
]
[z
H
z
X
z
Y
m
n
h
m
x
n
h
n
x
n
y
m
=
⇔
-
=
*
=∑∞
-∞
=
系统的转移函数为
N
N
M
M
z
a
z
a
a
z
b
z
b
b
z
A
z
B
z
H
-
-
-
-
+
+
+
+
+
+
=
=
...
...
)
(
)
(
)
(
1
1
1
1
分解因式
∏
∏
∑
∑
=
-
=
-
=
-
=
-
-
-
=
=
N
i
i
M
i
i
N
i
i
k
M
i
i
k
z
z
K
z
a
z
b
z
H
1
1
1
1
)
1(
)
1(
)
(
λ
ξ
，其中i
ξ
和i
λ
称为
零、极点。

在MATLAB中，可以用函数[z,p,K]=tf2zp（num,den）求得有理分式形式的系统转移函数的零、极点，用函数zplane（z，p）绘出零、极点分布图；也可以用函数zplane（num，den）直接绘出有理分式形式的系统转移函数的零、极
点分布图。

四、实验目的：
1、加深对离散系统转移函数、零极点概念的理解；
2、根据系统转移函数求系统零极点分布。

五、实验内容：
MATLAB 仿真
六、实验器材（设备、元器件）：
计算机、MATLAB 软件。

七、实验步骤：
对系统系统
2181.09.011
)(--+-=z z z H
1、 编程实现系统的参数输入，绘出幅度频率响应曲线和零、极点分布图。

2、 根据系统的零极点计算系统频率响应的幅值和相位。

定义
omega=[0:511]*pi/256和unitcirc=exp(j*omega)得到在单位圆上512个等分点，
在这些点上将要对频率响应
)(jw e H 求值。

（a ）定义polevectors1是一个2×512的矩阵，其中每一行包含这样一些复数，这些复数是由unitcirc 的相应列减去一个极点位置得到的。

如果ps1是一个列向量，它包含了极点的位置，可以用下列命令来完成
》polevectors1 = ones(2, 1) * unitcirc – ps1 * ones(1, 512)
用abs 和astan2定义polelength1和poleangle1作为polevectors1中每一元素的幅值和相位。

（b ）类似与polevectors1定义zerovectors1，使得它是2×512的矩阵，其
中包含从零点位置到unitcirc元素的向量。

定义zerolength1和zeroangle1分别是这些向量的幅值和相位。

（c）画出polelength1和zerolength1对于omega的图。

根据这些图，预
计
)
(
1
jw
e
H
在哪里有最大值和最小值？
（d）利用MATLAB命令H = freqz (b,a,512,’whole’) 导出幅值和相位，
与（c）中的结果相比较。

八、实验数据及结果分析：
(a)
b=[1];
a=[1 -0.9 0.81]
la=length(a);
lb=length(b);
if(la>lb)
b=[b zeros(1,la-lb)];
elseif(lb>la)
a=[a zeros(1,lb-la)];
end
ps=roots(a);
zs=roots(b);
mx=max(abs([ps+0.95zs+0.95]))+0.05; clg
axis('equal');
hold on
w=[0:0.01:2*pi];
plot(cos(w),sin(w),'.');
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(-0.1, -0.1, num2str(numz)); elseif nump>1
text(-0.1,-0.1, num2str(nump)); end
hold off
omega=[0:511]*pi/256;
un=exp(j*omega);
ps1=roots(a);
polev1=ones(2,1)*un-ps1*ones(1,512);
(b)
b=[1];
a=[1 -0.9 0.81];
la=length(a);
lb=length(b);
if(la>lb)
b=[b zeros(1,la-lb)];
elseif(lb>la)
a=[a zeros(1,lb-la)];
end
ps=roots(a);
zs=roots(b);
mx=max(abs([ps+0.95zs+0.95]))+0.05; clg
axis('equal');
w=[0:0.01:2*pi];
omega=[0:511]*pi/256;
un=exp(j*omega);
ps1=roots(a);
polev1=ones(2,1)*un-ps1*ones(1,512); polelength1=abs(polev1);
poleangle1=atan(polev1)
(c)
b=[1];
a=[1 -0.9 0.81];
la=length(a);
lb=length(b);
if(la>lb)
b=[b zeros(1,la-lb)];
elseif(lb>la)
a=[a zeros(1,lb-la)];
end
ps=roots(a);
zs=roots(b);
mx=max(abs([ps+0.95zs+0.95]))+0.05; clg
axis('equal');
w=[0:0.01:2*pi];
omega=[0:511]*pi/256;
un=exp(j*omega);
ps2=roots(b);
zerov=ones(2,1)*un-ps2*ones(1,512); zerolength1=abs(zerov)
zeroangle1=atan(zerov)
(d)
b=[1];
a=[1 -0.9 0.81];
la=length(a);
lb=length(b);
if(la>lb)
b=[b zeros(1,la-lb)];
elseif(lb>la)
a=[a zeros(1,lb-la)];
end
ps=roots(a);
zs=roots(b);
mx=max(abs([ps+0.95zs+0.95]))+0.05; clg
axis('equal');
w=[0:0.01:2*pi];
omega=[0:511]*pi/256;
un=exp(j*omega);
ps2=roots(b);
zerov=ones(2,1)*un-ps2*ones(1,512); zerolength1=abs(zerov);
zeroangle1=atan(zerov);
ps1=roots(a);
polev1=ones(2,1)*un-ps1*ones(1,512); polelength1=abs(polev1);
poleangle1=atan(polev1);
plot(omega,polelength1)
plot(omega,zerolength1)
(e)
b=[1];
a=[1 -0.9 0.81];
la=length(a);
lb=length(b);
if(la>lb)
b=[b zeros(1,la-lb)];
elseif(lb>la)
a=[a zeros(1,lb-la)];
end
ps=roots(a);
zs=roots(b);
mx=max(abs([ps+0.95zs+0.95]))+0.05; clg
axis('equal');
w=[0:0.01:2*pi];
omega=[0:511]*pi/256;
un=exp(j*omega);
ps2=roots(b);
zerov=ones(2,1)*un-ps2*ones(1,512);
zerolength1=abs(zerov);
zeroangle1=atan(zerov);
ps1=roots(a);
polev1=ones(2,1)*un-ps1*ones(1,512); polelength1=abs(polev1);
poleangle1=atan(polev1);
geonH1mag=1./prod(polelength1)
plot(omega,geonH1mag)
geomH1phase=sum(zeroangle1)-sum(poleangle1); plot(omega,geomH1phase)
H1=freqz(b,a,512,'whole')
plot(omega,H1)
对比发现H1和所画geonH1mag图象一致。