Implement range to moc function

healpix_alchemy/tests/test_types.py

@@ -0,0 +1,32 @@
+import numpy as np
+from mocpy import MOC
+import itertools
+from astropy_healpix import HEALPix
+from astropy.coordinates import Angle, Longitude, Latitude
+import astropy.units as u
+
+from .. import types
+
+
+def to_ranges(iterable):
+    iterable = sorted(set(iterable))
+    for key, group in itertools.groupby(enumerate(iterable),
+                                        lambda t: t[1] - t[0]):
+        group = list(group)
+        yield group[0][1], group[-1][1]+1
+
+
+def test_to_moc():
+    depth = 10
+    rng = np.random.default_rng()
+    lon = Longitude(360 * rng.random() * u.deg)
+    lat = Latitude((180 * rng.random() - 90) * u.deg)
+    radius = Angle(rng.normal(loc=10, scale=2.5) * u.deg)
+    hp = HEALPix(nside=2**depth, order='nested')
+    healpix_list = hp.cone_search_lonlat(lon, lat, radius)
+    depth_list = depth * np.ones(np.shape(healpix_list))
+    moc = MOC.from_healpix_cells(healpix_list, depth_list)
+    nested_hpx_ranges = [item for item in to_ranges(healpix_list)]
+
+    assert types.Point.to_moc(rangeSet=nested_hpx_ranges,
+                              nside=2**depth) == moc

healpix_alchemy/types.py

@@ -28,6 +28,13 @@ def process_bind_param(self, value, dialect):
             value = int(value)
         return value
 
+    @classmethod
+    def to_moc(cls, rangeSet, nside):
+        depth = int(np.log2(nside))
+        MOCstr = f'{depth}/' + ' '.join(map(lambda x: f'{x[0]}-{x[1]-1}',
+                                            rangeSet))
+        return MOC.from_string(MOCstr)
+
 
 class Tile(sa.TypeDecorator):