lilit.functions

  1import numpy as np
  2
  3__all__ = [
  4    "get_keys",
  5    "get_Gauss_keys",
  6    "cov_filling",
  7    "find_spectrum",
  8    "sigma",
  9    "inv_sigma",
 10    "CAMBres2dict",
 11    "txt2dict",
 12]
 13
 14
 15def get_keys(fields: list, *, debug: bool = False):
 16    """Extracts the keys that has to be used as a function of the requested fields. These will be the usual 2-points, e.g., tt, te, ee, etc.
 17
 18    Parameters:
 19        fields (list):
 20            List of fields
 21        debug (bool, optional):
 22            If True, print out the requested keys
 23    """
 24    n = len(fields)
 25    res = [fields[i] + fields[j] for i in range(n) for j in range(i, n)]
 26    if debug:
 27        print(f"\nThe requested keys are {res}")
 28    return res
 29
 30
 31def get_Gauss_keys(n: int, keys: list, *, debug: bool = False):
 32    """Find the proper dictionary keys for the requested fields.
 33
 34    Extracts the keys that has to be used as a function of the requested fields for the Gaussian likelihood. Indeed, the Gaussian likelihood is computed using 4-points, so the keys are different. E.g., there will be keys such as tttt, ttee, tete, etc.
 35
 36    Parameters:
 37        n (int):
 38            Number of fields.
 39        keys (list):
 40            List of keys to use for the computation.
 41        debug (bool, optional):
 42            If set, print the keys that are used, by default False.
 43    """
 44    n = int(n * (n + 1) / 2)
 45    res = np.zeros((n, n, 4), dtype=str)
 46    for i in range(n):
 47        for j in range(i, n):
 48            elem = keys[i] + keys[j]
 49            for k in range(4):
 50                res[i, j, k] = np.asarray(list(elem)[k])
 51                res[j, i, k] = res[i, j, k]
 52    if debug:
 53        print(f"\nThe requested keys are {res}")
 54    return res
 55
 56
 57def cov_filling(
 58    fields: list,
 59    absolute_lmin: int,
 60    absolute_lmax: int,
 61    cov_dict: dict,
 62    lmins: dict = {},
 63    lmaxs: dict = {},
 64):
 65    """Fill covariance matrix with appropriate spectra.
 66
 67    Computes the covariance matrix once given a dictionary. Returns the covariance matrix of the considered fields, in a shape equal to (num_fields x num_fields x lmax). Note that if more than one lmax, or lmin, is specified, there will be null values in the matrices, making them singular. This will be handled in another method.
 68
 69    Parameters:
 70        fields (list):
 71            The list of fields to consider.
 72        absolute_lmin (int):
 73            The minimum multipole to consider.
 74        absolute_lmax (int):
 75            The maximum multipole to consider.
 76        cov_dict (dict):
 77            The input dictionary of spectra.
 78        lmins (dict):
 79            The dictionary of minimum multipole to consider for each field pair.
 80        lmaxs (dict):
 81            The dictionary of maximum multipole to consider for each field pair.
 82    """
 83    n = len(fields)
 84    res = np.zeros((n, n, absolute_lmax + 1))
 85
 86    for i, field1 in enumerate(fields):
 87        for j, field2 in enumerate(fields[i:]):
 88            j += i
 89
 90            key = field1 + field2
 91
 92            lmin = lmins.get(key, absolute_lmin)
 93            lmax = lmaxs.get(key, absolute_lmax)
 94
 95            cov = cov_dict.get(key, np.zeros(lmax + 1))
 96
 97            res[i, j, lmin : lmax + 1] = cov[lmin : lmax + 1]
 98            res[j, i] = res[i, j]
 99
