remove twopoint dependency on mapper

heracles/maps/_healpix.py

@@ -186,6 +186,7 @@ def create(
             kernel="healpix",
             nside=self.__nside,
             lmax=self.__lmax,
+            deconv=self.__deconv,
             spin=spin,
         )
         return m
@@ -229,7 +230,7 @@ def transform(
         maps = np.asanyarray(maps)
         md: Mapping[str, Any] = maps.dtype.metadata or {}
         spin = md.get("spin", 0)
-        pw: NDArray[np.floating] | None = None
+        pw: NDArray[Any] | None = None
 
         if spin == 0:
             pol = False

heracles/twopoint.py

@@ -30,8 +30,7 @@
 import healpy as hp
 import numpy as np
 
-from .core import TocDict, items_with_suffix, toc_match, update_metadata
-from .maps import mapper_from_dict
+from .core import TocDict, toc_match, update_metadata
 
 if TYPE_CHECKING:
     from collections.abc import Mapping, MutableMapping
@@ -56,6 +55,9 @@ def _debias_cl(
 ) -> NDArray[Any]:
     """
     Remove additive bias from angular power spectrum.
+
+    This function special-cases the bias from HEALPix maps.
+
     """
 
     if md is None:
@@ -65,39 +67,46 @@ def _debias_cl(
         cl = cl.copy()
         update_metadata(cl, **md)
 
-    lmax = len(cl) - 1
+    # use explicit bias, if given, or bias value from metadata
+    if bias is None:
+        bias = md.get("bias")
+        # return early if there is no bias to be subtracted
+        if bias is None:
+            return cl
 
     # spins of the spectrum
     spin1, spin2 = md.get("spin_1", 0), md.get("spin_2", 0)
 
-    # minimum angular mode for bias correction
+    # minimum and maximum angular mode for bias correction
     lmin = max(abs(spin1), abs(spin2))
+    lmax = len(cl) - 1
 
-    # use explicit bias, if given, or bias value from metadata
-    b = np.zeros(lmax + 1)
-    b[lmin:] = md.get("bias", 0.0) if bias is None else bias
-
-    # apply biasing kernel from mappers
-    try:
-        mapper = mapper_from_dict(items_with_suffix(md, "_1"))
-    except ValueError:
-        pass
-    else:
-        if (bl := mapper.bl(lmax=lmax, spin=spin1)) is not None:
-            b *= bl
-    try:
-        mapper = mapper_from_dict(items_with_suffix(md, "_2"))
-    except ValueError:
-        pass
-    else:
-        if (bl := mapper.bl(lmax=lmax, spin=spin2)) is not None:
-            b *= bl
+    # this will be subtracted from the cl
+    # modes up to lmin are ignored
+    bl = np.full(lmax + 1, bias)
+    bl[:lmin] = 0.0
+
+    # handle HEALPix pseudo-convolution
+    for i, s in (1, spin1), (2, spin2):
+        if md.get(f"kernel_{i}") == "healpix":
+            nside: int | None = md.get(f"nside_{i}")
+            deconv: bool = md.get(f"deconv_{i}", True)
+            if nside is not None and deconv:
+                pw: NDArray[Any] | None
+                if s == 0:
+                    pw = hp.pixwin(nside, lmax=lmax, pol=False)
+                elif s == 2:
+                    pw = hp.pixwin(nside, lmax=lmax, pol=True)[1]
+                else:
+                    pw = None
+                if pw is not None:
+                    bl[lmin:] /= pw[lmin:]
 
     # remove bias
     if cl.dtype.names is None:
-        cl -= b
+        cl -= bl
     else:
-        cl["CL"] -= b
+        cl["CL"] -= bl
 
     return cl
 
@@ -221,23 +230,13 @@ def angular_power_spectra(
 def debias_cls(cls, bias=None, *, inplace=False):
     """remove bias from cls"""
 
-    logger.info("debiasing %d cl(s)%s", len(cls), " in place" if inplace else "")
-    t = time.monotonic()
-
     # the output toc dict
     out = cls if inplace else TocDict()
 
     # subtract bias of each cl in turn
     for key in cls:
-        logger.info("debiasing %s x %s cl for bins %s, %s", *key)
         out[key] = _debias_cl(cls[key], bias and bias.get(key), inplace=inplace)
 
-    logger.info(
-        "debiased %d cl(s) in %s",
-        len(out),
-        timedelta(seconds=(time.monotonic() - t)),
-    )
-
     # return the toc dict of debiased cls
     return out
 

tests/test_maps.py

@@ -34,6 +34,7 @@ def test_healpix_maps(rng):
         "kernel": "healpix",
         "nside": nside,
         "lmax": lmax,
+        "deconv": deconv,
         "spin": -3,
     }
 

tests/test_twopoint.py

@@ -116,16 +116,72 @@ def test_debias_cls():
     from heracles.twopoint import debias_cls
 
     cls = {
-        ("PP", 0, 0): np.zeros(100),
+        0: np.zeros(100),
+        2: np.zeros(100, dtype=np.dtype(float, metadata={"bias": 4.56, "spin_2": 2})),
     }
 
     nbs = {
-        ("PP", 0, 0): 1.23,
+        0: 1.23,
     }
 
     debias_cls(cls, nbs, inplace=True)
 
-    assert np.all(cls["PP", 0, 0] == -1.23)
+    assert np.all(cls[0] == -1.23)
+
+    assert np.all(cls[2][:2] == 0.0)
+    assert np.all(cls[2][2:] == -4.56)
+
+
+def test_debias_cls_healpix():
+    import healpy as hp
+
+    from heracles.twopoint import debias_cls
+
+    pw0, pw2 = hp.pixwin(64, lmax=99, pol=True)
+
+    md1 = {
+        "kernel_1": "healpix",
+        "nside_1": 64,
+        "kernel_2": "healpix",
+        "nside_2": 64,
+        "spin_2": 2,
+    }
+    md2 = {
+        "kernel_1": "healpix",
+        "nside_1": 64,
+        "spin_2": 2,
+    }
+    md3 = {
+        "kernel_1": "healpix",
+        "nside_1": 64,
+        "kernel_2": "healpix",
+        "nside_2": 64,
+        "spin_2": 2,
+        "deconv_2": False,
+    }
+
+    cls = {
+        1: np.zeros(100, dtype=np.dtype(float, metadata=md1)),
+        2: np.zeros(100, dtype=np.dtype(float, metadata=md2)),
+        3: np.zeros(100, dtype=np.dtype(float, metadata=md3)),
+    }
+
+    nbs = {
+        1: 1.23,
+        2: 4.56,
+        3: 7.89,
+    }
+
+    debias_cls(cls, nbs, inplace=True)
+
+    assert np.all(cls[1][:2] == 0.0)
+    assert np.all(cls[1][2:] == -1.23 / pw0[2:] / pw2[2:])
+
+    assert np.all(cls[2][:2] == 0.0)
+    assert np.all(cls[2][2:] == -4.56 / pw0[2:])
+
+    assert np.all(cls[3][:2] == 0.0)
+    assert np.all(cls[3][2:] == -7.89 / pw0[2:])
 
 
 @patch("convolvecl.mixmat_eb")