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]