gh-289: notebook for Gaussian covariance

heracles/healpy.py

@@ -130,15 +130,7 @@ def create(
         Create a new HEALPix map.
         """
         m = np.zeros((*dims, hp.nside2npix(self.__nside)), dtype=self.__dtype)
-        update_metadata(
-            m,
-            geometry="healpix",
-            kernel="healpix",
-            nside=self.__nside,
-            lmax=self.__lmax,
-            deconv=self.__deconv,
-            spin=spin,
-        )
+        self.update_metadata(m, spin=spin)
         return m
 
     def map_values(
@@ -202,6 +194,21 @@ def transform(self, data: NDArray[Any]) -> NDArray[Any]:
 
         return alm
 
+    def update_metadata(self, m: NDArray[Any], spin: int) -> NDArray[Any]:
+        """
+        Reads external data into the mapper's format.
+        """
+        update_metadata(
+            m,
+            geometry="healpix",
+            kernel="healpix",
+            nside=self.__nside,
+            lmax=self.__lmax,
+            deconv=self.__deconv,
+            spin=spin,
+        )
+        return m
+
     def resample(self, data: NDArray[Any]) -> NDArray[Any]:
         """
         Change resolution of HEALPix map.

heracles/mapper.py

@@ -67,6 +67,11 @@ def transform(self, data: NDArray[Any]) -> NDArray[Any]:
         The spherical harmonic transform for this mapper.
         """
 
+    def update_metadata(self, m: NDArray[Any], spin: int) -> NDArray[Any]:
+        """
+        Update metadata of external data to this mapper's format.
+        """
+
     def resample(self, data: NDArray[Any]) -> NDArray[Any]:
         """
         Change resolution of data, which must be in this mapper's format.

tests/test_mapping.py

@@ -50,6 +50,37 @@ def test_map_catalogs_match():
     assert set(maps.keys()) == {("b", "x"), ("b", "y"), ("c", "x"), ("c", "y")}
 
 
+def test_add_metadata_to_external_map():
+    import heracles.healpy
+    import numpy as np
+
+    nside = 16
+    lmax = 10
+    mapper = heracles.healpy.HealpixMapper(nside, lmax)
+
+    # spin = 0 case
+    m = np.ones(12 * nside**2)
+    m = mapper.update_metadata(m, spin=0)
+
+    assert m.dtype.metadata is not None
+    assert m.dtype.metadata["spin"] == 0
+    assert m.dtype.metadata["geometry"] == "healpix"
+    assert m.dtype.metadata["kernel"] == "healpix"
+    assert m.dtype.metadata["deconv"] is True
+    assert m.dtype.metadata["nside"] == nside
+    assert m.dtype.metadata["lmax"] == lmax
+
+    m2 = np.ones((2, 12 * nside**2))
+    m2 = mapper.update_metadata(m2, spin=2)
+    assert m2.dtype.metadata is not None
+    assert m2.dtype.metadata["spin"] == 2
+    assert m2.dtype.metadata["geometry"] == "healpix"
+    assert m2.dtype.metadata["kernel"] == "healpix"
+    assert m2.dtype.metadata["deconv"] is True
+    assert m2.dtype.metadata["nside"] == nside
+    assert m2.dtype.metadata["lmax"] == lmax
+
+
 def test_transform(rng):
     from heracles.mapping import transform