poliastro.core.angles

Module Contents

Functions

newton_factory(func, fprime)

D_to_nu(D)

True anomaly from parabolic anomaly.

nu_to_D(nu)

Parabolic anomaly from true anomaly.

nu_to_E(nu, ecc)

Eccentric anomaly from true anomaly.

nu_to_F(nu, ecc)

Hyperbolic anomaly from true anomaly.

E_to_nu(E, ecc)

True anomaly from eccentric anomaly.

F_to_nu(F, ecc)

True anomaly from hyperbolic anomaly.

M_to_E(M, ecc)

Eccentric anomaly from mean anomaly.

M_to_F(M, ecc)

Hyperbolic anomaly from mean anomaly.

M_to_D(M)

Parabolic anomaly from mean anomaly.

E_to_M(E, ecc)

Mean anomaly from eccentric anomaly.

F_to_M(F, ecc)

Mean anomaly from eccentric anomaly.

D_to_M(D)

Mean anomaly from parabolic anomaly.

fp_angle(nu, ecc)

Returns the flight path angle.

poliastro.core.angles.newton_factory(func, fprime)
poliastro.core.angles.D_to_nu(D)

True anomaly from parabolic anomaly.

Parameters

D (float) – Eccentric anomaly.

Returns

nu – True anomaly.

Return type

float

Notes

From 1:

\[\nu = 2 \arctan{D}\]
poliastro.core.angles.nu_to_D(nu)

Parabolic anomaly from true anomaly.

Parameters

nu (float) – True anomaly in radians.

Returns

D – Parabolic anomaly.

Return type

float

Warning

The parabolic anomaly will be continuous in (-∞, ∞) only if the true anomaly is in (-π, π]. No validation or wrapping is performed.

Notes

The treatment of the parabolic case is heterogeneous in the literature, and that includes the use of an equivalent quantity to the eccentric anomaly: 1 calls it “parabolic eccentric anomaly” D, 2 also uses the letter D but calls it just “parabolic anomaly”, 3 uses the letter B citing indirectly 4 (which however calls it “parabolic time argument”), and 5 does not bother to define it.

We use this definition:

\[B = \tan{\frac{\nu}{2}}\]

References

1(1,2)

Farnocchia, Davide, Davide Bracali Cioci, and Andrea Milani. “Robust resolution of Kepler’s equation in all eccentricity regimes.”

2

Bate, Muller, White.

3(1,2,3,4,5,6)

Vallado, David. “Fundamentals of Astrodynamics and Applications”, 2013.

4

IAU VIth General Assembly, 1938.

5(1,2,3)

Battin, Richard H. “An introduction to the Mathematics and Methods of Astrodynamics, Revised Edition”, 1999.

poliastro.core.angles.nu_to_E(nu, ecc)

Eccentric anomaly from true anomaly.

New in version 0.4.0.

Parameters
  • nu (float) – True anomaly in radians.

  • ecc (float) – Eccentricity.

Returns

E – Eccentric anomaly, between -π and π radians.

Return type

float

Warning

The eccentric anomaly will be between -π and π radians, no matter the value of the true anomaly.

Notes

The implementation uses the half-angle formula from 3:

\[E = 2 \arctan \left ( \sqrt{\frac{1 - e}{1 + e}} \tan{\frac{\nu}{2}} \right) \in (-\pi, \pi]\]
poliastro.core.angles.nu_to_F(nu, ecc)

Hyperbolic anomaly from true anomaly.

Parameters
  • nu (float) – True anomaly in radians.

  • ecc (float) – Eccentricity (>1).

Returns

F – Hyperbolic anomaly.

Return type

float

Warning

The hyperbolic anomaly will be continuous in (-∞, ∞) only if the true anomaly is in (-π, π], which should happen anyway because the true anomaly is limited for hyperbolic orbits. No validation or wrapping is performed.

Notes

The implementation uses the half-angle formula from 3:

\[F = 2 \operatorname{arctanh} \left( \sqrt{\frac{e-1}{e+1}} \tan{\frac{\nu}{2}} \right)\]
poliastro.core.angles.E_to_nu(E, ecc)

