====================================== Orbit Data Objects (`sbpy.data.Orbit`) ====================================== Orbit Queries ============= `~sbpy.data.Orbit.from_horizons` enables the query of Solar System body osculating elements from the `JPL Horizons service `_: .. doctest-remote-data:: >>> from sbpy.data import Orbit >>> from astropy.time import Time >>> epoch = Time('2018-05-14', scale='utc') >>> elem = Orbit.from_horizons('Ceres', epochs=epoch) >>> elem # doctest: +SKIP <...>/sbpy/data/orbit.py:113: TimeScaleWarning: converting utc epochs to tdb for use in astroquery.jplhorizons TimeScaleWarning) targetname H G ... P epoch mag ... d str7 float64 float64 ... float64 object ---------- ------- ------- ... ----------------- ----------------- 1 Ceres 3.34 0.12 ... 1681.218136185659 2458252.500800755 Please note the ``TimeScaleWarning``, which is raised when the time scale of the desired epoch is not supported by the query function and hence converted to a scale that is supported (``tdb`` in this case). The following fields are available in this object: >>> elem.field_names # doctest: +REMOTE_DATA ['targetname', 'H', 'G', 'e', 'q', 'incl', 'Omega', 'w', 'n', 'M', 'nu', 'a', 'Q', 'P', 'epoch', 'Tp'] If ``epochs`` is not set, the osculating elements for the current epoch (current time) are queried. Similar to `~sbpy.data.Ephem.from_horizons`, this function is a wrapper for `~astroquery.jplhorizons.HorizonsClass.elements` and passes optional parameter on to that function. Furthermore, it is possible to query orbital elements for a number of targets: .. doctest-remote-data:: >>> epoch = Time('2018-08-03 14:20', scale='tdb') >>> elem = Orbit.from_horizons(['3749', '2009 BR60'], ... epochs=epoch, ... refplane='earth') >>> elem # doctest: +SKIP targetname H G ... P epoch mag ... d str21 float64 float64 ... float64 object --------------------- ------- ------- ... ----------------- ----------------- 3749 Balam (1982 BG1) 13.3 0.15 ... 1221.865723743203 2458334.097222222 312497 (2009 BR60) 17.7 0.15 ... 1221.776912893334 2458334.097222222 Alternatively, orbital elements can also be queried from the `Minor Planet Center `_, although in this case only the most recent elements are accessible: .. doctest-remote-data:: >>> elem = Orbit.from_mpc(['3552', '12893']) >>> elem # doctest: +SKIP absmag Q arc w ... a Tj moid_uranus moid_venus mag AU d deg ... AU AU AU float64 float64 float64 float64 ... float64 float64 float64 float64 ------- ------- ------- --------- ... --------- ------- ----------- ---------- 12.9 7.278 12955.0 316.44802 ... 4.2589272 2.3 11.7518 0.56105 13.9 3.028 12990.0 184.31369 ... 2.8281991 3.3 15.1738 1.90535 Orbit Transformations ===================== An existing `~Orbit` instance can be transformed to a different orbital element definition system (e.g., Keplerian, cometary, cartesian) using `~sbpy.data.Orbit.oo_transform` or it can be propagated into the future or past using `~sbpy.data.Orbit.oo_propagate`. Both functions are implemented in `sbpy` to provide an interface to `pyoorb `_, a Python module using `OpenOrb `_. In order to transform some current orbits to a state vector in cartesian coordinates, one could use the following code: .. doctest-remote-data:: >>> elem = Orbit.from_horizons(['Ceres', 'Pallas', 'Vesta']) >>> statevec = elem.oo_transform('CART') # doctest: +SKIP >>> statevec # doctest: +SKIP id x y ... H G AU AU ... mag str8 float64 float64 ... float64 float64 -------- ------------------- ------------------- ... ------- ------- 1 Ceres -0.4867631007775121 -2.7702346649193696 ... 3.34 0.12 2 Pallas -1.7745931352186222 -1.7169356664520194 ... 4.13 0.11 4 Vesta 2.24552918427612 1.0169886872736296 ... 3.2 0.32 Orbits can currently be transformed to the following definitions: cartesian (``'CART'``), Keplerian (``'KEP'``), and cometary (``'COM'``). Orbit Propagations ================== Orbit propagation requires the epoch to which the orbit should be propagated to either as `~astropy.time.Time` object, or as float in terms of Julian date. The following example propagates the current orbit of Ceres back to year 2000: .. doctest-remote-data:: >>> elem = Orbit.from_horizons('Ceres') >>> epoch = Time('2000-01-01', scale='tdb') >>> newelem = elem.oo_propagate(epoch) >>> newelem # doctest: +SKIP id a e ... epoch H G AU ... mag str7 float64 float64 ... object float64 float64 ------- ----------------- ------------------- ... --------- ------- ------- 1 Ceres 2.766494220549446 0.07837504411299284 ... 2451544.5 3.34 0.12 Note that both functions require `pyoorb `_ to be installed. Calculate dynamical parameters ============================== The Tisserand parameter is a commonly used dynamic parameter to characterize the orbit of a small body, especially a comet, when its orbital evolution is dominated by the gravitational effect of a particular planet. The Tisserand parameter with respect to Jupiter is used in the dynamical classification of comets. The Tisserand parameter can be calculated by `~sbpy.Orbit.tisserand` as follows: .. doctest-remote-data:: >>> epoch = Time(2449400.5, format='jd', scale='tdb') >>> halley = Orbit.from_horizons('1P', id_type='designation', ... closest_apparition=True, epochs=epoch) >>> T = halley.tisserand() >>> print('{:.4f}'.format(T)) # doctest: +SKIP -0.6050 One can also specify the planet with respect to which the Tisserand parameter is calculated with optional parameter `planet`. It also allows multiple planet to be specified simultaneously: .. doctest-remote-data:: >>> import numpy as np >>> chariklo = Orbit.from_horizons('chariklo', id_type='name') >>> T = chariklo.tisserand(planet=['599', '699', '799', '899']) >>> with np.printoptions(precision=3): ... print(T) # doctest: +FLOAT_CMP [3.485 2.931 2.858 3.224] `~sbpy.Orbit` also provides a method to compare the orbits of two objects in terms of the "D-criterion" (`Jopek 1993 `_). The `~sbpy.Orbit.D_criterion` method implements all three versions of the D-criterion, including Southworth & Hawkins function (`Southworth and Hawkins 1963 `_), Drummond function (`Drummond 1991 `_), and the hybrid function (`Jopek 1993 `_). The code example below demonstrates the calculation of three versions of D_criterion: .. doctest-remote-data:: >>> comets = Orbit.from_horizons(['252P', 'P/2016 BA14'], ... id_type='designation', closest_apparition=True) >>> >>> # Southworth & Hawkins function >>> D_SH = comets[0].D_criterion(comets[1]) >>> # Drummond function >>> D_D = comets[0].D_criterion(comets[1], version='d') >>> # hybrid function >>> D_H = comets[0].D_criterion(comets[1], version='h') >>> print('D_SH = {:.4f}, D_D = {:.4f}, D_H = {:.4f}'. ... format(D_SH, D_D, D_H)) D_SH = 0.1560, D_D = 0.0502, D_H = 0.1556