# Coordinates

“Coordinates” are the numbers we use to name astronomical positions. For any given position, Skyfield can return any of at least a half-dozen different sets of coordinates naming that position.

Through the discussion that follows, remember that coordinates are merely names for a position. Asking for a different coordinate doesn’t change the position itself. As explained in the Positions chapter, there is a real physical difference between, say, an astrometric position and an apparent position: they are two distinct points in the sky. But given a single position — a single point in the sky — all the various coordinates you can produce mean exactly the same thing: the direction you would point a telescope to put that point at its center.

## Cartesian coordinates versus Spherical coordinates

A given coordinate system can be used to produce two kinds of coordinate: Cartesian and spherical.

Cartesian coordinates specify a position by providing three distances (x,y,z) that are measured along axes that stand at right angles to each other.

Spherical coordinates instead specify a location using two angles plus a distance. Just as a geographic latitude and longitude designate a position on the Earth, so the two angles of a spherical coordinate pick out a particular point in the sky — a particular point on the celestial sphere — and the distance says how far you would need to travel in that direction to reach the target.

Through the following sections, we will see how to ask for both Cartesian and spherical coordinates in a variety of coordinate systems.

## ICRS right ascension and declination

Equatorial coordinates, which measure angles from the Earth’s equator and poles, have already been described in the Positions chapter:

- Astrometric right ascension and declination

- Apparent right ascension and declination

Here’s how the axes and angles are defined:

> Equatorial Coordinates
> xy-plane: Earth’s equator
> x-axis: March equinox
> z-axis: North celestial pole
> ↔ Right ascension 0h–24h east around the equator
> ↕ Declination ±90° from equator toward poles

There’s a problem with this coordinate system: it slowly changes. Thanks to precession, the Earth’s poles are in constant — if gradual — motion. So equatorial coordinates always have to specify whether they are measuring against the Earth’s equator and equinox as of 1950, or 2000, or some other specific year. In the late twentieth century “J2000” coordinates became popular, measured against the Earth’s mean equator as of noon Terrestrial Time on 2000 January 1.

But in 1998, astronomers decided it was time to get off of the equatorial treadmill. They defined a new coordinate system: the ICRS (International Celestial Reference System). Its axes point almost exactly along the old J2000 axes, but are permanent and immobile. It can achieve unbounded precision because, instead of using the Earth, its axes are tied to the locations of very distant quasars.

In the Positions chapter you’ve already seen how to retrieve ICRS coordinates in both Cartesian and spherical form:

    # Building an example position.

    from skyfield.api import load

    ts = load.timescale()
    t = ts.utc(2021, 12, 2, 14, 7)
    eph = load('de421.bsp')
    earth, mars = eph['earth'], eph['mars']
    position = earth.at(t).observe(mars)

    # Producing coordinates.

    print('Cartesian ICRS:')

    x, y, z = position.xyz.au

    print('  x = {:.3f} au'.format(x))
    print('  y = {:.3f} au'.format(y))
    print('  z = {:.3f} au'.format(z))
    print()

    print('Spherical ICRS:')

    ra, dec, distance = position.radec()

    print(' ', ra, 'right ascension')
    print(' ', dec, 'declination')
    print(' ', distance, 'distance')

    Cartesian ICRS:
      x = -1.521 au
      y = -1.800 au
      z = -0.772 au

    Spherical ICRS:
      15h 19m 15.08s right ascension
      -18deg 08' 37.0" declination
      2.47983 au distance

Note that `.xyz` is a `Distance` object and supports not only `.au` but also several other units of measure.

## Equatorial right ascension and declination

If instead of using the permanent ICRS you want right ascension and declination measured against the Earth’s true equator and axes, then retrieving spherical coordinates is easy. Simply provide a date to `radec()` and it will return a right ascension and declination relative to the Earth’s equator and equinox as of that date. Or provide the string `'date'` and the date of the position itself is used. You can find example code in the ‘Positions’ chapter in the sections Astrometric right ascension and declination and Apparent right ascension and declination.

But maybe your situation is more complicated. You might want (x,y,z) coordinates instead of right ascension and declination. Or maybe you want coordinates measured against the mean equator and equinox instead of the true equator and equinox.

In that case you will want to import one of these two reference frames:

- `true_equator_and_equinox_of_date` — a reference frame that accounts for nearly every wobble of the Earth’s real axis. The effects of both precession and nutation are included. (The only effect omitted is the tiny daily offset 𝑥 and 𝑦 that you may have loaded from an IERS data file as described in the Polar Motion section.)

