如何在 Matlab 中循环表示和绘制数据?
How to respresent and plot data while on loop in Matlab?
我想做的是模拟两颗卫星的位置并确定它们是否会发生碰撞。为此,我使用了 Rk4 方法。如您所见,在循环中,我有 drawnow
命令在每次迭代后在 3d 图形上添加一个点。我感兴趣的是找到一种方法来表示 3d 图,直到满足碰撞条件(我有消息),并在 .txt
文件中获取值,而不是在循环之外.当我使用 return 命令并尝试绘制或添加文本时,出现错误 "vectors must be the same length" 并且不让我获取数据。发生这种情况是因为假设 sat1=1X200 而 tspan=1X201。我可以使用其他命令代替 return 或更好的解决方案吗?
这是我的代码:
clear; close all; clc;
figure(1);
orbitsat1 = animatedline('Color','r');
msat1 = 124;
msat2 = 234;
Asat1 = 2.2;
Asat2 = 3.2;
sat10(1) = 453322.3616;
sat10(2) = -2346021.219;
sat10(3) = -7894131.349;
sat10(4) = 2142.38067;
sat10(5) = 7487.44895;
sat10(6) = -9864.00872;
sat10 = sat10';
sat20(1) = 543322.3616;
sat20(2) = -3246021.219;
sat20(3) = -8794131.349;
sat20(4) = 1242.38067;
sat20(5) = 4787.44895;
sat20(6) = -8964.00872;
sat20 = sat20';
tspan = 0:200;
secpday = 60*60*24;
a1 = 2018;
la1 = 1;
z1 = 2;
o1 = 12;
m1 = 0;
s1 = 0;
n1 = datenum(a1,la1,z1,o1,m1,s1);
n1 = n1 + tspan/secpday;
[an,luna,zi,ora,min,sec] = datevec(n1);
h = 1;
sat1 = zeros(6, tspan(end)/h);
sat1(:,1) = sat10;
sat2 = zeros(6, tspan(end)/h);
sat2(:,1) = sat20;
for i = 1:tspan(end)/h
k_1 = simsat1(tspan(i), sat1(:,i), msat1, Asat1, sat1(4:6,i));
k_2 = simsat1(tspan(i)+0.5*h, sat1(:,i)+0.5*h*k_1,msat1, Asat1, sat1(4:6,i));
k_3 = simsat1(tspan(i)+0.5*h, sat1(:,i)+0.5*h*k_2, msat1, Asat1, sat1(4:6,i));
k_4 = simsat1(tspan(i)+h, sat1(:,i)+k_3*h, msat1, Asat1, sat1(4:6,i));
k_12 = simsat2(tspan(i), sat2(:,i), msat2, Asat2, sat2(4:6,i));
k_22 = simsat2(tspan(i)+0.5*h, sat2(:,i)+0.5*h*k_12, msat2, Asat2, sat2(4:6,i));
k_32 = simsat2(tspan(i)+0.5*h, sat2(:,i)+0.5*h*k_22, msat2, Asat2, sat2(4:6,i));
k_42 = simsat2(tspan(i)+h, sat2(:,i)+k_32*h, msat2, Asat2, sat2(4:6,i));
sat2(:,i+1) = sat2(:,i) + (1/6)*(k_12+2*k_22+2*k_32+k_42)*h;
sat1(:,i+1) = sat1(:,i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h;
addpoints(orbitsat1, sat1(1,i), sat1(2,i), sat1(3,i));
drawnow update;
hold off;
Rorb = sqrt(sum(sat1(1:3,i).^2));
Rcsc = sqrt(sum(sat2(1:3,i).^2));
dif = Rorb-Rcsc;
if (dif >0.5e10 && dif < 2e5) && i>350
Message = sprintf('Data:\nan:%d\nluna:%d\nzi:%d\nora:%d\nminut:%d\nsecunda:%d',an(i),luna(i),zi(i),ora(i),min(i),sec(i));
msgbox(Message)
msgbox('Coliziune!','Coliziune','warn',Message);
return
end
end
Pozx = sat1(1,:);
Pozy = sat1(2,:);
Pozz = sat1(3,:);
Tbody = [an luna zi ora min sec Pozx Pozy Pozz]';
twobody = fopen('nobody.txt','w');
fprintf(twobody,'Rezultate Twobody (metri) \n\n');
fprintf(twobody,' An luna zi ora min sec pozitia pe x pozitia pe y pozitia pe z viteza pe x viteza pe y viteza pe z\n\n');
fprintf(twobody,'---------------------------------------------------------------------------------------------------------------------------------------------\n%6d-%3d-%3d\t%3d:%3d:%3d\t\t\t\t%12.6f\t%12.6f\t\t%12.6f\t\t%12.6f\t%\n',Tbody);
fclose(twobody);
其中 simsat1
function sat1prim = simsat1(t,sat1,msat1,Asat1,vit)
miu = 398600.4418e9;
magn = sum(sat1(1:3).^2)^(3/2);
sat1prim = zeros(6,1);
sat1prim(1:3) = sat1(4:6);
sat1prim(4:6) = -miu.*sat1(1:3)./magn;
end
和 simsat2
function sat2prim = simsat2(t,sat2,msat2,Asat2,vit)
miu = 398600.4418e9;
magn = sum(sat2(1:3).^2)^(3/2);
sat2prim = zeros(6,1);
sat2prim(1:3) = sat2(4:6);
sat2prim(4:6) = -miu.*sat2(1:3)./magn;
end
好的;如果我没理解错的话这应该是你想要的。
我已经将 Tbody
变量移到“for”循环中,并且我正在将每个时间步写入文件中。另外,您搜索的关键字是 break
,而不是 return
。 break
只是让你脱离当前的 for 或 while 循环(只向上一级),而 return
returns 控制调用函数。
我稍微改变了初始化时间向量的方式,使其更清晰。 "sat1=1X200 while tspan=1X201" 的问题是因为你的 tspan = 0:200
是 201 点,而 tpsan(end)
值只有 200.这很好并且很有意义,因为您将初始位置包含在向量中。基本上确保他们所有人都有正确的点数。我通过定义时间步 dt=1
小时然后定义时间步数 nt=200
来完成此操作。这样我的时间向量就变成 tspan = 0:nt*dt
并且它将有 nt+1
个点。然后我可以像你一样迭代 for i = 1:(nt-1)
。
我添加了注释,我鼓励您在代码中也这样做。反正;这里是:
clear; close all; clc;
warning('on'); % Turn on warnings (Because they may be off)
figure(1);
orbitsat1 = animatedline('Color','r');
%% Initial conditions
msat1 = 124; msat2 = 234; Asat1 = 2.2; Asat2 = 3.2;
% Initialise satellites
sat10 = [453322.3616 -2346021.219 -7894131.349 2142.38067 7487.44895 -9864.00872];
sat20 = [543322.3616 -3246021.219 -8794131.349 1242.38067 4787.44895 -8964.00872];
dt = 1; % On hour timestep
nt = 200; % Number of timesteps
tspan = (0:nt)*dt; % Time from 0 to 200 timesteps in hours
% Make the time vectors
n1 = datenum(2018,1,2,12,0,0); % Starting point (serial date number in days)
n1 = n1 + tspan/24; % Timespan from hours -> days
[an,luna,zi,ora,min,sec] = datevec(n1);
%% Initialise containers
% Satelites
% nt+1 points because (:,1) is the initial condition
sat1 = zeros(6,nt+1); sat1(:,1) = sat10;
sat2 = zeros(6,nt+1); sat2(:,1) = sat20;
% Initialise Tbody here
Tbody = zeros(9,nt+1); Tbody(1:6,:) = [an; luna; zi; ora; min; sec];
%% Open the file
f_twobody = fopen('nobody.txt','w');
fprintf2(f_twobody,'Rezultate Twobody (metri)\n');
fprintf2(f_twobody,'Numar\tAn\tluna\tzi\tora\tmin\tsec\tpozitia_pe_x\tpozitia_pe_y\tpozitia_pe_z\n');
fprintf2(f_twobody,'%06d\t%6d\t%3d\t%3d\t%3d\t%3d\t%3d\t%12.6f\t%12.6f\t%12.6f\n',0,Tbody(:,1));
%% Go though all the timesteps
for i = 1:(nt-1)
% Calculate timestep
k_1 = simsat1(tspan(i), sat1(:,i), msat1, Asat1, sat1(4:6,i));
k_2 = simsat1(tspan(i)+0.5*dt, sat1(:,i)+0.5*dt*k_1,msat1, Asat1, sat1(4:6,i));
k_3 = simsat1(tspan(i)+0.5*dt, sat1(:,i)+0.5*dt*k_2, msat1, Asat1, sat1(4:6,i));
k_4 = simsat1(tspan(i)+dt, sat1(:,i)+k_3*dt, msat1, Asat1, sat1(4:6,i));
k_12 = simsat2(tspan(i), sat2(:,i), msat2, Asat2, sat2(4:6,i));
k_22 = simsat2(tspan(i)+0.5*dt, sat2(:,i)+0.5*dt*k_12, msat2, Asat2, sat2(4:6,i));
k_32 = simsat2(tspan(i)+0.5*dt, sat2(:,i)+0.5*dt*k_22, msat2, Asat2, sat2(4:6,i));
k_42 = simsat2(tspan(i)+dt, sat2(:,i)+k_32*dt, msat2, Asat2, sat2(4:6,i));
sat2(:,i+1) = sat2(:,i) + (1/6)*(k_12+2*k_22+2*k_32+k_42)*dt;
sat1(:,i+1) = sat1(:,i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*dt;
% Put the calculated positions into Tbody
Tbody(7:9,i+1) = [sat1(1,i+1) sat1(2,i+1) sat1(3,i+1)]';
% Write next line to the file
fprintf2(f_twobody,'%06d\t%6d\t%3d\t%3d\t%3d\t%3d\t%3d\t%12.6f\t%12.6f\t%12.6f\n',i,Tbody(:,i+1));
% Update plot
addpoints(orbitsat1, sat1(1,i), sat1(2,i), sat1(3,i));
drawnow update; hold off;
% Check for collision
Rorb = sqrt(sum(sat1(1:3,i+1).^2));
Rcsc = sqrt(sum(sat2(1:3,i+1).^2));
dif = abs(Rorb-Rcsc);
%if (dif >0.5e10 && dif < 2e5) && i>350
if sat1(1,i)>sat2(1,i) % Some kindof collision condition
warning('Coliziune! %d/%d/%d %d:%d:%d',Tbody(1:6,i+1));
fprintf(f_twobody,'Coliziune!\n');
break % Break out of the for loop (use this instead of 'return')
end
end
fclose(f_twobody); % Close the file
function sat1prim = simsat1(t,sat1,msat1,Asat1,vit)
miu = 398600.4418e9;
magn = sum(sat1(1:3).^2)^(3/2);
sat1prim = zeros(6,1);
sat1prim(1:3) = sat1(4:6);
sat1prim(4:6) = -miu.*sat1(1:3)./magn;
end
function sat2prim = simsat2(t,sat2,msat2,Asat2,vit)
miu = 398600.4418e9;
magn = sum(sat2(1:3).^2)^(3/2);
sat2prim = zeros(6,1);
sat2prim(1:3) = sat2(4:6);
sat2prim(4:6) = -miu.*sat2(1:3)./magn;
end
function fprintf2(f,varargin)
% Write both to stdout (console) and file
fprintf(varargin{:}); fprintf(f,varargin{:});
end
我想做的是模拟两颗卫星的位置并确定它们是否会发生碰撞。为此,我使用了 Rk4 方法。如您所见,在循环中,我有 drawnow
命令在每次迭代后在 3d 图形上添加一个点。我感兴趣的是找到一种方法来表示 3d 图,直到满足碰撞条件(我有消息),并在 .txt
文件中获取值,而不是在循环之外.当我使用 return 命令并尝试绘制或添加文本时,出现错误 "vectors must be the same length" 并且不让我获取数据。发生这种情况是因为假设 sat1=1X200 而 tspan=1X201。我可以使用其他命令代替 return 或更好的解决方案吗?
这是我的代码:
clear; close all; clc;
figure(1);
orbitsat1 = animatedline('Color','r');
msat1 = 124;
msat2 = 234;
Asat1 = 2.2;
Asat2 = 3.2;
sat10(1) = 453322.3616;
sat10(2) = -2346021.219;
sat10(3) = -7894131.349;
sat10(4) = 2142.38067;
sat10(5) = 7487.44895;
sat10(6) = -9864.00872;
sat10 = sat10';
sat20(1) = 543322.3616;
sat20(2) = -3246021.219;
sat20(3) = -8794131.349;
sat20(4) = 1242.38067;
sat20(5) = 4787.44895;
sat20(6) = -8964.00872;
sat20 = sat20';
tspan = 0:200;
secpday = 60*60*24;
a1 = 2018;
la1 = 1;
z1 = 2;
o1 = 12;
m1 = 0;
s1 = 0;
n1 = datenum(a1,la1,z1,o1,m1,s1);
n1 = n1 + tspan/secpday;
[an,luna,zi,ora,min,sec] = datevec(n1);
h = 1;
sat1 = zeros(6, tspan(end)/h);
sat1(:,1) = sat10;
sat2 = zeros(6, tspan(end)/h);
sat2(:,1) = sat20;
for i = 1:tspan(end)/h
k_1 = simsat1(tspan(i), sat1(:,i), msat1, Asat1, sat1(4:6,i));
k_2 = simsat1(tspan(i)+0.5*h, sat1(:,i)+0.5*h*k_1,msat1, Asat1, sat1(4:6,i));
k_3 = simsat1(tspan(i)+0.5*h, sat1(:,i)+0.5*h*k_2, msat1, Asat1, sat1(4:6,i));
k_4 = simsat1(tspan(i)+h, sat1(:,i)+k_3*h, msat1, Asat1, sat1(4:6,i));
k_12 = simsat2(tspan(i), sat2(:,i), msat2, Asat2, sat2(4:6,i));
k_22 = simsat2(tspan(i)+0.5*h, sat2(:,i)+0.5*h*k_12, msat2, Asat2, sat2(4:6,i));
k_32 = simsat2(tspan(i)+0.5*h, sat2(:,i)+0.5*h*k_22, msat2, Asat2, sat2(4:6,i));
k_42 = simsat2(tspan(i)+h, sat2(:,i)+k_32*h, msat2, Asat2, sat2(4:6,i));
sat2(:,i+1) = sat2(:,i) + (1/6)*(k_12+2*k_22+2*k_32+k_42)*h;
sat1(:,i+1) = sat1(:,i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h;
addpoints(orbitsat1, sat1(1,i), sat1(2,i), sat1(3,i));
drawnow update;
hold off;
Rorb = sqrt(sum(sat1(1:3,i).^2));
Rcsc = sqrt(sum(sat2(1:3,i).^2));
dif = Rorb-Rcsc;
if (dif >0.5e10 && dif < 2e5) && i>350
Message = sprintf('Data:\nan:%d\nluna:%d\nzi:%d\nora:%d\nminut:%d\nsecunda:%d',an(i),luna(i),zi(i),ora(i),min(i),sec(i));
msgbox(Message)
msgbox('Coliziune!','Coliziune','warn',Message);
return
end
end
Pozx = sat1(1,:);
Pozy = sat1(2,:);
Pozz = sat1(3,:);
Tbody = [an luna zi ora min sec Pozx Pozy Pozz]';
twobody = fopen('nobody.txt','w');
fprintf(twobody,'Rezultate Twobody (metri) \n\n');
fprintf(twobody,' An luna zi ora min sec pozitia pe x pozitia pe y pozitia pe z viteza pe x viteza pe y viteza pe z\n\n');
fprintf(twobody,'---------------------------------------------------------------------------------------------------------------------------------------------\n%6d-%3d-%3d\t%3d:%3d:%3d\t\t\t\t%12.6f\t%12.6f\t\t%12.6f\t\t%12.6f\t%\n',Tbody);
fclose(twobody);
其中 simsat1
function sat1prim = simsat1(t,sat1,msat1,Asat1,vit)
miu = 398600.4418e9;
magn = sum(sat1(1:3).^2)^(3/2);
sat1prim = zeros(6,1);
sat1prim(1:3) = sat1(4:6);
sat1prim(4:6) = -miu.*sat1(1:3)./magn;
end
和 simsat2
function sat2prim = simsat2(t,sat2,msat2,Asat2,vit)
miu = 398600.4418e9;
magn = sum(sat2(1:3).^2)^(3/2);
sat2prim = zeros(6,1);
sat2prim(1:3) = sat2(4:6);
sat2prim(4:6) = -miu.*sat2(1:3)./magn;
end
好的;如果我没理解错的话这应该是你想要的。
我已经将 Tbody
变量移到“for”循环中,并且我正在将每个时间步写入文件中。另外,您搜索的关键字是 break
,而不是 return
。 break
只是让你脱离当前的 for 或 while 循环(只向上一级),而 return
returns 控制调用函数。
我稍微改变了初始化时间向量的方式,使其更清晰。 "sat1=1X200 while tspan=1X201" 的问题是因为你的 tspan = 0:200
是 201 点,而 tpsan(end)
值只有 200.这很好并且很有意义,因为您将初始位置包含在向量中。基本上确保他们所有人都有正确的点数。我通过定义时间步 dt=1
小时然后定义时间步数 nt=200
来完成此操作。这样我的时间向量就变成 tspan = 0:nt*dt
并且它将有 nt+1
个点。然后我可以像你一样迭代 for i = 1:(nt-1)
。
我添加了注释,我鼓励您在代码中也这样做。反正;这里是:
clear; close all; clc;
warning('on'); % Turn on warnings (Because they may be off)
figure(1);
orbitsat1 = animatedline('Color','r');
%% Initial conditions
msat1 = 124; msat2 = 234; Asat1 = 2.2; Asat2 = 3.2;
% Initialise satellites
sat10 = [453322.3616 -2346021.219 -7894131.349 2142.38067 7487.44895 -9864.00872];
sat20 = [543322.3616 -3246021.219 -8794131.349 1242.38067 4787.44895 -8964.00872];
dt = 1; % On hour timestep
nt = 200; % Number of timesteps
tspan = (0:nt)*dt; % Time from 0 to 200 timesteps in hours
% Make the time vectors
n1 = datenum(2018,1,2,12,0,0); % Starting point (serial date number in days)
n1 = n1 + tspan/24; % Timespan from hours -> days
[an,luna,zi,ora,min,sec] = datevec(n1);
%% Initialise containers
% Satelites
% nt+1 points because (:,1) is the initial condition
sat1 = zeros(6,nt+1); sat1(:,1) = sat10;
sat2 = zeros(6,nt+1); sat2(:,1) = sat20;
% Initialise Tbody here
Tbody = zeros(9,nt+1); Tbody(1:6,:) = [an; luna; zi; ora; min; sec];
%% Open the file
f_twobody = fopen('nobody.txt','w');
fprintf2(f_twobody,'Rezultate Twobody (metri)\n');
fprintf2(f_twobody,'Numar\tAn\tluna\tzi\tora\tmin\tsec\tpozitia_pe_x\tpozitia_pe_y\tpozitia_pe_z\n');
fprintf2(f_twobody,'%06d\t%6d\t%3d\t%3d\t%3d\t%3d\t%3d\t%12.6f\t%12.6f\t%12.6f\n',0,Tbody(:,1));
%% Go though all the timesteps
for i = 1:(nt-1)
% Calculate timestep
k_1 = simsat1(tspan(i), sat1(:,i), msat1, Asat1, sat1(4:6,i));
k_2 = simsat1(tspan(i)+0.5*dt, sat1(:,i)+0.5*dt*k_1,msat1, Asat1, sat1(4:6,i));
k_3 = simsat1(tspan(i)+0.5*dt, sat1(:,i)+0.5*dt*k_2, msat1, Asat1, sat1(4:6,i));
k_4 = simsat1(tspan(i)+dt, sat1(:,i)+k_3*dt, msat1, Asat1, sat1(4:6,i));
k_12 = simsat2(tspan(i), sat2(:,i), msat2, Asat2, sat2(4:6,i));
k_22 = simsat2(tspan(i)+0.5*dt, sat2(:,i)+0.5*dt*k_12, msat2, Asat2, sat2(4:6,i));
k_32 = simsat2(tspan(i)+0.5*dt, sat2(:,i)+0.5*dt*k_22, msat2, Asat2, sat2(4:6,i));
k_42 = simsat2(tspan(i)+dt, sat2(:,i)+k_32*dt, msat2, Asat2, sat2(4:6,i));
sat2(:,i+1) = sat2(:,i) + (1/6)*(k_12+2*k_22+2*k_32+k_42)*dt;
sat1(:,i+1) = sat1(:,i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*dt;
% Put the calculated positions into Tbody
Tbody(7:9,i+1) = [sat1(1,i+1) sat1(2,i+1) sat1(3,i+1)]';
% Write next line to the file
fprintf2(f_twobody,'%06d\t%6d\t%3d\t%3d\t%3d\t%3d\t%3d\t%12.6f\t%12.6f\t%12.6f\n',i,Tbody(:,i+1));
% Update plot
addpoints(orbitsat1, sat1(1,i), sat1(2,i), sat1(3,i));
drawnow update; hold off;
% Check for collision
Rorb = sqrt(sum(sat1(1:3,i+1).^2));
Rcsc = sqrt(sum(sat2(1:3,i+1).^2));
dif = abs(Rorb-Rcsc);
%if (dif >0.5e10 && dif < 2e5) && i>350
if sat1(1,i)>sat2(1,i) % Some kindof collision condition
warning('Coliziune! %d/%d/%d %d:%d:%d',Tbody(1:6,i+1));
fprintf(f_twobody,'Coliziune!\n');
break % Break out of the for loop (use this instead of 'return')
end
end
fclose(f_twobody); % Close the file
function sat1prim = simsat1(t,sat1,msat1,Asat1,vit)
miu = 398600.4418e9;
magn = sum(sat1(1:3).^2)^(3/2);
sat1prim = zeros(6,1);
sat1prim(1:3) = sat1(4:6);
sat1prim(4:6) = -miu.*sat1(1:3)./magn;
end
function sat2prim = simsat2(t,sat2,msat2,Asat2,vit)
miu = 398600.4418e9;
magn = sum(sat2(1:3).^2)^(3/2);
sat2prim = zeros(6,1);
sat2prim(1:3) = sat2(4:6);
sat2prim(4:6) = -miu.*sat2(1:3)./magn;
end
function fprintf2(f,varargin)
% Write both to stdout (console) and file
fprintf(varargin{:}); fprintf(f,varargin{:});
end