lsst.afw.geom._hpxUtils Namespace Reference

Classes
class	_ZOrderCurve2DInt

Functions
def	makeHpxWcs (hips_order, pixel, shift_order=9)

Function Documentation

makeHpxWcs()

def lsst.afw.geom._hpxUtils.makeHpxWcs	(		hips_order,
			pixel,
			shift_order = 9
	)

Make a SkyWcs object with HEALPix grid projection (HPX).

The SkyWcs generated by this function is suitable to be used with a
Hierarchical Progressive Survey (HiPS) FITS image as described in
https://www.ivoa.net/documents/HiPS/20170519/REC-HIPS-1.0-20170519.pdf

A HiPS image covers one HEALPix cell, with the HEALPix nside equal to
2**hips_order. Each cell is 'shift_order' orders deeper than the HEALPix
cell, with 2**shift_order x 2**shift_order sub-pixels on a side, which
defines the target resolution of the HiPS image. The IVOA recommends
shift_order=9, for 2**9=512 pixels on a side.  See Notes below to
convert from hips_order to image resolution.

Parameters
----------
hips_order : `int`
    HiPS order, such that HEALPix nside=2**hips_order.
    Must be a positive integer.
pixel : `int`
    Pixel number in the nest ordering scheme.
shift_order : `int`, optional
    Shift order for subpixels, such that there are 2**shift_order
    sub-pixels on a side of the HiPS cell.
    Must be a positive integer.

Returns
-------
wcs : `lsst.geom.SkyWcs`

Raises
------
`ValueError`: Raise if hips_order is <=0, or if shift_order is <=0, or
    if pixel number is out of range for the given hips_order
    (0 <= pixel < 12*nside*nside).

Notes
-----
Table 5 from
https://www.ivoa.net/documents/HiPS/20170519/REC-HIPS-1.0-20170519.pdf
shows the relationship between hips_order, number of tiles (full
sky coverage), cell size, and sub-pixel size/image resolution (with
the default shift_order=9):

+------------+-----------------+--------------+------------------+
| hips_order | Number of Tiles | Cell Size    | Image Resolution |
+============+=================+==============+==================+
| 0          | 12              | 58.63 deg    | 6.871 arcmin     |
| 1          | 48              | 29.32 deg    | 3.435 arcmin     |
| 2          | 192             | 14.66 deg    | 1.718 arcmin     |
| 3          | 768             | 7.329 deg    | 51.53 arcsec     |
| 4          | 3072            | 3.665 deg    | 25.77 arcsec     |
| 5          | 12288           | 1.832 deg    | 12.88 arcsec     |
| 6          | 49152           | 54.97 arcmin | 6.442 arcsec     |
| 7          | 196608          | 27.48 arcmin | 3.221 arcsec     |
| 8          | 786432          | 13.74 arcmin | 1.61 arcsec      |
| 9          | 3145728         | 6.871 arcmin | 805.2mas         |
| 10         | 12582912        | 3.435 arcmin | 402.6mas         |
| 11         | 50331648        | 1.718 arcmin | 201.3mas         |
| 12         | 201326592       | 51.53 arcsec | 100.6mas         |
| 13         | 805306368       | 25.77 arcsec | 50.32mas         |
+------------+-----------------+--------------+------------------+

Definition at line 31 of file _hpxUtils.py.