- `mean_equator_and_equinox_of_date` — a reference frame that includes only the smooth long-term motion of precession, but omits all of the smaller wobbles in the Earth’s axis that we call nutation. The mean equator and equinox are sometimes used to tabulate coordinates in books like the Astronomical Almanac.

You can get both spherical and (x,y,z) coordinates for either of these reference frames by importing them and using them like this:

    from skyfield.framelib import true_equator_and_equinox_of_date

    print('Spherical equinox-of-date coordinates:')

    dec, ra, distance = position.frame_latlon(true_equator_and_equinox_of_date)

    ra.preference = 'hours'  # print hours instead of degrees

    print(' ', ra, 'right ascension')
    print(' ', dec, 'declination')
    print(' ', distance, 'distance')

    print()
    print('Cartesian equinox-of-date coordinates:')

    x, y, z = position.frame_xyz(true_equator_and_equinox_of_date).au

    print('  x = {:.6f} au'.format(x))
    print('  y = {:.6f} au'.format(y))
    print('  z = {:.6f} au'.format(z))

    Spherical equinox-of-date coordinates:
      15h 20m 28.70s right ascension
      -18deg 13' 18.5" declination
      2.47983 au distance

    Cartesian equinox-of-date coordinates:
      x = -1.510302 au
      y = -1.807555 au
      z = -0.775435 au

The code will look almost exactly the same for the mean equator and equinox, except for which frame is imported and used.

## Altitude and azimuth (‘horizonal’ coordinates)

Altitude and azimuth have already been explained in the Positions chapter, so you can start reading about them there:

- Azimuth and altitude from a geographic position

The coordinate system is called horizonal in the sense of “pertaining to the horizon.”

> Horizonal Coordinates
> xy-plane: Horizon
> x-axis: North point on the horizon
> y-axis: East point on the horizon (left-handed)
> z-axis: Zenith
> ↕ Altitude ±90° above or below horizon
> ↔ Azimuth 0°–360° measured clockwise from north

As with the equatorial system, the angles associated with horizontal coordinates are so popular that Skyfield provides a built-in method `altaz()` to retrieve them, while (x,y,z) coordinates require a call to `frame_xyz()` with the geographic location itself passed as the reference frame:

    # From the chapter on Positions:
    # computing altitude and azimuth.

    from skyfield.api import load, wgs84

    bluffton = wgs84.latlon(+40.8939, -83.8917)
    astrometric = (earth + bluffton).at(t).observe(mars)
    position = astrometric.apparent()

    print('Cartesian:')

    x, y, z = position.frame_xyz(bluffton).au

    print('  x = {:.3f} au north'.format(x))
    print('  y = {:.3f} au east'.format(y))
    print('  z = {:.3f} au up'.format(z))
    print()

    print('Spherical:')

    alt, az, distance = position.altaz()

    print('  Altitude:', alt)
    print('  Azimuth:', az)
    print('  Distance:', distance)

    Cartesian:
      x = -1.913 au north
      y = 1.200 au east
      z = 1.025 au up

    Spherical:
      Altitude: 24deg 24' 20.6"
      Azimuth: 147deg 54' 28.8"
      Distance: 2.47981 au

Note that some astronomers use the term “elevation” for what Skyfield calls “altitude”: the angle at which a target stands above the horizon. Obviously both words are ambiguous, since “elevation” can also mean a site’s vertical distance above sea level, and “altitude” can mean an airplane’s height above either sea level or the ground.

## Hour Angle and Declination

If you are pointing a telescope or other instrument, you might be interested in a variation on equatorial coordinates: replacing right ascension with hour angle, which measures ±180° from your own local meridian, whose hour angle is 0°.

    ha, dec, distance = position.hadec()

    print('Hour Angle:', ha)
    print('Declination:', dec, )
    print('Distance:', distance)

    Hour Angle: -02h 02m 28.88s
    Declination: -18deg 13' 16.4"
    Distance: 2.47981 au

To make the hour angle and declination even more useful for pointing real-world instruments, Skyfield includes the effect of polar motion if you have loaded a polar motion table. In that case the declination you get from `hadec()` will vary slightly from the declination returned by `radec()`, which doesn’t include polar motion.

## ECI versus ECEF coordinates

Here’s a quick explanation of two acronyms that you’re likely to run across in discussions about coordinates.

ECI stands for Earth-Centered Inertial and specifies coordinates that are (a) measured from the Earth’s center and (b) that don’t rotate with the Earth itself. The very first coordinates we computed in this chapter, for example, qualify as ECI coordinates, because the `position` used the Earth as its center and because the ICRS system of right ascension and declination stays fixed on the celestial sphere even as the Earth rotates beneath it.

