KeplerOrbit Class Reference
[GeoDynamics]

#include <KeplerOrbit.hpp>

List of all members.


Detailed Description

This class do some useful Keplerian orbit computation.

Reference: Satellite orbits models methods applications Montenbruck, E. Gill

Definition at line 49 of file KeplerOrbit.hpp.

Static Public Member Functions

Vector< double > State (double GM, const Vector< double > &Kep, double dt)
 Computes the satellite state vector from osculating Keplerian elements for elliptic orbits.
gpstk::Matrix< double > StatePartials (double GM, const Vector< double > &Kep, double dt)
 Computes the partial derivatives of the satellite state vector with respect to the orbital elements for elliptic, Keplerian orbits.
Vector< double > Elements (double GM, const Vector< double > &y)
 Computes the osculating Keplerian elements from the satellite state vector for elliptic orbits.
Vector< double > Elements (double GM, double Mjda, double Mjdb, const Vector< double > &ra, const Vector< double > &rb)
 Computes orbital elements from two given position vectors and associated times.
void TwoBody (double GM, const Vector< double > &Y0, double dt, Vector< double > &Y, Matrix< double > &dYdY0)
 Propagates a given state vector and computes the state transition matrix for elliptical Keplerian orbits.
double EccentricAnomaly (double M, double e)
 Computes the eccentric anomaly for elliptic orbits.
double TrueAnomaly (double M, double e)
 Computes the true anomaly for elliptic orbits.
double MeanAnomaly (double cta, double e)
 Computes the true anomaly for elliptic orbits.
double getPeriod (double GM, const Vector< double > &Kep)
 Get the period of the orbit.
double getApogee (double GM, const Vector< double > &Kep)
 Get the distance to the apogee point.
double getPerigee (double GM, const Vector< double > &Kep)
 Get the distance to the perigee point.
void test ()

Protected Member Functions

 KeplerOrbit ()
 Default constructor.
 ~KeplerOrbit ()
 Default destructor.

Static Protected Member Functions

double FindEta (const Vector< double > &r_a, const Vector< double > &r_b, double tau)
 Computes the sector-triangle ratio from two position vectors and the intermediate time.
double Frac (double x)
 Fractional part of a number (y=x-[x]).
double Modulo (double x, double y)
 x mod y
double F (double eta, double m, double l)
 local function for use by FindEta()
Vector< double > Stack (Vector< double > r, Vector< double > v)
 connect two vector


Constructor & Destructor Documentation

KeplerOrbit  )  [inline, protected]
 

Default constructor.

Definition at line 163 of file KeplerOrbit.hpp.

~KeplerOrbit  )  [inline, protected]
 

Default destructor.

Definition at line 166 of file KeplerOrbit.hpp.


Member Function Documentation

double EccentricAnomaly double  M,
double  e
[static]
 

Computes the eccentric anomaly for elliptic orbits.

Parameters:
M Mean anomaly in [rad]
e Eccentricity of the orbit [0,1]
Returns:
Eccentric anomaly in [rad]

Definition at line 48 of file KeplerOrbit.cpp.

References gpstk::cos(), KeplerOrbit::Modulo(), and gpstk::sin().

Referenced by KeplerOrbit::State(), KeplerOrbit::StatePartials(), and KeplerOrbit::TrueAnomaly().

Vector< double > Elements double  GM,
double  Mjda,
double  Mjdb,
const Vector< double > &  ra,
const Vector< double > &  rb
[static]
 

Computes orbital elements from two given position vectors and associated times.

Parameters:
GM Gravitational coefficient
Mjda Time ta (Modified Julian Date)
Mjdb Time tb (Modified Julian Date)
ra Position vector at time t_a
rb Position vector at time t_b
Returns:
> Keplerian elements (a,e,i,Omega,omega,M) at ta
Warning:
The function cannot be used with state vectors describing a circular or non-inclined orbit.

Definition at line 447 of file KeplerOrbit.cpp.

References gpstk::cross(), gpstk::dot(), KeplerOrbit::FindEta(), log, KeplerOrbit::Modulo(), gpstk::norm(), gpstk::pow(), and gpstk::sqrt().

Vector< double > Elements double  GM,
const Vector< double > &  y
[static]
 

Computes the osculating Keplerian elements from the satellite state vector for elliptic orbits.

@ GM Gravitational coefficient @ y State vector(position and velocity) @ return Keplerian elements(a e i OGM omg M) @ warning The state vector and GM must be given in consistent units

Definition at line 383 of file KeplerOrbit.cpp.

References gpstk::cross(), gpstk::dot(), KeplerOrbit::Modulo(), gpstk::norm(), and gpstk::sqrt().

Referenced by KeplerOrbit::test(), and KeplerOrbit::TwoBody().

double F double  eta,
double  m,
double  l
[static, protected]
 

local function for use by FindEta()

Definition at line 612 of file KeplerOrbit.cpp.

References gpstk::asin(), log, gpstk::pow(), gpstk::sin(), and gpstk::sinh().

Referenced by KeplerOrbit::FindEta().

double FindEta const Vector< double > &  r_a,
const Vector< double > &  r_b,
double  tau
[static, protected]
 

Computes the sector-triangle ratio from two position vectors and the intermediate time.

Parameters:
r_a Position at time t_a
r_a Position at time t_b
tau Normalized time (sqrt(GM)*(t_a-t_b))
Returns:
Sector-triangle ratio