100    return res
101
102
103def find_spectrum(lmin, lmax, input_dict, key):
104    """Find a spectrum in a given dictionary.
105
106    Returns the corresponding power sepctrum for a given key. If the key is not found, it will try to find the reverse key. Otherwise it will fill the array with zeros.
107
108    Parameters:
109        lmin (int):
110            The minimum multipole to consider.
111        lmax (int):
112            The maximum multipole to consider.
113        input_dict (dict):
114            Dictionary where you want to search for keys.
115        key (str):
116            Key to search for.
117    """
118    res = np.zeros(lmax + 1)
119
120    if key in input_dict:
121        cov = input_dict[key]
122    # if the key is not found, try the reverse key, otherwise fill with zeros
123    else:
124        cov = input_dict.get(key[::-1], np.zeros(lmax + 1))
125
126    res[lmin : lmax + 1] = cov[lmin : lmax + 1]
127
128    return res
129
130
131def sigma(n, lmin, lmax, keys, fiduDICT, noiseDICT, fsky=None, fskies=[]):
132    """Define the covariance matrix for the Gaussian case.
133
134    In case of Gaussian likelihood, this returns the covariance matrix needed for the computation of the chi2. Note that the inversion is done in a separate funciton.
135
136    Parameters:
137        n (int):
138            Number of fields.
139        lmin (int):
140            The minimum multipole to consider.
141        lmax (int):
142            The maximum multipole to consider.
143        keys (dict):
144            Keys for the covariance elements.
145        fiduDICT (dict):
146            Dictionary with the fiducial spectra.
147        noiseDICT (dict):
148            Dictionary with the noise spectra.
149    """
150    # The covariance matrix has to be symmetric.
151    # The number of parameters in the likelihood is n.
152    # The covariance matrix is a (n x n x lmax+1) ndarray.
153    # We will store the covariance matrix in a (n x n x lmax+1) ndarray,
154    # where n = int(n * (n + 1) / 2).
155    n = int(n * (n + 1) / 2)
156    res = np.zeros((n, n, lmax + 1))
157    for i in range(n):  # Loop over all combinations of pairs of spectra
158        for j in range(i, n):
159            C_AC = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 0] + keys[i, j, 2])
160            C_BD = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 1] + keys[i, j, 3])
161            C_AD = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 0] + keys[i, j, 3])
162            C_BC = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 1] + keys[i, j, 2])
163            N_AC = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 0] + keys[i, j, 2])
164            N_BD = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 1] + keys[i, j, 3])
165            N_AD = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 0] + keys[i, j, 3])
166            N_BC = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 1] + keys[i, j, 2])
167            ell = np.arange(len(C_AC))
168            if fsky is not None:
169                res[i, j] = (
170                    (C_AC + N_AC) * (C_BD + N_BD) + (C_AD + N_AD) * (C_BC + N_BC)
171                ) / fsky
172            else:
173                AC = keys[i, j, 0] + keys[i, j, 2]
174                BD = keys[i, j, 1] + keys[i, j, 3]
175                AD = keys[i, j, 0] + keys[i, j, 3]
176                BC = keys[i, j, 1] + keys[i, j, 2]
177                AB = keys[i, j, 0] + keys[i, j, 1]
178                CD = keys[i, j, 2] + keys[i, j, 3]
179                res[i, j] = (
180                    np.sqrt(fskies[AC] * fskies[BD]) * (C_AC + N_AC) * (C_BD + N_BD)
181                    + np.sqrt(fskies[AD] * fskies[BC]) * (C_AD + N_AD) * (C_BC + N_BC)
182                ) / (fskies[AB] * fskies[CD])
183            res[i, j, 2:] /= 2 * ell[2:] + 1
184            res[j, i] = res[i, j]
185    return res
186
187
188def inv_sigma(lmin, lmax, sigma):
189    """Invert the covariance matrix of the Gaussian case.
190
191    Inverts the previously calculated sigma ndarray. Note that some elements may be null, thus the covariance may be singular. If so, this also reduces the dimension of the matrix by deleting the corresponding row and column.
192
193    Parameters:
194        lmin (int):
195            The minimum multipole to consider.
196        lmax (int):
197            The maximum multipole to consider.
198        sigma (np.ndarray):
199            (n x n x lmax+1) ndarray with the previously computed sigma (not inverted).
200    """
201    res = np.zeros(sigma.shape)
202
203    for ell in range(lmin, lmax + 1):
204        # Check if matrix is singular
205        COV = sigma[:, :, ell]
206        if np.linalg.det(COV) == 0:
207            # Get indices of null diagonal elements
208            idx = np.where(np.diag(COV) == 0)[0]
209            # Remove corresponding rows and columns
210            COV = np.delete(COV, idx, axis=0)
211            COV = np.delete(COV, idx, axis=1)
212
213        res[:, :, ell] = np.linalg.inv(COV)
214    return res[:, :, lmin:]
215
216
217def CAMBres2dict(camb_results, probes):
218    """Takes the CAMB result product from get_cmb_power_spectra and convert it to a dictionary with the proper keys.
219
220    Parameters:
221        camb_results (CAMBdata):
222            CAMB result product from the method get_cmb_power_spectra.
223        probes (dict):
224            Keys for the spectra dictionary we want to save.
225    """
226    ls = np.arange(camb_results["total"].shape[0], dtype=np.int64)
227    # Mapping between the CAMB keys and the ones we want
228    mapping = {"tt": 0, "ee": 1, "bb": 2, "te": 3, "et": 3}
229    res = {"ell": ls}
230    for probe, i in mapping.items():
231        res[probe] = camb_results["total"][:, i]
232    if "pp" in probes:
233        cl_lens = camb_results.get("lens_potential")
234        if cl_lens is not None:
235            # Save it with the normalization to obtain phiphi
236            array = cl_lens[:, 0].copy()
237            array[2:] /= (res["ell"] * (res["ell"] + 1))[2:]
238            res["pp"] = array
239            if "pt" in probes and "pe" in probes:
240                # Loop over the cross terms correcting the normalization to phiX
241                for i, cross in enumerate(["pt", "pe"]):
242                    array = cl_lens[:, i + 1].copy()
243                    array[2:] /= np.sqrt(res["ell"] * (res["ell"] + 1))[2:]
244                    res[cross] = array
245                    res[cross[::-1]] = res[cross]
246    return res
247
248
249def txt2dict(txt, mapping_probe2colnum=None, apply_ellfactor=None):
250    """Takes a txt file and convert it to a dictionary. This requires a way to map the columns to the keys. Also, it is possible to apply an ell factor to the Cls.
251
252    Parameters:
253        txt (str):
254            Path to txt file containing the spectra as columns.
255        mapping (dict):
256            Dictionary containing the mapping. Keywords will become the new keywords and values represent the index of the corresponding column.
257    """
258    assert (
259        mapping_probe2colnum is not None
260    ), "You must provide a way to map the columns of your txt to the keys of a dictionary"
261    res = {}
262    for probe, i in mapping_probe2colnum.items():
263        ls = np.arange(len(txt[:, i]), dtype=np.int64)
264        res["ell"] = ls
265        if apply_ellfactor:
266            res[probe] = txt[:, i] * ls * (ls + 1) / 2 / np.pi
267        else:
268            res[probe] = txt[:, i]
269    return res
270
271
272__docformat__ = "google"
273__pdoc__ = {}
274__pdoc__[
275    "Likelihood"
276] = "Likelihood class from Cobaya, refer to Cobaya documentation for more information."
def get_keys(fields: list, *, debug: bool = False):
16def get_keys(fields: list, *, debug: bool = False):
17    """Extracts the keys that has to be used as a function of the requested fields. These will be the usual 2-points, e.g., tt, te, ee, etc.
18
19    Parameters:
20        fields (list):
21            List of fields
22        debug (bool, optional):
23            If True, print out the requested keys
24    """
25    n = len(fields)
26    res = [fields[i] + fields[j] for i in range(n) for j in range(i, n)]
27    if debug:
28        print(f"\nThe requested keys are {res}")
29    return res