ECEF stands for Earth-Centered Earth-Fixed and specifies coordinates that are (a) measured from the Earth’s center but (b) which rotate with the Earth instead of staying fixed in space. An example would be the latitude and longitude of the Lowell Observatory, which stays in place on the Earth’s surface but from the point of view of the rest of the Solar System is rotating with the Earth. A fixed on the Earth’s surface is a good example. We will learn about generating ECEF coordinates in the next section.

## Geographic ITRS latitude and longitude

Skyfield uses the standard ITRS reference frame to specify positions that are fixed relative to the Earth’s surface.

> ITRS Coordinates
> xy-plane: Earth’s equator
> x-axis: 0° longitude on the equator
> y-axis: 90° east longitude on the equator
> z-axis: North pole
> ↕ Latitude ±90° from equator toward poles
> ↔ Longitude ±180° from prime meridian with east positive

A location’s latitude will vary slightly depending on whether you model the Earth as a simple sphere or more realistically as a slightly flattened ellipsoid. The most popular choice today is to use the WGS84 ellipsoid, which is the one used by the GPS system.

    from skyfield.api import wgs84
    from skyfield.framelib import itrs

    # Important: must start with a position
    # measured from the Earth’s center.
    position = earth.at(t).observe(mars)

    print('Cartesian:')

    x, y, z = position.frame_xyz(itrs).au

    print('  x = {:.3f} au'.format(x))
    print('  y = {:.3f} au'.format(y))
    print('  z = {:.3f} au'.format(z))
    print()

    print('Geographic:')

    lat, lon = wgs84.latlon_of(position)
    height = wgs84.height_of(position)

    print(' {:.4f}° latitude'.format(lat.degrees))
    print(' {:.4f}° longitude'.format(lon.degrees))
    print(' {:.0f} km above sea level'.format(distance.km))

    Cartesian:
      x = 1.409 au
      y = -1.888 au
      z = -0.775 au

    Geographic:
     -18.2218° latitude
     -53.2662° longitude
     370974969 km above sea level

Note that height is measured from sea level, not from the center of the Earth.

The code above is slightly inefficient, because `height_of()` will wind up recomputing several values that were already computed in `latlon_of()`. If you need both, it’s more efficient to call `geographic_position_of()`.

There’s also a `subpoint_of()` method if you want Skyfield to compute the geographic position of the sea-level point beneath a given celestial object.

## Ecliptic coordinates

Ecliptic coordinates are measured from the plane of the Earth’s orbit. They are useful when making maps and diagrams of the Solar System and when exploring the properties of orbits around the Sun, because they place the orbits of the major planets nearly flat against the xy-plane — unlike right ascension and declination, which twist the Solar System up at a 23° angle because of the tilt of the Earth’s axis.

You might be tempted to ask why we measure against the plane of the Earth’s orbit, instead of averaging together all the planets to compute the “invariable plane” of the whole Solar System (to which the Earth’s orbit is inclined by something like 1.57°). The answer is: precision. We know the plane of the Earth’s orbit to many decimal places, because the Earth carries all of our highest-precision observatories along with it as it revolves around the Sun. Our estimate of the invariable plane, by contrast, is a mere average that changes — at least slightly — every time we discover a new asteroid, comet, or trans-Neptunian object. So the Earth’s own orbit winds up being the most pragmatic choice for a coordinate system oriented to the Solar System.

> Ecliptic Coordinates
> xy-plane: Ecliptic plane (plane of Earth’s orbit)
> x-axis: March equinox
> z-axis: North ecliptic pole
> ↕ Latitude ±90° above or below the ecliptic
> ↔ Longitude 0°–360° measured east from March equinox

Let’s take the `position` we already calculated above, which happens to measure the position of Mars as viewed from Earth’s center, and express that vector’s direction in ecliptic coordinates:

    from skyfield.framelib import ecliptic_frame

    print('Cartesian ecliptic coordinates:')

    x, y, z = position.frame_xyz(ecliptic_frame).au

    print('  x = {:.3f} au'.format(x))
    print('  y = {:.3f} au'.format(y))
    print('  z = {:.3f} au'.format(z))
    print()

    print('Spherical ecliptic coordinates:')

    lat, lon, distance = position.frame_latlon(ecliptic_frame)

    print(' {:.4f} latitude'.format(lat.degrees))
    print(' {:.4f} longitude'.format(lon.degrees))
    print(' {:.3f} au distant'.format(distance.au))

    Cartesian ecliptic coordinates:
      x = -1.510 au
      y = -1.967 au
      z = 0.007 au

    Spherical ecliptic coordinates:
     0.1732 latitude
     232.4801 longitude
     2.480 au distant

Note the very small values returned for the ecliptic z coordinate and for the ecliptic latitude, because Mars revolves around the Sun in very nearly the same plane as the Earth.

