%% Apollotemplate.m  
%  18.2.2018 
%  Please, include a narrative with intermediate headlines, comments, etc.
%  and publish! 
%%
close all
%%
mu=1/82.45;
lam=1-mu;
w0=[1.2;0;0;-1.04935751];
% w0=[1.2;0;0;-0.5];  % Another set of initial values, try your own! 
Tperiod=6.19217;
%tmax=4*pi;
Tspan=[0 Tperiod];
[T,Z]=ode45(@Apollosys,Tspan,w0);
%%  Plots of orbit
hold on
plot(Z(:,1),Z(:,3));shg
plot(w0(1),w0(3),'*')    % Initial point.
plot(Z(1,1),Z(1,3),'*')  % Same as the previous line.-
plot(-mu,0,'*k',lam,0,'or') % plot moon and earth
%% Test: Some orbit points from the start and to the end.
plot(Z(10:10:30,1),Z(10:10:30,3),'r*')  % Next 10,20,30 time steps
plot(Z(end:-1:end-3,1),Z(end:-1:end-3,3),'bo')      % Last points
shg
%% Next animate
x=Z(:,1);y=Z(:,3);
for j=1:length(T)
    plot([x(j)],[y(j)],'or')
    pause(0.02)   % adjust
end
%% For fun: animate the rays from earth to orbit.
for j=1:length(T)
    plot([0 x(j)],[0 y(j)],'y')
    pause(0.02)   % adjust
end
%% Look at stepsizes
figure
tmax=max(T);
t=linspace(0,tmax,length(T));
plot(t(2:end),diff(T),'.');shg
%% How close does the Apollo come to the surface of the earth?
% Real sizes in kilometers: 
%
distMoonEarth = 384400  
R_earth= 6371            
R_moon = 1737

%%
%  Distance to centre: 
%
% $$ d(t)=\sqrt{y(t)^2 + (x(t)+\mu)^2}$$
%
%% Find the minimum of $d2(t) = d(t)^2$
% 
X=Z(:,1);
Y=Z(:,3);
D2=Y.^2+(X+mu).^2;
[minval2,minind]=min(D2)  % study help find
T(minind)
figure(1)
hold on
plot(X(minind),Y(minind),'*r')
D=sqrt(minval2)*distMoonEarth
%% 
% Remember to test with different tolerances, and other options, like 
% 'stats','on
%