Merge 516f08d into 528292b

astropy_healpix/core.py

@@ -18,6 +18,10 @@
     'pixel_resolution_to_nside',
     'nside_to_npix',
     'npix_to_nside',
+    'level_to_nside',
+    'nside_to_level',
+    'levelipix_to_uniq',
+    'uniq_to_levelipix',
     'lonlat_to_healpix',
     'healpix_to_lonlat',
     'bilinear_interpolation_weights',
@@ -61,6 +65,9 @@ def _validate_nside(nside):
     if not np.all(2 ** log_2_nside == nside):
         raise ValueError('nside must be a power of two')
 
+def _validate_npix(level, ipix):
+    if not np.all(ipix < (3 << 2*(level + 1))):
+        raise ValueError('ipix for a specific level must be inferior to npix')
 
 def level_to_nside(level):
     """
@@ -77,10 +84,92 @@ def level_to_nside(level):
     nside : int
         The number of pixels on the side of one of the 12 'top-level' HEALPix tiles.
     """
+    level = np.asarray(level, dtype=np.int64)
+
     _validate_level(level)
     return 2 ** level
 
 
+def nside_to_level(nside):
+    """
+    Find the HEALPix level for a given nside. (this is given by log2(nside))
+
+    This function is the inverse of `level_to_nside`.
+
+    Parameters
+    ----------
+    nside : int
+        The number of pixels on the side of one of the 12 'top-level' HEALPix tiles.
+        Must be a power of two.
+
+    Returns
+    -------
+    level : int
+        The level of the HEALPix cells
+    """
+    nside = np.asarray(nside, dtype=np.int64)
+
+    _validate_nside(nside)
+    return np.log2(nside)
+
+
+def uniq_to_levelipix(uniq):
+    """
+	Convert a HEALPix cell uniq number to its (level, ipix) equivalent
+
+    A uniq number is a 64 bit integer equals to : ipix + 4*(4**level). Please read
+    this `paper <http://ivoa.net/documents/MOC/20140602/REC-MOC-1.0-20140602.pdf>`_
+    for more details about uniq numbers.
+
+    Parameters
+    ----------
+    uniq : int
+        The uniq number of a HEALPix cell.
+
+    Returns
+    -------
+    level, ipix: int, int
+        The level and index of the HEALPix cell computed from ``uniq``.
+    """
+    uniq = np.asarray(uniq, dtype=np.int64)
+
+    level = ((np.log2(uniq//4)) // 2)
+    level = level.astype(np.int64)
+    _validate_level(level)
+
+    ipix = uniq - (1 << 2*(level + 1))
+    _validate_npix(level, ipix)
+
+    return level, ipix
+
+
+def levelipix_to_uniq(level, ipix):
+    """
+    Convert a level and HEALPix index into a uniq number representing the cell
+
+    This function is the inverse of `uniq_to_levelipix`.
+
+    Parameters
+    ----------
+    level : int
+        The level of the HEALPix cell
+    ipix : int
+        The index of the HEALPix cell
+
+    Returns
+    -------
+    uniq : int
+        The uniq number representing the HEALPix cell.
+    """
+    level = np.asarray(level, dtype=np.int64)
+    ipix = np.asarray(ipix, dtype=np.int64)
+
+    _validate_level(level)
+    _validate_npix(level, ipix)
+
+    return ipix + (1 << 2*(level + 1))
+
+
 def nside_to_pixel_area(nside):
     """
     Find the area of HEALPix pixels given the pixel dimensions of one of

astropy_healpix/tests/test_core.py

@@ -17,7 +17,9 @@
                     nside_to_npix, npix_to_nside, healpix_to_lonlat,
                     lonlat_to_healpix, interpolate_bilinear_lonlat,
                     neighbours, healpix_cone_search, boundaries_lonlat,
-                    level_to_nside, nested_to_ring, ring_to_nested,
+                    level_to_nside, nside_to_level,
+                    nested_to_ring, ring_to_nested,
+                    levelipix_to_uniq, uniq_to_levelipix,
                     bilinear_interpolation_weights)
 
 
@@ -28,6 +30,35 @@ def test_level_to_nside():
     assert exc.value.args[0] == 'level must be positive'
 
 
+def test_nside_to_level():
+    assert nside_to_level(1024) == 10
+    with pytest.raises(ValueError) as exc:
+        nside_to_level(511)
+    assert exc.value.args[0] == 'nside must be a power of two'
+
+
+@pytest.mark.parametrize("level, ipix, expected_nuniq", [
+    (0, 11, 11+4),
+    (15, 62540, 62540 + 4*4**15),
+])
+def test_levelipix_to_uniq(level, ipix, expected_nuniq):
+    assert ipix + 4*4**level == levelipix_to_uniq(level, ipix)
+
+
+@pytest.mark.parametrize("level", [
+    0, 5, 10, 15, 20, 22, 25, 26, 27, 28, 29
+])
+def test_uniq_to_levelipix(level):
+    # Generate 10 HEALPix indexes for each level
+    size = 10
+    npix = 3 << 2*(level + 1)
+    ipix_arr = np.random.randint(npix, size=size)
+    level_arr = np.ones(size) * level
+
+    level_res_arr, ipix_res_arr = uniq_to_levelipix(levelipix_to_uniq(level_arr, ipix_arr))
+    assert np.all(level_res_arr == level_arr) & np.all(ipix_res_arr == ipix_arr)
+
+
 def test_nside_to_pixel_area():
     resolution = nside_to_pixel_area(256)
     assert_allclose(resolution.value, 1.5978966540475428e-05)