Extracts the keys that has to be used as a function of the requested fields. These will be the usual 2-points, e.g., tt, te, ee, etc.

Arguments:
  • fields (list): List of fields
  • debug (bool, optional): If True, print out the requested keys
def get_Gauss_keys(n: int, keys: list, *, debug: bool = False):
32def get_Gauss_keys(n: int, keys: list, *, debug: bool = False):
33    """Find the proper dictionary keys for the requested fields.
34
35    Extracts the keys that has to be used as a function of the requested fields for the Gaussian likelihood. Indeed, the Gaussian likelihood is computed using 4-points, so the keys are different. E.g., there will be keys such as tttt, ttee, tete, etc.
36
37    Parameters:
38        n (int):
39            Number of fields.
40        keys (list):
41            List of keys to use for the computation.
42        debug (bool, optional):
43            If set, print the keys that are used, by default False.
44    """
45    n = int(n * (n + 1) / 2)
46    res = np.zeros((n, n, 4), dtype=str)
47    for i in range(n):
48        for j in range(i, n):
49            elem = keys[i] + keys[j]
50            for k in range(4):
51                res[i, j, k] = np.asarray(list(elem)[k])
52                res[j, i, k] = res[i, j, k]
53    if debug:
54        print(f"\nThe requested keys are {res}")
55    return res

