astropy:docs
  • Index
  • Modules

Navigation

  • next »
  • « previous |
  • Astropy v1.2 »
  • Astronomical Coordinate Systems (astropy.coordinates) »
  • Separations, Catalog Matching, and Related Functionality

Separations, Catalog Matching, and Related Functionality¶

astropy.coordinates contains commonly-used tools for comparing or matching coordinate objects. Of particular importance are those for determining separations between coordinates and those for matching a coordinate (or coordinates) to a catalog. These are mainly implemented as methods on the coordinate objects.

Separations¶

The on-sky separation is easily computed with the astropy.coordinates.BaseCoordinateFrame.separation() or astropy.coordinates.SkyCoord.separation() methods, which computes the great-circle distance (not the small-angle approximation):

>>> import numpy as np
>>> from astropy import units as u
>>> from astropy.coordinates import SkyCoord
>>> c1 = SkyCoord('5h23m34.5s', '-69d45m22s', frame='icrs')
>>> c2 = SkyCoord('0h52m44.8s', '-72d49m43s', frame='fk5')
>>> sep = c1.separation(c2)
>>> sep  
<Angle 20.74611447604398 deg>

The returned object is an Angle instance, so it is straightforward to access the angle in any of several equivalent angular units:

>>> sep.radian  
0.36208800460262575
>>> sep.hour  
1.3830742984029323
>>> sep.arcminute  
1244.7668685626388
>>> sep.arcsecond  
74686.01211375833

Also note that the two input coordinates were not in the same frame - one is automatically converted to match the other, ensuring that even though they are in different frames, the separation is determined consistently.

In addition to the on-sky separation described above, astropy.coordinates.BaseCoordinateFrame.separation_3d() or astropy.coordinates.SkyCoord.separation_3d() methods will determine the 3D distance between two coordinates that have distance defined:

>>> from astropy.coordinates import SkyCoord
>>> c1 = SkyCoord('5h23m34.5s', '-69d45m22s', distance=70*u.kpc, frame='icrs')
>>> c2 = SkyCoord('0h52m44.8s', '-72d49m43s', distance=80*u.kpc, frame='icrs')
>>> sep = c1.separation_3d(c2)
>>> sep  
<Distance 28.743988157814094 kpc>

There is also a spherical_offsets_to() method for computing angular offsets (e.g., small shifts like you might give a telescope operator to move from a bright star to a fainter target.):

>>> from astropy.coordinates import SkyCoord
>>> bright_star = SkyCoord('8h50m59.75s', '+11d39m22.15s', frame='icrs')
>>> faint_galaxy = SkyCoord('8h50m47.92s', '+11d39m32.74s', frame='icrs')
>>> dra, ddec = bright_star.spherical_offsets_to(faint_galaxy)
>>> dra.to(u.arcsec)  
<Angle -173.78873354064126 arcsec>
>>> ddec.to(u.arcsec)  
<Angle 10.605103417374696 arcsec>

“Sky Offset” Frames¶

To extend the concept of spherical offsets, coordinates has a frame class SkyOffsetFrame which creates distinct frames that are centered on a specific point. These are known as “sky offset frames”, as they are a convenient way to create a frame centered on an arbitrary position on the sky, suitable for computing positional offsets (e.g., for astrometry):

>>> from astropy.coordinates import SkyOffsetFrame, ICRS
>>> center = ICRS(10*u.deg, 45*u.deg)
>>> center.transform_to(SkyOffsetFrame(origin=center))  
<SkyOffsetICRS Coordinate (rotation=0.0 deg, origin=<ICRS Coordinate: (ra, dec) in deg
    (10.0, 45.0)>): (lon, lat) in deg
    (0.0, 0.0)>
>>> target = ICRS(11*u.deg, 46*u.deg)
>>> target.transform_to(SkyOffsetFrame(origin=center))  
<SkyOffsetICRS Coordinate (rotation=0.0 deg, origin=<ICRS Coordinate: (ra, dec) in deg
    (10.0, 45.0)>): (lon, lat) in deg
    (0.69474685, 1.00428706)>

Alternatively, the convenience method skyoffset_frame() lets you create an skyoffset frame from an already-existing SkyCoord:

>>> center = SkyCoord(10*u.deg, 45*u.deg)
>>> aframe = center.skyoffset_frame()
>>> target.transform_to(aframe)  
<SkyOffsetICRS Coordinate (rotation=0.0 deg, origin=<ICRS Coordinate: (ra, dec) in deg
    (10.0, 45.0)>): (lon, lat) in deg
    (0.69474685, 1.00428706)>
>>> other = SkyCoord(9*u.deg, 44*u.deg, frame='fk5')
>>> other.transform_to(aframe)  
<SkyCoord (SkyOffsetICRS: rotation=0.0 deg, origin=<ICRS Coordinate: (ra, dec) in deg
    (10.0, 45.0)>): (lon, lat) in deg
    (-0.71943945, -0.99556216)>

Note

While sky offset frames appear to be all the same class, this not the case: the sky offset frame for each different type of frame for origin is actually a distinct class. E.g., SkyOffsetFrame(origin=ICRS(...)) yields an object of class SkyOffsetICRS, not SkyOffsetFrame. While this is not important for most uses of this class, it is important for things like type-checking, because something like SkyOffsetFrame(origin=ICRS(...)).__class__ is SkyOffsetFrame will not be True, as it would be for most classes.