## Galactic coordinates

Galactic coordinates are measured against the disc of our own Milky Way galaxy, as measured from our vantage point here inside the Orion Arm:

> Galactic Coordinates
> xy-plane: Galactic plane
> x-axis: Galactic center
> z-axis: North galactic pole
> ↕ Latitude ±90° above galactic plane
> ↔ Longitude 0°–360° east from galactic center

    from skyfield.framelib import galactic_frame

    print('Cartesian galactic coordinates:')

    x, y, z = position.frame_xyz(galactic_frame).au

    print('  x = {:.3f} au'.format(x))
    print('  y = {:.3f} au'.format(y))
    print('  z = {:.3f} au'.format(z))
    print()

    print('Spherical galactic coordinates:')

    lat, lon, distance = position.frame_latlon(galactic_frame)

    print(' {:.4f} latitude'.format(lat.degrees))
    print(' {:.4f} longitude'.format(lon.degrees))
    print(' {:.3f} au distant'.format(distance.au))

    Cartesian galactic coordinates:
      x = 2.029 au
      y = -0.527 au
      z = 1.324 au

    Spherical galactic coordinates:
     32.2664 latitude
     345.4330 longitude
     2.480 au distant

Astronomers have generated a series of more and more precise estimates of our galaxy’s orientation over the past hundred years. Skyfield uses the IAU 1958 Galactic System II, which is believed to be accurate to within ±0.1°.

## Velocity

If you need not only a position vector but a velocity vector relative to a particular reference frame, then switch from the `frame_xyz()` method to the `frame_xyz_and_velocity()` method. As explained in its documentation, its return values include the velocity vector.

## Turning coordinates into a position

All of the above examples take a Skyfield position and return coordinates, but sometimes you start with coordinates and want to produce a position.

If you happen to start with ICRS (x,y,z) coordinates, then you can create a position with a function call:

    from skyfield.positionlib import build_position

    icrs_xyz_au = [-1.521, -1.800, -0.772]
    position = build_position(icrs_xyz_au, t=t)

But it’s probably more common for you to have been given a target’s right ascension and declination. One common approach is to create a `Star` as described in the Stars and Distant Objects chapter, which will give you an object that you can pass to `.observe()` like any other Skyfield body. But you can also create a position directly by using the `from_radec()` method carried by each position class. To create an apparent position, for example:

    from skyfield.positionlib import Apparent

    position = Apparent.from_radec(ra_hours=5.59, dec_degrees=5.45)

A final common situation is that you measured an altitude and azimuth relative to your horizon and want to learn the right ascension and declination of that position. The solution is to use the `from_altaz()` method, but there’s a catch: because the true coordinates behind any particular altitude and azimuth are changing every moment as the Earth spins, you first need to compute the position of your observatory `.at()` the moment you measured the altitude and azimuth, and only then call the method.

    # What are the coordinates of the zenith?

    b = bluffton.at(t)
    apparent = b.from_altaz(alt_degrees=90.0, az_degrees=0.0)

    ra, dec, distance = apparent.radec()
    print('Right ascension:', ra)
    print('Declination:', dec)

    Right ascension: 13h 17m 00.26s
    Declination: +41deg 00' 27.7"

If you find yourself in an even less common situation, like needing to build a position from ecliptic or galactic coordinates, then — while there aren’t yet any documented examples for you to follow — you might be able to assemble a solution together from these pieces:

- The position constructor method `from_time_and_frame_vectors()`.

- The `from_spherical(r, theta, phi)` method in `skyfield/functions.py`.

## Rotation Matrices

If you are doing some of your own mathematics, you might want access to the low-level 3×3 rotation matrices that define the relationship between each coordinate reference frame and the ICRS. To compute a rotation matrix, simply pass a time to the frame’s `rotation_at()` method:

    # 3×3 rotation matrix: ICRS → frame

    R = ecliptic_frame.rotation_at(t)
    print(R)

    [[ 0.99998613 -0.00482998 -0.00209857]
     [ 0.00526619  0.9174892   0.39772583]
     [ 0.00000441 -0.39773137  0.91750191]]

You should find the matrix easy to work with if `t` is a single time, but if `t` is a whole array of times then the resulting matrix will have a third dimension with the same number of elements as the time vector. NumPy provides no direct support for rotation matrices with an extra dimension, so avoid using NumPy’s multiplication operators. Instead, use Skyfield utility functions:

    from skyfield.functions import T, mxm, mxv

    T(R)        # reverse rotation matrix: frame → ICRS
    mxm(R, R2)  # matrix × matrix: combines rotations
    mxv(R, v)   # matrix × vector: rotates a vector