Find the proper dictionary keys for the requested fields.

Extracts the keys that has to be used as a function of the requested fields for the Gaussian likelihood. Indeed, the Gaussian likelihood is computed using 4-points, so the keys are different. E.g., there will be keys such as tttt, ttee, tete, etc.

Arguments:
  • n (int): Number of fields.
  • keys (list): List of keys to use for the computation.
  • debug (bool, optional): If set, print the keys that are used, by default False.
def cov_filling( fields: list, absolute_lmin: int, absolute_lmax: int, cov_dict: dict, lmins: dict = {}, lmaxs: dict = {}):
 58def cov_filling(
 59    fields: list,
 60    absolute_lmin: int,
 61    absolute_lmax: int,
 62    cov_dict: dict,
 63    lmins: dict = {},
 64    lmaxs: dict = {},
 65):
 66    """Fill covariance matrix with appropriate spectra.
 67
 68    Computes the covariance matrix once given a dictionary. Returns the covariance matrix of the considered fields, in a shape equal to (num_fields x num_fields x lmax). Note that if more than one lmax, or lmin, is specified, there will be null values in the matrices, making them singular. This will be handled in another method.
 69
 70    Parameters:
 71        fields (list):
 72            The list of fields to consider.
 73        absolute_lmin (int):
 74            The minimum multipole to consider.
 75        absolute_lmax (int):
 76            The maximum multipole to consider.
 77        cov_dict (dict):
 78            The input dictionary of spectra.
 79        lmins (dict):
 80            The dictionary of minimum multipole to consider for each field pair.
 81        lmaxs (dict):
 82            The dictionary of maximum multipole to consider for each field pair.
 83    """
 84    n = len(fields)
 85    res = np.zeros((n, n, absolute_lmax + 1))
 86
 87    for i, field1 in enumerate(fields):
 88        for j, field2 in enumerate(fields[i:]):
 89            j += i
 90
 91            key = field1 + field2
 92
 93            lmin = lmins.get(key, absolute_lmin)
 94            lmax = lmaxs.get(key, absolute_lmax)
 95
 96            cov = cov_dict.get(key, np.zeros(lmax + 1))
 97
 98            res[i, j, lmin : lmax + 1] = cov[lmin : lmax + 1]
 99            res[j, i] = res[i, j]
100
101    return res

Fill covariance matrix with appropriate spectra.

Computes the covariance matrix once given a dictionary. Returns the covariance matrix of the considered fields, in a shape equal to (num_fields x num_fields x lmax). Note that if more than one lmax, or lmin, is specified, there will be null values in the matrices, making them singular. This will be handled in another method.

Arguments:
  • fields (list): The list of fields to consider.
  • absolute_lmin (int): The minimum multipole to consider.
  • absolute_lmax (int): The maximum multipole to consider.
  • cov_dict (dict): The input dictionary of spectra.
  • lmins (dict): The dictionary of minimum multipole to consider for each field pair.
  • lmaxs (dict): The dictionary of maximum multipole to consider for each field pair.
def find_spectrum(lmin, lmax, input_dict, key):
104def find_spectrum(lmin, lmax, input_dict, key):
105    """Find a spectrum in a given dictionary.
106
107    Returns the corresponding power sepctrum for a given key. If the key is not found, it will try to find the reverse key. Otherwise it will fill the array with zeros.
108
109    Parameters:
110        lmin (int):
111            The minimum multipole to consider.
112        lmax (int):
113            The maximum multipole to consider.
114        input_dict (dict):
115            Dictionary where you want to search for keys.
116        key (str):
117            Key to search for.
118    """
119    res = np.zeros(lmax + 1)
120
121    if key in input_dict:
122        cov = input_dict[key]
123    # if the key is not found, try the reverse key, otherwise fill with zeros
124    else:
125        cov = input_dict.get(key[::-1], np.zeros(lmax + 1))
126
127    res[lmin : lmax + 1] = cov[lmin : lmax + 1]
128
129    return res

Find a spectrum in a given dictionary.

Returns the corresponding power sepctrum for a given key. If the key is not found, it will try to find the reverse key. Otherwise it will fill the array with zeros.

