通信仿真.docx
- 文档编号:748799
- 上传时间:2023-04-30
- 格式:DOCX
- 页数:20
- 大小:408.75KB
通信仿真.docx
《通信仿真.docx》由会员分享,可在线阅读,更多相关《通信仿真.docx(20页珍藏版)》请在冰点文库上搜索。
通信仿真
%ch1example1prg1.m
g=9.8;%重力加速度
v=0;%设定初始速度条件
s=0;%设定初始位移条件
t=0;%设定起始时间
dt=0.1;%设置计算步长
N=60;%设置仿真递推次数.仿真时间等于N与dt的乘积
fork=1:
N
v=v+g*dt;%计算新时刻的速度
s(k+1)=s(k)+v*dt;%新位移
t(k+1)=t(k)+dt;%时间更新
end
%理论计算,以便与仿真结果对照
t_theory=0:
0.01:
N*dt;%设置解析计算的时间点
v_theory=g*t_theory;%解析计算的瞬时速度
s_theory=1/2*g*t_theory.^2;%解析计算的瞬时位移
%作图:
仿真结果与解析结果对比
t=0:
dt:
N*dt;
plot(t,s,'o',t_theory,s_theory,'-');
xlabel('时间t');ylabel('位移s');
legend('仿真结果','理论结果');
N=20
%ch1example1prg2.m
g=9.8;%重力加速度
forL=1:
5%仿真重复5次以便于观察
v=0;%初始速度
s=0;%初始位置
t=0;
dt=0.01;%计算步长
fork=1:
200
v=v+g*dt;%速度
s=s+v*dt;%位移
t=t+dt;%时间
plot(0,-s,'o');
axis([-22-200]);%坐标范围固定
text(0.5,-1,['当前时间:
t=',num2str(t)]);
text(0.5,-2,['当前速度:
v=',num2str(v)]);
text(0.5,-3,['当前位置:
s=',num2str(s)]);
set(gcf,'DoubleBuffer','on');%双缓冲避免作图闪烁
drawnow;%立即作图
end
end
k=1:
100
%ch1example2prg1.m
g=9.8;%重力加速度
v0=0;%初始速度
y0=1;%初始位置
m=1;%小球质量
t0=0;%起始时间
K=0.85;%弹跳的损耗系数
N=5000;%仿真的总步进数
dt=0.001;%仿真步长
v=v0;%初状态
y=y0;
fork=1:
N
if(y>0)|(v>0)%小球在空中的(含刚刚弹起瞬间)动力方程计算
v=v-g*dt;
y=y+v*dt;
else%碰击瞬间的计算
y=y-K.*v*dt;
v=-K.*v-g*dt;
end
s(k)=y;%将当前位移记录到s数组中以便作图
end
t=t0:
dt:
dt*(N-1);%仿真时间长度
plot(t,s);
xlabel('时间t');
ylabel('位移y(t)');
axis([0501.1]);
holdon
g=9.8;%重力加速度
v0=0;%初始速度
y0=1;%初始位置
m=1;%小球质量
t0=0;%起始时间
K=1;%弹跳的损耗系数
N=5000;%仿真的总步进数
dt=0.001;%仿真步长
v=v0;%初状态
y=y0;
fork=1:
N
if(y>0)|(v>0)%小球在空中的(含刚刚弹起瞬间)动力方程计算
v=v-g*dt;
y=y+v*dt;
else%碰击瞬间的计算
y=y-K.*v*dt;
v=-K.*v-g*dt;
end
s(k)=y;%将当前位移记录到s数组中以便作图
end
t=t0:
dt:
dt*(N-1);%仿真时间长度
plot(t,s);
xlabel('时间t');
ylabel('位移y(t)');
axis([0501.1]);
%ch1example3prg1.m
sita=0:
0.01:
2*pi;
x=sin(sita);
y=cos(sita);%计算半径为1的圆周上的点,以便作出圆周观察
m=0;%在圆内在落点计数器
x1=2*rand(1000,1)-1;%产生均匀分布于[-1,+1]直接的两个独立随机数x1,y1
y1=2*rand(1000,1)-1;
N=1000;%设置试验次数
forn=1:
N%循环进行重复试验并统计
p1=x1(1:
n);
q1=y1(1:
n);
if(x1(n)*x1(n)+y1(n)*y1(n))<1%计算落点到坐标原点的距离,判别落点是否在圆内
m=m+1;%如果落入圆中,计数器加1
end
plot(p1,q1,'.',x,y,'-k',[-1-111-1],[-111-1-1],'-k');
axisequal;%坐标纵横比例相同
axis([-22-22]);%固定坐标范围
text(-1,-1.2,['试验总次数n=',num2str(n)]);%显示试验结果
text(-1,-1.4,['落入圆中数m=',num2str(m)]);
text(-1,-1.6,['近似圆面积S_c=',num2str(m/n*4)]);
set(gcf,'DoubleBuffer','on');%双缓冲避免作图闪烁
drawnow;%显示结果
end
%ch1example3prg2.m
tic%计时器启动
n=10000;%每次随机落点10000个
fork=1:
1000%重复试验1000次
x1=2*rand(n,1)-1;
y1=2*rand(n,1)-1;%随机落点产生
m(k)=sum((x1.*x1+y1.*y1)<1);%求落入圆中的点数和
end
S_c=mean(m).*4./n%计算并显示结果
time=toc%显示耗时
S_c=3.1417
time=0.5460
%ch1example4prg1.m
g=9.8;%重力加速度
v0=0;%初始速度
y0=1;%初始位置
m=0.1;%小球质量
t0=0;%起始时间
K=0.8;%弹跳的损耗系数
N=500;%仿真的总步进数
dt=0.005;%仿真步长
v=v0;%初状态
y=y0;
vx=0;
vz=0;
sx=0;
sz=0;
fork=1:
N
ify>0%小球在空中的动力方程计算
v=v-g*dt;
y=y+v*dt;
else%碰击瞬间的计算
y=-K.*v*dt;
v=-K.*v-g*dt;
end
Fx=randn;%x水平方向的随机力,方差为1
ax=Fx./m;%Fx导致的x水平方向的加速度
vx=vx+ax*dt;%小球在x水平方向的瞬时速度
sx=sx+vx*dt;%小球在x水平方向的位移
Fz=randn;%z水平方向的随机力,方差为1
az=Fz./m;%Fz导致的z水平方向的加速度
vz=vz+az*dt;%小球在z水平方向的瞬时速度
sz=sz+vz*dt;%小球在z水平方向的位移
plot3(sx,sz,y,'r.');gridon;holdon;
axis([-22-2201]);%坐标范围固定
set(gcf,'DoubleBuffer','on');%双缓冲避免作图闪烁
xlabel('水平方向x');ylabel('水平方向z');zlabel('垂直方向y');
drawnow;%作图
end
%ch1problem2.m
g=9.8;%重力加速度
k=-1;%空气阻力系数
dt=0.1;%设置计算步长
N=30;%设置仿真递推次数.仿真时间等于N与dt的乘积
form=[1,2,10]%三种落体质量
v=0;%设定初始速度条件
s=0;%设定初始位移条件
t=0;%设定起始时间
fori=1:
N
a=g+k/m*v;%计算加速度
v=v+a*dt;%计算新时刻的速度
s(i+1)=s(i)+v*dt;%新位移
t(i+1)=t(i)+dt;%时间更新
end
plot(t,s,'ro');
holdon;
end
%理论计算,以便与仿真结果对照
t_theory=0:
0.01:
N*dt;%设置解析计算的时间点
v_theory=g*t_theory;%解析计算的瞬时速度
s_theory=1/2*g*t_theory.^2;%解析计算的瞬时位移
%作图:
仿真结果与解析结果对比
t=0:
dt:
N*dt;
plot(t_theory,s_theory,'b-');
xlabel('时间t');ylabel('位移s');
%ch1problem5.m
M=10000;%每次随机落点10000个
fork=1:
1000%重复试验1000次
x1=sqrt
(2)*(2*rand(M,1)-1);
y1=2*rand(M,1);%随机落点产生
m(k)=sum(2-x1.^2>y1);%求落入抛物线所围区域的点数和
end
S_c=mean(m)/M*2*2*sqrt
(2)%计算并显示结果
S_c=3.7713
%ch2example1prg1.m
dt=1e-5;%仿真采样间隔
T=3*1e-3;%仿真终止时间
t=0:
dt:
T;
input=2*cos(2*pi*1000*t);%输入被调信号
carrier=5*cos(2*pi*1e4*t);%载波
output=(2+0.5*input).*carrier;%调制输出
%作图:
观察输入信号,载波,以及调制输出
subplot(3,1,1);plot(t,input);xlabel('时间t');ylabel('被调信号');
subplot(3,1,2);plot(t,carrier);xlabel('时间t');ylabel('载波');
subplot(3,1,3);plot(t,output);xlabel('时间t');ylabel('调幅输出');
%ch2example1prg2.m
clear;%清空内存变量,以避免以往运行的结果影响本程序
dt=1e-5;%仿真采样间隔
T=3*1e-3;%仿真终止时间
t=0:
dt:
T;
fork=1:
length(t)%基于时间流的仿真计算
input(k)=2*cos(2*pi*1000*t(k));%第k个仿真步进时的输入被调信号
carrier(k)=5*cos(2*pi*1e4*t(k));%第k个仿真步进时的载波
output(k)=(2+0.5*input(k)).*carrier(k);%第k个仿真步进时的调制输出
end
%作图:
观察输入信号,载波,以及调制输出
subplot(3,1,1);plot(t,input);xlabel('时间t');ylabel('被调信号');
subplot(3,1,2);plot(t,carrier);xlabel('时间t');ylabel('载波');
subplot(3,1,3);plot(t,output);xlabel('时间t');ylabel('调幅输出');
%ch2example1prg3.m
dt=1e-6;%仿真采样间隔
T=2*1e-3;%仿真的帧周期
forN=0:
500%总共仿真的帧数
t=N*T+(0:
dt:
T);%帧中的取样时刻
input=2*cos(2*pi*1005*t);%输入被调信号
carrier=5*cos(2*pi*(1e4)*t+0.1*randn);%载波
output=(2+0.5*input).*carrier;%调制输出
noise=randn(size(t));%噪声
r=output+noise;%调制信号通过加性噪声信道
%作图:
观察输入信号,载波,以及调制输出
subplot(3,1,1);plot([0:
dt:
T],input);xlabel('时间t');
ylabel('被调信号');text(T*2/3,1.5,['当前帧数:
N=',num2str(N)]);
subplot(3,1,2);plot([0:
dt:
T],carrier);
xlabel('时间t');ylabel('载波');
subplot(3,1,3);plot([0:
dt:
T],r);
xlabel('时间t');ylabel('调幅输出');
set(gcf,'DoubleBuffer','on');%双缓冲避免作图闪烁
drawnow;
end
%ch2example2prg1.m
dt=1e-5;%仿真采样间隔
R=1e3;%电阻值
C=1e-6;%电容量
T=5*1e-3;%仿真区间从-T到+T
t=-T:
dt:
T;%计算的离散时刻序列
y
(1)=0;%电容电压初始值,在时间小于零区间将保持不变
%如果要仿真零输入响应,可设置y
(1)=1等非零值.
%----输入信号设定:
可选择:
零输入,阶跃输入,正弦输入,方波输入等----
x=zeros(size(t));%初始化输入信号存储矩阵
x=1*(t>=0);%在0时刻的输入信号跃变为1,即输入为阶跃信号.
%如果要仿真零输入响应,这里可设x=0即可
%x=sin(2*pi*1000*t).*(t>=0);%这是从0时刻开始的1000Hz的正弦信号
%x=square(2*pi*500*t).*(t>=0);%这是从0时刻开始的500Hz的方波信号
%仿真开始,注意:
设零时刻之前电路不工作,系统状态保持不变
fork=1:
length(t)
time=-T+k*dt;
iftime>=0
y(k+1)=y(k)+1./(R*C)*(x(k)-y(k))*dt;%递推求解下一个仿真时刻的状态值
else
y(k+1)=y(k);%在时间小于零时设电路断开,系统不工作
end
end
subplot(2,1,1);plot(t,x(1:
length(t)));axis([-TT-1.11.1]);
xlabel('t');ylabel('input');
subplot(2,1,2);plot(t,y(1:
length(t)));axis([-TT-1.11.1]);
xlabel('t');ylabel('output');
%ch2example3prg1.m
dt=0.0001;%仿真步进
T=15;%仿真时间长度
t=0:
dt:
T;%仿真计算时间序列
g=9.8;%重力加速度
L=1;%摆线长度
m=10;%摆锤质量
k=5;%空气阻力比例系数
theta0=3.1;%初始摆角设置
v0=0;%初始摆速设置
v=zeros(size(t));%程序存储变量预先初始化,可提高执行速度
theta=zeros(size(t));
v
(1)=v0;%初始值赋值
theta
(1)=theta0;
forn=1:
length(t)%仿真求解开始
v(n+1)=v(n)+(g*sin(theta(n))-k./m.*v(n)).*dt;
theta(n+1)=theta(n)-1./L.*v(n).*dt;
end
%使用双坐标系统来作图,注意作图和图形标注的技巧
[AX,H1,H2]=plotyy(t,v(1:
length(t)),t,theta(1:
length(t)),'plot');
set(H1,'LineStyle','-');%设置作图线型
set(H2,'LineStyle','-.');
set(get(AX
(1),'Ylabel'),'String','线速度v(t)m/s');%坐标标注
set(get(AX
(2),'Ylabel'),'String','角位移\theta(t)rad');
xlabel('时间t');
legend(H1,'线速度v(t)',2);
legend(H2,'角位移\theta(t)',1);
%ch2example5prg1.m
dt=0.01;%仿真步进
T=15;%仿真时间长度
t=0:
dt:
T;%仿真计算时间序列
g=9.8;%重力加速度
L=1;%摆线长度
m=10;%摆锤质量
k=5;%空气阻力比例系数
theta0=3.1;%初始摆角设置
v0=0;%初始摆速设置
x0=[v0;theta0];%初始状态赋值
par=[g;k;m;L];%系统参数赋值
%以欧拉算法计算并作图
[t_out,x_out]=eulerode('pendulumstateeq',t',x0,[],par);
plot(t_out,x_out(:
1),'b-.');holdon;
%以ode45算法计算并作图对比
[t_out,x_out]=ode45('pendulumstateeq',t',x0,[],par);
plot(t_out,x_out(:
1),'r-');
xlabel('时间t');ylabel('线速度m/s');
legend('欧拉算法','ode45算法');
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 通信 仿真