True anomaly from eccentric anomaly.

New in version 0.4.0.

Parameters
  • E (float) – Eccentric anomaly in radians.

  • ecc (float) – Eccentricity.

Returns

nu – True anomaly, between -π and π radians.

Return type

float

Warning

The true anomaly will be between -π and π radians, no matter the value of the eccentric anomaly.

Notes

The implementation uses the half-angle formula from 3:

\[\nu = 2 \arctan \left( \sqrt{\frac{1 + e}{1 - e}} \tan{\frac{E}{2}} \right) \in (-\pi, \pi]\]
poliastro.core.angles.F_to_nu(F, ecc)

True anomaly from hyperbolic anomaly.

Parameters
  • F (float) – Hyperbolic anomaly.

  • ecc (float) – Eccentricity (>1).

Returns

nu – True anomaly.

Return type

float

Notes

The implementation uses the half-angle formula from 3:

\[\nu = 2 \arctan \left( \sqrt{\frac{e + 1}{e - 1}} \tanh{\frac{F}{2}} \right) \in (-\pi, \pi]\]
poliastro.core.angles.M_to_E(M, ecc)

Eccentric anomaly from mean anomaly.

New in version 0.4.0.

Parameters
  • M (float) – Mean anomaly in radians.

  • ecc (float) – Eccentricity.

Returns

E – Eccentric anomaly.

Return type

float

Notes

This uses a Newton iteration on the Kepler equation.

poliastro.core.angles.M_to_F(M, ecc)

Hyperbolic anomaly from mean anomaly.

Parameters
  • M (float) – Mean anomaly in radians.

  • ecc (float) – Eccentricity (>1).

Returns

F – Hyperbolic anomaly.

Return type

float

Notes

This uses a Newton iteration on the hyperbolic Kepler equation.

poliastro.core.angles.M_to_D(M)

Parabolic anomaly from mean anomaly.

Parameters

M (float) – Mean anomaly in radians.

Returns

D – Parabolic anomaly.

Return type

float

Notes

This uses the analytical solution of Barker’s equation from 5.

poliastro.core.angles.E_to_M(E, ecc)

Mean anomaly from eccentric anomaly.

New in version 0.4.0.

Parameters
  • E (float) – Eccentric anomaly in radians.

  • ecc (float) – Eccentricity.

Returns

M – Mean anomaly.

Return type

float

Warning

The mean anomaly will be outside of (-π, π] if the eccentric anomaly is. No validation or wrapping is performed.

Notes

The implementation uses the plain original Kepler equation:

\[M = E - e \sin{E}\]
poliastro.core.angles.F_to_M(F, ecc)

Mean anomaly from eccentric anomaly.

Parameters
  • F (float) – Hyperbolic anomaly.

  • ecc (float) – Eccentricity (>1).

Returns

M – Mean anomaly.

Return type

float

Notes

As noted in 5, by manipulating the parametric equations of the hyperbola we can derive a quantity that is equivalent to the eccentric anomaly in the elliptic case:

\[M = e \sinh{F} - F\]
poliastro.core.angles.D_to_M(D)

Mean anomaly from parabolic anomaly.

Parameters

D (float) – Parabolic anomaly.

Returns

M – Mean anomaly.

Return type

float

Notes

We use this definition:

\[M = B + \frac{B^3}{3}\]

Notice that M < ν until ν ~ 100 degrees, then it reaches π when ν ~ 120 degrees, and grows without bounds after that. Therefore, it can hardly be called an “anomaly” since it is by no means an angle.

poliastro.core.angles.fp_angle(nu, ecc)

Returns the flight path angle.

Parameters
  • nu (float) – True anomaly in radians.

  • ecc (float) – Eccentricity.

Returns

fp_angle – Flight path angle

Return type

float

Notes

From 3, pp. 113:

\[\phi = \arctan(\frac {e \sin{\nu}}{1 + e \cos{\nu}})\]