boinor.core.elements ==================== .. py:module:: boinor.core.elements .. autoapi-nested-parse:: This module contains a set of functions that can be used to convert between different elements that define the orbit of a body. Functions --------- .. autoapisummary:: boinor.core.elements.eccentricity_vector boinor.core.elements.circular_velocity boinor.core.elements.rv_pqw boinor.core.elements.coe_rotation_matrix boinor.core.elements.coe2rv boinor.core.elements.coe2rv_many boinor.core.elements.coe2mee boinor.core.elements.rv2coe boinor.core.elements.mee2coe boinor.core.elements.mee2rv Module Contents --------------- .. py:function:: eccentricity_vector(k, r, v) Eccentricity vector. .. math:: \vec{e} = \frac{1}{\mu} \left( \left( v^{2} - \frac{\mu}{r}\right ) \vec{r} - (\vec{r} \cdot \vec{v})\vec{v} \right) :param k: Standard gravitational parameter (km^3 / s^2). :type k: float :param r: Position vector (km) :type r: numpy.ndarray :param v: Velocity vector (km / s) :type v: numpy.ndarray .. py:function:: circular_velocity(k, a) Compute circular velocity for a given body given thegravitational parameter and the semimajor axis. .. math:: v = \sqrt{\frac{\mu}{a}} :param k: Gravitational Parameter :type k: float :param a: Semimajor Axis :type a: float .. py:function:: rv_pqw(k, p, ecc, nu) Returns r and v vectors in perifocal frame. :param k: Standard gravitational parameter (km^3 / s^2). :type k: float :param p: Semi-latus rectum or parameter (km). :type p: float :param ecc: Eccentricity. :type ecc: float :param nu: True anomaly (rad). :type nu: float :returns: * **r** (*numpy.ndarray*) -- Position. Dimension 3 vector * **v** (*numpy.ndarray*) -- Velocity. Dimension 3 vector .. rubric:: Notes These formulas can be checked at :cite:t:`Curtis2013`, 3rd Edition, page 110. Also the example proposed is 2.11 of :cite:t:`Curtis2013`, 3rd Edition book. .. math:: \vec{r} = \frac{h^2}{\mu}\frac{1}{1 + e\cos(\theta)}\begin{bmatrix} \cos(\theta)\\ \sin(\theta)\\ 0 \end{bmatrix} \\\\\\ \vec{v} = \frac{h^2}{\mu}\begin{bmatrix} -\sin(\theta)\\ e+\cos(\theta)\\ 0 \end{bmatrix} .. rubric:: Examples >>> from boinor.constants import GM_earth >>> k = GM_earth.value # Earth gravitational parameter >>> ecc = 0.3 # Eccentricity >>> h = 60000e6 # Angular momentum of the orbit (m**2 / s) >>> nu = np.deg2rad(120) # True Anomaly (rad) >>> p = h**2 / k # Parameter of the orbit >>> r, v = rv_pqw(k, p, ecc, nu) >>> # Printing the results r = [-5312706.25105345 9201877.15251336 0] [m] v = [-5753.30180931 -1328.66813933 0] [m]/[s] .. py:function:: coe_rotation_matrix(inc, raan, argp) Create a rotation matrix for coe transformation. .. py:function:: coe2rv(k, p, ecc, inc, raan, argp, nu) Converts from classical orbital to state vectors. Classical orbital elements are converted into position and velocity vectors by `rv_pqw` algorithm. A rotation matrix is applied to position and velocity vectors to get them expressed in terms of an IJK basis. :param k: Standard gravitational parameter (km^3 / s^2). :type k: float :param p: Semi-latus rectum or parameter (km). :type p: float :param ecc: Eccentricity. :type ecc: float :param inc: Inclination (rad). :type inc: float :param raan: Longitude of ascending node, omega (rad). :type raan: float :param argp: Argument of perigee (rad). :type argp: float :param nu: True anomaly (rad). :type nu: float :returns: * **r_ijk** (*numpy.ndarray*) -- Position vector in basis ijk. * **v_ijk** (*numpy.ndarray*) -- Velocity vector in basis ijk. .. rubric:: Notes .. math:: \begin{aligned} \vec{r}_{IJK} &= [ROT3(-\Omega)][ROT1(-i)][ROT3(-\omega)]\vec{r}_{PQW} = \left [ \frac{IJK}{PQW} \right ]\vec{r}_{PQW}\\ \vec{v}_{IJK} &= [ROT3(-\Omega)][ROT1(-i)][ROT3(-\omega)]\vec{v}_{PQW} = \left [ \frac{IJK}{PQW} \right ]\vec{v}_{PQW}\\ \end{aligned} Previous rotations (3-1-3) can be expressed in terms of a single rotation matrix: .. math:: \left [ \frac{IJK}{PQW} \right ] .. math:: \begin{bmatrix} \cos(\Omega)\cos(\omega) - \sin(\Omega)\sin(\omega)\cos(i) & -\cos(\Omega)\sin(\omega) - \sin(\Omega)\cos(\omega)\cos(i) & \sin(\Omega)\sin(i)\\ \sin(\Omega)\cos(\omega) + \cos(\Omega)\sin(\omega)\cos(i) & -\sin(\Omega)\sin(\omega) + \cos(\Omega)\cos(\omega)\cos(i) & -\cos(\Omega)\sin(i)\\ \sin(\omega)\sin(i) & \cos(\omega)\sin(i) & \cos(i) \end{bmatrix} .. py:function:: coe2rv_many(k, p, ecc, inc, raan, argp, nu) Parallel version of coe2rv. .. py:function:: coe2mee(p, ecc, inc, raan, argp, nu) Converts from classical orbital elements to modified equinoctial orbital elements. The definition of the modified equinoctial orbital elements is taken from :cite:t:`Walker1985`. The modified equinoctial orbital elements are a set of orbital elements that are useful for trajectory analysis and optimization. They are valid for circular, elliptic, and hyperbolic orbits. These direct modified equinoctial equations exhibit no singularity for zero eccentricity and orbital inclinations equal to 0 and 90 degrees. However, the h and k components are singular for an orbital inclination of 180 degrees since the retrograde factor is not yet supported. :param p: Semi-latus rectum or parameter (km). :type p: float :param ecc: Eccentricity. :type ecc: float :param inc: Inclination (rad). :type inc: float :param raan: Longitude of ascending node (rad). :type raan: float :param argp: Argument of perigee (rad). :type argp: float :param nu: True anomaly (rad). :type nu: float :returns: * **p** (*float*) -- Semi-latus rectum or parameter * **f** (*float*) -- Equinoctial parameter f * **g** (*float*) -- Equinoctial parameter g * **h** (*float*) -- Equinoctial parameter h * **k** (*float*) -- Equinoctial parameter k * **L** (*float*) -- Longitude .. rubric:: Notes The conversion equations are taken directly from the original paper: .. math:: \begin{aligned} p &= a(1-e^2) \\ f &= e\cos(\omega + \Omega) \\ g &= e\sin(\omega + \Omega) \\ h &= \tan(\frac{i}{2})\cos(\Omega) \\ k &= \tan(\frac{i}{2})\sin(\Omega) \\ L &= \Omega + \omega + \theta \\ \end{aligned} .. py:function:: rv2coe(k, r, v, tol=1e-08) Converts from vectors to classical orbital elements. :param k: Standard gravitational parameter (km^3 / s^2) :type k: float :param r: Position vector (km) :type r: numpy.ndarray :param v: Velocity vector (km / s) :type v: numpy.ndarray :param tol: Tolerance for eccentricity and inclination checks, default to 1e-8 :type tol: float, optional :returns: * **p** (*float*) -- Semi-latus rectum of parameter (km) * **ecc** (*float*) -- Eccentricity * **inc** (*float*) -- Inclination (rad) * **raan** (*float*) -- Right ascension of the ascending nod (rad) * **argp** (*float*) -- Argument of Perigee (rad) * **nu** (*float*) -- True Anomaly (rad) .. rubric:: Notes This example is a real exercise from :cite:t:`Curtis2013`, 3rd Edition, exercise is 4.3, page 200. 1. First the angular momentum is computed: .. math:: \vec{h} = \vec{r} \times \vec{v} 2. With it the eccentricity can be solved: .. math:: \begin{aligned} \vec{e} &= \frac{1}{\mu}\left [ \left ( v^{2} - \frac{\mu}{r}\right ) \vec{r} - (\vec{r} \cdot \vec{v})\vec{v} \right ] \\ e &= \sqrt{\vec{e}\cdot\vec{e}} \\ \end{aligned} 3. The node vector line is solved: .. math:: \begin{aligned} \vec{N} &= \vec{k} \times \vec{h} \\ N &= \sqrt{\vec{N}\cdot\vec{N}} \end{aligned} 4. The rigth ascension node is computed: .. math:: \Omega = \left\{ \begin{array}{lcc} cos^{-1}{\left ( \frac{N_{x}}{N} \right )} & if & N_{y} \geq 0 \\ \\ 360^{o} -cos^{-1}{\left ( \frac{N_{x}}{N} \right )} & if & N_{y} < 0 \\ \end{array} \right. 5. The argument of perigee: .. math:: \omega = \left\{ \begin{array}{lcc} cos^{-1}{\left ( \frac{\vec{N}\vec{e}}{Ne} \right )} & if & e_{z} \geq 0 \\ \\ 360^{o} -cos^{-1}{\left ( \frac{\vec{N}\vec{e}}{Ne} \right )} & if & e_{z} < 0 \\ \end{array} \right. 6. And finally the true anomaly: .. math:: \nu = \left\{ \begin{array}{lcc} cos^{-1}{\left ( \frac{\vec{e}\vec{r}}{er} \right )} & if & v_{r} \geq 0 \\ \\ 360^{o} -cos^{-1}{\left ( \frac{\vec{e}\vec{r}}{er} \right )} & if & v_{r} < 0 \\ \end{array} \right. .. rubric:: Examples >>> from boinor.bodies import Earth >>> from astropy import units as u >>> k = Earth.k.to_value(u.km ** 3 / u.s ** 2) >>> r = np.array([-6045., -3490., 2500.]) >>> v = np.array([-3.457, 6.618, 2.533]) >>> p, ecc, inc, raan, argp, nu = rv2coe(k, r, v) >>> print("p:", p, "[km]") # doctest: +FLOAT_CMP p: 8530.47436396927 [km] >>> print("ecc:", ecc) # doctest: +FLOAT_CMP ecc: 0.17121118195416898 >>> print("inc:", np.rad2deg(inc), "[deg]") # doctest: +FLOAT_CMP inc: 153.2492285182475 [deg] >>> print("raan:", np.rad2deg(raan), "[deg]") # doctest: +FLOAT_CMP raan: 255.27928533439618 [deg] >>> print("argp:", np.rad2deg(argp), "[deg]") # doctest: +FLOAT_CMP argp: 20.068139973005362 [deg] >>> print("nu:", np.rad2deg(nu), "[deg]") # doctest: +FLOAT_CMP nu: 28.445804984192122 [deg] .. py:function:: mee2coe(p, f, g, h, k, L) Converts from modified equinoctial orbital elements to classical orbital elements. The definition of the modified equinoctial orbital elements is taken from :cite:t:`Walker1985`. .. math:: \begin{aligned} p &= a(1 - e^{2})\\ e &= \sqrt{f^{2} + g^{2}}\\ i &= 2\arctan{(\sqrt{h^{2} + k^{2}})}\\ raan &= atan2(k, h) \pmod{2\pi}\\ argp &= (atan2(g, f) - raan) \pmod{2\pi}\\ nu &= (L - atan2(g, f)) \pmod{2\pi}\\ \end{aligned} :param p: Semi-latus rectum :type p: float :param f: Equinoctial parameter f :type f: float :param g: Equinoctial parameter g :type g: float :param h: Equinoctial parameter h :type h: float :param k: Equinoctial parameter k :type k: float :param L: Longitude :type L: float :returns: * **p** (*float*) -- Semi-latus rectum * **ecc** (*float*) -- Eccentricity of the orbit * **inc** (*float*) -- Inclination of the orbit * **raan** (*float*) -- RAAN of orbit * **argp** (*float*) -- Argument of the periapsis * **nu** (*float*) -- True anomaly .. rubric:: Notes The conversion is always safe because arctan2 works also for 0, 0 arguments. .. py:function:: mee2rv(p, f, g, h, k, L) Calculates position and velocity vector from modified equinoctial elements. :param p: Semi-latus rectum :type p: float :param f: Equinoctial parameter f :type f: float :param g: Equinoctial parameter g :type g: float :param h: Equinoctial parameter h :type h: float :param k: Equinoctial parameter k :type k: float :param L: Longitude :type L: float :returns: * **r** (*numpy.ndarray*) -- Position vector. * **v** (*numpy.ndarray*) -- Velocity vector. .. note:: The definition of `r` and `v` is taken from https://spsweb.fltops.jpl.nasa.gov/portaldataops/mpg/MPG_Docs/Source%20Docs/EquinoctalElements-modified.pdf, Equation 3a and 3b.