get_horizons_coord#

sunpy.coordinates.get_horizons_coord(body, time='now', id_type=None, *, include_velocity=False)[source]#

Queries JPL HORIZONS and returns a SkyCoord for the location of a solar-system body at a specified time. This location is the instantaneous or “true” location, and is not corrected for light travel time or observer motion.

Parameters:
  • body (str) – The solar-system body for which to calculate positions. One can also use the search form linked below to find valid names or ID numbers.

  • id_type (None, str) – Defaults to None, which searches major bodies first, and then searches small bodies (comets and asteroids) if no major body is found. If 'smallbody', the search is limited to only small bodies. If 'designation', the search is limited to only small-body designations. If 'name', the search is limited to only small-body names. If 'asteroid_name' or 'comet_name', the search is limited to only asteroid names or only comet names, respectively.

  • time (tuple, list, str, pandas.Timestamp, pandas.Series, pandas.DatetimeIndex, datetime.datetime, datetime.date, numpy.datetime64, numpy.ndarray, astropy.time.Time, dict) – Time to use in a parse_time-compatible format.

    Alternatively, this can be a dictionary defining a range of times and dates; the range dictionary has to be of the form {‘start’: start_time, ‘stop’: stop_time, ‘step’:’n[y|d|m|s]’}. start_time and stop_time must be in a parse_time-compatible format, and are interpreted as UTC time. step must be a string with either a number and interval length (e.g. for every 10 seconds, '10s'), or a plain number for a number of evenly spaced intervals.

  • include_velocity (bool, optional) – If True, include the body’s velocity in the output coordinate. Defaults to False.

Returns:

SkyCoord – Location of the solar-system body

Notes

Be aware that there can be discrepancies between the coordinates returned by JPL HORIZONS, the coordinates reported in mission data files, and the coordinates returned by get_body_heliographic_stonyhurst.

References

Examples

Obtaining a spacecraft trajectory from JPL Horizons

Obtaining a spacecraft trajectory from JPL Horizons

Setting the correct position for SOHO in a LASCO C3 Map

Setting the correct position for SOHO in a LASCO C3 Map

Reproducing the “Where is STEREO Today?” plot

Reproducing the "Where is STEREO Today?" plot
>>> from sunpy.coordinates import get_horizons_coord

Query the location of Venus

>>> get_horizons_coord('Venus barycenter', '2001-02-03 04:05:06')  
INFO: Obtained JPL HORIZONS location for Venus Barycenter (2) [sunpy.coordinates.ephemeris]
<SkyCoord (HeliographicStonyhurst: obstime=2001-02-03T04:05:06.000, rsun=695700.0 km): (lon, lat, radius) in (deg, deg, AU)
    (-33.93155836, -1.64998443, 0.71915147)>

Query the location of the SDO spacecraft

>>> get_horizons_coord('SDO', '2011-11-11 11:11:11')  
INFO: Obtained JPL HORIZONS location for Solar Dynamics Observatory (spacecraft) (-136395) [sunpy.coordinates.ephemeris]
<SkyCoord (HeliographicStonyhurst: obstime=2011-11-11T11:11:11.000, rsun=695700.0 km): (lon, lat, radius) in (deg, deg, AU)
    (0.01019118, 3.29640728, 0.99011042)>

Query the location of the SOHO spacecraft via its ID number (-21)

>>> get_horizons_coord(-21, '2004-05-06 11:22:33')  
INFO: Obtained JPL HORIZONS location for SOHO (spacecraft) (-21) [sunpy.coordinates.ephemeris]
<SkyCoord (HeliographicStonyhurst: obstime=2004-05-06T11:22:33.000, rsun=695700.0 km): (lon, lat, radius) in (deg, deg, AU)
    (0.25234902, -3.55863633, 0.99923086)>

Query the location and velocity of the asteroid Juno

>>> get_horizons_coord('Juno', '1995-07-18 07:17', 'smallbody', include_velocity=True)  
INFO: Obtained JPL HORIZONS location for 3 Juno (A804 RA) [sunpy.coordinates.ephemeris]
<SkyCoord (HeliographicStonyhurst: obstime=1995-07-18T07:17:00.000, rsun=695700.0 km): (lon, lat, radius) in (deg, deg, AU)
    (-25.16107532, 14.59098438, 3.17667664)
 (d_lon, d_lat, d_radius) in (arcsec / s, arcsec / s, km / s)
    (-0.03306548, 0.00052415, -2.66709222)>

Query the location of Solar Orbiter at a set of 12 regularly sampled times

>>> get_horizons_coord('Solar Orbiter',
...                    time={'start': '2020-12-01',
...                           'stop': '2020-12-02',
...                           'step': '12'})  
INFO: Obtained JPL HORIZONS location for Solar Orbiter (spacecraft) (-144) [sunpy.coordinates.ephemeris]
...