测控专业计算机仿真实验指导书
实验一
123，，(1)已知矩阵，b=[7 8;9 10;11 12]，利用MATLAB进行a,,,456，，
矩阵计算:a*b a/b a\b a^3 a.*b a./b a.\b a.^3
解:
>> a=[1 2 3;4 5 6]
a =
1 2 3
4 5 6
>> b=[7 8;9 10;11 12]
b =
7 8
9 10
11 12
>> a*b
ans =
58 64
139 154
>> a=[1 2 3;4 5 6]
a =
1 2 3
4 5 6
>> b=[7 8 9;10 11 12]
b =
1
7 8 9
10 11 12
>> a/b
ans =
3.0000 -2.0000
2.0000 -1.0000
>>
(2)多项式计算:求[12+2×(7-4)]?32的算数计算结果。
解:>> a=12;
>> b=2*(7-4);
>> c=(a+b)/32
c =
0.5625
2(编写M文件绘制函数的图形。 【例1-1】:绘制下面函数在[0 3]区间的图形。
sinx, x,0,
,y(x),x, 0,x,3 ,
,,x，6, x,3,
MATLAB程序:
x=-6:0.1:6;
leng=length(x);
for m=1:leng
if x(m)<=0
y(m)=sin(x(m));
elseif x(m)<=3
y(m)=x(m);
else
2
y(m)=-x(m)+6;
end
end
plot(x,y,'*'),grid;
练习:绘制下面的图形:(1)sin(1/t),-1> t=-0.1:0.0001:0.1; >> leng=length(t);
>> for m=1:leng
if t(m)<=-1
y(m)=0;
elseif t(m)>=1
y(m)=0;
else
y(m)=sin(1/t(m));
end
end
% Warning: Divide by zero. >> plot(t,y,'*'),grid;
3
(2)1-cos3(7t)
>> t=-1:0.0001:1;
>> leng=length(t);
>> for m=1:leng
if t(m)<=0
y(m)=1-cos(3*(7*t(m)));
else
y(m)=1-cos(3*(7*t(m))); end
end
>> plot(t,y,'*'),grid; >>
4
3(建立Simulink模型，进行系统仿真。
【例1-2】:信号发生器发出幅值为1，频率为0.2Hz(下图实为0.5Hz)的正弦信号。信号分两路输出:一路直接送到示波器，一路扩大5倍送到另外的示波器。
启动Simulink，建立如图模型，并运行，观察运行结果。
5
解:
练习:(1)应用Simulink建立系统模型，并对系统的阶跃响应进行仿真。
1G(s), 2s，4s，8
解:
6
(2)建立一个简单的模型，用信号发生器产生一个幅值为2V，频率为0.5Hz的正弦波，并叠加一个0.1V的噪声信号，将叠加后的信号显示在示波器上。
7
8
实验2
f(t),sin(,t，,)(1)正弦交流信号
实现 ，MATLAB程序如下:
t=-0:0.001:1;
y=sin(2*pi*t);
plot(t,y,'k');
xlabel('时间(t)');ylabel('幅值(f)');title('正弦交流信号');
f(t),,(t)(2)单位阶跃信号
MATLAB程序如下:
t=-2:0.01:6;
u=(t>=0);
plot(t,u);
axis([-2,6,0,1.2])
xlabel('时间(t)');ylabel('幅值(f)');title('单位阶跃信号');
,3t，4jtf(t),e(3)复指数信号
MATLAB程序如下:
t=0:0.01:3;
a=-3;b=4;
z=exp((a+i*b)*t);
subplot(221);
plot(t,imag(z)),title('实部');xlabel('时间');ylabel('幅值');
f(t),sin(2,t)subplot(222);
plot(t,real(z)),title('虚部');xlabel('时间');ylabel('幅值'); subplot(223);
plot(t,abs(z)),title('模');xlabel('时间');ylabel('幅值'); subplot(224);
plot(t,angle(z)),title('相角');xlabel('时间');ylabel('幅值');
9
2.连续时间信号的基本运算
(1)加(减)、乘运算
实现两个信号的加法和乘法运算:要求两个信号运算的时间序列长度相同。
MATLAB程序如下:
t=0:0.01:2;
f1=exp(-3*t);
f2=0.2*sin(4*pi*t);
f3=f1+f2;
f4=f1.*f2;
subplot(221);plot(t,f1);title('f1(t)'); subplot(222);plot(t,f2);title('f2(t)'); subplot(223);plot(t,f3);title('f1+f2'); subplot(224);plot(t,f4);title('f1*f2');
(2)信号的反褶、移位、尺度变换
f(t)f(,at，b)(a,0)由到步骤:
移位尺度反褶f(t),,,,f(t，b),,,,f(at，b),,,,f(,at，b)
f(t),sin(t)/tf(t)【例2-1】已知。试通过反褶、移位、尺度变换由
f(,2t，3)的波形得到的波形。
MATLAB程序如下
syms t;
f=sym(‘sin(t)/t’);
f1=subs(f,t,t+3);
f2=subs(f1,t,2*t);
f3=sub(f2,t,-t);
subplot(2,2,1);ezplot(f,[-8,8]);grid on; subplot(2,2,2);ezplot(f1,[-8,8]);grid on; subplot(2,2,3);ezplot(f2,[-8,8]);grid on; subplot(2,2,4);ezplot(f3,[-8,8]);grid on;
10
(3)连续信号的卷积运算
y,conv(x,h)用于实现x、h两个序列的卷积，假定都是从n=0开始。Y序列的长度为x、h序列的长度之和再减1.
【例2-2】求两个方波信号的卷积。
MATLAB程序如下:
y1=[ones(1,20),zeros(1,20)];
y2=[ones(1,10),zeros(1,20)];
y=conv(y1,y2);
n1=1:length(y1);
n2=1:length(y2);
L=length(y);
subplot(3,1,1);plot(n1,y1);axis([1 L 0 2]); subplot(3,1,2);plot(n2,y2);axis([1 L 0 2]); n=1:L;
subplot(3,1,3);plot(n,y),axis([1 L 0 20]);
【例2-3】求两个指数信号的卷积。
MATLAB程序如下:
11
t=0:0.01:1;
y1=exp(-6*t);
y2=exp(-3*t);
y=conv(y1,y2);
L1=length(y1);
L2=length(y2);
L=length(y);
subplot(3,1,1);plot(t,y1); subplot(3,1,2);plot(t,y2); t1=0:0.01:2;
subplot(3,1,3);plot(t1,y);
,(t,3),cos3t，sin2t1. 编制程序，生成如下连续信号:
解:>> t=-2:0.01:6;
>> u=(t>=3);
>> plot(t,u);
12
>> t=-2:0.01:6; >> u=(t>=3);
>> plot(t,u); >> t=-0:0.001:1; >> y=sin(2*t); >> z=cos(3*t); >> w=y+z;
>> plot(w);
13
14
实验3
''''的冲激 【例3-1】求系统(t)，(t)，(t),(t)，(t)y6y3y3x9x响应和阶跃响应。
(1)系统的冲激响应的MATLAB程序如下:
b=[3,9];a=[1,6,8];
sys=tf(b,a)
t=0:0.1:10;
y=impulse(sys,t);
plot(t,y);
xlabel('time(t)');ylabel('y(t)');title('单位冲激响应');
(2)系统的阶跃响应的MATLAB程序如下:
b=[3,9];a=[1,6,8];
sys=tf(b,a)
15
t=0:0.1:10;
y=step(sys,t);
plot(t,y);
xlabel('time(t)');ylabel('y(t)');title('单位阶跃响应');
'''，，【例3-2】求系统的全响(t)，(t),cos(t),y(),(),0yyy00应。
MATLAB程序如下:
b=[1];a=[1,0,1];
sys=tf(b,a)
t=0:0.1:10;
x=cos(t);
y=lsim(sys,x,t);
plot(t,y);
xlabel('time(t)');ylabel('y(t)');title('零状态响应');
16
利用函数直接求的程序:
b=[1];a=[1,0,1];
[A B C D]=tf2ss(b,a);
sys=ss(A,B,C,D);
t=0:0.1:10;
x=cos(t);
y=lsim(sys,x,t);
plot(t,y);
xlabel('time(t)');ylabel('y(t)');title('系统响应');
17
?琠????????氠湥?敬杮桴琨??‾潦???敬杮植?琠洨??礍?氠湥?敬杮桴琨??‾潦???敬杮植?琠洨??礍洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩????? ?‾潦???敬杮植?琠洨??礍洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩?????
‾潦???敬杮植?琠洨??礍洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩?????
植?琠洨??礍洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩?????
?琠洨??礍洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩?????
礍洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩?????
洨??潣???琪洨??攠獬????挭獯???????湥?湥??†汰瑯琨礬???牧摩?????
????
???
??
18
?
19
20
21
0
0週
0週?
0週?
0週?倲
0週?倲
0週?倲 ?
0週?倲 ??
0週?倲 ??
0週?倲 ???
0週?倲 ????
22
0週?倲 ????
0週?倲 ????
0週?倲 ????
0週?倲 ????????
0週?倲 ????????
0週?倲 ????????連
0週?倲 ????????連?
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤
0週?倲 ????????連?逤?? ????
0週?倲 ????????連?逤?? ????
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌
0週?倲 ????????連?逤?? ????逌0週?倲 ????????連?逤?? ????逌週?倲 ????????連?逤?? ????逌?倲 ????????連?逤?? ????逌倲 ????????連?逤?? ????逌? ? ? ? ? ú
????????連?逤?? ????逌? ? ? ? ? ú ,????????連?逤?? ????逌? ? ? ? ? ú ,???????連?逤?? ????逌? ? ? ? ? ú ,??????連?逤?? ????逌? ? ? ? ? ú ,
23
?????連?逤?? ????逌? ? ? ? ? ú ,????連?逤?? ????逌? ? ? ? ? ú ,???連?逤?? ????逌? ? ? ? ? ú ,??連?逤?? ????逌? ? ? ? ? ú ,?連?逤?? ????逌? ? ? ? ? ú ,連?逤?? ????逌? ? ? ? ? ú ,?逤?? ????逌? ? ? ? ? ú , ú 逤?? ????逌? ? ? ? ? ú , ú ,?? ????逌? ? ? ? ? ú , ú ,? ????逌? ? ? ? ? ú , ú ,
????逌? ? ? ? ? ú , ú ,????逌? ? ? ? ? ú , ú ,???逌? ? ? ? ? ú , ú ,??逌? ? ? ? ? ú , ú ,?逌? ? ? ? ? ú , ú ,逌? ? ? ? ? ú , ú ,? ? ? ? ? ú , ú ,
? ? ? ? ú , ú , ú ? ? ? ? ú , ú , ú ,
? ? ? ú , ú , ú ,? ? ? ú , ú , ú ,
? ? ú , ú , ú ,? ? ú , ú , ú ,
? ú , ú , ú ,? ú , ú , ú ,
ú , ú , ú ,ú , ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú , ú
, ú , ú , ú ,, ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
24
ú , ú , ú ,
ú , ú , ú ,ú , ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú , ú
, ú , ú , ú ,, ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,ú , ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú , ú
, ú , ú , ú ,, ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,ú , ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú , ú
, ú , ú , ú ,, ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
25
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,
ú , ú , ú ,ú , ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,
, ú , ú ,, ú , ú ,
ú , ú ,
ú , ú ,
ú , ú ,
ú , ú ,
ú , ú ,
ú , ú ,ú , ú ,
, ú ,
, ú ,
, ú ,
, ú ,
, ú ,, ú ,
ú ,
ú ,
ú ,
ú ,
ú ,
ú ,ú ,
,
,
,
,
,,
26
ú
27
ú
28
ú
29
ú
30
ú
31
ú
32
ú
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
实验4
1(周期信号的分解
【例4-1】用正弦信号的叠加近似合成一个频率为50Hz，幅度为3的方波。
MATLAB程序如下:
clear all
fs=10000;
t=[0:1/fs:0.1];
f0=50;
sum=0;
subplot(211)
for n=1:2:9
plot(t,4/pi*1/n*sin(2*pi*n*f0*t),'k');
hold on;
end
title('信号叠加前');
subplot(212)
for n=1:2:9;
sum=sum+4/pi*1/n*sin(2*pi*n*f0*t);
end
plot(t,sum,'k');
title('信号叠加后');
69
2(傅里叶变换和逆变换的实现
求傅里叶变换，可以调用forier函数，调用格式为F=fourier(f,u,v),是关于u的函数f的傅里叶变换，返回函数F是关于v的函数。
求傅里叶逆变换，可以调用iforier函数，调用格式为f=ifourier(F,u,v),是关于v的函数F的傅里叶逆变换，返回函数f是关于u的函数。
,2|t|【例4-2】已知连续信号，通过程序完成其傅里叶变换。 f(t),e
MATLAB程序如下:
syms t
f=fourier(exp(-2*abs(t)));
ezplot(f);
70
1【例4-3】已知连续信号，通过程序完成其傅里叶F(j),,21，,变换。
MATLAB程序如下:
syms t w
ifourier(1/(1+w^2),t)
3(傅里叶变换的性质
举例练习傅里叶变换的时移特性和频移特性。
1,2t【例4-4】分别绘出信号和的频谱，求f(t,1)f(t),,(t)e2
1,2t的频谱。 f(t),,(t)e2
MATLAB程序如下:
r=0.02;t=-5:r:5;
N=200;W=2*pi;
k=-N:N;w=k*W/N;
f1=1/2*exp(-2*t).*stepfun(t,0);
71
F=r*f1*exp(-j*t'*w);
F1=abs(F);P1=angle(F);
subplot(3,1,1);plot(t,f1);grid xlabel('t'),ylabel('f(t)');title('f(t)'); subplot(3,1,2);
plot(w,F1);xlabel('w');grid;ylabel('F(jw)的模');
subplot(3,1,3);
(度)'); plot(w,P1*180/pi);grid;xlabel('w');ylabel('相位
再求信号f(t-1)的频谱:
r=0.02;t=-5:r:5;N=200;W=2*pi;k=-N:N;w=k*W/N;
f1=1/2*exp(-2*(t-1)).*stepfun(t,1); F=r*f1*exp(-j*t'*w);
F1=abs(F);P1=angle(F);
subplot(3,1,1);plot(t,f1);grid xlabel('t'),ylabel('f(t)');title('f(t)'); subplot(3,1,2);
plot(w,F1);xlabel('w');grid;ylabel('F(jw)的模');
subplot(3,1,3);
plot(w,P1*180/pi);grid;xlabel('w');ylabel('相位(度)');
72
【例4-5】傅里叶变换的频移特性:信号为门信号，绘出f(t),g(t)
,j10tj10t(t),f(t)(t),f(t)ff信号和信号的频谱，并与原信号的ee12
频谱图进行比较。
(1)，求其频谱可以采用数值计算的f(t),g(t),,(t,1),,(t，1)
方法。MATLAB程序如下:
R=0.02;t=-2:R:2;
f=stepfun(t,-1)-stepfun(t,1); w1=2*pi*5;
N=500;
k=0:N;
W=k*w1/N;
F=f*exp(-j*t'*W)*R;
F=real(F);W=[-fliplr(W),W(2:501)]; F=[fliplr(F),F(2:501)]; subplot(211);plot(t,f);
73
xlabel('t'),ylabel('f(t)');axis([-2,2,-0.5,2]);
title('f(t)=u(t+1)-u(t-1)'); subplot(212);plot(W,F);
xlabel('w'),ylabel('F(w)');title('f(t)的傅里叶变换');
,j10tj10t(t),f(t)(t),f(t)ff (2)得到，的频谱的ee12
MATLAB程序如下:
R=0.02;t=-2:R:2;
f=stepfun(t,-1)-stepfun(t,1); f1=f.*exp(-j*10*t);f2=f.*exp(j*10*t); w1=2*pi*5;
N=500;
k=-N:N;
W=k*w1/N;
F1=f1*exp(-j*t'*W)*R;F2=f2*exp(-j*t'*W)*R; F1=real(F1);F2=real(F2);
subplot(211);plot(W,F1);
xlabel('w'),ylabel('F1(w)');title('F1(jw)的频谱');
74
subplot(212);plot(W,F2); xlabel('w'),ylabel('F2(w)');title('F2(jw)的频谱');
实验5
【例5-1】信号，绘制MATLAB曲面观察拉普f(t),,(t),,(t,2)
拉斯变换与傅里叶变换的关系，可以求得信号的拉普拉斯变换与傅里叶变
,s21,,j,e换:F(s),,F(j,),2Sa(,)es
MATLAB程序如下:
%绘制矩形信号拉氏变换曲面图
clf;
a=-0:0.1:5;
75
b=-20:0.1:20;
[a,b]=meshgrid(a,b); c=a+i*b;
c=(1-exp(-2*c))./c; c=abs(c);
mesh(a,b,c);
surf(a,b,c);
view(-60,20);
axis([-0,5,-20,20,0,2]); title('拉氏变换(s域像函数)');
colormap(hsv);
MATLAB程序如下:
%绘制矩形信号傅氏变换曲线
w=-20:0.1:20; %确定频率范围 Fw=(2*sin(w).*exp(i*w))./w; %计算傅里叶变换 plot(w,abs(Fw)); %绘制信号振幅频谱曲线 title('傅里叶变换(振幅频谱曲线)'); xlabel('频率');
76
2(应用拉普拉斯变换求解系统的响应
【例5-2】一个非时变电路的转移函数为:
4(s，6000)u100， H(s),,26，875s，88，ugs10
,12.5cos(8000t)V，求的稳态响应。 若uu0g
应用拉普拉斯变换法求解，MATLAB程序如下:
syms s t;
Hs=sym('(10^4*(s+6000))/(s^2+875*s+88*10^6)');
Vs=laplace(12.5*cos(8000*t));Vos=Hs*Vs; Vo=ilaplace(Vos);
Vo=vpa(Vo,4);
ezplot(Vo,[1,1+5e-3]);hold on; ezplot('12.5*cos(8000*t)',[1,1+5e-3]);
50,50]); axis([1,1+2e-3,-
text(0.25,0.85,'输出电压','sc');
text(0.07,0.35,'输入电压','sc');
title('稳态响应');xlabel('电压');ylabel('时间');
77
3(零极点分布与系统稳定性的关系
2，3，2ss【例4-3】已知连续系统的函数为，H(s),432，2，3，s，5sss编写程序求出系统的零极点，并画出系统的零极点图，判断其稳定性。
系统零极点的位置可以用MATLAB的多项式求根的函数roots，来求得，调用的格式为p=roots(A),其中A为待求根的多项式的行向量，返回的p为包含该多项式所有根的列向量。
系统零极点求取的MATLAB程序如下:
%求零极点
clear;
b=[1,3,2];a=[1,2,3,1,5];
zs=roots(b)
ps=roots(a)
plot(real(zs),imag(zs),'go',real(ps),imag(ps),'mx','markersize',12)
;grid on;
legend('零点','极点');
78
4.求系统的冲激响应
已知系统的零极点分布，可以利用MATLAB绘制其连续系统冲激响应的函数impluse()，来绘制系统的冲激响应。
16-4】连续系统的系统函数为【例2-，求系统的冲激响H(s),s，5
应。
MATLAB程序如下:
clear;
a=[1,2];b=[1];
impulse(b,a)
79
80
本文档为【测控专业计算机仿真实验指导书】，请使用软件OFFICE或WPS软件打开。作品中的文字与图均可以修改和编辑，
图片更改请在作品中右键图片并更换，文字修改请直接点击文字进行修改，也可以新增和删除文档中的内容。
该文档来自用户分享，如有侵权行为请发邮件ishare@vip.sina.com联系网站客服，我们会及时删除。
[版权声明] 本站所有资料为用户分享产生，若发现您的权利被侵害，请联系客服邮件isharekefu@iask.cn，我们尽快处理。
本作品所展示的图片、画像、字体、音乐的版权可能需版权方额外授权，请谨慎使用。
网站提供的党政主题相关内容(国旗、国徽、党徽..)目的在于配合国家政策宣传，仅限个人学习分享使用，禁止用于任何广告和商用目的。
浙公网安备 33021202002483