This same frame is also useful as a tool for defining frames that are relative to a specific known object, useful for hierarchical physical systems like galaxy groups. For example, objects around M31 are sometimes shown in a coordinate frame aligned with standard ICRA RA/Dec, but on M31:

>>> m31 = SkyCoord(10.6847083*u.deg, 41.26875*u.deg, frame='icrs')
>>> ngc147 = SkyCoord(8.3005*u.deg, 48.5087389*u.deg, frame='icrs')
>>> ngc147_inm31 = ngc147.transform_to(m31.skyoffset_frame())
>>> xi, eta = ngc147_inm31.lon, ngc147_inm31.lat
>>> xi  
<Longitude -1.5920694752086249 deg>
>>> eta  
<Latitude 7.261837574183891 deg>

Matching Catalogs¶

coordinates supports leverages the coordinate framework to make it straightforward to find the closest coordinates in a catalog to a desired set of other coordinates. For example, assuming ra1/dec1 and ra2/dec2 are numpy arrays loaded from some file:

>>> from astropy.coordinates import SkyCoord
>>> from astropy import units as u
>>> c = SkyCoord(ra=ra1*u.degree, dec=dec1*u.degree)  
>>> catalog = SkyCoord(ra=ra2*u.degree, dec=dec2*u.degree)  
>>> idx, d2d, d3d = c.match_to_catalog_sky(catalog)  

You can also find the nearest 3d matches, different from the on-sky separation shown above only when the coordinates were initialized with a distance:

>>> c = SkyCoord(ra=ra1*u.degree, dec=dec1*u.degree, distance=distance1*u.kpc)  
>>> catalog = SkyCoord(ra=ra2*u.degree, dec=dec2*u.degree, distance=distance2*u.kpc)  
>>> idx, d2d, d3d = c.match_to_catalog_3d(catalog)  

Now idx are indices into catalog that are the closest objects to each of the coordinates in c, d2d are the on-sky distances between them, and d3d are the 3-dimensional distances. Because coordinate objects support indexing, idx enables easy access to the matched set of coordinates in the catalog:

>>> matches = catalog[idx]  
>>> (matches.separation_3d(c) == d3d).all()  
True
>>> dra, ddec = c.spherical_offsets_to(matches)  

This functionality can also be accessed from the match_coordinates_sky() and match_coordinates_3d() functions. These will work on either SkyCoord objects or the lower-level frame classes:

>>> from astropy.coordinates import match_coordinates_sky
>>> idx, d2d, d3d = match_coordinates_sky(c, catalog)  
>>> idx, d2d, d3d = match_coordinates_sky(c.frame, catalog.frame)  

Searching Around Coordinates¶

Closely-related functionality can be used to search for all coordinates within a certain distance (either 3D distance or on-sky) of another set of coordinates. The search_around_* methods (and functions) provide this functionality, with an interface very similar to match_coordinates_*:

>>> idxc, idxcatalog, d2d, d3d = catalog.search_around_sky(c, 1*u.deg)  
>>> np.all(d2d < 1*u.deg)  
True
>>> idxc, idxcatalog, d2d, d3d = catalog.search_around_3d(c, 1*u.kpc)  
>>> np.all(d3d < 1*u.kpc)  
True

The key difference for these methods is that there can be multiple (or no) matches in catalog around any locations in c. Hence, indices into both c and catalog are returned instead of just indices into catalog. These can then be indexed back into the two SkyCoord objects, or, for that matter, any array with the same order:

>>> np.all(c[idxc].separation(catalog[idxcatalog]) == d2d)  
True
>>> np.all(c[idxc].separation_3d(catalog[idxcatalog]) == d3d)  
True
>>> print catalog_objectnames[idxcatalog]  
['NGC 1234' 'NGC 4567' ...]

Note, though, that this dual-indexing means that search_around_* does not work well if one of the coordinates is a scalar, because the returned index would not make sense for a scalar:

>>> scalarc = SkyCoord(1*u.deg, 2*u.deg)  
>>> idxscalarc, idxcatalog, d2d, d3d = catalog.search_around_sky(scalarc, 1*u.deg)  # THIS DOESN'T ACTUALLY WORK  
>>> scalarc[idxscalarc]  
IndexError: 0-d arrays can't be indexed

As a result (and because the search_around_* algorithm is inefficient in the scalar case, anyway), the best approach for this scenario is to instead use the separation* methods:

>>> d2d = scalarc.separation(catalog)  
>>> catalogmsk = d2d < 1*u.deg  
>>> d3d = scalarc.separation_3d(catalog)  
>>> catalog3dmsk = d3d < 1*u.kpc  

The resulting catalogmsk or catalog3dmsk variables are boolean arrays rather than arrays of indices, but in practice they usually can be used in the same way as idxcatalog from the above examples. If you definitely do need indices instead of boolean masks, you can do:

>>> idxcatalog = np.where(catalogmsk)[0]  
>>> idxcatalog3d = np.where(catalog3dmsk)[0]  

Page Contents

  • Separations, Catalog Matching, and Related Functionality
    • Separations
    • “Sky Offset” Frames
    • Matching Catalogs
    • Searching Around Coordinates

Edit This Page on Github   Page Source   Back to Top

© Copyright 2011-2016, The Astropy Developers.
Created using Sphinx 1.3.6.   Last built 20 Jun 2016.