Adding missing utils.py

pygdsm/utils.py

@@ -0,0 +1,33 @@
+import numpy as np
+import healpy as hp
+from astropy.coordinates import SkyCoord
+
+
+def hpix2sky(nside: int, pix_ids: np.ndarray) -> SkyCoord:
+    """ Convert a healpix pixel_id into a SkyCoord
+
+    Args:
+        nside (int): Healpix NSIDE parameter
+        pix_ids (np.array): Array of pixel IDs
+
+    Returns:
+        sc (SkyCoord): Corresponding SkyCoordinates
+    """
+    gl, gb = hp.pix2ang(nside, pix_ids, lonlat=True)
+    sc = SkyCoord(gl, gb, frame='galactic', unit=('deg', 'deg'))
+    return sc
+
+
+def sky2hpix(nside: int, sc: SkyCoord) -> np.ndarray:
+    """ Convert a SkyCoord into a healpix pixel_id
+
+    Args:
+        nside (int): Healpix NSIDE parameter
+        sc (SkyCoord): Astropy sky coordinates array
+
+    Returns:
+        pix (np.array): Array of healpix pixel IDs
+    """
+    gl, gb = sc.galactic.l.to('deg').value, sc.galactic.b.to('deg').value
+    pix = hp.ang2pix(nside, gl, gb, lonlat=True)
+    return pix

tests/test_init.py

@@ -0,0 +1,28 @@
+"""
+test_init.py
+===========
+
+Tests for GSM init commands
+"""
+
+from pygdsm import GlobalSkyModel, GlobalSkyModel08, GlobalSkyModel16, LowFrequencySkyModel, HaslamSkyModel
+from pygdsm import GSMObserver, GSMObserver16, LFSMObserver, HaslamObserver
+from pygdsm import init_gsm, init_observer
+
+
+def test_init():
+    assert isinstance(init_gsm('gsm'), GlobalSkyModel)
+    assert isinstance(init_gsm('gsm08'), GlobalSkyModel)
+    assert isinstance(init_gsm('gsm16'), GlobalSkyModel16)
+    assert isinstance(init_gsm('lfsm'), LowFrequencySkyModel)
+    assert isinstance(init_gsm('haslam'), HaslamSkyModel)
+
+    assert isinstance(init_observer('gsm'), GSMObserver)
+    assert isinstance(init_observer('gsm08'), GSMObserver)
+    assert isinstance(init_observer('gsm16'), GSMObserver16)
+    assert isinstance(init_observer('lfsm'), LFSMObserver)
+    assert isinstance(init_observer('haslam'), HaslamObserver)
+
+
+if __name__ == "__main__":
+    test_init()

tests/test_utils.py

@@ -0,0 +1,16 @@
+import numpy as np
+import healpy as hp
+from pygdsm.utils import hpix2sky, sky2hpix
+
+
+def test_pix2sky():
+    """ Small test routine for converting healpix pixel_id to and from SkyCoords """
+    NSIDE = 32
+    pix = np.arange(hp.nside2npix(NSIDE))
+    sc  = hpix2sky(NSIDE, pix)
+    pix_roundtrip = sky2hpix(NSIDE, sc)
+    assert np.allclose(pix, pix_roundtrip)
+
+
+if __name__ == "__main__":
+    test_pix2sky()