Investigate orbit coordinate systems.

Since version 1.
------------------------------------------------------------------------
See also FOrb, Constant, NewFig, TimeGUI, TitleS, XLabelS, YLabelS, RK4,
Period
------------------------------------------------------------------------

Contents

%-------------------------------------------------------------------------------
%   Copyright 1994 Princeton Satellite Systems, Inc. All rights reserved.
%-------------------------------------------------------------------------------

Global for the time GUI

%------------------------
global simulationAction
simulationAction = ' ';

Constants

%----------
mu     = Constant('mu');
rE     = Constant('Earth Radius Mean');

Simulate a circular orbit with an altitude of 200 km

%-----------------------------------------------------
rCirc  = rE+200;
r      = [rCirc;0;0];
v      = [0;sqrt(mu/rCirc);0];
x      = [r;v];
y      = [[rCirc;0;0];[0;sqrt(mu/rCirc^3);0]];

nSim   = 500;
xPlot  = zeros(2,nSim);
yPlot  = zeros(2,nSim);
aExt   = [0;0;0];
period = 2*pi*sqrt(rCirc^3/mu);
dTSim  = period/nSim;
t      = 0;

Initialize the time display

%----------------------------
tToGoMem.lastJD        = 0;
tToGoMem.lastStepsDone = 0;
tToGoMem.kAve          = 0;
ratioRealTime          = 0;
[ ratioRealTime, tToGoMem ] =  TimeGUI( nSim, 0, tToGoMem, 0, dTSim, 'TOrbit Sim' );

for k = 1:nSim;

  % Display the status message
  %---------------------------
  [ ratioRealTime, tToGoMem ] = TimeGUI( nSim, k, tToGoMem, ratioRealTime, dTSim );

  xPlot(:,k) = x(1:2);
  yPlot(:,k) = y(1)*[cos(y(2));sin(y(2))];
  x          = RK4('FOrb',x,dTSim,t,'car',aExt,mu);
  y          = RK4('FOrb',y,dTSim,t,'cyl',aExt,mu);
  t          = t + dTSim;

  % Time control
  %-------------
  switch simulationAction
    case 'pause'
      pause
      simulationAction = ' ';
    case 'stop'
      return;
    case 'plot'
      break;
  end

end

theta = linspace(0,2*pi*(1-1/k),k);
j     = 1:k;

NewFig('Orbit Integration Error');
subplot(211)
plot(xPlot(1,j)-rCirc*cos(theta),xPlot(2,j)-rCirc*sin(theta));
XLabelS('x Error (km)');
YLabelS('y Error (km)');
TitleS('Cartesian Orbit Integration Errors');
grid
subplot(212)
plot(yPlot(1,j)-rCirc*cos(theta),yPlot(2,j)-rCirc*sin(theta));
XLabelS('x Error (km)');
YLabelS('y Error (km)');
TitleS('Cylindrical Orbit Integration Errors');
grid


%--------------------------------------
% PSS internal file version information
%--------------------------------------