Merge branch 'main' into ana_cov_comp

pspipe_utils/io.py

@@ -132,3 +132,54 @@ def port2sacc(
 
     log.info(f"Writing {sacc_file_name} \n")
     s.save_fits(sacc_file_name, overwrite=True)
+
+
+def extract_sacc_spectra(likelihood_name, input_file, cov_Bbl_file=None):
+    """
+    This function extracts spectra from a sacc file through an "mflike"-like likelihood.
+    It returns spectra and covariance block as python dictionnary the same way
+    the likelihood reads and parses the sacc file content.
+
+    Parameters
+    ----------
+    likelihood_name: str
+      the likelihood name. Must equivalent to the name set in the cobaya yaml file
+      i.e. "mflike.MFLike" for default SO likelihood or "act_dr6_mflike.ACTDR6MFLike"
+      for ACT DR6 likelihood
+    input_file: path
+      the path to the input sacc file
+    cov_Bbl_file: path
+      the path to the covariance file if not inside the input file (default: None).
+      The dirname **must** be the same as the input file.
+    """
+    if cov_Bbl_file and os.path.dirname(cov_Bbl_file) != os.path.dirname(input_file):
+        raise ValueError(
+            "The directory path of the covariance file is different from the input file!"
+        )
+
+    likelihood_module, likelihood_class = likelihood_name.rsplit(".", 1)
+    likelihood_module = importlib.import_module(likelihood_module)
+    likelihood_class = getattr(likelihood_module, likelihood_class)
+
+    # Do not check installation of likelihood
+    likelihood_class.install_options = None
+
+    my_like = likelihood_class(
+        {
+            "path": os.path.realpath(input_file).replace(input_file, ""),
+            "data_folder": os.path.dirname(input_file),
+            "input_file": os.path.basename(input_file),
+            "cov_Bbl_file": os.path.basename(cov_Bbl_file) if cov_Bbl_file else None,
+        }
+    )
+
+    spectra = {}
+    for data in my_like.spec_meta:
+        lb, db = data.get("leff"), data.get("cl_data")
+        cross = (data.get("t1"), data.get("t2"))
+        mode = data.get("pol") if not data.get("hasYX_xsp") else "et"
+        ids = data.get("ids")
+        cov = my_like.cov[np.ix_(ids, ids)]
+        spectra.setdefault((mode, *cross), []).append(dict(lb=lb, db=db, cov=cov))
+
+    return spectra

pspipe_utils/leakage.py

@@ -97,7 +97,7 @@ def apply_leakage_model_to_alm(alms, gamma_TE, gamma_TB):
     alms[2] = alms[2] + curvedsky.almxfl(alms[0], gamma_TB)
     return alms
 
-def read_leakage_model(leakage_file_dir, file_name, lmax, lmin=0, include_error_modes=True):
+def read_leakage_model_old(leakage_file_dir, file_name, lmax, lmin=0, include_error_modes=True):
     """
     This routine serves to read the leakage model in the ACT format, both for
     the mean value of the leakage and its covariance.
@@ -133,6 +133,49 @@ def read_leakage_model(leakage_file_dir, file_name, lmax, lmin=0, include_error_
 
     return l, gamma_TE, error_modes_gTE,  gamma_TB, error_modes_gTB
 
+
+def read_leakage_model(leakage_file_dir,
+                       file_name_TE,
+                       file_name_TB,
+                       lmax,
+                       lmin=0,
+                       pol_eff=1):
+    """
+    This routine serves to read the leakage model in the ACT format, both for
+    the mean value of the leakage and its covariance.
+    not that the error modes file is expected to be of the form
+    l, err_modegE1, err_modegE2, err_modegE3, err_modegB1, err_modegE2, err_modegE3
+
+    Parameters
+    ----------
+    leakage_file_dir : str
+        location of the files describing the leakage
+    file_name_TE : str
+        name of the  file  that contain gamma_TE and error_modes (e.g pa4_f150_gamma_t2e.txt)
+    file_name_TB : str
+        name of the  file  that contain gamma_TB and error_modes (e.g pa4_f150_gamma_t2b.txt)
+    lmin : integer
+        minimum multipole to consider
+    lmax : integer
+        maximum multipole to consider
+    pol_eff: float
+        the polarisation efficiency of the data
+    """
+
+    def extract_beam_leakage_and_error_modes(file_name):
+        data = np.loadtxt(f"{leakage_file_dir}/{file_name}")
+        l, gamma = data[lmin: lmax, 0], data[lmin: lmax, 1]
+        error_modes = data[lmin: lmax, 2:]
+        return l, gamma, error_modes
+
+    l, gamma_TE, error_modes_gTE = extract_beam_leakage_and_error_modes(file_name_TE)
+    l, gamma_TB, error_modes_gTB = extract_beam_leakage_and_error_modes(file_name_TB)
+
+    return l, gamma_TE/pol_eff, error_modes_gTE/pol_eff,  gamma_TB/pol_eff, error_modes_gTB/pol_eff
+
+
+
+
 def error_modes_to_cov(error_modes):
     """
     Use the beam leakage error modes to reconstruct the beam leakage covarance matrix

setup.cfg

@@ -1,5 +1,5 @@
 [bumpversion]
-current_version = 0.1.5
+current_version = 0.1.6
 commit = True
 tag = True