astropy:docs

Astronomical Coordinate Systems (astropy.coordinates)

Introduction

The coordinates package provides classes for representing celestial coordinates, as well as tools for converting between standard systems in a uniform way.

Warning

coordinates is currently a work-in-progress, and thus it is possible there will be significant API changes in later versions of Astropy. If you have specific ideas for how it might be improved, feel free to let us know on the astropy-dev mailing list or at http://feedback.astropy.org

Getting Started

Coordinate objects are instantiated with a flexible and natural approach that supports both numeric angle values, (limited) string parsing, and can optionally include lists of multiple coordinates in one object:

>>> from astropy.coordinates import ICRS, Galactic
>>> from astropy import units as u
>>> ICRS(ra=10.68458, dec=41.26917, unit=(u.degree, u.degree))
<ICRS RA=10.68458 deg, Dec=41.26917 deg>
>>> ICRS('00h42m44.3s +41d16m9s')
<ICRS RA=10.68458 deg, Dec=41.26917 deg>
>>> ICRS('00h42m44.3s +41d16m9s')
<ICRS RA=10.68458 deg, Dec=41.26917 deg>
>>> ICRS(ra=[10.68458, 83.82208], dec=[41.26917, -5.39111], unit=(u.degree, u.degree))
<ICRS RA=[ 10.68458  83.82208] deg, Dec=[ 41.26917  -5.39111] deg>

The individual components of a coordinate are Longitude or Latitude objects, which are specialized versions of the general Angle class. The component values are accessed using aptly named attributes:

>>> c = ICRS(ra=10.68458, dec=41.26917,
...          unit=(u.degree, u.degree))
>>> c.ra
<Longitude 10.684579999999983 deg>
>>> c.ra.hour
0.7123053333333323
>>> c.ra.hms
(0.0, 42.0, 44.299199999996262)
>>> c.dec
<Latitude 41.26917 deg>
>>> c.dec.radian
0.7202828960652683

Coordinates can easily be converted to strings using the to_string method:

>>> c.to_string()
'0h42m44.2992s 41d16m09.012s'
>>> c.to_string(precision=1)
'0h42m44.3s 41d16m09.0s'
>>> c.to_string(precision=1, sep=' ')
'0 42 44.3 41 16 09.0'

To convert to some other coordinate system, the easiest method is to use attribute-style access with short names for the built-in systems, but explicit transformations via the transform_to method are also available:

>>> c.galactic
<Galactic l=121.17430 deg, b=-21.57280 deg>
>>> c.transform_to(Galactic)
<Galactic l=121.17430 deg, b=-21.57280 deg>

Distance from the origin (which is system-dependent, but often the Earth center) can also be assigned to a coordinate. This specifies a unique point in 3D space, which also allows conversion to Cartesian coordinates:

>>> from astropy.coordinates import Distance
>>> c = ICRS(ra=10.68458, dec=41.26917,
...          unit=(u.degree, u.degree),
...          distance=Distance(770, u.kpc))
>>> c.x
<Quantity 568.7128654235232 kpc>
>>> c.y
<Quantity 107.30089740420232 kpc>
>>> c.z
<Quantity 507.88994291875713 kpc>

Coordinate objects can also store arrays of coordinates instead of a single coordinate. This has a major performance advantage over transforming many individual coordinate objects separtely. It also allows coordinate objects to be used to find matches between two sets of coordinates:

>>> #assume ra1/dec1 and ra2/dec2 are arrays loaded from some file
>>> c = ICRS(ra1, dec1, unit=(u.degree, u.degree))  
>>> catalog = ICRS(ra2, dec2, unit=(u.degree, u.degree))  
>>> idx, d2d, d3d = c1.match_to_catalog_sky(catalog)  

These array coordinates can also be indexed in the same way as numpy arrays:

>>> len(c[0].ra) 
TypeError: 'Longitude' object with a scalar value has no len()
>>> len(c[1:5].ra) 
4
>>> matches = catalog[idx]  
>>> len(matches) == len(c)  
True

The astropy.coordinates subpackage also provides a quick way to get coordinates for named objects (if you have an active internet connection). All subclasses of SphericalCoordinatesBase have a special class method, from_name(), that accepts a string and queries Sesame to retrieve coordinates for that object:

>>> c = ICRS.from_name("M42")
>>> c.ra, c.dec
(<Longitude 83.82208... deg>, <Latitude -5.39111... deg>)

This works for any subclass of SphericalCoordinatesBase:

>>> c = Galactic.from_name("M42")
>>> c.l, c.b
(<Longitude 3.64797... rad>, <Latitude -0.33827... rad>)

Note

This is intended to be a convenience, and is very simple. If you need precise coordinates for an object you should find the appropriate reference for that measurement and input the coordinates manually.

Using astropy.coordinates

More details of using astropy.coordinates are provided in the following sections:

In addition, another resource for the capabilities of this package is the astropy.coordinates.tests.test_api testing file. It showcases most of the major capabilities of the package, and hence is a useful supplement to this document. You can see it by either looking at it directly if you downloaded a copy of the astropy source code, or typing the following in an IPython session:

In [1]: from astropy.coordinates.tests import test_api
In [2]: test_api??

See Also

Some references particularly useful in understanding subtleties of the coordinate systems implemented here include:

  • Standards Of Fundamental Astronomy

    The definitive implementation of IAU-defined algorithms. The “SOFA Tools for Earth Attitude” document is particularly valuable for understanding the latest IAU standards in detail.

  • USNO Circular 179

    A useful guide to the IAU 2000/2003 work surrounding ICRS/IERS/CIRS and related problems in precision coordinate system work.

  • Meeus, J. “Astronomical Algorithms”

    A valuable text describing details of a wide range of coordinate-related problems and concepts.

Reference/API

astropy.coordinates Module

This subpackage contains classes and functions for celestial coordinates of astronomical objects. It also contains a framework for conversions between coordinate systems.

The diagram below shows all of the coordinate systems built into the coordinates package, their aliases (usable for converting other coordinates to them using attribute-style access) and the pre-defined transformations between them. The user is free to override any of these transformations by defining new trasnformation between these systems, but the pre-defined transformations should be sufficient for typical usage.

The graph also indicates the priority for each transformation as a number next to the arrow. These priorities are used to decide the preferred order when two trasnformation paths have the same number of steps. These priorities are defined such that path with a smaller total priority are favored over one with larger.

digraph AstropyCoordinateTransformGraph {
FK5 [shape=oval label="FK5\n`fk5`"]; FK4NoETerms [shape=oval label="FK4NoETerms\n`fk4_no_e`"]; Galactic [shape=oval label="Galactic\n`galactic`"]; ICRS [shape=oval label="ICRS\n`icrs`"]; FK4 [shape=oval label="FK4\n`fk4`"]; AltAz [shape=oval label="AltAz\n`altaz`"];
FK5 -> FK4NoETerms[ label = "1" ];
FK5 -> Galactic[ label = "1" ];
FK5 -> ICRS[ label = "1" ];
FK4NoETerms -> FK5[ label = "1" ];
FK4NoETerms -> FK4[ label = "1" ];
FK4NoETerms -> Galactic[ label = "1" ];
FK4 -> FK4NoETerms[ label = "1" ];
ICRS -> FK5[ label = "1" ];
Galactic -> FK5[ label = "1" ];
Galactic -> FK4NoETerms[ label = "1" ];

overlap=false
}

Functions

cartesian_to_spherical(x, y, z) Converts 3D rectangular cartesian coordinates to spherical polar coordinates.
coordinate_alias(name[, coordcls]) Gives a short name to this coordinate system, allowing other coordinate objects to convert to this one using attribute-style access.
dynamic_transform_matrix(fromsys, tosys[, ...]) A function decorator for defining transformations between coordinate systems using a function that yields a matrix.
get_icrs_coordinates(name) Retrieve an ICRS object by using an online name resolving service to retrieve coordinates for the specified name.
match_coordinates_3d(matchcoord, catalogcoord) Finds the nearest 3-dimensional matches of a coordinate or coordinates in a set of catalog coordinates.
match_coordinates_sky(matchcoord, catalogcoord) Finds the nearest on-sky matches of a coordinate or coordinates in a set of catalog coordinates.
spherical_to_cartesian(r, lat, lon) Converts spherical polar coordinates to rectangular cartesian coordinates.
static_transform_matrix(fromsys, tosys[, ...]) A function decorator for defining transformations between coordinate systems using a matrix.
transform_function(fromsys, tosys[, ...]) A function decorator for defining transformations between coordinate systems.

