boinor.core.angles ================== .. py:module:: boinor.core.angles .. autoapi-nested-parse:: module containing functions related to angles in the core sub-package Functions --------- .. autoapisummary:: boinor.core.angles.newton_factory boinor.core.angles.D_to_nu boinor.core.angles.nu_to_D boinor.core.angles.nu_to_E boinor.core.angles.nu_to_F boinor.core.angles.E_to_nu boinor.core.angles.F_to_nu boinor.core.angles.M_to_E_scalar boinor.core.angles.M_to_E_vector boinor.core.angles.M_to_E boinor.core.angles.M_to_F_vector boinor.core.angles.M_to_F_scalar boinor.core.angles.M_to_F boinor.core.angles.M_to_D boinor.core.angles.E_to_M boinor.core.angles.F_to_M boinor.core.angles.D_to_M boinor.core.angles.fp_angle Module Contents --------------- .. py:function:: newton_factory(func, fprime) core function for the Newton method of finding zeroes .. py:function:: D_to_nu(D) True anomaly from parabolic anomaly. :param D: Eccentric anomaly. :type D: float :returns: **nu** -- True anomaly. :rtype: float .. rubric:: Notes From :cite:t:`Farnocchia2013`: .. math:: \nu = 2 \arctan{D} .. py:function:: nu_to_D(nu) Parabolic anomaly from true anomaly. :param nu: True anomaly in radians. :type nu: float :returns: **D** -- Parabolic anomaly. :rtype: float .. warning:: The parabolic anomaly will be continuous in (-∞, ∞) only if the true anomaly is in (-π, π]. No validation or wrapping is performed. .. rubric:: 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: :cite:t:`Farnocchia2013` calls it "parabolic eccentric anomaly" D, :cite:t:`Bate_et_al2020` also uses the letter D but calls it just "parabolic anomaly", :cite:p:`Vallado2013` uses the letter B citing indirectly :cite:t:`IAU_IV_GA_1938` (which however calls it "parabolic time argument"), and :cite:t:`Battin1999` does not bother to define it. We use this definition: .. math:: B = \tan{\frac{\nu}{2}} .. py:function:: nu_to_E(nu, ecc) Eccentric anomaly from true anomaly. .. versionadded:: 0.4.0 :param nu: True anomaly in radians. :type nu: float :param ecc: Eccentricity. :type ecc: float :returns: **E** -- Eccentric anomaly, between -π and π radians. :rtype: float .. warning:: The eccentric anomaly will be between -π and π radians, no matter the value of the true anomaly. .. rubric:: Notes The implementation uses the half-angle formula from :cite:t:`Vallado2013`: .. math:: E = 2 \arctan \left ( \sqrt{\frac{1 - e}{1 + e}} \tan{\frac{\nu}{2}} \right) \in (-\pi, \pi] .. py:function:: nu_to_F(nu, ecc) Hyperbolic anomaly from true anomaly. :param nu: True anomaly in radians. :type nu: float :param ecc: Eccentricity (>1). :type ecc: float :returns: **F** -- Hyperbolic anomaly. :rtype: 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. .. rubric:: Notes The implementation uses the half-angle formula from :cite:t:`Vallado2013`: .. math:: F = 2 \operatorname{arctanh} \left( \sqrt{\frac{e-1}{e+1}} \tan{\frac{\nu}{2}} \right) .. py:function:: E_to_nu(E, ecc) True anomaly from eccentric anomaly. .. versionadded:: 0.4.0 :param E: Eccentric anomaly in radians. :type E: float :param ecc: Eccentricity. :type ecc: float :returns: **nu** -- True anomaly, between -π and π radians. :rtype: float .. warning:: The true anomaly will be between -π and π radians, no matter the value of the eccentric anomaly. .. rubric:: Notes The implementation uses the half-angle formula from :cite:t:`Vallado2013`: .. math:: \nu = 2 \arctan \left( \sqrt{\frac{1 + e}{1 - e}} \tan{\frac{E}{2}} \right) \in (-\pi, \pi] .. py:function:: F_to_nu(F, ecc) True anomaly from hyperbolic anomaly. :param F: Hyperbolic anomaly. :type F: float :param ecc: Eccentricity (>1). :type ecc: float :returns: **nu** -- True anomaly. :rtype: float .. rubric:: Notes The implementation uses the half-angle formula from :cite:t:`Vallado2013`: .. math:: \nu = 2 \arctan \left( \sqrt{\frac{e + 1}{e - 1}} \tanh{\frac{F}{2}} \right) \in (-\pi, \pi] .. py:function:: M_to_E_scalar(M, ecc) Eccentric anomaly from mean anomaly (parameter is scalar). .. versionadded:: 0.4.0 :param M: Mean anomaly in radians. :type M: float :param ecc: Eccentricity. :type ecc: float :returns: **E** -- Eccentric anomaly. :rtype: float .. rubric:: Notes This uses a Newton iteration on the Kepler equation. .. py:function:: M_to_E_vector(M, ecc) Eccentric anomaly from mean anomaly (parameter is vector). .. versionadded:: 0.4.0 :param M: Mean anomaly in radians. :type M: np.array(float) :param ecc: Eccentricity. :type ecc: np.array(float) :returns: **E** -- Eccentric anomaly. :rtype: np.array(float) .. rubric:: Notes This uses a Newton iteration on the Kepler equation. .. py:function:: M_to_E(M, ecc) Eccentric anomaly from mean anomaly. .. versionadded:: 0.4.0 :param M: Mean anomaly in radians. :type M: float :param ecc: Eccentricity. :type ecc: float :returns: **E** -- Eccentric anomaly. :rtype: float .. rubric:: Notes This uses a Newton iteration on the Kepler equation. .. py:function:: M_to_F_vector(M, ecc) Hyperbolic anomaly from mean anomaly (parameter is vector). :param M: Mean anomaly in radians. :type M: float :param ecc: Eccentricity (>1). :type ecc: float :returns: **F** -- Hyperbolic anomaly. :rtype: float .. rubric:: Notes This uses a Newton iteration on the hyperbolic Kepler equation. .. py:function:: M_to_F_scalar(M, ecc) Hyperbolic anomaly from mean anomaly (parameter is scalar). :param M: Mean anomaly in radians. :type M: float :param ecc: Eccentricity (>1). :type ecc: float :returns: **F** -- Hyperbolic anomaly. :rtype: float .. rubric:: Notes This uses a Newton iteration on the hyperbolic Kepler equation. .. py:function:: M_to_F(M, ecc) Hyperbolic anomaly from mean anomaly. :param M: Mean anomaly in radians. :type M: float :param ecc: Eccentricity (>1). :type ecc: float :returns: **F** -- Hyperbolic anomaly. :rtype: float .. rubric:: Notes This uses a Newton iteration on the hyperbolic Kepler equation. .. py:function:: M_to_D(M) Parabolic anomaly from mean anomaly. :param M: Mean anomaly in radians. :type M: float :returns: **D** -- Parabolic anomaly. :rtype: float .. rubric:: Notes This uses the analytical solution of Barker's equation from :cite:t:`Battin1999`. .. py:function:: E_to_M(E, ecc) Mean anomaly from eccentric anomaly. .. versionadded:: 0.4.0 :param E: Eccentric anomaly in radians. :type E: float :param ecc: Eccentricity. :type ecc: float :returns: **M** -- Mean anomaly. :rtype: float .. warning:: The mean anomaly will be outside of (-π, π] if the eccentric anomaly is. No validation or wrapping is performed. .. rubric:: Notes The implementation uses the plain original Kepler equation: .. math:: M = E - e \sin{E} .. py:function:: F_to_M(F, ecc) Mean anomaly from hyperbolic anomaly. :param F: Hyperbolic anomaly. :type F: float :param ecc: Eccentricity (>1). :type ecc: float :returns: **M** -- Mean anomaly. :rtype: float .. rubric:: Notes As noted in :cite:t:`Battin1999`, by manipulating the parametric equations of the hyperbola we can derive a quantity that is equivalent to the mean anomaly in the elliptic case: .. math:: M = e \sinh{F} - F .. py:function:: D_to_M(D) Mean anomaly from parabolic anomaly. :param D: Parabolic anomaly. :type D: float :returns: **M** -- Mean anomaly. :rtype: float .. rubric:: Notes We use this definition: .. math:: 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. .. py:function:: fp_angle(nu, ecc) Returns the flight path angle. :param nu: True anomaly in radians. :type nu: float :param ecc: Eccentricity. :type ecc: float :returns: **fp_angle** -- Flight path angle :rtype: float .. rubric:: Notes From :cite:t:`Vallado2013`, pp. 113: .. math:: \phi = \arctan(\frac {e \sin{\nu}}{1 + e \cos{\nu}})