def makeHpxWcs(hips_order, pixel, shift_order=9):
    """
    Make a SkyWcs object with HEALPix grid projection (HPX).
 
    The SkyWcs generated by this function is suitable to be used with a
    Hierarchical Progressive Survey (HiPS) FITS image as described in
    https://www.ivoa.net/documents/HiPS/20170519/REC-HIPS-1.0-20170519.pdf
 
    A HiPS image covers one HEALPix cell, with the HEALPix nside equal to
    2**hips_order. Each cell is 'shift_order' orders deeper than the HEALPix
    cell, with 2**shift_order x 2**shift_order sub-pixels on a side, which
    defines the target resolution of the HiPS image. The IVOA recommends
    shift_order=9, for 2**9=512 pixels on a side.  See Notes below to
    convert from hips_order to image resolution.
 
    Parameters
    ----------
    hips_order : `int`
        HiPS order, such that HEALPix nside=2**hips_order.
        Must be a positive integer.
    pixel : `int`
        Pixel number in the nest ordering scheme.
    shift_order : `int`, optional
        Shift order for subpixels, such that there are 2**shift_order
        sub-pixels on a side of the HiPS cell.
        Must be a positive integer.
 
    Returns
    -------
    wcs : `lsst.geom.SkyWcs`
 
    Raises
    ------
    `ValueError`: Raise if hips_order is <=0, or if shift_order is <=0, or
        if pixel number is out of range for the given hips_order
        (0 <= pixel < 12*nside*nside).
 
    Notes
    -----
    Table 5 from
    https://www.ivoa.net/documents/HiPS/20170519/REC-HIPS-1.0-20170519.pdf
    shows the relationship between hips_order, number of tiles (full
    sky coverage), cell size, and sub-pixel size/image resolution (with
    the default shift_order=9):
 
    +------------+-----------------+--------------+------------------+
    | hips_order | Number of Tiles | Cell Size    | Image Resolution |
    +============+=================+==============+==================+
    | 0          | 12              | 58.63 deg    | 6.871 arcmin     |
    | 1          | 48              | 29.32 deg    | 3.435 arcmin     |
    | 2          | 192             | 14.66 deg    | 1.718 arcmin     |
    | 3          | 768             | 7.329 deg    | 51.53 arcsec     |
    | 4          | 3072            | 3.665 deg    | 25.77 arcsec     |
    | 5          | 12288           | 1.832 deg    | 12.88 arcsec     |
    | 6          | 49152           | 54.97 arcmin | 6.442 arcsec     |
    | 7          | 196608          | 27.48 arcmin | 3.221 arcsec     |
    | 8          | 786432          | 13.74 arcmin | 1.61 arcsec      |
    | 9          | 3145728         | 6.871 arcmin | 805.2mas         |
    | 10         | 12582912        | 3.435 arcmin | 402.6mas         |
    | 11         | 50331648        | 1.718 arcmin | 201.3mas         |
    | 12         | 201326592       | 51.53 arcsec | 100.6mas         |
    | 13         | 805306368       | 25.77 arcsec | 50.32mas         |
    +------------+-----------------+--------------+------------------+
    """
    if shift_order <= 0:
        raise ValueError(f"shift_order {shift_order} must be positive.")
    hips_tilepix = 2**shift_order
 
    if hips_order <= 0:
        raise ValueError(f"order {hips_order} must be positive.")
    nside_cell = 2**hips_order
 
    if pixel < 0 or pixel >= 12*nside_cell*nside_cell:
        raise ValueError(f"pixel value {pixel} out of range.")
 
    # The HEALPix grid projection (HPX) is defined in the FITS standard
    # https://fits.gsfc.nasa.gov/standard40/fits_standard40aa-le.pdf
    # from Calabretta & Roukema (2007)
    # https://ui.adsabs.harvard.edu/abs/2007MNRAS.381..865C/abstract
    # which defines the standard H = 4, K = 3 pixelization parameters
    # encoded in PV2_1 = H, PV2_2 = K.
    # The CRVAL1, CRVAL2 values should always be 0, 0 according to
    # the FITS standard.
    # The CD matrix is defined in wcslib HPXcvt.c.
    # The Calabretta & Roukema (2007) paper and wcslib HPXcvt.c only
    # define full-sky HPX projections.  For single pixels we
    # use the code from AladinSrc.jar Tile2HPX.java to compute
    # CRPIX1, CRPIX2.
 
    # The nside of the sub-pixels is the product of the tile nside
    # and the number of sub-pixels on a side.
    nside_pix = nside_cell*hips_tilepix
    # All tiles are rotated 45 degrees.
    cos45 = np.sqrt(2.0)/2.0
    # This defines the pixel scale.
    scale = 90.0/nside_pix/np.sqrt(2.0)
    cos45_scale = cos45*scale
    # The projected center of the pixel used for the HPX header is
    # a non-trivial computation, and typically is outside of the
    # tile pixel itself.  Therefore, these values are not the same
    # as the values computed from healpy.pix2ang().
    cent_ra_proj, cent_dec_proj = _hpx_projected_center(hips_order, pixel)
 
    md = PropertySet()
    md['CD1_1'] = -cos45_scale
    md['CD1_2'] = -cos45_scale
    md['CD2_1'] = cos45_scale
    md['CD2_2'] = -cos45_scale
    md['CTYPE1'] = 'RA---HPX'
    md['CTYPE2'] = 'DEC--HPX'
    md['CRVAL1'] = 0.0
    md['CRVAL2'] = 0.0
    md['PV2_1'] = 4
    md['PV2_2'] = 3
    md['CRPIX1'] = ((hips_tilepix + 1)/2.0) - 0.5*(-cent_ra_proj/cos45_scale + cent_dec_proj/cos45_scale)
    md['CRPIX2'] = ((hips_tilepix + 1)/2.0) - 0.5*(-cent_ra_proj/cos45_scale - cent_dec_proj/cos45_scale)
 
    return makeSkyWcs(md)