% RT2

constants
mu=G*M_E;
%---
h = 300e3;

a1 = R_E+h;  % parking orbit
a2 = k3_SMA(12*3600); % final orbit for T=12h

m0 = 3000;
ve = 2.5e3;
f = 50; % mass flow.
%--- a)

a = (a1+a2)/2

v1 = sqrt(mu/(a1)) % parking circ. orbit
vGTO1 = sqrt(mu*(2/a1 - 1/a)) % GTO at periapsis.
vGTO2 = sqrt(mu*(2/a2 - 1/a)) % GTO at apoapsis
v2 = sqrt(mu/(a2)) % Geo circ. orbit.

% Burn 1
deltaV1 = vGTO1 - v1
m1 = m0/exp(deltaV1/ve) 
deltam1 = m0-m1
tb1 = deltam1/f

% Burn 2
deltaV2 = v2 - vGTO2
m2 = m1/exp(deltaV2/ve)
deltam2 = m1 - m2
tb2 = deltam2/f
%---b)

mfinal = m0-deltam1-deltam2

%--- c)
deltaVia = sqrt(vGTO2^2 + visviva_circ(a2)^2 - 2*vGTO2*visviva_circ(a2)*cos(d2r(28.5)))
deltam2i = m0*(1-exp(-deltaVia/ve)) - deltam1*(1-exp(-deltaVia/ve)); % mass change in Burn 2 (with inclination change)
mfinali = m0-deltam1-deltam2i

%m2i = m1/exp(deltaVia/ve)
%deltam2i = m1 - m2i
tb2i = deltam2i/f

%--- d)
GTO_time = k3_period(a)/2/3600 % 

