MATLAB Answers

How to plot the point where two complicated graphs cross?

1 次查看(过去 30 天)
Rebekah Ericson
Rebekah Ericson2021-4-20
回答: Star Strider ,2021-4-20
Hello! I want to plot blue points where the red circle graph meets the dark blue graph, but I'm having trouble. I'd like to get this figured out before my presentation tomorrow afternoon.
Here's what my code plots out:
Here's my current code:
clear; clc; close all;
% variables used
jd = juliandate(2020,10,13); %launch date; first opposition
jd_1 = juliandate(2022, 12,8); % next opposition
jd_2 = juliandate(2025,1,24); %end of time date
t = linspace (jd, jd_2, 4*366+10);
syms tt ttt
transfer = double(solve(149.6+160*(tt-jd)/687-227.9,tt)); %moment spaceship makes it to Mars
i=3;
transfer4 = transfer - jd
transfer2 = transfer + 883; %Spaceship starts return to earth
transfer3 = transfer2 + transfer4 %time spaceship completed return to earth
% creating arrays to plot the earth, moon, mars, and my spaceship
x_earth=149.6*cos(2*pi*(t-jd)/365.2);
y_earth=149.6*sin(2*pi*(t-jd)/365.2);
x_moon=20*cos(2*pi*(t-jd)/29.53) + x_earth ;
y_moon=20*sin(2*pi*(t-jd)/29.53) + y_earth;
x_mars=227.9*cos(2*pi*(t-jd)/687);
y_mars=227.9*sin(2*pi*(t-jd)/687);
x_earthtoMars=(149.6+160*(t-jd)/687).*cos(2*pi*(t-jd)/687);
y_earthtoMars=(149.6+160*(t-jd)/687).*sin(2*pi*(t-jd)/687);
x_MarsSatellite = 10*cos(2*pi*(t-transfer+2)/10.23) + x_mars;
y_MarsSatellite = 10*sin(2*pi*(t-transfer+2)/10.23) + y_mars;
x_MarstoEarth=(227.9-160*(t-transfer2)/687).*cos(2*pi*(t-jd)/687);
y_MarstoEarth = (227.9-160*(t-transfer2)/687).*sin(2*pi*(t-jd)/687);
%while ((double(x_MarstoEarth(ttt)) == double(x_earth(ttt)))) && ((double(y_MarstoEarth(ttt)) == double(y_earth(ttt))))
% transfer3 = double(solve(227.9-160*(ttt-transfer2)/687 - 149.6, ttt));
%end
msinitial_velocity = 9.9814681210623879384561120582467; %miles per second
initial_velocity = 1.3878963905307975340027267293406; %million km / day
audience_intial_velocity_ms = input ("Pick a number between 8 and 11.") %miles per second
audience_intial_velocity_kms = audience_intial_velocity_ms/0.6213712 %km per second
audience_intial_velocity_mkms = (audience_intial_velocity_kms/1000000) %million km per second
audience_initial_velocity_mkmh = audience_intial_velocity_mkms*3600 %million km per hour
audience_intial_velocity_mkmday = audience_initial_velocity_mkmh*24 %million km per day
x_earthtoMarsalt=(149.6+(230*audience_intial_velocity_mkmday)*(t-jd)/687).*cos(2*pi*(t-jd)/687);
y_earthtoMarsalt=(149.6+(230*audience_intial_velocity_mkmday)*(t-jd)/687).*sin(2*pi*(t-jd)/687);
x_MarstoEarthalt=(227.9-(230*audience_intial_velocity_mkmday)*(t-transfer2)/687).*cos(2*pi*(t-jd)/687);
y_MarstoEarthalt = (227.9-(230*audience_intial_velocity_mkmday)*(t-transfer2)/687).*sin(2*pi*(t-jd)/687);
%while ((double(x_MarstoEarth(ttt)) == double(x_earth(ttt)))) && ((double(y_MarstoEarth(ttt)) == double(y_earth(ttt))))
% transfer3 = double(solve(227.9-160*(ttt-transfer2)/687 - 149.6, ttt));
%end
%{
i3=0
for i=1:length(t)
if (t(i)> transfer3)
i3 = i
i = length(t)
end
end
%}
% plotting the sun
fhandle = figure
fhandle.Position(3) = 2 * fhandle.Position(3)
% Plot both versions simultaneously
subplot(1,2,1)
plot (0,0, '*', 'color', 'yellow', 'linewidth', 8)
hold on
xlabel('horizontal movement in space')
ylabel('vertical movement in space')
title('Planets around the Sun')
subplot(1,2,2)
plot (0,0, '*', 'color', 'yellow', 'linewidth', 8)
xlabel('horizontal movement in space')
ylabel('vertical movement in space')
title('Planets around the Sun')
hold on
while i<=4*368
subplot(1,2,1)
hold on
xlim([-300 300])
ylim([-300 300])
plot([x_moon(i) x_moon(i+2)], [y_moon(i) y_moon(i + 2)], 'k')
plot([x_moon(i-2) x_moon(i)], [y_moon(i-2) y_moon(i)], 'Color', 'white')
if t(i) <= transfer
plot([x_earth(i) x_earth(i+2)], [y_earth(i) y_earth(i + 2)], 'g')
plot([x_earthtoMars(i) x_earthtoMars(i+2)], [y_earthtoMars(i) y_earthtoMars(i + 2)], 'b')
plot([x_mars(i) x_mars(i+2)], [y_mars(i) y_mars(i + 2)], 'r')
else if ((transfer < t(i))) && ((t(i) <= transfer2))
plot([x_MarsSatellite(i) x_MarsSatellite(i+2)], [y_MarsSatellite(i) y_MarsSatellite(i + 2)], 'b')
plot([x_MarsSatellite(i-2) x_MarsSatellite(i)], [y_MarsSatellite(i-2) y_MarsSatellite(i)], 'Color', 'white')
plot([x_earth(i) x_earth(i+2)], [y_earth(i) y_earth(i + 2)], 'Color', 'cyan')
plot([x_mars(i) x_mars(i+2)], [y_mars(i) y_mars(i + 2)], 'Color', 'magenta')
else if ((transfer2 < t(i))) && ((t(i) <= transfer3))
plot([x_MarstoEarth(i) x_MarstoEarth(i+2)], [y_MarstoEarth(i) y_MarstoEarth(i + 2)], 'b')
plot([x_earth(i) x_earth(i+2)], [y_earth(i) y_earth(i + 2)], 'g')
plot([x_mars(i) x_mars(i+2)], [y_mars(i) y_mars(i + 2)], 'r')
end
end
end
subplot(1,2,2)
plot([x_moon(i) x_moon(i+2)], [y_moon(i) y_moon(i + 2)], 'k')
plot([x_moon(i-2) x_moon(i)], [y_moon(i-2) y_moon(i)], 'Color', 'white')
if t(i) <= transfer3
plot([x_earth(i) x_earth(i+2)], [y_earth(i) y_earth(i + 2)], 'g')
plot([x_earthtoMarsalt(i) x_earthtoMarsalt(i+2)], [y_earthtoMarsalt(i) y_earthtoMarsalt(i + 2)], 'b')
plot([x_mars(i) x_mars(i+2)], [y_mars(i) y_mars(i + 2)], 'r')
end
hold on
xlim([-300 300])
ylim([-300 300])
i = i + 2;
pause(0.0000000001)
end
subplot(1,2,1)
plot (x_earthtoMars(fix(transfer-jd)),y_earthtoMars(fix(transfer-jd)), '*', 'color', 'blue', 'linewidth', 5)
plot (x_MarstoEarth(fix(transfer2-jd)),y_MarstoEarth(fix(transfer2-jd)), '*', 'color', 'blue', 'linewidth', 5)
subplot(1,2,2)
plot (x_earthtoMarsalt(fix(transfer-jd)),y_earthtoMarsalt(fix(transfer-jd)), '*', 'color', 'blue', 'linewidth', 5)
% position is in millions of kilometers
% time is in days

回答(2 个)


Star Strider
Star Strider 2021-4-20
Clever animation!
Since I assume you wnat to know where the ‘x_eartoMars’ and ‘x_mars’ intersect, the easiest way is to simply find where those two curves cross:
int_idx = find(diff(sign(x_earthtoMars-x_mars)))
Plotting them as:
figure
plot(t,x_earthtoMars,'-b', t,x_mars,'-r')
hold on
plot(t(int_idx),x_earthtoMars(int_idx),'sb', t(int_idx),x_mars(int_idx),'dr')
hold off
grid
demonstrates that the indices are likely precise enough that it is not necessary to interpolate them to get more precise values. (In any event, when I attempted that it was possible to get accurate x-intercept values, however not y-intercept values, so using the indices is likely all that is necessary.)
To get other intersections, use this idea with other vector differences.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by