Arguments:
  • lmin (int): The minimum multipole to consider.
  • lmax (int): The maximum multipole to consider.
  • input_dict (dict): Dictionary where you want to search for keys.
  • key (str): Key to search for.
def sigma(n, lmin, lmax, keys, fiduDICT, noiseDICT, fsky=None, fskies=[]):
132def sigma(n, lmin, lmax, keys, fiduDICT, noiseDICT, fsky=None, fskies=[]):
133    """Define the covariance matrix for the Gaussian case.
134
135    In case of Gaussian likelihood, this returns the covariance matrix needed for the computation of the chi2. Note that the inversion is done in a separate funciton.
136
137    Parameters:
138        n (int):
139            Number of fields.
140        lmin (int):
141            The minimum multipole to consider.
142        lmax (int):
143            The maximum multipole to consider.
144        keys (dict):
145            Keys for the covariance elements.
146        fiduDICT (dict):
147            Dictionary with the fiducial spectra.
148        noiseDICT (dict):
149            Dictionary with the noise spectra.
150    """
151    # The covariance matrix has to be symmetric.
152    # The number of parameters in the likelihood is n.
153    # The covariance matrix is a (n x n x lmax+1) ndarray.
154    # We will store the covariance matrix in a (n x n x lmax+1) ndarray,
155    # where n = int(n * (n + 1) / 2).
156    n = int(n * (n + 1) / 2)
157    res = np.zeros((n, n, lmax + 1))
158    for i in range(n):  # Loop over all combinations of pairs of spectra
159        for j in range(i, n):
160            C_AC = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 0] + keys[i, j, 2])
161            C_BD = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 1] + keys[i, j, 3])
162            C_AD = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 0] + keys[i, j, 3])
163            C_BC = find_spectrum(lmin, lmax, fiduDICT, keys[i, j, 1] + keys[i, j, 2])
164            N_AC = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 0] + keys[i, j, 2])
165            N_BD = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 1] + keys[i, j, 3])
166            N_AD = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 0] + keys[i, j, 3])
167            N_BC = find_spectrum(lmin, lmax, noiseDICT, keys[i, j, 1] + keys[i, j, 2])
168            ell = np.arange(len(C_AC))
169            if fsky is not None:
170                res[i, j] = (
171                    (C_AC + N_AC) * (C_BD + N_BD) + (C_AD + N_AD) * (C_BC + N_BC)
172                ) / fsky
173            else:
174                AC = keys[i, j, 0] + keys[i, j, 2]
175                BD = keys[i, j, 1] + keys[i, j, 3]
176                AD = keys[i, j, 0] + keys[i, j, 3]
177                BC = keys[i, j, 1] + keys[i, j, 2]
178                AB = keys[i, j, 0] + keys[i, j, 1]
179                CD = keys[i, j, 2] + keys[i, j, 3]
180                res[i, j] = (
181                    np.sqrt(fskies[AC] * fskies[BD]) * (C_AC + N_AC) * (C_BD + N_BD)
182                    + np.sqrt(fskies[AD] * fskies[BC]) * (C_AD + N_AD) * (C_BC + N_BC)
183                ) / (fskies[AB] * fskies[CD])
184            res[i, j, 2:] /= 2 * ell[2:] + 1
185            res[j, i] = res[i, j]
186    return res

Define the covariance matrix for the Gaussian case.

In case of Gaussian likelihood, this returns the covariance matrix needed for the computation of the chi2. Note that the inversion is done in a separate funciton.

Arguments:
  • n (int): Number of fields.
  • lmin (int): The minimum multipole to consider.
  • lmax (int): The maximum multipole to consider.
  • keys (dict): Keys for the covariance elements.
  • fiduDICT (dict): Dictionary with the fiducial spectra.
  • noiseDICT (dict): Dictionary with the noise spectra.