Classes

AltAz(*args, **kwargs) A coordinate in the altitude/azimuth or “horizontal” system.
Angle One or more angular value(s) with units equivalent to radians or degrees.
BoundsError Raised when an angle is outside of its user-specified bounds.
CartesianPoints A cartesian representation of a point in three-dimensional space.
CompositeStaticMatrixTransform(fromsys, ...) A MatrixTransform constructed by combining a sequence of matricies together.
ConvertError Raised if a coordinate system cannot be converted to another
Distance A one-dimensional distance.
DynamicMatrixTransform(fromsys, tosys, ...) A coordinate transformation specified as a function that yields a 3 x 3 cartesian transformation matrix.
FK4(*args, **kwargs) A coordinate in the FK4 system.
FK4NoETerms(*args, **kwargs) A coordinate in the FK4 system.
FK5(*args, **kwargs) A coordinate in the FK5 system.
FunctionTransform(fromsys, tosys, func[, ...]) A coordinate transformation defined by a function that simply accepts a coordinate object and returns the transformed coordinate object.
Galactic(*args, **kwargs) A coordinate in Galactic Coordinates.
ICRS(*args, **kwargs) A coordinate in the ICRS.
IllegalHourError(hour) Raised when an hour value is not in the range [0,24).
IllegalHourWarning(hour[, alternativeactionstr]) Raised when an hour value is 24.
IllegalMinuteError(minute) Raised when an minute value is not in the range [0,60].
IllegalMinuteWarning(minute[, ...]) Raised when a minute value is 60.
IllegalSecondError(second) Raised when an second value (time) is not in the range [0,60].
IllegalSecondWarning(second[, ...]) Raised when a second value is 60.
Latitude Latitude-like angle(s) which must be in the range -90 to +90 deg.
Longitude Longitude-like angle(s) which are wrapped within a contiguous 360 degree range.
RangeError Raised when some part of an angle is out of its valid range.
SphericalCoordinatesBase(*args, **kwargs) Abstract superclass for all coordinate classes representing points in three dimensions.
StaticMatrixTransform(fromsys, tosys, matrix) A coordinate transformation defined as a 3 x 3 cartesian transformation matrix.

Class Inheritance Diagram

Inheritance diagram of astropy.coordinates.builtin_systems.AltAz, astropy.coordinates.angles.Angle, astropy.coordinates.errors.BoundsError, astropy.coordinates.distances.CartesianPoints, astropy.coordinates.transformations.CompositeStaticMatrixTransform, astropy.coordinates.errors.ConvertError, astropy.coordinates.distances.Distance, astropy.coordinates.transformations.DynamicMatrixTransform, astropy.coordinates.builtin_systems.FK4, astropy.coordinates.old_builtin_systems_names.FK4Coordinates, astropy.coordinates.old_builtin_systems_names.FK4NoETermCoordinates, astropy.coordinates.builtin_systems.FK4NoETerms, astropy.coordinates.builtin_systems.FK5, astropy.coordinates.old_builtin_systems_names.FK5Coordinates, astropy.coordinates.transformations.FunctionTransform, astropy.coordinates.builtin_systems.Galactic, astropy.coordinates.old_builtin_systems_names.GalacticCoordinates, astropy.coordinates.old_builtin_systems_names.HorizontalCoordinates, astropy.coordinates.builtin_systems.ICRS, astropy.coordinates.old_builtin_systems_names.ICRSCoordinates, astropy.coordinates.errors.IllegalHourError, astropy.coordinates.errors.IllegalHourWarning, astropy.coordinates.errors.IllegalMinuteError, astropy.coordinates.errors.IllegalMinuteWarning, astropy.coordinates.errors.IllegalSecondError, astropy.coordinates.errors.IllegalSecondWarning, astropy.coordinates.angles.Latitude, astropy.coordinates.angles.Longitude, astropy.coordinates.errors.RangeError, astropy.coordinates.coordsystems.SphericalCoordinatesBase, astropy.coordinates.transformations.StaticMatrixTransform