Definition at line 118 of file KeplerOrbit.cpp.

References gpstk::dot(), KeplerOrbit::F(), F1, F2, gpstk::norm(), gpstk::pow(), and gpstk::sqrt().

Referenced by KeplerOrbit::Elements().

double Frac double  x  )  [inline, static, protected]
 

Fractional part of a number (y=x-[x]).

Definition at line 182 of file KeplerOrbit.hpp.

double getApogee double  GM,
const Vector< double > &  Kep
[static]
 

Get the distance to the apogee point.

Definition at line 654 of file KeplerOrbit.cpp.

double getPerigee double  GM,
const Vector< double > &  Kep
[static]
 

Get the distance to the perigee point.

Definition at line 661 of file KeplerOrbit.cpp.

double getPeriod double  GM,
const Vector< double > &  Kep
[static]
 

Get the period of the orbit.

Definition at line 648 of file KeplerOrbit.cpp.

References gpstk::sqrt().

double MeanAnomaly double  cta,
double  e
[static]
 

Computes the true anomaly for elliptic orbits.

Parameters:
cta True anomaly in [rad]
e Eccentricity of the orbit [0,1]
Returns:
Mean anomaly in [rad]

Definition at line 102 of file KeplerOrbit.cpp.

References gpstk::cos(), gpstk::sin(), and gpstk::sqrt().

double Modulo double  x,
double  y
[inline, static, protected]
 

x mod y

Definition at line 186 of file KeplerOrbit.hpp.

Referenced by KeplerOrbit::EccentricAnomaly(), and KeplerOrbit::Elements().

Vector< double > Stack Vector< double >  r,
Vector< double >  v
[static, protected]
 

connect two vector

Definition at line 235 of file KeplerOrbit.cpp.

References Vector::size().

Referenced by KeplerOrbit::State(), and KeplerOrbit::StatePartials().

Vector< double > State double  GM,
const Vector< double > &  Kep,
double  dt
[static]
 

Computes the satellite state vector from osculating Keplerian elements for elliptic orbits.

Parameters:
GM Gravitational coefficient
Kep Keplerian elements(a e i OMG omg M)
dt Time since epoch
Returns:
State vector(position and velocity)
Warning:
The semimajor axis a=Kep(0), dt and GM must be given in consistent units

Definition at line 184 of file KeplerOrbit.cpp.

References gpstk::cos(), KeplerOrbit::EccentricAnomaly(), gpstk::Rx(), gpstk::Rz(), gpstk::sin(), gpstk::sqrt(), and KeplerOrbit::Stack().

Referenced by KeplerOrbit::test(), and KeplerOrbit::TwoBody().

Matrix< double > StatePartials double  GM,
const Vector< double > &  Kep,
double  dt
[static]
 

Computes the partial derivatives of the satellite state vector with respect to the orbital elements for elliptic, Keplerian orbits.

Parameters:
GM Gravitational coefficient
Kep Keplerian elements (a,e,i,Omega,omega,M)
dt Time since epoch
Returns:
Partials derivatives of the state vector (x,y,z,vx,vy,vz) at time dt with respect to the epoch orbital elements warning The semimajor axis a=Kep(0), dt and GM must be given in consistent units, The function cannot be used with circular or non-inclined orbit.

Definition at line 267 of file KeplerOrbit.cpp.

References gpstk::cos(), gpstk::cross(), KeplerOrbit::EccentricAnomaly(), gpstk::norm(), gpstk::Rx(), gpstk::Rz(), gpstk::sin(), gpstk::sqrt(), and KeplerOrbit::Stack().

Referenced by KeplerOrbit::TwoBody().

void test  )  [static]
 

Definition at line 671 of file KeplerOrbit.cpp.

References KeplerOrbit::Elements(), and KeplerOrbit::State().

double TrueAnomaly double  M,
double  e
[static]
 

Computes the true anomaly for elliptic orbits.

Parameters:
M Mean anomaly in [rad]
e Eccentricity of the orbit [0,1]
Returns:
True anomaly in [rad]

Definition at line 86 of file KeplerOrbit.cpp.

References gpstk::cos(), KeplerOrbit::EccentricAnomaly(), gpstk::sin(), and gpstk::sqrt().

void TwoBody double  GM,
const Vector< double > &  Y0,
double  dt,
Vector< double > &  Y,
Matrix< double > &  dYdY0
[static]
 

Propagates a given state vector and computes the state transition matrix for elliptical Keplerian orbits.

Parameters:
GM Gravitational coefficient
Y0 Epoch state vector (position and velocity)
dt Time since epoch
Y State vector (position and velocity)
dYdY0 State transition matrix d(x,y,z,vx,vy,vz)/d(x,y,z,vx,vy,vz)_0
Warning:
The state vector, dt and GM must be given in consistent units. Due to the internal use of Keplerian elements, the function cannot be used with epoch state vectors describing a circular or non-inclined orbit.

Definition at line 543 of file KeplerOrbit.cpp.

References KeplerOrbit::Elements(), gpstk::sin(), gpstk::sqrt(), KeplerOrbit::State(), KeplerOrbit::StatePartials(), and gpstk::tan().


The documentation for this class was generated from the following files:
Generated on Thu May 23 03:31:38 2013 for GPS ToolKit Software Library by  doxygen 1.3.9.1