def inv_sigma(lmin, lmax, sigma):
189def inv_sigma(lmin, lmax, sigma):
190    """Invert the covariance matrix of the Gaussian case.
191
192    Inverts the previously calculated sigma ndarray. Note that some elements may be null, thus the covariance may be singular. If so, this also reduces the dimension of the matrix by deleting the corresponding row and column.
193
194    Parameters:
195        lmin (int):
196            The minimum multipole to consider.
197        lmax (int):
198            The maximum multipole to consider.
199        sigma (np.ndarray):
200            (n x n x lmax+1) ndarray with the previously computed sigma (not inverted).
201    """
202    res = np.zeros(sigma.shape)
203
204    for ell in range(lmin, lmax + 1):
205        # Check if matrix is singular
206        COV = sigma[:, :, ell]
207        if np.linalg.det(COV) == 0:
208            # Get indices of null diagonal elements
209            idx = np.where(np.diag(COV) == 0)[0]
210            # Remove corresponding rows and columns
211            COV = np.delete(COV, idx, axis=0)
212            COV = np.delete(COV, idx, axis=1)
213
214        res[:, :, ell] = np.linalg.inv(COV)
215    return res[:, :, lmin:]

Invert the covariance matrix of the Gaussian case.

Inverts the previously calculated sigma ndarray. Note that some elements may be null, thus the covariance may be singular. If so, this also reduces the dimension of the matrix by deleting the corresponding row and column.

Arguments:
  • lmin (int): The minimum multipole to consider.
  • lmax (int): The maximum multipole to consider.
  • sigma (np.ndarray): (n x n x lmax+1) ndarray with the previously computed sigma (not inverted).
def CAMBres2dict(camb_results, probes):
218def CAMBres2dict(camb_results, probes):
219    """Takes the CAMB result product from get_cmb_power_spectra and convert it to a dictionary with the proper keys.
220
221    Parameters:
222        camb_results (CAMBdata):
223            CAMB result product from the method get_cmb_power_spectra.
224        probes (dict):
225            Keys for the spectra dictionary we want to save.
226    """
227    ls = np.arange(camb_results["total"].shape[0], dtype=np.int64)
228    # Mapping between the CAMB keys and the ones we want
229    mapping = {"tt": 0, "ee": 1, "bb": 2, "te": 3, "et": 3}
230    res = {"ell": ls}
231    for probe, i in mapping.items():
232        res[probe] = camb_results["total"][:, i]
233    if "pp" in probes:
234        cl_lens = camb_results.get("lens_potential")
235        if cl_lens is not None:
236            # Save it with the normalization to obtain phiphi
237            array = cl_lens[:, 0].copy()
238            array[2:] /= (res["ell"] * (res["ell"] + 1))[2:]
239            res["pp"] = array
240            if "pt" in probes and "pe" in probes:
241                # Loop over the cross terms correcting the normalization to phiX
242                for i, cross in enumerate(["pt", "pe"]):
243                    array = cl_lens[:, i + 1].copy()
244                    array[2:] /= np.sqrt(res["ell"] * (res["ell"] + 1))[2:]
245                    res[cross] = array
246                    res[cross[::-1]] = res[cross]
247    return res

Takes the CAMB result product from get_cmb_power_spectra and convert it to a dictionary with the proper keys.

Arguments:
  • camb_results (CAMBdata): CAMB result product from the method get_cmb_power_spectra.
  • probes (dict): Keys for the spectra dictionary we want to save.
def txt2dict(txt, mapping_probe2colnum=None, apply_ellfactor=None):
250def txt2dict(txt, mapping_probe2colnum=None, apply_ellfactor=None):
251    """Takes a txt file and convert it to a dictionary. This requires a way to map the columns to the keys. Also, it is possible to apply an ell factor to the Cls.
252
253    Parameters:
254        txt (str):
255            Path to txt file containing the spectra as columns.
256        mapping (dict):
257            Dictionary containing the mapping. Keywords will become the new keywords and values represent the index of the corresponding column.
258    """
259    assert (
260        mapping_probe2colnum is not None
261    ), "You must provide a way to map the columns of your txt to the keys of a dictionary"
262    res = {}
263    for probe, i in mapping_probe2colnum.items():
264        ls = np.arange(len(txt[:, i]), dtype=np.int64)
265        res["ell"] = ls
266        if apply_ellfactor:
267            res[probe] = txt[:, i] * ls * (ls + 1) / 2 / np.pi
268        else:
269            res[probe] = txt[:, i]
270    return res

Takes a txt file and convert it to a dictionary. This requires a way to map the columns to the keys. Also, it is possible to apply an ell factor to the Cls.

Arguments:
  • txt (str): Path to txt file containing the spectra as columns.
  • mapping (dict): Dictionary containing the mapping. Keywords will become the new keywords and values represent the index of the corresponding column.