diff --git a/paramfiles/paramfile_SAT.yaml b/paramfiles/paramfile_SAT.yaml index b09fe79..7c800f9 100644 --- a/paramfiles/paramfile_SAT.yaml +++ b/paramfiles/paramfile_SAT.yaml @@ -77,8 +77,8 @@ general_pars: lmax: 600 deltal: 10 binning_file: "binning.npz" - pure_B: True - tf_est_pure_B: True + pure_B: False # Compute data and covariance sims with purified B-modes + tf_est_pure_B: True # Compute transfer function with purified B-modes ## Simulation related sim_pars: @@ -99,9 +99,9 @@ sim_pars: ## Filtering transfer function related parameters tf_settings: - filtering_type: "m_filterer" # "m_filterer" or "toast" - m_cut: 30 # necessary only if `filtering_type` is set to "m_filterer" - toast: # necessary only if `filtering_type` is set to "toast" + filtering_type: "m_filterer" # "m_filterer" or "toast" + m_cut: 30 # necessary only if `filtering_type` is set to "m_filterer" + toast: # necessary only if `filtering_type` is set to "toast" schedule: "../outputs/schedules/schedule_sat.txt" instrument: "SAT1" band: "SAT_f090" diff --git a/pipeline/mcmer.py b/pipeline/mcmer.py index f9e2478..9e125ee 100644 --- a/pipeline/mcmer.py +++ b/pipeline/mcmer.py @@ -29,7 +29,9 @@ def mcmer(args): # Create dummy NaMaster fields, spin-2 is optionally purified. field_spin0 = nmt.NmtField(mask, None, spin=0) field_spin2 = nmt.NmtField(mask, None, spin=2) - if meta.pure_B: + # Whether we want to B-purify data or transfer function estimations, + # in either case we will need a purified mode coupling matrix. + if meta.pure_B or meta.tf_est_pure_B: field_spin2_pure = nmt.NmtField(mask, None, spin=2, purify_b=meta.pure_B) @@ -45,7 +47,7 @@ def mcmer(args): print("Computing MCM") w = nmt.NmtWorkspace() w.compute_coupling_matrix(field_spin0, field_spin2, nmt_bins, is_teb=True) - if meta.pure_B: + if meta.pure_B or meta.tf_est_pure_B: w_pure = nmt.NmtWorkspace() w_pure.compute_coupling_matrix(field_spin0, field_spin2_pure, nmt_bins, is_teb=True) @@ -53,7 +55,7 @@ def mcmer(args): mcm = np.transpose(w.get_coupling_matrix().reshape([nl, nspec, nl, nspec]), axes=[1, 0, 3, 2]) mcm_binned = np.einsum('ij,kjlm->kilm', binner, mcm) - if meta.pure_B: + if meta.pure_B or meta.tf_est_pure_B: mcm_pure = np.transpose( w_pure.get_coupling_matrix().reshape([nl, nspec, nl, nspec]), axes=[1, 0, 3, 2] @@ -87,8 +89,9 @@ def mcmer(args): spin2xspin2_binned=mcm_binned[3:, :, 3:, :] ) - # If we B-purify, store also the un-beamed purified MCM - if meta.pure_B: + # If we B-purify for transfer function estimation, store the + # un-beamed purified mode coupling matrix + if meta.tf_est_pure_B: fname = f"{coupling_dir}/mcm_pure.npz" np.savez( fname, diff --git a/pipeline/pcler.py b/pipeline/pcler.py index 5161752..7642baf 100644 --- a/pipeline/pcler.py +++ b/pipeline/pcler.py @@ -375,7 +375,8 @@ def pcler(args): if args.tf_val: inv_couplings = {} for filter_flag in ["filtered", "unfiltered"]: - couplings = np.load(f"{meta.coupling_directory}/couplings_{filter_flag}_pure.npz") # noqa + pure_str = "_pure" if meta.tf_est_pure_B else "" + couplings = np.load(f"{meta.coupling_directory}/couplings_{filter_flag}{pure_str}.npz") # noqa inv_couplings[filter_flag] = { k1: couplings[f"inv_coupling_{k2}"].reshape([ncl*n_bins, ncl*n_bins]) diff --git a/pipeline/transfer.py b/pipeline/transfer.py index c0df2e6..3f5c29c 100644 --- a/pipeline/transfer.py +++ b/pipeline/transfer.py @@ -61,11 +61,13 @@ def transfer(args): spin_pair: mcm[f"{spin_pair}_binned"] for spin_pair in spin_pairs } - # If we B-purify, load also the un-beamed purified mode coupling matrix - mcm_pure = np.load(f"{coupling_dir}/mcm_pure.npz") - mcms_dict_nobeam_pure = { - spin_pair: mcm_pure[f"{spin_pair}_binned"] for spin_pair in spin_pairs - } + # If we B-purify for transfer function estimation, load the un-beamed + # purified mode coupling matrix + if meta.tf_est_pure_B: + mcm_pure = np.load(f"{coupling_dir}/mcm_pure.npz") + mcms_dict_nobeam_pure = { + spin_pair: mcm_pure[f"{spin_pair}_binned"] for spin_pair in spin_pairs + } meta.timer.stop("Load mode-coupling matrices", verbose=True) @@ -156,11 +158,15 @@ def transfer(args): for filter_tag in ["filtered", "unfiltered"]: couplings_nobeam = {} - couplings_nobeam_pure = {} - for couplings, mcms_dict in zip([couplings_nobeam, - couplings_nobeam_pure], - [mcms_dict_nobeam, - mcms_dict_nobeam_pure]): + couplings_list = [couplings_nobeam] + mcms_list = [mcms_dict_nobeam] + + if meta.tf_est_pure_B: + couplings_nobeam_pure = {} + couplings_list.append(couplings_nobeam_pure) + mcms_list.append(mcms_dict_nobeam_pure) + + for couplings, mcms_dict in zip(couplings_list, mcms_list): for spin_pair in spin_pairs: if filter_tag == "filtered": tbmcm = np.einsum('ijk,jklm->iklm', trans[spin_pair], @@ -194,10 +200,11 @@ def transfer(args): f"{coupling_dir}/couplings_{filter_tag}.npz", **couplings_nobeam ) - np.savez( - f"{coupling_dir}/couplings_{filter_tag}_pure.npz", - **couplings_nobeam_pure - ) + if meta.tf_est_pure_B: + np.savez( + f"{coupling_dir}/couplings_{filter_tag}_pure.npz", + **couplings_nobeam_pure + ) meta.timer.stop("Compute full coupling", verbose=True)