Coordinate conversions

Converting between pixel indices and spherical coordinates

As described in Getting started, coordinates in a HEALPix pixellization can follow either the ‘ring’ or ‘nested’ convention. Let’s start by setting up an example pixellization:

>>> from astropy_healpix import HEALPix
>>> hp = HEALPix(nside=16, order='nested')

The healpix_to_lonlat() method can be used to convert HEALPix indices to Longitude and Latitude objects:

>>> lon, lat = hp.healpix_to_lonlat([1, 442, 2200])
>>> lon  
<Longitude [ 0.83448555, 1.63624617, 0.4712389 ] rad>
>>> lat  
<Latitude [ 0.08343009, 0.94842784,-0.78529135] rad>

The Longitude and Latitude objects are fully-fledged Quantity objects and also include shortcuts to get the values in various units:

>>> lon.hourangle  
array([ 3.1875,  6.25  ,  1.8   ])
>>> lat.degree  
array([  4.78019185,  54.3409123 , -44.99388015])

Conversely, given longitudes and latitudes as Quantity objects, it is possible to recover HEALPix pixel indices:

>>> from astropy import units as u
>>> print(hp.lonlat_to_healpix([1, 3, 4] * u.deg, [5, 6, 9] * u.deg))
[1217 1217 1222]

In these examples, what is being converted is the position of the center of each pixel. In fact, the lonlat_to_healpix() method can also take or give the fractional position inside each HEALPix pixel, e.g.:

>>> index, dx, dy = hp.lonlat_to_healpix([1, 3, 4] * u.deg, [5, 6, 9] * u.deg,
...                                      return_offsets=True)
>>> print(index)
[1217 1217 1222]
>>> dx  
array([ 0.22364669,  0.78767489,  0.58832469])
>>> dy  
array([ 0.86809114,  0.72100823,  0.16610247])

and the healpix_to_lonlat() method can take offset positions - for example we can use this to find the position of the corners of a given pixel:

>>> dx = [0., 1., 1., 0.]
>>> dy = [0., 0., 1., 1.]
>>> lon, lat = hp.healpix_to_lonlat([133, 133, 133, 133], dx=dx, dy=dy)
>>> lon  
<Longitude [ 0.53996124, 0.58904862, 0.53996124, 0.49087385] rad>
>>> lat  
<Latitude [ 0.47611906, 0.52359878, 0.57241857, 0.52359878] rad>

Celestial coordinates

For cases where the HEALPix pixellization is of the celestial sphere, a frame argument can be passed to HEALPix. This argument should specify the celestial frame (using an astropy.coordinates frame) in which the HEALPix pixellization is defined:

>>> from astropy_healpix import HEALPix
>>> from astropy.coordinates import Galactic
>>> hp = HEALPix(nside=16, order='nested', frame=Galactic())

Each method defined in HEALPix and ending in lonlat has an equivalent method ending in skycoord which can be used if the frame is set. For example, to convert from HEALPix indices to celestial coordinates, you can use the healpix_to_skycoord() method:

>>> hp.healpix_to_skycoord([144, 231])  
<SkyCoord (Galactic): (l, b) in deg
    [( 33.75      ,  32.7971683 ), ( 32.14285714,  69.42254649)]>

and to convert from celestial coordinates to HEALPix indices you can use the skycoord_to_healpix() method, e.g:

>>> from astropy.coordinates import SkyCoord
>>> coord = SkyCoord('00h42m44.3503s +41d16m08.634s')
>>> hp.skycoord_to_healpix(coord)
2537

Converting between ring and nested conventions

The HEALPix class has methods that can be used to convert HEALPix pixel indices between the ring and nested convention. These are nested_to_ring():

>>> print(hp.nested_to_ring([30]))
[873]

and ring_to_nested():

>>> print(hp.ring_to_nested([1, 2, 3]))
[ 511  767 1023]