cosmology

cosmology

Author: Corentin Ravoux

Description : Classes to return a Dachshund input based on data that can be treated via picca software.

DeltaModifier

DeltaModifier(pwd, delta_path)

Bases: object

Transform delta files: shuffle, redshift-cut and subsample lines of sight.

Used both to build null (shuffled) mocks and to thin the line-of-sight density, optionally matching a reference density/separation profile.

Store the output and input delta directories.

Parameters:

Name	Type	Description	Default
pwd	str	Output directory for the modified deltas.	required
delta_path	str	Input delta directory.	required

Source code in lelantos/cosmology.py

def __init__(self, pwd, delta_path):
    """Store the output and input delta directories.

    Args:
        pwd (str): Output directory for the modified deltas.
        delta_path (str): Input delta directory.
    """
    self.pwd = pwd
    self.delta_path = delta_path

shuffle_deltas

shuffle_deltas(other_delta_path=None, other_path_out=None, seed=None)

Randomly permute all delta/ivar pixels to build a null field.

Parameters:

Name	Type	Description	Default
other_delta_path	str	Second delta set shuffled in lock-step (same permutation).	None
other_path_out	str	Output directory for the second set.	None
seed	int	RNG seed (random if None).	None

Source code in lelantos/cosmology.py

def shuffle_deltas(self, other_delta_path=None, other_path_out=None, seed=None):
    """Randomly permute all delta/ivar pixels to build a null field.

    Args:
        other_delta_path (str, optional): Second delta set shuffled in
            lock-step (same permutation).
        other_path_out (str, optional): Output directory for the second set.
        seed (int, optional): RNG seed (random if None).
    """
    namefile = get_delta_list(self.delta_path)
    namefile_other = None
    if other_delta_path is not None:
        namefile_other = get_delta_list(other_delta_path)

    (delta, ivar, delta_other, weight_other) = self.get_delta_sigma_array(
        namefile, namefile_other=namefile_other
    )
    if seed is None:
        seed = np.random.randint(10000000)
    np.random.seed(seed)
    ivar_rand = np.random.permutation(ivar)
    delta_rand = np.random.permutation(delta)
    if other_delta_path is not None:
        np.random.seed(seed)
        weight_other_rand = np.random.permutation(weight_other)
        delta_other_rand = np.random.permutation(delta_other)

    self.write_delta_sigma_array(
        delta_rand,
        ivar_rand,
        namefile,
        other_delta_path=other_delta_path,
        namefile_other=namefile_other,
        delta_other_rand=delta_other_rand,
        weight_other_rand=weight_other_rand,
        other_path_out=other_path_out,
    )

shuffle_deltas_cut_z

shuffle_deltas_cut_z(n_cut, zmin, zmax, other_delta_path=None, other_path_out=None, seed=None)

Shuffle delta/ivar pixels independently within redshift bins.

Parameters:

Name	Type	Description	Default
n_cut	int	Number of redshift bins the shuffling is confined to.	required
zmin	float	Minimum redshift of the binning.	required
zmax	float	Maximum redshift of the binning.	required
other_delta_path	str	Second delta set (lock-step).	None
other_path_out	str	Output directory for the second set.	None
seed	int	RNG seed (random if None).	None

Source code in lelantos/cosmology.py

def shuffle_deltas_cut_z(
    self, n_cut, zmin, zmax, other_delta_path=None, other_path_out=None, seed=None
):
    """Shuffle delta/ivar pixels independently within redshift bins.

    Args:
        n_cut (int): Number of redshift bins the shuffling is confined to.
        zmin (float): Minimum redshift of the binning.
        zmax (float): Maximum redshift of the binning.
        other_delta_path (str, optional): Second delta set (lock-step).
        other_path_out (str, optional): Output directory for the second set.
        seed (int, optional): RNG seed (random if None).
    """
    namefile = get_delta_list(self.delta_path)
    namefile_other = None
    if other_delta_path is not None:
        namefile_other = get_delta_list(other_delta_path)
    redshift_cut = np.linspace(zmin, zmax, n_cut + 1)
    (delta, ivar, delta_other, weight_other) = self.get_delta_sigma_array_cut_z(
        namefile, redshift_cut, namefile_other=namefile_other
    )

    if seed is None:
        seed = np.random.randint(10000000)
    ivar_rand, delta_rand = [], []
    np.random.seed(seed)
    for k in range(n_cut):
        ivar_rand.append(np.random.permutation(ivar[k]))
    np.random.seed(seed)
    for k in range(n_cut):
        delta_rand.append(np.random.permutation(delta[k]))
    if other_delta_path is not None:
        weight_other_rand, delta_other_rand = [], []
        np.random.seed(seed)
        for k in range(n_cut):
            weight_other_rand.append(np.random.permutation(weight_other[k]))
        np.random.seed(seed)
        for k in range(n_cut):
            delta_other_rand.append(np.random.permutation(delta_other[k]))

    self.write_delta_sigma_array_cut_z(
        delta_rand,
        ivar_rand,
        namefile,
        n_cut,
        redshift_cut,
        other_delta_path=other_delta_path,
        namefile_other=namefile_other,
        delta_other_rand=delta_other_rand,
        weight_other_rand=weight_other_rand,
        other_path_out=other_path_out,
    )

get_delta_sigma_array

get_delta_sigma_array(namefile, namefile_other=None)

Concatenate the delta and ivar pixels of all delta files.

Parameters:

Name	Type	Description	Default
namefile	list[str]	Delta file paths.	required
namefile_other	list[str]	Second delta set.	None

Returns:

Name	Type	Description
tuple		(delta, ivar, delta_other, weight_other) flat arrays
		(the *_other entries are None if no second set is given).

Source code in lelantos/cosmology.py

def get_delta_sigma_array(self, namefile, namefile_other=None):
    """Concatenate the delta and ivar pixels of all delta files.

    Args:
        namefile (list[str]): Delta file paths.
        namefile_other (list[str], optional): Second delta set.

    Returns:
        tuple: ``(delta, ivar, delta_other, weight_other)`` flat arrays
        (the ``*_other`` entries are None if no second set is given).
    """
    weight_other, delta_other = None, None
    if namefile_other is not None:
        weight_other, delta_other = [], []
    ivar, delta = [], []
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=True)
        delta_tomo.read()
        for j in range(len(delta_tomo.delta_array)):
            ivar.append(delta_tomo.delta_array[j].ivar)
            delta.append(delta_tomo.delta_array[j].delta)
        if namefile_other is not None:
            delta_tomo_other = tomographic_objects.Delta(
                name=namefile_other[i], pk1d_type=False
            )
            delta_tomo_other.read()
            for j in range(len(delta_tomo_other.delta_array)):
                weight_other.append(delta_tomo_other.delta_array[j].weights)
                delta_other.append(delta_tomo_other.delta_array[j].delta)
    ivar = np.concatenate(ivar, axis=0)
    delta = np.concatenate(delta, axis=0)
    if namefile_other is not None:
        weight_other = np.concatenate(weight_other, axis=0)
        delta_other = np.concatenate(delta_other, axis=0)
    return (delta, ivar, delta_other, weight_other)

write_delta_sigma_array

write_delta_sigma_array(delta_rand, ivar_rand, namefile, other_delta_path=None, namefile_other=None, delta_other_rand=None, weight_other_rand=None, other_path_out=None)

Write shuffled delta/ivar pixels back into per-file delta objects.

Parameters:

Name	Type	Description	Default
delta_rand	ndarray	Shuffled delta pixels.	required
ivar_rand	ndarray	Shuffled ivar pixels.	required
namefile	list[str]	Input delta file paths (define the LOS lengths).	required
other_delta_path	str	Second delta set input directory.	None
namefile_other	list[str]	Second delta file paths.	None
delta_other_rand	ndarray	Shuffled second deltas.	None
weight_other_rand	ndarray	Shuffled second weights.	None
other_path_out	str	Second delta set output directory.	None

Source code in lelantos/cosmology.py

def write_delta_sigma_array(
    self,
    delta_rand,
    ivar_rand,
    namefile,
    other_delta_path=None,
    namefile_other=None,
    delta_other_rand=None,
    weight_other_rand=None,
    other_path_out=None,
):
    """Write shuffled delta/ivar pixels back into per-file delta objects.

    Args:
        delta_rand (numpy.ndarray): Shuffled delta pixels.
        ivar_rand (numpy.ndarray): Shuffled ivar pixels.
        namefile (list[str]): Input delta file paths (define the LOS lengths).
        other_delta_path (str, optional): Second delta set input directory.
        namefile_other (list[str], optional): Second delta file paths.
        delta_other_rand (numpy.ndarray, optional): Shuffled second deltas.
        weight_other_rand (numpy.ndarray, optional): Shuffled second weights.
        other_path_out (str, optional): Second delta set output directory.
    """
    ibegin = 0
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=True)
        delta_tomo.read()
        delta_obj_list = []
        for j in range(len(delta_tomo.delta_array)):
            delta_tomo.delta_array[j].de = delta_rand[
                ibegin : ibegin + len(delta_tomo.delta_array[j].de)
            ]
            delta_tomo.delta_array[j].iv = ivar_rand[
                ibegin : ibegin + len(delta_tomo.delta_array[j].iv)
            ]
            ibegin = ibegin + len(delta_tomo.delta_array[j].de)
            delta_obj_list.append(delta_tomo.delta_array[j])
        name_delta = os.path.join(
            self.pwd, namefile[i].split(self.delta_path)[-1].split("/")[-1]
        )
        new_delta_tomo = tomographic_objects.Delta(name=name_delta, pk1d_type=True)
        new_delta_tomo.delta_array = delta_obj_list
        new_delta_tomo.write()

    if other_delta_path is not None:
        ibegin = 0
        for i in range(len(namefile_other)):
            delta_tomo_other = tomographic_objects.Delta(
                name=namefile_other[i], pk1d_type=False
            )
            delta_tomo_other.read()
            delta_obj_list = []
            for j in range(len(delta_tomo_other.delta_array)):
                delta_tomo_other.delta_array[j].de = delta_other_rand[
                    ibegin : ibegin + len(delta_tomo_other.delta_array[j].de)
                ]
                delta_tomo_other.delta_array[j].we = weight_other_rand[
                    ibegin : ibegin + len(delta_tomo_other.delta_array[j].we)
                ]
                ibegin = ibegin + len(delta_tomo_other.delta_array[j].de)
                delta_obj_list.append(delta_tomo_other.delta_array[j])
            name_delta = os.path.join(
                other_path_out,
                namefile_other[i].split(other_delta_path)[-1].split("/")[-1],
            )
            new_delta_tomo = tomographic_objects.Delta(
                name=name_delta, pk1d_type=False
            )
            new_delta_tomo.delta_array = delta_obj_list
            new_delta_tomo.write()

get_delta_sigma_array_cut_z

get_delta_sigma_array_cut_z(namefile, redshift_cut, namefile_other=None)

Concatenate delta/ivar pixels grouped into redshift bins.

Parameters:

Name	Type	Description	Default
namefile	list[str]	Delta file paths.	required
redshift_cut	array - like	Bin edges (n_cut + 1 values).	required
namefile_other	list[str]	Second delta set.	None

Returns:

Name	Type	Description
tuple		(delta, ivar, delta_other, weight_other) as lists of
		per-bin flat arrays.

Source code in lelantos/cosmology.py

def get_delta_sigma_array_cut_z(self, namefile, redshift_cut, namefile_other=None):
    """Concatenate delta/ivar pixels grouped into redshift bins.

    Args:
        namefile (list[str]): Delta file paths.
        redshift_cut (array-like): Bin edges (``n_cut + 1`` values).
        namefile_other (list[str], optional): Second delta set.

    Returns:
        tuple: ``(delta, ivar, delta_other, weight_other)`` as lists of
        per-bin flat arrays.
    """
    n_cut = len(redshift_cut) - 1
    weight_other, delta_other = None, None
    if namefile_other is not None:
        weight_other, delta_other = [[] for i in range(n_cut)], [
            [] for i in range(n_cut)
        ]
    ivar, delta = [[] for i in range(n_cut)], [[] for i in range(n_cut)]
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=True)
        delta_tomo.read()
        for j in range(len(delta_tomo.delta_array)):
            redshift = (
                10 ** delta_tomo.delta_array[j].log_lambda / utils.lambdaLy
            ) - 1
            for k in range(n_cut):
                mask = (redshift >= redshift_cut[k]) & (
                    redshift < redshift_cut[k + 1]
                )
                ivar[k].append(delta_tomo.delta_array[j].ivar[mask])
                delta[k].append(delta_tomo.delta_array[j].delta[mask])
        if namefile_other is not None:
            delta_tomo_other = tomographic_objects.Delta(
                name=namefile_other[i], pk1d_type=False
            )
            delta_tomo_other.read()
            for j in range(len(delta_tomo_other.delta_array)):
                redshift = (
                    10 ** delta_tomo_other.delta_array[j].log_lambda
                    / utils.lambdaLy
                ) - 1
                for k in range(n_cut):
                    mask = (redshift >= redshift_cut[k]) & (
                        redshift < redshift_cut[k + 1]
                    )
                    weight_other[k].append(
                        delta_tomo_other.delta_array[j].weights[mask]
                    )
                    delta_other[k].append(
                        delta_tomo_other.delta_array[j].delta[mask]
                    )
    for k in range(n_cut):
        ivar[k] = np.concatenate(ivar[k], axis=0)
        delta[k] = np.concatenate(delta[k], axis=0)
    if namefile_other is not None:
        for k in range(n_cut):
            weight_other[k] = np.concatenate(weight_other[k], axis=0)
            delta_other[k] = np.concatenate(delta_other[k], axis=0)
    return (delta, ivar, delta_other, weight_other)

write_delta_sigma_array_cut_z

write_delta_sigma_array_cut_z(delta_rand, ivar_rand, namefile, n_cut, redshift_cut, other_delta_path=None, namefile_other=None, delta_other_rand=None, weight_other_rand=None, other_path_out=None)

Write per-redshift-bin shuffled pixels back into delta objects.

Parameters:

Name	Type	Description	Default
delta_rand	list[array]	Per-bin shuffled delta pixels.	required
ivar_rand	list[array]	Per-bin shuffled ivar pixels.	required
namefile	list[str]	Input delta file paths.	required
n_cut	int	Number of redshift bins.	required
redshift_cut	array - like	Bin edges.	required
other_delta_path	str	Second delta set input directory.	None
namefile_other	list[str]	Second delta file paths.	None
delta_other_rand	list[array]	Second shuffled deltas.	None
weight_other_rand	list[array]	Second shuffled weights.	None
other_path_out	str	Second delta set output directory.	None

Source code in lelantos/cosmology.py

def write_delta_sigma_array_cut_z(
    self,
    delta_rand,
    ivar_rand,
    namefile,
    n_cut,
    redshift_cut,
    other_delta_path=None,
    namefile_other=None,
    delta_other_rand=None,
    weight_other_rand=None,
    other_path_out=None,
):
    """Write per-redshift-bin shuffled pixels back into delta objects.

    Args:
        delta_rand (list[array]): Per-bin shuffled delta pixels.
        ivar_rand (list[array]): Per-bin shuffled ivar pixels.
        namefile (list[str]): Input delta file paths.
        n_cut (int): Number of redshift bins.
        redshift_cut (array-like): Bin edges.
        other_delta_path (str, optional): Second delta set input directory.
        namefile_other (list[str], optional): Second delta file paths.
        delta_other_rand (list[array], optional): Second shuffled deltas.
        weight_other_rand (list[array], optional): Second shuffled weights.
        other_path_out (str, optional): Second delta set output directory.
    """
    ibegin = [0 for i in range(n_cut)]
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=True)
        delta_tomo.read()
        delta_obj_list = []
        for j in range(len(delta_tomo.delta_array)):
            redshift = (
                10 ** delta_tomo.delta_array[j].log_lambda / utils.lambdaLy
            ) - 1
            for k in range(n_cut):
                mask = (redshift >= redshift_cut[k]) & (
                    redshift < redshift_cut[k + 1]
                )
                delta_tomo.delta_array[j].de[mask] = delta_rand[k][
                    ibegin[k] : ibegin[k] + len(delta_tomo.delta_array[j].de[mask])
                ]
                delta_tomo.delta_array[j].iv[mask] = ivar_rand[k][
                    ibegin[k] : ibegin[k] + len(delta_tomo.delta_array[j].iv[mask])
                ]
                ibegin[k] = ibegin[k] + len(delta_tomo.delta_array[j].de[mask])
            delta_obj_list.append(delta_tomo.delta_array[j])
        name_delta = os.path.join(
            self.pwd, namefile[i].split(self.delta_path)[-1].split("/")[-1]
        )
        new_delta_tomo = tomographic_objects.Delta(name=name_delta, pk1d_type=True)
        new_delta_tomo.delta_array = delta_obj_list
        new_delta_tomo.write()

    if other_delta_path is not None:
        ibegin = [0 for i in range(n_cut)]
        for i in range(len(namefile_other)):
            delta_tomo_other = tomographic_objects.Delta(
                name=namefile_other[i], pk1d_type=False
            )
            delta_tomo_other.read()
            delta_obj_list = []
            for j in range(len(delta_tomo_other.delta_array)):
                redshift = (
                    10 ** delta_tomo_other.delta_array[j].log_lambda
                    / utils.lambdaLy
                ) - 1
                for k in range(n_cut):
                    mask = (redshift >= redshift_cut[k]) & (
                        redshift < redshift_cut[k + 1]
                    )
                    delta_tomo_other.delta_array[j].de[mask] = delta_other_rand[k][
                        ibegin[k] : ibegin[k]
                        + len(delta_tomo_other.delta_array[j].de[mask])
                    ]
                    delta_tomo_other.delta_array[j].we[mask] = weight_other_rand[k][
                        ibegin[k] : ibegin[k]
                        + len(delta_tomo_other.delta_array[j].we[mask])
                    ]
                    ibegin[k] = ibegin[k] + len(
                        delta_tomo_other.delta_array[j].de[mask]
                    )
                delta_obj_list.append(delta_tomo_other.delta_array[j])
            name_delta = os.path.join(
                other_path_out,
                namefile_other[i].split(other_delta_path)[-1].split("/")[-1],
            )
            new_delta_tomo = tomographic_objects.Delta(
                name=name_delta, pk1d_type=False
            )
            new_delta_tomo.delta_array = delta_obj_list
            new_delta_tomo.write()

get_new_healpix

get_new_healpix(number_cut, random_density_parameter=None, number_repeat=1, iterative_selection_parameters=None, ra_cut_min=None, ra_cut_max=None, dec_cut_min=None, dec_cut_max=None, z_cut_min=None, z_cut_max=None, center_ra=True)

Build a (optionally thinned) set of lines of sight split by RA.

Dispatches to a plain, randomly-thinned, or iteratively density-matched selection depending on the arguments.

Parameters:

Name	Type	Description	Default
number_cut	int	Number of RA sub-regions.	required
random_density_parameter	float	Target LOS density; None keeps all lines of sight.	None
number_repeat	int	Number of iterative selection trials.	1
iterative_selection_parameters	dict	Reference density / separation config for the iterative selection.	None
ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max	float	Sky window.	required
z_cut_min, z_cut_max	float	Redshift window.	required
center_ra	bool	Recenter RA around 0.	True

Returns:

Name	Type	Description
dict		Per-RA-region lists of selected delta objects.

Source code in lelantos/cosmology.py

def get_new_healpix(
    self,
    number_cut,
    random_density_parameter=None,
    number_repeat=1,
    iterative_selection_parameters=None,
    ra_cut_min=None,
    ra_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    z_cut_min=None,
    z_cut_max=None,
    center_ra=True,
):
    """Build a (optionally thinned) set of lines of sight split by RA.

    Dispatches to a plain, randomly-thinned, or iteratively density-matched
    selection depending on the arguments.

    Args:
        number_cut (int): Number of RA sub-regions.
        random_density_parameter (float, optional): Target LOS density; None
            keeps all lines of sight.
        number_repeat (int, optional): Number of iterative selection trials.
        iterative_selection_parameters (dict, optional): Reference density /
            separation config for the iterative selection.
        ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max (float, optional):
            Sky window.
        z_cut_min, z_cut_max (float, optional): Redshift window.
        center_ra (bool, optional): Recenter RA around 0.

    Returns:
        dict: Per-RA-region lists of selected delta objects.
    """
    if random_density_parameter is None:
        deltas = self.create_healpix(
            number_cut,
            random=None,
            return_len_ra=False,
            ra_cut_min=ra_cut_min,
            ra_cut_max=ra_cut_max,
            dec_cut_min=dec_cut_min,
            dec_cut_max=dec_cut_max,
            center_ra=center_ra,
        )
    else:
        if number_repeat == 1:
            deltas = self.create_healpix(
                number_cut,
                random=random_density_parameter,
                return_len_ra=False,
                ra_cut_min=ra_cut_min,
                ra_cut_max=ra_cut_max,
                dec_cut_min=dec_cut_min,
                dec_cut_max=dec_cut_max,
                center_ra=center_ra,
            )
        else:
            if iterative_selection_parameters is None:
                return KeyError(
                    "Please dictionary parameter for the iterative selection"
                )
            deltas = self.iterate_healpix_creation(
                number_cut,
                random_density_parameter,
                number_repeat,
                iterative_selection_parameters,
                ra_cut_min=ra_cut_min,
                ra_cut_max=ra_cut_max,
                dec_cut_min=dec_cut_min,
                dec_cut_max=dec_cut_max,
                z_cut_min=z_cut_min,
                z_cut_max=z_cut_max,
            )
    return deltas

create_healpix

create_healpix(number_cut, random=None, return_len_ra=False, ra_cut_min=None, ra_cut_max=None, dec_cut_min=None, dec_cut_max=None, center_ra=True)

Split the lines of sight into RA sub-regions, optionally thinned.

Parameters:

Name	Type	Description	Default
number_cut	int	Number of RA sub-regions.	required
random	float	Target LOS density for random thinning.	None
return_len_ra	bool	Also return the number of LOS.	False
ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max	float	Sky window.	required
center_ra	bool	Recenter RA around 0.	True

Returns:

Type	Description
	dict | tuple: The per-region delta dict, or (deltas, n_los) when
	return_len_ra is True.

Source code in lelantos/cosmology.py

def create_healpix(
    self,
    number_cut,
    random=None,
    return_len_ra=False,
    ra_cut_min=None,
    ra_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    center_ra=True,
):
    """Split the lines of sight into RA sub-regions, optionally thinned.

    Args:
        number_cut (int): Number of RA sub-regions.
        random (float, optional): Target LOS density for random thinning.
        return_len_ra (bool, optional): Also return the number of LOS.
        ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max (float, optional):
            Sky window.
        center_ra (bool, optional): Recenter RA around 0.

    Returns:
        dict | tuple: The per-region delta dict, or ``(deltas, n_los)`` when
        ``return_len_ra`` is True.
    """
    namefile = get_delta_list(self.delta_path)
    deltas = {}
    ra_array = []
    dec_array = []
    for cut in range(number_cut):
        deltas[cut] = []
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=True)
        delta_tomo.read()
        for j in range(len(delta_tomo.delta_array)):
            if center_ra:
                if delta_tomo.delta_array[j].ra * 180 / np.pi > 180:
                    ra = (delta_tomo.delta_array[j].ra * 180 / np.pi) - 360
                else:
                    ra = delta_tomo.delta_array[j].ra * 180 / np.pi
            else:
                ra = delta_tomo.delta_array[j].ra * 180 / np.pi
            dec = delta_tomo.delta_array[j].dec * 180 / np.pi
            if (
                (ra > ra_cut_min)
                & (ra < ra_cut_max)
                & (dec > dec_cut_min)
                & (dec < dec_cut_max)
            ):
                ra_array.append(ra)
                dec_array.append(dec)
                for cut in range(number_cut):
                    interval_ra = ((cut) / (number_cut)) * (
                        ra_cut_max - ra_cut_min
                    ) + ra_cut_min, ((cut + 1) / (number_cut)) * (
                        ra_cut_max - ra_cut_min
                    ) + ra_cut_min
                    if (ra > interval_ra[0]) & (ra <= interval_ra[1]):
                        if center_ra:
                            delta_tomo.delta_array[j].ra = (
                                delta_tomo.delta_array[j].ra - 2 * np.pi
                            )
                        deltas[cut].append(delta_tomo.delta_array[j])
    if random is not None:
        deltas = self.randomize_choice_of_los(
            deltas,
            random,
            number_cut,
            len(ra_array),
            ra_cut_min=ra_cut_min,
            ra_cut_max=ra_cut_max,
            dec_cut_min=dec_cut_min,
            dec_cut_max=dec_cut_max,
        )
    if return_len_ra:
        return (deltas, len(ra_array))
    else:
        return deltas

randomize_choice_of_los

randomize_choice_of_los(deltas_dict, random, number_cut, number_ra, ra_cut_min=None, ra_cut_max=None, dec_cut_min=None, dec_cut_max=None)

Randomly downsample the lines of sight to a target density.

Parameters:

Name	Type	Description	Default
deltas_dict	dict	Per-region delta lists.	required
random	float	Target LOS density (per deg^2).	required
number_cut	int	Number of RA sub-regions.	required
number_ra	int	Current total number of lines of sight.	required
ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max	float	Sky window (defines the area).	required

Returns:

Name	Type	Description
dict		The thinned per-region delta lists.

Source code in lelantos/cosmology.py

def randomize_choice_of_los(
    self,
    deltas_dict,
    random,
    number_cut,
    number_ra,
    ra_cut_min=None,
    ra_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
):
    """Randomly downsample the lines of sight to a target density.

    Args:
        deltas_dict (dict): Per-region delta lists.
        random (float): Target LOS density (per deg^2).
        number_cut (int): Number of RA sub-regions.
        number_ra (int): Current total number of lines of sight.
        ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max (float, optional):
            Sky window (defines the area).

    Returns:
        dict: The thinned per-region delta lists.
    """
    deltas = deltas_dict.copy()
    density = number_ra / ((ra_cut_max - ra_cut_min) * (dec_cut_max - dec_cut_min))
    utils.Logger.add("density before random choice =" + str(density))
    random_cut = random / density
    ra_random = []
    dec_random = []
    for cut in range(number_cut):
        number_of_delta_to_select = int(round(len(deltas[cut]) * random_cut, 0))
        deltas[cut] = sample(deltas[cut], number_of_delta_to_select)
        for i in range(len(deltas[cut])):
            ra_random.append(deltas[cut][i].ra * 180 / np.pi)
            dec_random.append(deltas[cut][i].dec * 180 / np.pi)
    density = len(ra_random) / (
        (ra_cut_max - ra_cut_min) * (dec_cut_max - dec_cut_min)
    )
    utils.Logger.add("density after random choice =" + str(density))
    return deltas

iterate_healpix_creation

iterate_healpix_creation(number_cut, random_density_parameter, number_repeat, iterative_selection_parameters, property_file_name, ra_cut_min=None, ra_cut_max=None, dec_cut_min=None, dec_cut_max=None, z_cut_min=None, z_cut_max=None, center_ra=True)

Pick the random LOS thinning closest to a reference n(z)/separation.

Repeats the random selection number_repeat times and keeps the draw that best matches the reference density and mean-separation profiles.

Parameters:

Name	Type	Description	Default
number_cut	int	Number of RA sub-regions.	required
random_density_parameter	float	Target LOS density.	required
number_repeat	int	Number of random trials.	required
iterative_selection_parameters	dict	Reference density/separation file names, Om and coordinate transform.	required
property_file_name	str	Pixel-property file for the cartesian map.	required
ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max	float	Sky window.	required
z_cut_min, z_cut_max	float	Redshift window.	required
center_ra	bool	Recenter RA around 0.	True

Returns:

Name	Type	Description
dict		The best-matching per-region delta selection.

Source code in lelantos/cosmology.py

def iterate_healpix_creation(
    self,
    number_cut,
    random_density_parameter,
    number_repeat,
    iterative_selection_parameters,
    property_file_name,
    ra_cut_min=None,
    ra_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    z_cut_min=None,
    z_cut_max=None,
    center_ra=True,
):
    """Pick the random LOS thinning closest to a reference n(z)/separation.

    Repeats the random selection ``number_repeat`` times and keeps the draw
    that best matches the reference density and mean-separation profiles.

    Args:
        number_cut (int): Number of RA sub-regions.
        random_density_parameter (float): Target LOS density.
        number_repeat (int): Number of random trials.
        iterative_selection_parameters (dict): Reference density/separation
            file names, ``Om`` and coordinate transform.
        property_file_name (str): Pixel-property file for the cartesian map.
        ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max (float, optional):
            Sky window.
        z_cut_min, z_cut_max (float, optional): Redshift window.
        center_ra (bool, optional): Recenter RA around 0.

    Returns:
        dict: The best-matching per-region delta selection.
    """
    density_names = iterative_selection_parameters["density_names"]
    dperp_names = iterative_selection_parameters["separation_names"]
    Om = iterative_selection_parameters["Om"]
    (rcomov, distang, inv_rcomov, inv_distang) = utils.get_cosmo_function(Om)
    coordinate_transform = iterative_selection_parameters["coordinate_transform"]
    suplementary_parameters = utils.return_suplementary_parameters(
        coordinate_transform, zmin=z_cut_min, zmax=z_cut_max
    )
    density_ref = {}
    dperp_ref = {}
    for cut in range(number_cut):
        density_ref[cut] = PixelAnalizer.read_density_file(density_names[cut])[1]
        dperp_ref[cut] = PixelAnalizer.read_dperp_file(dperp_names[cut])[1]
    min_density, min_dperp = np.inf, np.inf
    delta_to_keep = {}
    utils.Logger.add("Beginning of the random iteration selection")
    deltas_dict, number_ra = self.create_healpix(
        number_cut, random=False, return_len_ra=True, center_ra=center_ra
    )
    deltas_random = {}
    for i in range(number_repeat):
        utils.Logger.add("repeat " + str(i))
        utils.Logger.add("deltas " + str(i) + " computed")
        diff_dperp, diff_density = [], []
        deltas_random = self.randomize_choice_of_los(
            deltas_dict,
            random_density_parameter,
            number_cut,
            number_ra,
            ra_cut_min=ra_cut_min,
            ra_cut_max=ra_cut_max,
            dec_cut_min=dec_cut_min,
            dec_cut_max=dec_cut_max,
        )
        for cut in range(number_cut):
            delta_file = tomographic_objects.Delta(
                delta_file=None,
                delta_array=deltas_random[cut],
                name="delta_" + str(cut) + "_to_test.pickle",
            )
            delta_file.write()
            namefile = "delta_" + str(cut) + "_to_test.pickle"
            (ra, dec, z, zqso, ids, sigmas, deltas) = get_deltas(namefile)
            sky_deltas = np.array(
                [
                    [ra[i], dec[i], z[i][j], sigmas[i][j], deltas[i][j]]
                    for i in range(len(ra))
                    for j in range(len(z[i]))
                ]
            )
            sky_deltas = sky_deltas[
                utils.cut_sky_catalog(
                    sky_deltas[:, 0],
                    sky_deltas[:, 1],
                    sky_deltas[:, 2],
                    ramin=ra_cut_min,
                    ramax=ra_cut_max,
                    decmin=dec_cut_min,
                    decmax=dec_cut_max,
                    zmin=z_cut_min,
                    zmax=z_cut_max,
                )
            ]
            cartesian_deltas = np.zeros(sky_deltas.shape)
            (
                cartesian_deltas[:, 0],
                cartesian_deltas[:, 1],
                cartesian_deltas[:, 2],
            ) = utils.convert_sky_to_cartesian(
                sky_deltas[:, 0],
                sky_deltas[:, 1],
                sky_deltas[:, 2],
                coordinate_transform,
                rcomov=rcomov,
                distang=distang,
                suplementary_parameters=suplementary_parameters,
            )
            pixel = tomographic_objects.Pixel.init_from_property_files(
                property_file_name, pixel_array=cartesian_deltas, name=None
            )
            pixel_analyzer = PixelAnalizer(pixel=pixel)
            (
                zpar,
                dperpz,
                densityz,
            ) = pixel_analyzer.compute_plot_mean_distance_density("", plot=False)
            diff_dperp.append(
                np.mean(abs(np.array(dperpz) - np.array(dperp_ref[cut])))
            )
            diff_density.append(
                np.mean(abs(np.array(densityz) - np.array(density_ref[cut])))
            )
        utils.Logger.add(
            "Mean difference in term LOS density : {}".format(np.mean(diff_density))
        )
        utils.Logger.add(
            "Mean difference in term of Mean LOS separation : {}".format(
                np.mean(diff_dperp)
            )
        )
        if (np.mean(diff_density) < min_density) & (
            np.mean(diff_dperp) < min_dperp
        ):
            utils.Logger.add("Better at the repeat " + str(i))
            min_density = np.mean(diff_density)
            min_dperp = np.mean(diff_dperp)
            delta_to_keep = deltas_random
    for cut in range(number_cut):
        os.remove("delta_" + str(cut) + "_to_test.pickle")
    utils.Logger.add("End of the random iteration selection")
    return delta_to_keep

save_deltas

save_deltas(deltas, name_out, number_cut)

Write each RA sub-region of selected deltas to its own FITS file.

Parameters:

Name	Type	Description	Default
deltas	dict	Per-region delta lists.	required
name_out	str	Output file base name.	required
number_cut	int	Number of RA sub-regions.	required

Source code in lelantos/cosmology.py

def save_deltas(self, deltas, name_out, number_cut):
    """Write each RA sub-region of selected deltas to its own FITS file.

    Args:
        deltas (dict): Per-region delta lists.
        name_out (str): Output file base name.
        number_cut (int): Number of RA sub-regions.
    """
    for cut in range(number_cut):
        delta = tomographic_objects.Delta(
            name=os.path.join(self.pwd, f"{name_out}_{cut}.fits"),
            delta_array=deltas[cut],
        )
        delta.write()

subsample_deltas

subsample_deltas(name_out, number_cut, random_density_parameter=None, number_repeat=1, iterative_selection_parameters=None, ra_cut_min=None, ra_cut_max=None, dec_cut_min=None, dec_cut_max=None, z_cut_min=None, z_cut_max=None, center_ra=True)

Build a (thinned) LOS selection and write it to disk.

Parameters:

Name	Type	Description	Default
name_out	str	Output file base name.	required
number_cut	int	Number of RA sub-regions.	required
random_density_parameter	float	Target LOS density.	None
number_repeat	int	Number of iterative trials.	1
iterative_selection_parameters	dict	Reference config.	None
ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max	float	Sky window.	required
z_cut_min, z_cut_max	float	Redshift window.	required
center_ra	bool	Recenter RA around 0.	True

Source code in lelantos/cosmology.py

def subsample_deltas(
    self,
    name_out,
    number_cut,
    random_density_parameter=None,
    number_repeat=1,
    iterative_selection_parameters=None,
    ra_cut_min=None,
    ra_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    z_cut_min=None,
    z_cut_max=None,
    center_ra=True,
):
    """Build a (thinned) LOS selection and write it to disk.

    Args:
        name_out (str): Output file base name.
        number_cut (int): Number of RA sub-regions.
        random_density_parameter (float, optional): Target LOS density.
        number_repeat (int, optional): Number of iterative trials.
        iterative_selection_parameters (dict, optional): Reference config.
        ra_cut_min, ra_cut_max, dec_cut_min, dec_cut_max (float, optional):
            Sky window.
        z_cut_min, z_cut_max (float, optional): Redshift window.
        center_ra (bool, optional): Recenter RA around 0.
    """
    deltas = self.get_new_healpix(
        number_cut,
        random_density_parameter=random_density_parameter,
        number_repeat=number_repeat,
        iterative_selection_parameters=iterative_selection_parameters,
        ra_cut_min=ra_cut_min,
        ra_cut_max=ra_cut_max,
        dec_cut_min=dec_cut_min,
        dec_cut_max=dec_cut_max,
        z_cut_min=z_cut_min,
        z_cut_max=z_cut_max,
        center_ra=center_ra,
    )
    self.save_deltas(deltas, name_out, number_cut)

DeltaConverter

DeltaConverter(pwd, Omega_m, delta_path, coordinate_transform, plot_pixel_properties, software, return_qso_catalog=None, return_dla_catalog=None, dla_catalog=None, return_sky_catalogs=False, repeat=False, center_ra=True)

Convert picca delta files into a tomographic solver input.

Applies the footprint / redshift / sigma cuts, transforms sky coordinates to cartesian (Mpc.h^-1), writes the pixel binary (serial or tiled parallel) and the map property file, and emits the QSO/DLA sky and cartesian catalogs.

Store the conversion configuration.

Parameters:

Name	Type	Description	Default
pwd	str	Output directory.	required
Omega_m	float	Fiducial matter density for the comoving transform.	required
delta_path	str	Input delta directory.	required
coordinate_transform	str	Sky<->cartesian transform mode.	required
plot_pixel_properties	bool	Produce LOS density/separation plots.	required
software	str	Solver backend (e.g. "dachshund").	required
return_qso_catalog	str	Output QSO catalog file name.	None
return_dla_catalog	str	Output DLA catalog file name.	None
dla_catalog	str	Input DLA catalog file.	None
return_sky_catalogs	bool	Also write sky-coordinate catalogs.	False
repeat	bool	Merge multiple exposures of each LOS.	False
center_ra	bool	Recenter RA around 0.	True

Source code in lelantos/cosmology.py

def __init__(
    self,
    pwd,
    Omega_m,
    delta_path,
    coordinate_transform,
    plot_pixel_properties,
    software,
    return_qso_catalog=None,
    return_dla_catalog=None,
    dla_catalog=None,
    return_sky_catalogs=False,
    repeat=False,
    center_ra=True,
):
    """Store the conversion configuration.

    Args:
        pwd (str): Output directory.
        Omega_m (float): Fiducial matter density for the comoving transform.
        delta_path (str): Input delta directory.
        coordinate_transform (str): Sky<->cartesian transform mode.
        plot_pixel_properties (bool): Produce LOS density/separation plots.
        software (str): Solver backend (e.g. ``"dachshund"``).
        return_qso_catalog (str, optional): Output QSO catalog file name.
        return_dla_catalog (str, optional): Output DLA catalog file name.
        dla_catalog (str, optional): Input DLA catalog file.
        return_sky_catalogs (bool, optional): Also write sky-coordinate
            catalogs.
        repeat (bool, optional): Merge multiple exposures of each LOS.
        center_ra (bool, optional): Recenter RA around 0.
    """
    self.pwd = pwd
    self.delta_path = delta_path
    self.Omega_m = Omega_m
    self.coordinate_transform = coordinate_transform
    self.plot_pixel_properties = plot_pixel_properties
    self.return_qso_catalog = return_qso_catalog
    self.return_dla_catalog = return_dla_catalog
    self.dla_catalog = dla_catalog
    self.return_sky_catalogs = return_sky_catalogs
    self.repeat = repeat
    self.software = software
    self.center_ra = center_ra

transform_delta_to_pixel_file

transform_delta_to_pixel_file(rebin=None, sigma_min=None, sigma_max=None, z_cut_min=None, z_cut_max=None, dec_cut_min=None, dec_cut_max=None, ra_cut_min=None, ra_cut_max=None)

Read the deltas and produce cartesian/sky pixel and catalog arrays.

Reads (optionally exposure-merged, rebinned) deltas, applies the sky / redshift cuts, converts to cartesian coordinates shifted to the box origin, clips sigma, and builds the matching QSO and DLA arrays.

Parameters:

Name	Type	Description	Default
rebin	int	Rebin factor along each LOS.	None
sigma_min, sigma_max	float	Clip the pixel noise.	required
z_cut_min, z_cut_max	float	Redshift window.	required
dec_cut_min, dec_cut_max, ra_cut_min, ra_cut_max	float	Sky window.	required

Returns:

Name	Type	Description
tuple		``(cartesian_deltas, cartesian_qso_catalog,
		cartesian_dla_catalog, sky_deltas, sky_qso_catalog,
		sky_dla_catalog, properties_map_pixels)``.

Source code in lelantos/cosmology.py

def transform_delta_to_pixel_file(
    self,
    rebin=None,
    sigma_min=None,
    sigma_max=None,
    z_cut_min=None,
    z_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    ra_cut_min=None,
    ra_cut_max=None,
):
    """Read the deltas and produce cartesian/sky pixel and catalog arrays.

    Reads (optionally exposure-merged, rebinned) deltas, applies the sky /
    redshift cuts, converts to cartesian coordinates shifted to the box
    origin, clips sigma, and builds the matching QSO and DLA arrays.

    Args:
        rebin (int, optional): Rebin factor along each LOS.
        sigma_min, sigma_max (float, optional): Clip the pixel noise.
        z_cut_min, z_cut_max (float, optional): Redshift window.
        dec_cut_min, dec_cut_max, ra_cut_min, ra_cut_max (float, optional):
            Sky window.

    Returns:
        tuple: ``(cartesian_deltas, cartesian_qso_catalog,
        cartesian_dla_catalog, sky_deltas, sky_qso_catalog,
        sky_dla_catalog, properties_map_pixels)``.
    """
    namefile = get_delta_list(self.delta_path)
    # namefile = preselect_deltas(namefile,
    #                             ramin=ra_cut_min,
    #                             ramax=ra_cut_max,
    #                             decmin=dec_cut_min,
    #                             decmax=dec_cut_max)
    properties_map_pixels = {}
    (rcomov, distang, inv_rcomov, inv_distang) = utils.get_cosmo_function(
        self.Omega_m
    )
    if self.repeat:
        (
            ra,
            dec,
            z,
            zqso,
            ids,
            sigmas,
            deltas,
        ) = get_merged_multiple_exposure_deltas(namefile)
    else:
        (ra, dec, z, zqso, ids, sigmas, deltas) = get_deltas(
            namefile, center_ra=self.center_ra
        )

    if rebin is not None:
        (z, deltas, sigmas) = self.rebin_data(z, deltas, sigmas, rebin)

    if self.return_dla_catalog is not None:
        zdlas, z_qso_dlas = [], []
        if self.dla_catalog is None:
            raise KeyError(
                "Please give a DLA catalog name or turn off the return_dla_catalog option"
            )
        dla_catalog = tomographic_objects.DLACatalog.init_from_fits(
            self.dla_catalog
        )
        for i in range(len(ids)):
            mask = dla_catalog.primary_key == ids[i]
            if z_cut_min is not None:
                mask &= dla_catalog.coord_z > z_cut_min
            if z_cut_max is not None:
                mask &= dla_catalog.coord_z < z_cut_max
            zdlas.append(dla_catalog.coord_z[mask])
            z_qso_dlas.append(dla_catalog.z_qso[mask])

    sky_deltas = np.array(
        [
            [ra[i], dec[i], z[i][j], sigmas[i][j], deltas[i][j]]
            for i in range(len(ra))
            for j in range(len(z[i]))
        ]
    )
    sky_deltas = sky_deltas[
        utils.cut_sky_catalog(
            sky_deltas[:, 0],
            sky_deltas[:, 1],
            sky_deltas[:, 2],
            ramin=ra_cut_min,
            ramax=ra_cut_max,
            decmin=dec_cut_min,
            decmax=dec_cut_max,
            zmin=z_cut_min,
            zmax=z_cut_max,
        )
    ]
    suplementary_parameters = utils.return_suplementary_parameters(
        self.coordinate_transform,
        zmin=np.min(sky_deltas[:, 2]),
        zmax=np.max(sky_deltas[:, 2]),
    )
    cartesian_deltas = np.zeros(sky_deltas.shape)
    (
        cartesian_deltas[:, 0],
        cartesian_deltas[:, 1],
        cartesian_deltas[:, 2],
    ) = utils.convert_sky_to_cartesian(
        sky_deltas[:, 0],
        sky_deltas[:, 1],
        sky_deltas[:, 2],
        self.coordinate_transform,
        rcomov=rcomov,
        distang=distang,
        suplementary_parameters=suplementary_parameters,
    )
    cartesian_deltas[:, 3], cartesian_deltas[:, 4] = (
        sky_deltas[:, 3],
        sky_deltas[:, 4],
    )
    (
        properties_map_pixels["minx"],
        properties_map_pixels["miny"],
        properties_map_pixels["minz"],
    ) = (
        np.min(cartesian_deltas[:, 0]),
        np.min(cartesian_deltas[:, 1]),
        np.min(cartesian_deltas[:, 2]),
    )
    (
        properties_map_pixels["maxx"],
        properties_map_pixels["maxy"],
        properties_map_pixels["maxz"],
    ) = (
        np.max(cartesian_deltas[:, 0]),
        np.max(cartesian_deltas[:, 1]),
        np.max(cartesian_deltas[:, 2]),
    )
    (
        properties_map_pixels["minra"],
        properties_map_pixels["mindec"],
        properties_map_pixels["minredshift"],
    ) = (
        np.min(sky_deltas[:, 0]),
        np.min(sky_deltas[:, 1]),
        np.min(sky_deltas[:, 2]),
    )
    (
        properties_map_pixels["maxra"],
        properties_map_pixels["maxdec"],
        properties_map_pixels["maxredshift"],
    ) = (
        np.max(sky_deltas[:, 0]),
        np.max(sky_deltas[:, 1]),
        np.max(sky_deltas[:, 2]),
    )
    cartesian_deltas = cartesian_deltas - np.array(
        [
            properties_map_pixels["minx"],
            properties_map_pixels["miny"],
            properties_map_pixels["minz"],
            0,
            0,
        ]
    )

    if sigma_min is not None:
        cartesian_deltas[:, 3][cartesian_deltas[:, 3] < sigma_min] = sigma_min
    if sigma_max is not None:
        cartesian_deltas[:, 3][cartesian_deltas[:, 3] > sigma_max] = sigma_max

    if self.return_qso_catalog is not None:
        sky_qso_catalog = np.array(
            [[ra[i], dec[i], zqso[i], ids[i]] for i in range(len(ra))]
        )
        sky_qso_catalog = sky_qso_catalog[
            utils.cut_sky_catalog(
                sky_qso_catalog[:, 0],
                sky_qso_catalog[:, 1],
                sky_qso_catalog[:, 2],
                ramin=ra_cut_min,
                ramax=ra_cut_max,
                decmin=dec_cut_min,
                decmax=dec_cut_max,
                zmin=z_cut_min,
                zmax=z_cut_max,
            ),
            :,
        ]
        cartesian_qso_catalog = np.zeros(sky_qso_catalog.shape)
        (
            cartesian_qso_catalog[:, 0],
            cartesian_qso_catalog[:, 1],
            cartesian_qso_catalog[:, 2],
        ) = utils.convert_sky_to_cartesian(
            sky_qso_catalog[:, 0],
            sky_qso_catalog[:, 1],
            sky_qso_catalog[:, 2],
            self.coordinate_transform,
            rcomov=rcomov,
            distang=distang,
            suplementary_parameters=suplementary_parameters,
        )
        cartesian_qso_catalog[:, 3] = sky_qso_catalog[:, 3]
        cartesian_qso_catalog = cartesian_qso_catalog - np.array(
            [
                properties_map_pixels["minx"],
                properties_map_pixels["miny"],
                properties_map_pixels["minz"],
                0,
            ]
        )

    else:
        sky_qso_catalog, cartesian_qso_catalog = None, None

    if self.return_dla_catalog is not None:
        sky_dla_catalog = np.array(
            [
                [ra[i], dec[i], zdlas[i][j], z_qso_dlas[i][j]]
                for i in range(len(ra))
                for j in range(len(zdlas[i]))
            ]
        )
        sky_dla_catalog = sky_dla_catalog[
            utils.cut_sky_catalog(
                sky_dla_catalog[:, 0],
                sky_dla_catalog[:, 1],
                sky_dla_catalog[:, 2],
                ramin=ra_cut_min,
                ramax=ra_cut_max,
                decmin=dec_cut_min,
                decmax=dec_cut_max,
                zmin=z_cut_min,
                zmax=z_cut_max,
            )
        ]
        cartesian_dla_catalog = np.zeros(sky_dla_catalog.shape)
        (
            cartesian_dla_catalog[:, 0],
            cartesian_dla_catalog[:, 1],
            cartesian_dla_catalog[:, 2],
        ) = utils.convert_sky_to_cartesian(
            sky_dla_catalog[:, 0],
            sky_dla_catalog[:, 1],
            sky_dla_catalog[:, 2],
            self.coordinate_transform,
            rcomov=rcomov,
            distang=distang,
            suplementary_parameters=suplementary_parameters,
        )
        cartesian_qso_catalog[:, 3] = sky_qso_catalog[:, 3]
        cartesian_dla_catalog = cartesian_dla_catalog - np.array(
            [
                properties_map_pixels["minx"],
                properties_map_pixels["miny"],
                properties_map_pixels["minz"],
                0,
            ]
        )
    else:
        sky_dla_catalog, cartesian_dla_catalog = None, None

    return (
        cartesian_deltas,
        cartesian_qso_catalog,
        cartesian_dla_catalog,
        sky_deltas,
        sky_qso_catalog,
        sky_dla_catalog,
        properties_map_pixels,
    )

rebin_data

rebin_data(z, deltas, sigmas, bin_pixel, method='gauss')

Rebin each line of sight by grouping bin_pixel pixels together.

Parameters:

Name	Type	Description	Default
z	list[array]	Per-LOS redshift arrays.	required
deltas	list[array]	Per-LOS delta arrays.	required
sigmas	list[array]	Per-LOS sigma arrays.	required
bin_pixel	int	Number of pixels per output bin.	required
method	str	Reduction ("gauss", "mean" ...).	'gauss'

Returns:

Name	Type	Description
tuple		(z, deltas, sigmas) rebinned in place.

Source code in lelantos/cosmology.py

def rebin_data(self, z, deltas, sigmas, bin_pixel, method="gauss"):
    """Rebin each line of sight by grouping ``bin_pixel`` pixels together.

    Args:
        z (list[array]): Per-LOS redshift arrays.
        deltas (list[array]): Per-LOS delta arrays.
        sigmas (list[array]): Per-LOS sigma arrays.
        bin_pixel (int): Number of pixels per output bin.
        method (str, optional): Reduction (``"gauss"``, ``"mean"`` ...).

    Returns:
        tuple: ``(z, deltas, sigmas)`` rebinned in place.
    """
    for i in range(len(z)):
        if len(z[i]) > 1:
            if len(z[i]) <= bin_pixel:
                z[i] = [np.mean(z[i])]
                deltas[i] = [np.mean(deltas[i])]
                sigmas[i] = [np.mean(sigmas[i])]
            else:
                new_shape = len(z[i]) // bin_pixel
                first_coord = len(z[i]) - (len(z[i]) // bin_pixel) * bin_pixel
                if first_coord == 0:
                    z[i] = utils.bin_ndarray(
                        np.array(z[i])[:], [new_shape], operation=method
                    )
                    deltas[i] = utils.bin_ndarray(
                        np.array(deltas[i])[:], [new_shape], operation=method
                    )
                    sigmas[i] = utils.bin_ndarray(
                        np.array(sigmas[i])[:], [new_shape], operation=method
                    )
                else:
                    z[i] = np.concatenate(
                        [
                            [np.mean(np.array(z[i])[:first_coord])],
                            utils.bin_ndarray(
                                np.array(z[i])[first_coord:],
                                [new_shape],
                                operation=method,
                            ),
                        ]
                    )
                    deltas[i] = np.concatenate(
                        [
                            [np.mean(np.array(deltas[i])[:first_coord])],
                            utils.bin_ndarray(
                                np.array(deltas[i])[first_coord:],
                                [new_shape],
                                operation=method,
                            ),
                        ]
                    )
                    sigmas[i] = np.concatenate(
                        [
                            [np.mean(np.array(sigmas[i])[:first_coord])],
                            utils.bin_ndarray(
                                np.array(sigmas[i])[first_coord:],
                                [new_shape],
                                operation=method,
                            ),
                        ]
                    )
    return (z, deltas, sigmas)

create_input_files

create_input_files(coordinates_to_write, properties, name_pixel, create_launcher=None)

Write the solver pixel input, dispatching on the solver backend.

Parameters:

Name	Type	Description	Default
coordinates_to_write	ndarray	Pixel coordinates + values.	required
properties	dict	Solver/geometry parameters.	required
name_pixel	str	Output pixel file name.	required
create_launcher	str	If set, also write the launcher.	None

Source code in lelantos/cosmology.py

def create_input_files(
    self, coordinates_to_write, properties, name_pixel, create_launcher=None
):
    """Write the solver pixel input, dispatching on the solver backend.

    Args:
        coordinates_to_write (numpy.ndarray): Pixel coordinates + values.
        properties (dict): Solver/geometry parameters.
        name_pixel (str): Output pixel file name.
        create_launcher (str, optional): If set, also write the launcher.
    """
    if self.software.lower() == "dachshund":
        self.create_dachshund_input_files(
            coordinates_to_write,
            properties,
            name_pixel,
            create_launcher=create_launcher,
        )

create_dachshund_input_files

create_dachshund_input_files(coordinates_to_write, properties, name_pixel, create_launcher=None)

Write the Dachshund pixel binary and optionally its .cfg launcher.

Parameters:

Name	Type	Description	Default
coordinates_to_write	ndarray	Pixel coordinates + values.	required
properties	dict	Solver/geometry parameters.	required
name_pixel	str	Output pixel file name.	required
create_launcher	str	If set, launcher base name to write.	None

Source code in lelantos/cosmology.py

def create_dachshund_input_files(
    self, coordinates_to_write, properties, name_pixel, create_launcher=None
):
    """Write the Dachshund pixel binary and optionally its ``.cfg`` launcher.

    Args:
        coordinates_to_write (numpy.ndarray): Pixel coordinates + values.
        properties (dict): Solver/geometry parameters.
        name_pixel (str): Output pixel file name.
        create_launcher (str, optional): If set, launcher base name to write.
    """
    pixel = tomographic_objects.Pixel(
        name=os.path.join(self.pwd, name_pixel), pixel_array=coordinates_to_write
    )
    pixel.write()
    if create_launcher is not None:
        self.create_dachshund_launcher(
            np.max(coordinates_to_write[:, 0]),
            np.max(coordinates_to_write[:, 1]),
            np.max(coordinates_to_write[:, 2]),
            len(coordinates_to_write),
            properties["shape"][0],
            properties["shape"][1],
            properties["shape"][2],
            properties["sigma_f"],
            properties["lperp"],
            properties["lpar"],
            properties["name_pixel"],
            properties["name_map"],
            create_launcher,
        )

create_dachshund_launcher

create_dachshund_launcher(lx, ly, lz, npix, nx, ny, nz, sigmaf, lperp, lpar, namepixel, namemap, nameinput)

Write a Dachshund .cfg launcher describing one (sub-)map.

Parameters:

Name	Type	Description	Default
lx, ly, lz	float	Domain size in each direction (Mpc.h^-1).	required
npix	int	Total number of pixels.	required
nx, ny, nz	int	Map grid sizes.	required
sigmaf	float	Signal covariance amplitude sigma_f^2.	required
lperp	float	Transverse correlation length.	required
lpar	float	Radial correlation length.	required
namepixel	str	Input pixel binary path.	required
namemap	str	Output map binary path.	required
nameinput	str	Launcher base name (.cfg appended).	required

Source code in lelantos/cosmology.py

def create_dachshund_launcher(
    self,
    lx,
    ly,
    lz,
    npix,
    nx,
    ny,
    nz,
    sigmaf,
    lperp,
    lpar,
    namepixel,
    namemap,
    nameinput,
):
    """Write a Dachshund ``.cfg`` launcher describing one (sub-)map.

    Args:
        lx, ly, lz (float): Domain size in each direction (Mpc.h^-1).
        npix (int): Total number of pixels.
        nx, ny, nz (int): Map grid sizes.
        sigmaf (float): Signal covariance amplitude sigma_f^2.
        lperp (float): Transverse correlation length.
        lpar (float): Radial correlation length.
        namepixel (str): Input pixel binary path.
        namemap (str): Output map binary path.
        nameinput (str): Launcher base name (``.cfg`` appended).
    """
    f = open(os.path.join(self.pwd, f"{nameinput}.cfg"), "w")
    f.write("#lx, ly, lz: the domain size in each direction.\n")
    f.write("#num_pixels: the *total* number of pixels.\n")
    f.write(
        "#map_nx, map_ny, map_nz: the number of map points. The map points are arbitrary but for now these n's are used to setup a uniform grid across the domain given above.\n"
    )
    f.write("#corr_var_s: the signal cov prefactor sigma_f^2\n")
    f.write("#corr_l_perp: the signal cov perp scale.\n")
    f.write("#corr_l_para: the signal cov para scale.\n")
    f.write(
        "#pcg_max_iter: the PCG max number of iterations. 100 should be good.\n"
    )
    f.write(
        "#pcg_tol: the PCG stopping tolerance. I found 1.0e-3 is good enough. Set it very small if you want the most accurate map.\n"
    )
    f.write("lx = {}\n".format(lx))
    f.write("ly = {}\n".format(ly))
    f.write("lz = {}\n".format(lz))
    f.write("\n")
    f.write("# From output of GEN_DACH_INPUT.PRO\n")
    f.write("num_pixels = {}\n".format(npix))
    f.write("\n")
    f.write("map_nx = {}\n".format(nx))
    f.write("map_ny = {}\n".format(ny))
    f.write("map_nz = {}\n".format(nz))
    f.write("\n")
    f.write("corr_var_s = {}\n".format(sigmaf))
    f.write("corr_l_perp = {}\n".format(lperp))
    f.write("corr_l_para = {}\n".format(lpar))
    f.write("\n")
    f.write("pcg_max_iter = 500\n")
    f.write("pcg_tol = 1.0e-3\n")
    f.write("#pcg_step_r = 1\n")
    f.write("\n")
    f.write("option_map_covar = 0\n")
    f.write("option_noise_covar = 0\n")
    f.write("pixel_data_path = {}\n".format(namepixel))
    f.write("map_path = {}\n".format(namemap))
    f.close()

create_dachshund_map_pixel_property_file

create_dachshund_map_pixel_property_file(name_out, cartesian_coordinates, sky_coordinates, shape, properties_map_pixels)

Build the map pixel-property object from the conversion metadata.

Parameters:

Name	Type	Description	Default
name_out	str	Output property file name.	required
cartesian_coordinates	ndarray	Cartesian pixel coordinates.	required
sky_coordinates	ndarray	Sky pixel coordinates.	required
shape	tuple[int]	Map pixel shape.	required
properties_map_pixels	dict	Cartesian/sky bounds of the box.	required

Returns:

Type	Description
	tomographic_objects.MapPixelProperty: The property object (unwritten).

Source code in lelantos/cosmology.py

def create_dachshund_map_pixel_property_file(
    self,
    name_out,
    cartesian_coordinates,
    sky_coordinates,
    shape,
    properties_map_pixels,
):
    """Build the map pixel-property object from the conversion metadata.

    Args:
        name_out (str): Output property file name.
        cartesian_coordinates (numpy.ndarray): Cartesian pixel coordinates.
        sky_coordinates (numpy.ndarray): Sky pixel coordinates.
        shape (tuple[int]): Map pixel shape.
        properties_map_pixels (dict): Cartesian/sky bounds of the box.

    Returns:
        tomographic_objects.MapPixelProperty: The property object (unwritten).
    """
    size = (
        np.max(cartesian_coordinates[:, 0]),
        np.max(cartesian_coordinates[:, 1]),
        np.max(cartesian_coordinates[:, 2]),
    )
    coordinate_transform = self.coordinate_transform
    boundary_cartesian_coord = (
        (
            properties_map_pixels["minx"],
            properties_map_pixels["miny"],
            properties_map_pixels["minz"],
        ),
        (
            properties_map_pixels["maxx"],
            properties_map_pixels["maxy"],
            properties_map_pixels["maxz"],
        ),
    )
    boundary_sky_coord = (
        (
            properties_map_pixels["minra"],
            properties_map_pixels["mindec"],
            properties_map_pixels["minredshift"],
        ),
        (
            properties_map_pixels["maxra"],
            properties_map_pixels["maxdec"],
            properties_map_pixels["maxredshift"],
        ),
    )
    property_file = tomographic_objects.MapPixelProperty(
        name=os.path.join(self.pwd, name_out),
        size=size,
        shape=shape,
        boundary_cartesian_coord=boundary_cartesian_coord,
        boundary_sky_coord=boundary_sky_coord,
        coordinate_transform=coordinate_transform,
        Omega_m=self.Omega_m,
    )
    return property_file

create_serial_input

create_serial_input(nameout, properties, cartesian_deltas, sky_deltas)

Write the single-box (serial) solver input and optional sky pixels.

Parameters:

Name	Type	Description	Default
nameout	str	Launcher base name.	required
properties	dict	Solver/geometry parameters.	required
cartesian_deltas	ndarray	Cartesian pixels.	required
sky_deltas	ndarray	Sky pixels (written if configured).	required

Returns:

Type	Description
	tuple[int]: The map pixel shape.

Source code in lelantos/cosmology.py

def create_serial_input(self, nameout, properties, cartesian_deltas, sky_deltas):
    """Write the single-box (serial) solver input and optional sky pixels.

    Args:
        nameout (str): Launcher base name.
        properties (dict): Solver/geometry parameters.
        cartesian_deltas (numpy.ndarray): Cartesian pixels.
        sky_deltas (numpy.ndarray): Sky pixels (written if configured).

    Returns:
        tuple[int]: The map pixel shape.
    """
    self.create_input_files(
        cartesian_deltas,
        properties,
        properties["name_pixel"],
        create_launcher=nameout,
    )
    if self.return_sky_catalogs:
        self.create_input_files(
            sky_deltas,
            properties,
            "{}_sky_coordinates".format(properties["name_pixel"]),
            create_launcher=None,
        )
    return properties["shape"]

cut_in_chunks

cut_in_chunks(cartesian_deltas, number_chunks, overlaping, shape_sub_map)

Split the cartesian pixels into overlapping transverse chunks.

Parameters:

Name	Type	Description	Default
cartesian_deltas	ndarray	Cartesian pixels + values.	required
number_chunks	tuple[int]	Chunk grid (nx, ny).	required
overlaping	float	Chunk overlap (Mpc.h^-1).	required
shape_sub_map	tuple[int]	Pixel shape of one sub-map.	required

Returns:

Name	Type	Description
tuple		(Chunks, (shape_x, shape_y)) — per-chunk pixel arrays and
		limits, and the merged transverse pixel shape.

Source code in lelantos/cosmology.py

def cut_in_chunks(self, cartesian_deltas, number_chunks, overlaping, shape_sub_map):
    """Split the cartesian pixels into overlapping transverse chunks.

    Args:
        cartesian_deltas (numpy.ndarray): Cartesian pixels + values.
        number_chunks (tuple[int]): Chunk grid ``(nx, ny)``.
        overlaping (float): Chunk overlap (Mpc.h^-1).
        shape_sub_map (tuple[int]): Pixel shape of one sub-map.

    Returns:
        tuple: ``(Chunks, (shape_x, shape_y))`` — per-chunk pixel arrays and
        limits, and the merged transverse pixel shape.
    """
    if overlaping is None:
        overlaping = 0.0
    minx, maxx = np.min(cartesian_deltas[:, 0]), np.max(cartesian_deltas[:, 0])
    miny, maxy = np.min(cartesian_deltas[:, 1]), np.max(cartesian_deltas[:, 1])
    minz, maxz = np.min(cartesian_deltas[:, 2]), np.max(cartesian_deltas[:, 2])
    intervalx = maxx - minx
    intervaly = maxy - miny
    intervalz = maxz - minz
    subIntervalx = intervalx / number_chunks[0]
    subIntervaly = intervaly / number_chunks[1]
    Chunks = {}
    shape_x = number_chunks[0] * shape_sub_map[0]
    shape_y = number_chunks[1] * shape_sub_map[1]
    remove_shape_x, remove_shape_y = 0, 0
    for i in range(number_chunks[0]):
        for j in range(number_chunks[1]):
            filename = f"{i:03d}" + f"{j:03d}"
            Chunks[filename] = {}
            if (i == number_chunks[0] - 1) & (i == 0):
                intervalxChunk = [i * subIntervalx, (i + 1) * subIntervalx]
            elif i == 0:
                intervalxChunk = [
                    i * subIntervalx,
                    (i + 1) * subIntervalx + overlaping,
                ]
            elif i == number_chunks[0] - 1:
                intervalxChunk = [i * subIntervalx - overlaping, intervalx]
            else:
                intervalxChunk = [
                    i * subIntervalx - overlaping,
                    (i + 1) * subIntervalx + overlaping,
                ]
            if (j == number_chunks[1] - 1) & (j == 0):
                intervalyChunk = [j * subIntervaly, (j + 1) * subIntervaly]
            elif j == 0:
                intervalyChunk = [
                    j * subIntervaly,
                    (j + 1) * subIntervaly + overlaping,
                ]
            elif j == number_chunks[1] - 1:
                intervalyChunk = [j * subIntervaly - overlaping, intervaly]
            else:
                intervalyChunk = [
                    j * subIntervaly - overlaping,
                    (j + 1) * subIntervaly + overlaping,
                ]
            mask = (cartesian_deltas[:, 0] < intervalxChunk[1]) & (
                cartesian_deltas[:, 0] >= intervalxChunk[0]
            )
            mask &= (cartesian_deltas[:, 1] < intervalyChunk[1]) & (
                cartesian_deltas[:, 1] >= intervalyChunk[0]
            )
            chunks_deltas = []
            chunks_deltas.append(cartesian_deltas[:, 0][mask] - intervalxChunk[0])
            chunks_deltas.append(cartesian_deltas[:, 1][mask] - intervalyChunk[0])
            chunks_deltas.append(cartesian_deltas[:, 2][mask])
            chunks_deltas.append(cartesian_deltas[:, 3][mask])
            chunks_deltas.append(cartesian_deltas[:, 4][mask])
            chunks_deltas = np.transpose(np.stack(chunks_deltas))
            Chunks[filename]["coord"] = chunks_deltas
            Chunks[filename]["limits"] = [
                intervalxChunk[0],
                intervalxChunk[1],
                intervalyChunk[0],
                intervalyChunk[1],
                np.min(cartesian_deltas[:, 2]),
                np.max(cartesian_deltas[:, 2]),
            ]
            size = (
                intervalxChunk[1] - intervalxChunk[0],
                intervalyChunk[1] - intervalyChunk[0],
                intervalz,
            )
            pixel_to_remove = np.around(
                utils.pixel_per_mpc(size, shape_sub_map) * overlaping, 0
            ).astype(int)
            if number_chunks[0] != 1:
                if (i == 0) | (i == number_chunks[0] - 1):
                    remove_shape_x = remove_shape_x + pixel_to_remove[0]
                else:
                    remove_shape_x = remove_shape_x + 2 * pixel_to_remove[0]
            if number_chunks[1] != 1:
                if (j == 0) | (j == number_chunks[1] - 1):
                    remove_shape_y = remove_shape_y + pixel_to_remove[1]
                else:
                    remove_shape_y = remove_shape_y + 2 * pixel_to_remove[1]
    shape_x = shape_x - remove_shape_x // number_chunks[1]
    shape_y = shape_y - remove_shape_y // number_chunks[0]
    Chunks["overlaping"] = overlaping
    return (Chunks, (shape_x, shape_y))

create_parallel_input

create_parallel_input(properties, cartesian_deltas, number_chunks, overlaping, shape_sub_map)

Build the per-chunk launcher parameters for a tiled solver run.

Parameters:

Name	Type	Description	Default
properties	dict	Solver parameters (sigma_f, l_perp, l_par ...).	required
cartesian_deltas	ndarray	Cartesian pixels + values.	required
number_chunks	tuple[int]	Chunk grid.	required
overlaping	float	Chunk overlap (Mpc.h^-1).	required
shape_sub_map	tuple[int]	Pixel shape of one sub-map.	required

Returns:

Name	Type	Description
tuple		(parallel_launcher_params, filename, chunks, shape).

Source code in lelantos/cosmology.py

def create_parallel_input(
    self, properties, cartesian_deltas, number_chunks, overlaping, shape_sub_map
):
    """Build the per-chunk launcher parameters for a tiled solver run.

    Args:
        properties (dict): Solver parameters (sigma_f, l_perp, l_par ...).
        cartesian_deltas (numpy.ndarray): Cartesian pixels + values.
        number_chunks (tuple[int]): Chunk grid.
        overlaping (float): Chunk overlap (Mpc.h^-1).
        shape_sub_map (tuple[int]): Pixel shape of one sub-map.

    Returns:
        tuple: ``(parallel_launcher_params, filename, chunks, shape)``.
    """
    chunks, shape = self.cut_in_chunks(
        cartesian_deltas, number_chunks, overlaping, shape_sub_map
    )
    shape = (shape[0], shape[1], shape_sub_map[2])
    filename = []
    parallel_launcher_params = []
    for i in range(len(list(chunks.keys()))):
        key = list(chunks.keys())[i]
        if key != "overlaping":
            filename.append(key)
            parallel_launcher_params.append({})
            parallel_launcher_params[i]["maxx"] = chunks[key]["limits"][1]
            parallel_launcher_params[i]["minx"] = chunks[key]["limits"][0]
            parallel_launcher_params[i]["maxy"] = chunks[key]["limits"][3]
            parallel_launcher_params[i]["miny"] = chunks[key]["limits"][2]
            parallel_launcher_params[i]["maxz"] = chunks[key]["limits"][5]
            parallel_launcher_params[i]["minz"] = chunks[key]["limits"][4]
            parallel_launcher_params[i]["lx"] = (
                chunks[key]["limits"][1] - chunks[key]["limits"][0]
            )
            parallel_launcher_params[i]["ly"] = (
                chunks[key]["limits"][3] - chunks[key]["limits"][2]
            )
            parallel_launcher_params[i]["lz"] = (
                chunks[key]["limits"][5] - chunks[key]["limits"][4]
            )
            parallel_launcher_params[i]["npix"] = len(chunks[key]["coord"])
            parallel_launcher_params[i]["nx"] = shape_sub_map[0]
            parallel_launcher_params[i]["ny"] = shape_sub_map[1]
            parallel_launcher_params[i]["nz"] = shape_sub_map[2]
            parallel_launcher_params[i]["sigmaf"] = properties["sigma_f"]
            parallel_launcher_params[i]["lperp"] = properties["lperp"]
            parallel_launcher_params[i]["lpar"] = properties["lpar"]
            parallel_launcher_params[i]["namepixel"] = "{}_{}".format(
                properties["name_pixel"], key
            )
            parallel_launcher_params[i]["namemap"] = "map_{}_{}".format(
                properties["name_pixel"], key
            )
            parallel_launcher_params[i]["nameinput"] = "input_{}.cfg".format(key)
    return (parallel_launcher_params, filename, chunks, shape)

write_parallel_input

write_parallel_input(cartesian_deltas, parallel_launcher_params, filename, chunks, properties, nameout, number_chunks, overlaping)

Write the full and per-chunk pixel files and the launch pickle.

Parameters:

Name	Type	Description	Default
cartesian_deltas	ndarray	Full cartesian pixel set.	required
parallel_launcher_params	list[dict]	Per-chunk solver params.	required
filename	list[str]	Per-chunk keys.	required
chunks	dict	Per-chunk pixel arrays.	required
properties	dict	Solver parameters.	required
nameout	str	Base name for the launch pickle.	required
number_chunks	tuple[int]	Chunk grid.	required
overlaping	float	Chunk overlap (Mpc.h^-1).	required

Source code in lelantos/cosmology.py

def write_parallel_input(
    self,
    cartesian_deltas,
    parallel_launcher_params,
    filename,
    chunks,
    properties,
    nameout,
    number_chunks,
    overlaping,
):
    """Write the full and per-chunk pixel files and the launch pickle.

    Args:
        cartesian_deltas (numpy.ndarray): Full cartesian pixel set.
        parallel_launcher_params (list[dict]): Per-chunk solver params.
        filename (list[str]): Per-chunk keys.
        chunks (dict): Per-chunk pixel arrays.
        properties (dict): Solver parameters.
        nameout (str): Base name for the launch pickle.
        number_chunks (tuple[int]): Chunk grid.
        overlaping (float): Chunk overlap (Mpc.h^-1).
    """
    self.create_input_files(
        cartesian_deltas, properties, properties["name_pixel"], create_launcher=None
    )
    for i in range(len(list(chunks.keys()))):
        key = list(chunks.keys())[i]
        if key != "overlaping":
            self.create_input_files(
                chunks[key]["coord"],
                properties,
                f"{properties['name_pixel']}_{key}",
                create_launcher=None,
            )
    pickle.dump(
        [filename, parallel_launcher_params, number_chunks, overlaping],
        open(os.path.join(self.pwd, f"{nameout}_launch_data.pickle"), "wb"),
    )

create_additional_catalogs

create_additional_catalogs(cartesian_qso_catalog, cartesian_dla_catalog, sky_qso_catalog, sky_dla_catalog, properties_map_pixels)

Build and write the QSO and DLA catalogs (cartesian and sky).

Parameters:

Name	Type	Description	Default
cartesian_qso_catalog	ndarray | None	Cartesian QSO array.	required
cartesian_dla_catalog	ndarray | None	Cartesian DLA array.	required
sky_qso_catalog	ndarray | None	Sky QSO array.	required
sky_dla_catalog	ndarray | None	Sky DLA array.	required
properties_map_pixels	dict	Cartesian/sky bounds of the box.	required

Source code in lelantos/cosmology.py

def create_additional_catalogs(
    self,
    cartesian_qso_catalog,
    cartesian_dla_catalog,
    sky_qso_catalog,
    sky_dla_catalog,
    properties_map_pixels,
):
    """Build and write the QSO and DLA catalogs (cartesian and sky).

    Args:
        cartesian_qso_catalog (numpy.ndarray | None): Cartesian QSO array.
        cartesian_dla_catalog (numpy.ndarray | None): Cartesian DLA array.
        sky_qso_catalog (numpy.ndarray | None): Sky QSO array.
        sky_dla_catalog (numpy.ndarray | None): Sky DLA array.
        properties_map_pixels (dict): Cartesian/sky bounds of the box.
    """
    boundary_cartesian_coord = (
        (
            properties_map_pixels["minx"],
            properties_map_pixels["miny"],
            properties_map_pixels["minz"],
        ),
        (
            properties_map_pixels["maxx"],
            properties_map_pixels["maxy"],
            properties_map_pixels["maxz"],
        ),
    )
    boundary_sky_coord = (
        (
            properties_map_pixels["minra"],
            properties_map_pixels["mindec"],
            properties_map_pixels["minredshift"],
        ),
        (
            properties_map_pixels["maxra"],
            properties_map_pixels["maxdec"],
            properties_map_pixels["maxredshift"],
        ),
    )
    if self.return_dla_catalog is not None:
        dla_catalog_cartesian = (
            tomographic_objects.DLACatalog.init_from_pixel_catalog(
                cartesian_dla_catalog,
                name=os.path.join(self.pwd, self.return_dla_catalog),
                coordinate_transform=self.coordinate_transform,
                Omega_m=self.Omega_m,
                boundary_cartesian_coord=boundary_cartesian_coord,
                boundary_sky_coord=boundary_sky_coord,
            )
        )
        dla_catalog_cartesian.write()
        if self.return_sky_catalogs:
            dla_catalog_sky = (
                tomographic_objects.DLACatalog.init_from_pixel_catalog(
                    sky_dla_catalog,
                    name=os.path.join(
                        self.pwd, f"{self.return_dla_catalog}_sky_coordinates"
                    ),
                    coordinate_transform=self.coordinate_transform,
                    Omega_m=self.Omega_m,
                    boundary_cartesian_coord=boundary_cartesian_coord,
                    boundary_sky_coord=boundary_sky_coord,
                )
            )
            dla_catalog_sky.write()
    if self.return_qso_catalog is not None:
        quasar_catalog_cartesian = (
            tomographic_objects.QSOCatalog.init_from_pixel_catalog(
                cartesian_qso_catalog,
                name=os.path.join(self.pwd, self.return_qso_catalog),
                coordinate_transform=self.coordinate_transform,
                Omega_m=self.Omega_m,
                boundary_cartesian_coord=boundary_cartesian_coord,
                boundary_sky_coord=boundary_sky_coord,
            )
        )
        quasar_catalog_cartesian.write()
        if self.return_sky_catalogs:
            quasar_catalog_sky = (
                tomographic_objects.QSOCatalog.init_from_pixel_catalog(
                    sky_qso_catalog,
                    name=os.path.join(
                        self.pwd, f"{self.return_qso_catalog}_sky_coordinates"
                    ),
                    coordinate_transform=self.coordinate_transform,
                    Omega_m=self.Omega_m,
                    catalog_type="sky",
                    boundary_cartesian_coord=boundary_cartesian_coord,
                    boundary_sky_coord=boundary_sky_coord,
                )
            )
            quasar_catalog_sky.write()

write_additional_catalogs

write_additional_catalogs(dla_catalog_sky, dla_catalog_cartesian, quasar_catalog_sky, quasar_catalog_cartesian)

Write any provided QSO/DLA catalog objects to disk.

Parameters:

Name	Type	Description	Default
dla_catalog_sky		Sky DLA catalog object or None.	required
dla_catalog_cartesian		Cartesian DLA catalog object or None.	required
quasar_catalog_sky		Sky QSO catalog object or None.	required
quasar_catalog_cartesian		Cartesian QSO catalog object or None.	required

Source code in lelantos/cosmology.py

def write_additional_catalogs(
    self,
    dla_catalog_sky,
    dla_catalog_cartesian,
    quasar_catalog_sky,
    quasar_catalog_cartesian,
):
    """Write any provided QSO/DLA catalog objects to disk.

    Args:
        dla_catalog_sky: Sky DLA catalog object or None.
        dla_catalog_cartesian: Cartesian DLA catalog object or None.
        quasar_catalog_sky: Sky QSO catalog object or None.
        quasar_catalog_cartesian: Cartesian QSO catalog object or None.
    """
    if dla_catalog_sky is not None:
        dla_catalog_sky.write()
    if dla_catalog_cartesian is not None:
        dla_catalog_cartesian.write()
    if quasar_catalog_sky is not None:
        quasar_catalog_sky.write()
    if quasar_catalog_cartesian is not None:
        quasar_catalog_cartesian.write()

transform_delta

transform_delta(mode, nameout, properties, property_file_name, rebin=False, sigma_min=None, sigma_max=None, z_cut_min=None, z_cut_max=None, dec_cut_min=None, dec_cut_max=None, ra_cut_min=None, ra_cut_max=None, number_chunks=None, overlaping=None, shape_sub_map=None)

Full delta->solver-input conversion (main entry point).

Reads and cuts the deltas, writes the serial or parallel solver input, writes the map property file and QSO/DLA catalogs, and optionally produces the LOS density/separation diagnostic plots.

Parameters:

Name	Type	Description	Default
mode	str	"serial" or "parallel".	required
nameout	str	Run base name (launcher / launch pickle).	required
properties	dict	Solver parameters.	required
property_file_name	str	Output map property file name.	required
rebin	int	Rebin factor along each LOS.	False
sigma_min, sigma_max	float	Clip the pixel noise.	required
z_cut_min, z_cut_max	float	Redshift window.	required
dec_cut_min, dec_cut_max, ra_cut_min, ra_cut_max	float	Sky window.	required
number_chunks	tuple[int]	Chunk grid (parallel mode).	None
overlaping	float	Chunk overlap (parallel mode).	None
shape_sub_map	tuple[int]	Sub-map shape (parallel mode).	None

Raises:

Type	Description
KeyError	If mode is neither "serial" nor "parallel".

Source code in lelantos/cosmology.py

def transform_delta(
    self,
    mode,
    nameout,
    properties,
    property_file_name,
    rebin=False,
    sigma_min=None,
    sigma_max=None,
    z_cut_min=None,
    z_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    ra_cut_min=None,
    ra_cut_max=None,
    number_chunks=None,
    overlaping=None,
    shape_sub_map=None,
):
    """Full delta->solver-input conversion (main entry point).

    Reads and cuts the deltas, writes the serial or parallel solver input,
    writes the map property file and QSO/DLA catalogs, and optionally
    produces the LOS density/separation diagnostic plots.

    Args:
        mode (str): ``"serial"`` or ``"parallel"``.
        nameout (str): Run base name (launcher / launch pickle).
        properties (dict): Solver parameters.
        property_file_name (str): Output map property file name.
        rebin (int, optional): Rebin factor along each LOS.
        sigma_min, sigma_max (float, optional): Clip the pixel noise.
        z_cut_min, z_cut_max (float, optional): Redshift window.
        dec_cut_min, dec_cut_max, ra_cut_min, ra_cut_max (float, optional):
            Sky window.
        number_chunks (tuple[int], optional): Chunk grid (parallel mode).
        overlaping (float, optional): Chunk overlap (parallel mode).
        shape_sub_map (tuple[int], optional): Sub-map shape (parallel mode).

    Raises:
        KeyError: If ``mode`` is neither ``"serial"`` nor ``"parallel"``.
    """
    (
        cartesian_deltas,
        cartesian_qso_catalog,
        cartesian_dla_catalog,
        sky_deltas,
        sky_qso_catalog,
        sky_dla_catalog,
        properties_map_pixels,
    ) = self.transform_delta_to_pixel_file(
        rebin=rebin,
        sigma_min=sigma_min,
        sigma_max=sigma_max,
        z_cut_min=z_cut_min,
        z_cut_max=z_cut_max,
        dec_cut_min=dec_cut_min,
        dec_cut_max=dec_cut_max,
        ra_cut_min=ra_cut_min,
        ra_cut_max=ra_cut_max,
    )
    if mode.lower() == "serial":
        shape = self.create_serial_input(
            nameout, properties, cartesian_deltas, sky_deltas
        )
    elif mode.lower() == "parallel":
        (
            parallel_launcher_params,
            filename,
            chunks,
            shape,
        ) = self.create_parallel_input(
            properties, cartesian_deltas, number_chunks, overlaping, shape_sub_map
        )
        self.write_parallel_input(
            cartesian_deltas,
            parallel_launcher_params,
            filename,
            chunks,
            properties,
            nameout,
            number_chunks,
            overlaping,
        )
    else:
        raise KeyError("Please choose a mode between serial and parallel")
    property_file = self.create_dachshund_map_pixel_property_file(
        property_file_name,
        cartesian_deltas,
        sky_deltas,
        shape,
        properties_map_pixels,
    )
    property_file.write()
    self.create_additional_catalogs(
        cartesian_qso_catalog,
        cartesian_dla_catalog,
        sky_qso_catalog,
        sky_dla_catalog,
        properties_map_pixels,
    )
    if self.plot_pixel_properties:
        pixel_analyzer = PixelAnalizer(
            pixel=os.path.join(self.pwd, properties["name_pixel"]),
            property_file=os.path.join(self.pwd, property_file_name),
        )
        pixel_analyzer.analyze_pixels(
            False,
            True,
            name_dperp=os.path.join(self.pwd, nameout),
            coupled_plot=True,
        )

PixelAnalizer

PixelAnalizer(pixel=None, property_file=None)

Bases: object

Analyse a pixel file: line-of-sight density and mean separation.

Computes and plots the redshift dependence of the mean nearest-LOS separation and the LOS surface density, and can compare two datasets.

Load the pixel object to analyse.

Parameters:

Name	Type	Description	Default
pixel	str | Pixel	Pixel file path or a Pixel object.	None
property_file	str	Property file (if pixel is a path).	None

Source code in lelantos/cosmology.py

def __init__(self, pixel=None, property_file=None):
    """Load the pixel object to analyse.

    Args:
        pixel (str | Pixel, optional): Pixel file path or a Pixel object.
        property_file (str, optional): Property file (if ``pixel`` is a path).
    """
    if type(pixel) == str:
        pixel_class = tomographic_objects.Pixel.init_from_property_files(
            property_file, name=pixel
        )
        pixel_class.read()
    else:
        pixel_class = pixel
    self.pixel = pixel_class

write_density_file staticmethod

write_density_file(z, dperp, name)

Pickle a (z, density) profile to disk.

Parameters:

Name	Type	Description	Default
z	array - like	Redshift bin centres.	required
dperp	array - like	Profile values (density here).	required
name	str	Output pickle file.	required

Source code in lelantos/cosmology.py

@staticmethod
def write_density_file(z, dperp, name):
    """Pickle a ``(z, density)`` profile to disk.

    Args:
        z (array-like): Redshift bin centres.
        dperp (array-like): Profile values (density here).
        name (str): Output pickle file.
    """
    pickle.dump([z, dperp], open(name, "wb"))

write_dperp_file staticmethod

write_dperp_file(z, dperp, name)

Pickle a (z, mean-separation) profile to disk.

Parameters:

Name	Type	Description	Default
z	array - like	Redshift bin centres.	required
dperp	array - like	Mean nearest-LOS separation per bin.	required
name	str	Output pickle file.	required

Source code in lelantos/cosmology.py

@staticmethod
def write_dperp_file(z, dperp, name):
    """Pickle a ``(z, mean-separation)`` profile to disk.

    Args:
        z (array-like): Redshift bin centres.
        dperp (array-like): Mean nearest-LOS separation per bin.
        name (str): Output pickle file.
    """
    pickle.dump([z, dperp], open(name, "wb"))

read_density_file staticmethod

read_density_file(name)

Read a pickled (z, density) profile.

Parameters:

Name	Type	Description	Default
name	str	Pickle file.	required

Returns:

Name	Type	Description
tuple		(z, density).

Source code in lelantos/cosmology.py

@staticmethod
def read_density_file(name):
    """Read a pickled ``(z, density)`` profile.

    Args:
        name (str): Pickle file.

    Returns:
        tuple: ``(z, density)``.
    """
    a = pickle.load(open(name, "rb"))
    return (a[0], a[1])

read_dperp_file staticmethod

read_dperp_file(name)

Read a pickled (z, mean-separation) profile.

Parameters:

Name	Type	Description	Default
name	str	Pickle file.	required

Returns:

Name	Type	Description
tuple		(z, dperp).

Source code in lelantos/cosmology.py

@staticmethod
def read_dperp_file(name):
    """Read a pickled ``(z, mean-separation)`` profile.

    Args:
        name (str): Pickle file.

    Returns:
        tuple: ``(z, dperp)``.
    """
    a = pickle.load(open(name, "rb"))
    return (a[0], a[1])

plot_histogram_mean_distance staticmethod

plot_histogram_mean_distance(zpar, dmin, name_histo, nb_bins=50)

Plot the histogram of nearest-LOS distances at a fixed redshift.

Parameters:

Name	Type	Description	Default
zpar	float	Redshift of the slab.	required
dmin	array - like	Nearest-LOS distances.	required
name_histo	str	Output figure base name.	required
nb_bins	int	Number of histogram bins.	50

Source code in lelantos/cosmology.py

@staticmethod
def plot_histogram_mean_distance(zpar, dmin, name_histo, nb_bins=50):
    """Plot the histogram of nearest-LOS distances at a fixed redshift.

    Args:
        zpar (float): Redshift of the slab.
        dmin (array-like): Nearest-LOS distances.
        name_histo (str): Output figure base name.
        nb_bins (int, optional): Number of histogram bins.
    """
    plt.figure()
    plt.hist(dmin, nb_bins)
    plt.xlabel("minimal distance histogram at Z={}".format(zpar))
    plt.savefig(f"{name_histo}_at_Z{zpar}.pdf", format="pdf")

plot_mean_distance_density staticmethod

plot_mean_distance_density(zpar, dperpz, densityz, nameout, coupled_plot=False, comparison=False, dperp_comparison=None, density_comparison=None, zpar_comparison=None, legend=None, dperp_other=None, density_other=None, **kwargs)

Plot mean LOS separation and density vs redshift (coupled or not).

Parameters:

Name	Type	Description	Default
zpar	array - like	Redshift bin centres.	required
dperpz	array - like	Mean nearest-LOS separation per bin.	required
densityz	array - like	LOS density per bin.	required
nameout	str	Output figure base name.	required
coupled_plot	bool	Draw both curves on twin axes.	False
comparison	bool	Overlay a comparison dataset.	False
dperp_comparison, density_comparison, zpar_comparison	optional	Comparison profiles.	required
legend	list[str]	Legend labels.	None
dperp_other, density_other	optional	Extra comparison profiles.	required
**kwargs		Styling options.	{}

Source code in lelantos/cosmology.py

@staticmethod
def plot_mean_distance_density(
    zpar,
    dperpz,
    densityz,
    nameout,
    coupled_plot=False,
    comparison=False,
    dperp_comparison=None,
    density_comparison=None,
    zpar_comparison=None,
    legend=None,
    dperp_other=None,
    density_other=None,
    **kwargs,
):
    """Plot mean LOS separation and density vs redshift (coupled or not).

    Args:
        zpar (array-like): Redshift bin centres.
        dperpz (array-like): Mean nearest-LOS separation per bin.
        densityz (array-like): LOS density per bin.
        nameout (str): Output figure base name.
        coupled_plot (bool, optional): Draw both curves on twin axes.
        comparison (bool, optional): Overlay a comparison dataset.
        dperp_comparison, density_comparison, zpar_comparison (optional):
            Comparison profiles.
        legend (list[str], optional): Legend labels.
        dperp_other, density_other (optional): Extra comparison profiles.
        **kwargs: Styling options.
    """
    if coupled_plot:
        PixelAnalizer.plot_mean_distance_density_coupled(
            zpar,
            dperpz,
            densityz,
            nameout,
            comparison=comparison,
            dperp_comparison=dperp_comparison,
            density_comparison=density_comparison,
            zpar_comparison=zpar_comparison,
            legend=legend,
            dperp_other=dperp_other,
            density_other=density_other,
            **kwargs,
        )
    else:
        PixelAnalizer.plot_mean_distance_density_not_coupled(
            zpar,
            dperpz,
            densityz,
            nameout,
            comparison=comparison,
            dperp_comparison=dperp_comparison,
            density_comparison=density_comparison,
            zpar_comparison=zpar_comparison,
            legend=legend,
            dperp_other=dperp_other,
            density_other=density_other,
            **kwargs,
        )

plot_mean_distance_density_not_coupled staticmethod

plot_mean_distance_density_not_coupled(zpar, dperpz, densityz, nameout, comparison=False, dperp_comparison=None, density_comparison=None, zpar_comparison=None, legend=None, dperp_other=None, density_other=None, **kwargs)

Plot LOS separation and density vs redshift on separate figures.

Parameters:

Name	Type	Description	Default
zpar	array - like	Redshift bin centres.	required
dperpz	array - like	Mean nearest-LOS separation per bin.	required
densityz	array - like	LOS density per bin.	required
nameout	str	Output figure base name.	required
comparison	bool	Overlay a comparison dataset.	False
dperp_comparison, density_comparison, zpar_comparison	optional	Comparison profiles.	required
legend	list[str]	Legend labels.	None
dperp_other, density_other	optional	Extra comparison profiles.	required
**kwargs		Styling options.	{}

Source code in lelantos/cosmology.py

@staticmethod
def plot_mean_distance_density_not_coupled(
    zpar,
    dperpz,
    densityz,
    nameout,
    comparison=False,
    dperp_comparison=None,
    density_comparison=None,
    zpar_comparison=None,
    legend=None,
    dperp_other=None,
    density_other=None,
    **kwargs,
):
    """Plot LOS separation and density vs redshift on separate figures.

    Args:
        zpar (array-like): Redshift bin centres.
        dperpz (array-like): Mean nearest-LOS separation per bin.
        densityz (array-like): LOS density per bin.
        nameout (str): Output figure base name.
        comparison (bool, optional): Overlay a comparison dataset.
        dperp_comparison, density_comparison, zpar_comparison (optional):
            Comparison profiles.
        legend (list[str], optional): Legend labels.
        dperp_other, density_other (optional): Extra comparison profiles.
        **kwargs: Styling options.
    """
    plt.figure()
    plt.xlabel("Redshift")
    plt.ylabel("Mean los separation [" + r"$\mathrm{Mpc\cdot h^{-1}}$" + "]")
    plt.grid()
    plt.plot(zpar, dperpz)
    if comparison:
        plt.plot(zpar_comparison, dperp_comparison)
        plt.legend(legend)
        if dperp_other is not None:
            for i in range(len(dperp_other)):
                plt.plot(dperp_other[i][0], dperp_other[i][1])
    plt.savefig(f"{nameout}_separation.pdf", format="pdf")

    plt.figure()
    plt.xlabel("Redshift")
    plt.ylabel("Density [" + r"$\mathrm{deg^{-2}}$" + "]")
    plt.grid()
    plt.plot(zpar, densityz)
    if comparison:
        plt.plot(zpar_comparison, density_comparison)
        plt.legend(legend)
        if density_other is not None:
            for i in range(len(density_other)):
                plt.plot(density_other[i][0], density_other[i][1])
    plt.savefig(f"{nameout}_density.pdf", format="pdf")

plot_mean_distance_density_coupled staticmethod

plot_mean_distance_density_coupled(zpar, dperpz, densityz, nameout, comparison=False, dperp_comparison=None, density_comparison=None, zpar_comparison=None, legend=None, dperp_other=None, density_other=None, **kwargs)

Plot LOS separation and density vs redshift on shared twin axes.

Parameters:

Name	Type	Description	Default
zpar	array - like	Redshift bin centres.	required
dperpz	array - like	Mean nearest-LOS separation per bin.	required
densityz	array - like	LOS density per bin.	required
nameout	str	Output figure base name.	required
comparison	bool	Overlay a comparison dataset.	False
dperp_comparison, density_comparison, zpar_comparison	optional	Comparison profiles.	required
legend	list[str]	Legend labels.	None
dperp_other, density_other	optional	Extra comparison profiles.	required
**kwargs		Styling options.	{}

Source code in lelantos/cosmology.py

@staticmethod
def plot_mean_distance_density_coupled(
    zpar,
    dperpz,
    densityz,
    nameout,
    comparison=False,
    dperp_comparison=None,
    density_comparison=None,
    zpar_comparison=None,
    legend=None,
    dperp_other=None,
    density_other=None,
    **kwargs,
):
    """Plot LOS separation and density vs redshift on shared twin axes.

    Args:
        zpar (array-like): Redshift bin centres.
        dperpz (array-like): Mean nearest-LOS separation per bin.
        densityz (array-like): LOS density per bin.
        nameout (str): Output figure base name.
        comparison (bool, optional): Overlay a comparison dataset.
        dperp_comparison, density_comparison, zpar_comparison (optional):
            Comparison profiles.
        legend (list[str], optional): Legend labels.
        dperp_other, density_other (optional): Extra comparison profiles.
        **kwargs: Styling options.
    """
    figsize = utils.return_key(kwargs, "figsize", (8, 6))
    grid = utils.return_key(kwargs, "grid", True)
    fontsize = utils.return_key(kwargs, "fontsize", 13)
    fontsize_scale = utils.return_key(kwargs, "fontscalesize", 13)
    ylabel1 = utils.return_key(
        kwargs,
        "ylabel1",
        "Mean separation between nearest los ["
        + r"$h^{-1}$"
        + r"$\cdot$"
        + "Mpc"
        + "]",
    )
    ylabel2 = utils.return_key(
        kwargs, "ylabel2", "Density [" + r"$\mathrm{deg}^{-2}$" + "]"
    )
    xlabel = utils.return_key(kwargs, "xlabel", r"Redshift $z$")

    fig, ax1 = plt.subplots(1, 1, figsize=figsize)
    if grid:
        ax1.grid()
    line = ["dotted", "dashdot", "densely dashdotdotted"]

    color = "C0"
    ax1.set_xlabel(xlabel, fontsize=fontsize)
    ax1.set_ylabel(ylabel1, color=color, fontsize=fontsize)
    ax1.tick_params(axis="x", labelsize=fontsize_scale)
    ax1.tick_params(axis="y", labelsize=fontsize_scale)
    ax1.plot(zpar, dperpz, color=color)
    if comparison:
        ax1.plot(zpar_comparison, dperp_comparison, color=color, linestyle="--")
        if dperp_other is not None:
            for i in range(len(dperp_other)):
                ax1.plot(
                    dperp_other[i][0],
                    dperp_other[i][1],
                    color=color,
                    linestyle=line[i],
                )
    ax1.tick_params(axis="y", labelcolor=color)

    ax2 = ax1.twinx()
    if grid:
        ax2.grid()
    color = "C1"
    ax2.set_ylabel(ylabel2, color=color, fontsize=fontsize)
    ax2.tick_params(axis="y", labelsize=fontsize_scale)
    ax2.plot(zpar, densityz, color=color)
    if comparison:
        ax2.plot(zpar_comparison, density_comparison, color=color, linestyle="--")
        if density_other is not None:
            for i in range(len(density_other)):
                ax2.plot(
                    density_other[i][0],
                    density_other[i][1],
                    color=color,
                    linestyle=line[i],
                )
    ax2.tick_params(axis="y", labelcolor=color)

    fig.tight_layout()
    if comparison:
        legend_elements = [
            Line2D([0], [0], color="k", lw=1, label=legend[0]),
            Line2D([0], [0], color="k", linestyle="--", lw=1, label=legend[1]),
        ]
        if dperp_other is not None:
            for i in range(len(dperp_other)):
                legend_elements.append(
                    Line2D(
                        [0],
                        [0],
                        color="k",
                        linestyle=line[i],
                        lw=1,
                        label=legend[i + 2],
                    )
                )

        ax1.legend(handles=legend_elements, loc="upper center")
    plt.savefig(f"{nameout}_separation_density.pdf", format="pdf")

plot_histo_mean_distance_comparison staticmethod

plot_histo_mean_distance_comparison(dpername1, dpername2, densityname1, densityname2, nameout, legend, coupled_plot=False, **kwargs)

Compare the separation/density profiles of two datasets from files.

Parameters:

Name	Type	Description	Default
dpername1, dpername2	str | list[str]	Separation profile file(s) for datasets 1 and 2 (averaged if lists).	required
densityname1, densityname2	str | list[str]	Density profile file(s).	required
nameout	str	Output figure base name.	required
legend	list[str]	Legend labels.	required
coupled_plot	bool	Draw both curves on twin axes.	False
**kwargs		Styling options.	{}

Source code in lelantos/cosmology.py

@staticmethod
def plot_histo_mean_distance_comparison(
    dpername1,
    dpername2,
    densityname1,
    densityname2,
    nameout,
    legend,
    coupled_plot=False,
    **kwargs,
):
    """Compare the separation/density profiles of two datasets from files.

    Args:
        dpername1, dpername2 (str | list[str]): Separation profile file(s)
            for datasets 1 and 2 (averaged if lists).
        densityname1, densityname2 (str | list[str]): Density profile file(s).
        nameout (str): Output figure base name.
        legend (list[str]): Legend labels.
        coupled_plot (bool, optional): Draw both curves on twin axes.
        **kwargs: Styling options.
    """
    if type(dpername1) is str:
        zpar, dperp = PixelAnalizer.read_dperp_file(dpername1)
        zpar_comparison, dperp_comparison = PixelAnalizer.read_dperp_file(dpername2)
        zpar, density = PixelAnalizer.read_density_file(densityname1)
        zpar_comparison, density_comparison = PixelAnalizer.read_density_file(
            densityname2
        )
    else:
        zpar, dperp = np.mean(
            [
                PixelAnalizer.read_dperp_file(dpername1[i])
                for i in range(len(dpername1))
            ],
            axis=0,
        )
        zpar_comparison, dperp_comparison = np.mean(
            [
                PixelAnalizer.read_dperp_file(dpername2[i])
                for i in range(len(dpername2))
            ],
            axis=0,
        )
        zpar, density = np.mean(
            [
                PixelAnalizer.read_density_file(densityname1[i])
                for i in range(len(densityname1))
            ],
            axis=0,
        )
        zpar_comparison, density_comparison = np.mean(
            [
                PixelAnalizer.read_density_file(densityname2[i])
                for i in range(len(densityname2))
            ],
            axis=0,
        )
    PixelAnalizer.plot_mean_distance_density(
        zpar,
        dperp,
        density,
        nameout,
        coupled_plot=coupled_plot,
        comparison=True,
        dperp_comparison=dperp_comparison,
        density_comparison=density_comparison,
        zpar_comparison=zpar_comparison,
        legend=legend,
        **kwargs,
    )

compute_plot_histo_mean_distance

compute_plot_histo_mean_distance(zpar, name_histo)

Compute and plot the nearest-LOS distance histogram at a redshift.

Parameters:

Name	Type	Description	Default
zpar	float	Redshift of the slab.	required
name_histo	str	Output figure base name.	required

Source code in lelantos/cosmology.py

def compute_plot_histo_mean_distance(self, zpar, name_histo):
    """Compute and plot the nearest-LOS distance histogram at a redshift.

    Args:
        zpar (float): Redshift of the slab.
        name_histo (str): Output figure base name.
    """
    dmin = self.pixel.compute_mean_distance_histogram(zpar)
    PixelAnalizer.plot_histogram_mean_distance(zpar, dmin, name_histo, nb_bins=50)

compute_plot_mean_distance_density

compute_plot_mean_distance_density(nameout, coupled=False, plot=True)

Compute (and optionally plot/save) the separation/density profiles.

Parameters:

Name	Type	Description	Default
nameout	str	Output base name for the figures/pickles.	required
coupled	bool	Coupled twin-axes plot.	False
plot	bool	Write the profiles and figures.	True

Returns:

Name	Type	Description
tuple		(zpar, dperpz, densityz).

Source code in lelantos/cosmology.py

def compute_plot_mean_distance_density(self, nameout, coupled=False, plot=True):
    """Compute (and optionally plot/save) the separation/density profiles.

    Args:
        nameout (str): Output base name for the figures/pickles.
        coupled (bool, optional): Coupled twin-axes plot.
        plot (bool, optional): Write the profiles and figures.

    Returns:
        tuple: ``(zpar, dperpz, densityz)``.
    """
    (zpar, dperpz, densityz) = self.pixel.compute_mean_distance_density()
    if plot:
        PixelAnalizer.write_dperp_file(zpar, dperpz, f"{nameout}_dperp_file.pickle")
        PixelAnalizer.write_density_file(
            zpar, densityz, f"{nameout}_density_file.pickle"
        )
        PixelAnalizer.plot_mean_distance_density(
            zpar, dperpz, densityz, nameout, coupled_plot=coupled
        )
    return (zpar, dperpz, densityz)

analyze_pixels

analyze_pixels(compute_histo, compute_mean_distance_density, histo_zpar=None, name_histo='histogram_dperp', name_dperp='density_mean_distance', coupled_plot=False)

Run the requested pixel analyses (histogram and/or profiles).

Parameters:

Name	Type	Description	Default
compute_histo	bool	Compute the nearest-LOS distance histogram.	required
compute_mean_distance_density	bool	Compute the profiles.	required
histo_zpar	float	Redshift for the histogram.	None
name_histo	str	Histogram figure base name.	'histogram_dperp'
name_dperp	str	Profile figure base name.	'density_mean_distance'
coupled_plot	bool	Coupled twin-axes profile plot.	False

Source code in lelantos/cosmology.py

def analyze_pixels(
    self,
    compute_histo,
    compute_mean_distance_density,
    histo_zpar=None,
    name_histo="histogram_dperp",
    name_dperp="density_mean_distance",
    coupled_plot=False,
):
    """Run the requested pixel analyses (histogram and/or profiles).

    Args:
        compute_histo (bool): Compute the nearest-LOS distance histogram.
        compute_mean_distance_density (bool): Compute the profiles.
        histo_zpar (float, optional): Redshift for the histogram.
        name_histo (str, optional): Histogram figure base name.
        name_dperp (str, optional): Profile figure base name.
        coupled_plot (bool, optional): Coupled twin-axes profile plot.
    """
    if compute_histo:
        self.compute_plot_histo_mean_distance(histo_zpar, name_histo)
    if compute_mean_distance_density:
        self.compute_plot_mean_distance_density(name_dperp, coupled=coupled_plot)

DeltaAnalyzer

DeltaAnalyzer(pwd, delta_path, center_ra=True, z_cut_min=None, z_cut_max=None, dec_cut_min=None, dec_cut_max=None, ra_cut_min=None, ra_cut_max=None, degree=True, pk1d_type=True)

Bases: object

Diagnostic plots of picca delta files (histograms, sky maps, trends).

Reads the deltas within a sky/redshift window and plots the distribution and redshift dependence of the delta, sigma, SNR and redshift quantities, plus RA/Dec and LOS-density diagrams, with optional dataset comparisons.

Store the delta path, sky/redshift window and plotting options.

Parameters:

Name	Type	Description	Default
pwd	str	Output directory for the figures.	required
delta_path	str	Input delta directory.	required
center_ra	bool	Recenter RA around 0.	True
z_cut_min, z_cut_max	float	Redshift window.	required
dec_cut_min, dec_cut_max, ra_cut_min, ra_cut_max	float	Sky window.	required
degree	bool	Report/plot angles in degrees.	True
pk1d_type	bool	Whether the deltas are in pk1d format.	True

Source code in lelantos/cosmology.py

def __init__(
    self,
    pwd,
    delta_path,
    center_ra=True,
    z_cut_min=None,
    z_cut_max=None,
    dec_cut_min=None,
    dec_cut_max=None,
    ra_cut_min=None,
    ra_cut_max=None,
    degree=True,
    pk1d_type=True,
):
    """Store the delta path, sky/redshift window and plotting options.

    Args:
        pwd (str): Output directory for the figures.
        delta_path (str): Input delta directory.
        center_ra (bool, optional): Recenter RA around 0.
        z_cut_min, z_cut_max (float, optional): Redshift window.
        dec_cut_min, dec_cut_max, ra_cut_min, ra_cut_max (float, optional):
            Sky window.
        degree (bool, optional): Report/plot angles in degrees.
        pk1d_type (bool, optional): Whether the deltas are in pk1d format.
    """
    self.pwd = pwd
    self.delta_path = delta_path
    self.center_ra = center_ra
    self.z_cut_min = z_cut_min
    self.z_cut_max = z_cut_max
    self.dec_cut_min = dec_cut_min
    self.dec_cut_max = dec_cut_max
    self.ra_cut_min = ra_cut_min
    self.ra_cut_max = ra_cut_max
    self.degree = degree
    self.pk1d_type = pk1d_type

get_ra_dec

get_ra_dec(delta_path)

Obtain arrays of RA and DEC coordinates from a list or a name of a delta file in pickle, fits or ascii format

Source code in lelantos/cosmology.py

def get_ra_dec(self, delta_path):
    """Obtain arrays of RA and DEC coordinates from a list or a name of a delta file in pickle, fits or ascii format"""
    namefile = get_delta_list(delta_path)
    # namefile = preselect_deltas(namefile,
    #                             ramin=self.ra_cut_min,
    #                             ramax=self.ra_cut_max,
    #                             decmin=self.dec_cut_min,
    #                             decmax=self.dec_cut_max)
    (ra, dec, z, zqso, ids, sigmas, deltas) = get_deltas(
        namefile, center_ra=self.center_ra, pk1d_type=self.pk1d_type
    )
    pixel_coord = np.array(
        [
            [ra[i], dec[i], z[i][j], sigmas[i][j], deltas[i][j], zqso[i], ids[i]]
            for i in range(len(ra))
            for j in range(len(z[i]))
        ]
    )
    pixel_coord = pixel_coord[
        utils.cut_sky_catalog(
            pixel_coord[:, 0],
            pixel_coord[:, 1],
            pixel_coord[:, 2],
            ramin=self.ra_cut_min,
            ramax=self.ra_cut_max,
            decmin=self.dec_cut_min,
            decmax=self.dec_cut_max,
            zmin=self.z_cut_min,
            zmax=self.z_cut_max,
        )
    ]
    (redshift, redshift_qso, id, sigma, delta) = (
        pixel_coord[:, 2],
        pixel_coord[:, 5],
        pixel_coord[:, 6],
        pixel_coord[:, 3],
        pixel_coord[:, 4],
    )
    unique_coord = np.unique(pixel_coord[:, 0:2], axis=0)
    ra = unique_coord[:, 0]
    dec = unique_coord[:, 1]
    if self.degree:
        ra, dec = np.degrees(ra), np.degrees(dec)
    snr = np.abs((delta + 1) / sigma)
    return (ra, dec, redshift, redshift_qso, id, sigma, delta, snr)

get_ra_dec_comparison

get_ra_dec_comparison(comparison)

Read the plotting quantities for each comparison delta set.

Parameters:

Name	Type	Description	Default
comparison	list[str] | None	Comparison delta paths.	required

Returns:

Name	Type	Description
tuple		Per-comparison lists of ``(ra, dec, redshift, redshift_qso,
		id, sigma, delta, snr)(all None ifcomparison`` is None).

Source code in lelantos/cosmology.py

def get_ra_dec_comparison(self, comparison):
    """Read the plotting quantities for each comparison delta set.

    Args:
        comparison (list[str] | None): Comparison delta paths.

    Returns:
        tuple: Per-comparison lists of ``(ra, dec, redshift, redshift_qso,
        id, sigma, delta, snr)`` (all None if ``comparison`` is None).
    """
    if comparison is None:
        return (None, None, None, None, None, None, None, None)
    (
        ra_comp,
        dec_comp,
        redshift_comp,
        redshift_qso_comp,
        id_comp,
        sigma_comp,
        delta_comp,
        snr_comp,
    ) = ([], [], [], [], [], [], [], [])
    for i in range(len(comparison)):
        (ra, dec, redshift, redshift_qso, id, sigma, delta, snr) = self.get_ra_dec(
            comparison[i]
        )
        ra_comp.append(ra)
        dec_comp.append(dec)
        redshift_comp.append(redshift)
        redshift_qso_comp.append(redshift_qso)
        id_comp.append(id)
        sigma_comp.append(sigma)
        delta_comp.append(delta)
        snr_comp.append(snr)
    return (
        ra_comp,
        dec_comp,
        redshift_comp,
        redshift_qso_comp,
        id_comp,
        sigma_comp,
        delta_comp,
        snr_comp,
    )

plot

plot(value_names, name, comparison=None, comparison_legend=None, histo=True, mean_z_dependence=True, z_dependence=True, ra_dec_plots=True, print_stats=False, **kwargs)

Produce the requested delta diagnostic plots for several quantities.

Parameters:

Name	Type	Description	Default
value_names	list[str]	Quantities to plot (delta, sigma, snr, redshift ...).	required
name	str	Output figure base name.	required
comparison	list[str]	Comparison delta paths.	None
comparison_legend	list[str]	Legend labels.	None
histo	bool	Draw histograms.	True
mean_z_dependence	bool	Draw mean-vs-redshift plots.	True
z_dependence	bool	Draw value-vs-redshift plots.	True
ra_dec_plots	bool	Draw the RA/Dec and LOS-density plots.	True
print_stats	bool	Print comparison statistics.	False
**kwargs		Styling options.	{}

Source code in lelantos/cosmology.py

def plot(
    self,
    value_names,
    name,
    comparison=None,
    comparison_legend=None,
    histo=True,
    mean_z_dependence=True,
    z_dependence=True,
    ra_dec_plots=True,
    print_stats=False,
    **kwargs,
):
    """Produce the requested delta diagnostic plots for several quantities.

    Args:
        value_names (list[str]): Quantities to plot (``delta``, ``sigma``,
            ``snr``, ``redshift`` ...).
        name (str): Output figure base name.
        comparison (list[str], optional): Comparison delta paths.
        comparison_legend (list[str], optional): Legend labels.
        histo (bool, optional): Draw histograms.
        mean_z_dependence (bool, optional): Draw mean-vs-redshift plots.
        z_dependence (bool, optional): Draw value-vs-redshift plots.
        ra_dec_plots (bool, optional): Draw the RA/Dec and LOS-density plots.
        print_stats (bool, optional): Print comparison statistics.
        **kwargs: Styling options.
    """
    style = utils.return_key(kwargs, "style", None)
    if style is not None:
        plt.style.use(style)

    (ra, dec, redshift, redshift_qso, id, sigma, delta, snr) = self.get_ra_dec(
        self.delta_path
    )
    (
        ra_comp,
        dec_comp,
        redshift_comp,
        redshift_qso_comp,
        id_comp,
        sigma_comp,
        delta_comp,
        snr_comp,
    ) = self.get_ra_dec_comparison(comparison)

    for value_name in value_names:
        value = locals()[value_name]
        comparison_value = locals()[value_name + "_comp"]
        lambda_rest = utils.return_key(kwargs, f"{value_name}_lambda_rest", False)
        if lambda_rest:
            kwargs[f"{value_name}_redshift_qso"] = redshift_qso
        if histo:
            utils.save_histo(
                self.pwd,
                value,
                value_name,
                name,
                comparison=comparison_value,
                comparison_legend=comparison_legend,
                **kwargs,
            )
        if (mean_z_dependence) & (value_name not in ["redshift", "ra", "dec"]):
            utils.save_mean_redshift_dependence(
                self.pwd,
                value,
                redshift,
                value_name,
                name,
                comparison=comparison_value,
                comparison_redshift=redshift_comp,
                comparison_legend=None,
                **kwargs,
            )

        if (z_dependence) & (value_name not in ["redshift", "ra", "dec"]):
            utils.save_redshift_dependence(
                self.pwd,
                value,
                redshift,
                value_name,
                name,
                comparison=comparison_value,
                comparison_redshift=redshift_comp,
                comparison_legend=None,
                **kwargs,
            )
    if ra_dec_plots:
        self.plot_ra_dec(
            ra, dec, name, comparison=None, comparison_legend=None, **kwargs
        )
    if (comparison is not None) & (print_stats):
        self.stat_comparison(
            redshift, redshift_comp, sigma, sigma_comp, delta, delta_comp
        )

plot_ra_dec

plot_ra_dec(ra, dec, name, comparison_ra=None, comparison_dec=None, comparison_legend=None, **kwargs)

Plot the RA/Dec diagram and the LOS-density-vs-RA curve.

Parameters:

Name	Type	Description	Default
ra, dec	array - like	Object sky angles.	required
name	str	Output figure base name.	required
comparison_ra, comparison_dec	list[array]	Comparison sky coordinates.	required
comparison_legend	list[str]	Legend labels.	None
**kwargs		Styling options.	{}

Source code in lelantos/cosmology.py

def plot_ra_dec(
    self,
    ra,
    dec,
    name,
    comparison_ra=None,
    comparison_dec=None,
    comparison_legend=None,
    **kwargs,
):
    """Plot the RA/Dec diagram and the LOS-density-vs-RA curve.

    Args:
        ra, dec (array-like): Object sky angles.
        name (str): Output figure base name.
        comparison_ra, comparison_dec (list[array], optional): Comparison
            sky coordinates.
        comparison_legend (list[str], optional): Legend labels.
        **kwargs: Styling options.
    """
    utils.save_ra_dec(
        self.pwd,
        ra,
        dec,
        name,
        comparison_ra=comparison_ra,
        comparison_dec=comparison_dec,
        comparison_legend=comparison_legend,
        **kwargs,
    )

    DeltaAnalyzer.plot_los_density(self.pwd, ra, dec, name, **kwargs)

stat_comparison

stat_comparison(redshift, redshift_comp, sigma, sigma_comp, delta, delta_comp)

Print redshift-range, sigma and delta statistics for main vs comps.

Parameters:

Name	Type	Description	Default
redshift, sigma, delta	array - like	Main-dataset quantities.	required
redshift_comp, sigma_comp, delta_comp	list[array]	Per-comparison quantities.	required

Source code in lelantos/cosmology.py

def stat_comparison(
    self, redshift, redshift_comp, sigma, sigma_comp, delta, delta_comp
):
    """Print redshift-range, sigma and delta statistics for main vs comps.

    Args:
        redshift, sigma, delta (array-like): Main-dataset quantities.
        redshift_comp, sigma_comp, delta_comp (list[array]): Per-comparison
            quantities.
    """
    print("Redshift interval =", np.max(redshift), np.min(redshift))
    for i in range(len(redshift_comp)):
        print(
            f"Redshift interval for comp {i} =",
            np.max(redshift_comp[i]),
            np.min(redshift_comp[i]),
        )
    print("Maximal sigma =", np.max(sigma))
    for i in range(len(sigma_comp)):
        print(f"Maximal sigma for comp {i} =", np.max(sigma_comp[i]))
    print("Mean sigma =", np.mean(sigma))
    for i in range(len(sigma_comp)):
        print(f"Mean sigma for comp {i} =", np.mean(sigma_comp[i]))
    print("Median sigma =", np.median(sigma))
    for i in range(len(sigma_comp)):
        print(f"Median sigma for comp {i} =", np.median(sigma_comp[i]))
    print("Mean delta =", np.mean(delta))
    for i in range(len(delta_comp)):
        print(f"Mean delta for comp {i} =", np.mean(delta_comp[i]))

plot_los_density staticmethod

plot_los_density(pwd, ra, dec, plot_name, **kwargs)

Plot the line-of-sight surface density as a function of RA.

Parameters:

Name	Type	Description	Default
pwd	str	Output directory.	required
ra, dec	array - like	Object sky angles.	required
plot_name	str	Output figure base name.	required
**kwargs		nb_interval (RA bins) and different_sign_region (split by Dec sign) options.	{}

Source code in lelantos/cosmology.py

@staticmethod
def plot_los_density(pwd, ra, dec, plot_name, **kwargs):
    """Plot the line-of-sight surface density as a function of RA.

    Args:
        pwd (str): Output directory.
        ra, dec (array-like): Object sky angles.
        plot_name (str): Output figure base name.
        **kwargs: ``nb_interval`` (RA bins) and ``different_sign_region``
            (split by Dec sign) options.
    """
    nb_interval = utils.return_key(kwargs, "nb_interval", 20)
    different_sign_region = utils.return_key(kwargs, "different_sign_region", False)

    ra_interval = np.linspace(np.min(ra), np.max(ra), nb_interval)
    ra_size = abs((np.max(ra) - np.min(ra)) / nb_interval)
    ra_array = []
    for i in range(len(ra_interval) - 1):
        ra_array.append((ra_interval[i] + ra_interval[i + 1]) / 2)
    ra_array = np.asarray(ra_array)
    density_array = np.zeros(ra_array.shape)
    density_array_plus = np.zeros(ra_array.shape)
    density_array_minus = np.zeros(ra_array.shape)
    for i in range(len(ra_interval) - 1):
        mask = (ra > ra_interval[i]) & (ra < ra_interval[i + 1])
        density_array[i] = len(ra[mask])
        density_array_plus[i] = len(ra[mask & (dec >= 0)])
        density_array_minus[i] = len(ra[mask & (dec < 0)])
    maxdec = np.max(dec)
    mindec = np.min(dec)
    plt.figure()
    plt.plot(ra_array, density_array / abs(ra_size * (maxdec - mindec)))
    plt.title("LOS density in function of RA")
    plt.grid()
    plt.savefig(os.path.join(pwd, f"{plot_name}_los_density.pdf"), format="pdf")
    if different_sign_region:
        plt.figure()
        plt.plot(ra_array, density_array_plus / abs(ra_size * maxdec))
        plt.title("LOS density in function of RA for DEC >= 0")
        plt.grid()
        plt.savefig(
            os.path.join(pwd, f"{plot_name}_los_density_dec_positive.pdf"),
            format="pdf",
        )
        plt.figure()
        plt.plot(ra_array, density_array_minus / abs(ra_size * mindec))
        plt.title("LOS density in function of RA for DEC < 0")
        plt.grid()
        plt.savefig(
            os.path.join(pwd, f"{plot_name}_los_density_dec_negative.pdf"),
            format="pdf",
        )

plot_delta_gaussian_fit staticmethod

plot_delta_gaussian_fit(pwd, delta, plot_name, **kwargs)

Plot the delta histogram with an overlaid Gaussian fit.

Parameters:

Name	Type	Description	Default
pwd	str	Output directory.	required
delta	array - like	Delta values.	required
plot_name	str	Output figure base name.	required
**kwargs		Histogram styling options.	{}

Source code in lelantos/cosmology.py

@staticmethod
def plot_delta_gaussian_fit(pwd, delta, plot_name, **kwargs):
    """Plot the delta histogram with an overlaid Gaussian fit.

    Args:
        pwd (str): Output directory.
        delta (array-like): Delta values.
        plot_name (str): Output figure base name.
        **kwargs: Histogram styling options.
    """
    (name, data, bins, patches) = utils.plot_histo(delta, "delta", "", **kwargs)
    bin_centers = np.array(
        [0.5 * (bins[i] + bins[i + 1]) for i in range(len(bins) - 1)]
    )
    fit_function = lambda x, A, mu, sigma: A * np.exp(
        -1.0 * (x - mu) ** 2 / (2 * sigma**2)
    )
    popt, pcov = curve_fit(
        fit_function, xdata=bin_centers, ydata=data, p0=[1, 0.0, 0.1]
    )
    x = np.linspace(min(bins), max(bins), 1000)
    y = fit_function(x, *popt)
    plt.plot(x, y, "r--", linewidth=2)
    mu, sigma = popt[1], popt[2]
    plt.text(0.8, 2 * np.max(data) / 6, "mu = " + str(round(mu, 8)))
    plt.text(0.8, 1.5 * np.max(data) / 6, "sigma = " + str(round(sigma, 8)))

    plt.savefig(
        os.path.join(pwd, f"{plot_name}_histo_delta_gaussian_fit.pdf"), format="pdf"
    )

get_delta_list

get_delta_list(delta_path)

List all delta-*.fits* files under one or several directories.

Parameters:

Name	Type	Description	Default
delta_path	str | list[str]	Directory, or list of directories, holding picca delta files.	required

Returns:

Type	Description
	numpy.ndarray | list: Sorted list of delta file paths.

Raises:

Type	Description
KeyError	If no delta file is found.

Source code in lelantos/cosmology.py

def get_delta_list(delta_path):
    """List all ``delta-*.fits*`` files under one or several directories.

    Args:
        delta_path (str | list[str]): Directory, or list of directories, holding
            picca delta files.

    Returns:
        numpy.ndarray | list: Sorted list of delta file paths.

    Raises:
        KeyError: If no delta file is found.
    """
    if type(delta_path) == str:
        delta_list = np.sort(glob.glob(os.path.join(delta_path, "delta-*.fits*")))
    elif type(delta_path) == list:
        delta_list = []
        for i in range(len(delta_path)):
            delta_list = delta_list + list(
                np.sort(glob.glob(os.path.join(delta_path[i], "delta-*.fits*")))
            )
    if len(delta_list) == 0:
        raise KeyError("No delta file was found")
    return delta_list

preselect_deltas

preselect_deltas(namefile, ramin=None, ramax=None, decmin=None, decmax=None, center_ra=True, pk1d_type=True)

Keep only the delta files overlapping a sky window.

Parameters:

Name	Type	Description	Default
namefile	list[str]	Delta file paths.	required
ramin, ramax, decmin, decmax	float	Sky window bounds.	required
center_ra	bool	Recenter RA around 0 before comparing.	True
pk1d_type	bool	Whether the deltas are in pk1d format.	True

Returns:

Type	Description
	list[str]: The subset of files with at least one pixel in the window.

Raises:

Type	Description
KeyError	If no file falls inside the window.

Source code in lelantos/cosmology.py

def preselect_deltas(
    namefile,
    ramin=None,
    ramax=None,
    decmin=None,
    decmax=None,
    center_ra=True,
    pk1d_type=True,
):
    """Keep only the delta files overlapping a sky window.

    Args:
        namefile (list[str]): Delta file paths.
        ramin, ramax, decmin, decmax (float, optional): Sky window bounds.
        center_ra (bool, optional): Recenter RA around 0 before comparing.
        pk1d_type (bool, optional): Whether the deltas are in pk1d format.

    Returns:
        list[str]: The subset of files with at least one pixel in the window.

    Raises:
        KeyError: If no file falls inside the window.
    """
    subset_namefile = []
    if (ramin is None) & (ramax is None) & (decmin is None) & (decmax is None):
        return namefile
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=pk1d_type)
        delta_tomo.read()
        ra, dec, redshift, redshift_qso, id, sigma, delta = delta_tomo.return_params(
            center_ra=center_ra
        )
        mask = np.full(ra.shape, True)
        if ramin is not None:
            mask &= ra > ramin
        if ramax is not None:
            mask &= ra < ramax
        if decmin is not None:
            mask &= dec > decmin
        if ramin is not None:
            mask &= dec < decmax
        if len(mask[mask]) != 0:
            subset_namefile.append(namefile[i])
    print(len(subset_namefile))
    if len(subset_namefile) == 0:
        raise KeyError("Select window does not contain any delta files")
    return subset_namefile

get_deltas

get_deltas(namefile, center_ra=True, pk1d_type=True)

Extract delta properties

Source code in lelantos/cosmology.py

def get_deltas(namefile, center_ra=True, pk1d_type=True):
    """Extract delta properties"""
    ras, decs, redshifts, redshift_qsos, ids, sigmas, deltas = (
        [],
        [],
        [],
        [],
        [],
        [],
        [],
    )
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=pk1d_type)
        delta_tomo.read()
        ra, dec, redshift, redshift_qso, id, sigma, delta = delta_tomo.return_params(
            center_ra=center_ra
        )
        ras.append(ra)
        decs.append(dec)
        redshift_qsos.append(redshift_qso)
        ids.append(id)
        redshifts = redshifts + redshift
        sigmas = sigmas + sigma
        deltas = deltas + delta
    return (
        np.concatenate(ras),
        np.concatenate(decs),
        redshifts,
        np.concatenate(redshift_qsos),
        np.concatenate(ids),
        sigmas,
        deltas,
    )

get_merged_multiple_exposure_deltas

get_merged_multiple_exposure_deltas(namefile)

Merge deltas with repeated observation

Source code in lelantos/cosmology.py

def get_merged_multiple_exposure_deltas(namefile):
    """Merge deltas with repeated observation"""
    # Pack LOS by Id in the dict Deltas
    ra, dec, z, deltas, sigmas, zqso = [], [], [], [], [], []
    (Deltas, ids) = get_id_list(namefile)

    # For each pack of LOS
    for i in range(len(ids)):
        # Get the data
        zqso.append(tomographic_objects.Delta.z_qso(Deltas[ids[i]][0]))
        ra.append(tomographic_objects.Delta.ra(Deltas[ids[i]][0]))
        dec.append(tomographic_objects.Delta.dec(Deltas[ids[i]][0]))
        listsigmas, listz, listdelta = [], [], []
        for j in range(len(Deltas[ids[i]])):
            listsigmas.append(1 / np.sqrt(np.asarray(Deltas[ids[i]][j].ivar)))
            listdelta.append(Deltas[ids[i]][j].delta)
            listz.append(
                ((10 ** np.asarray(Deltas[ids[i]][j].log_lambda) / utils.lambdaLy) - 1)
            )

        (listz, listsigmas, listdelta) = delete_los_extrema(
            listz, listsigmas, listdelta
        )

        (listz, listsigmas, listdelta, lenlists) = delete_missing_pixels(
            listz, listsigmas, listdelta
        )

        # Weighted merging of the LOSs
        zmerged, sigmamerged, deltamerged = listz[0], [], []
        for k in range(len(listz[0])):
            sigma = 0
            delta = 0
            sumsigma = 0
            for m in range(len(listz)):
                sumsigma = sumsigma + 1 / (listsigmas[m][k] ** 2)
            for m in range(len(listz)):
                delta = delta + listdelta[m][k] / ((listsigmas[m][k] ** 2) * sumsigma)

                sigma = sigma + 1 / ((listsigmas[m][k] ** 2) * (sumsigma**2))
            sigma = np.sqrt(sigma / len(listz))
            sigmamerged.append(sigma)
            deltamerged.append(delta)
        deltas.append(deltamerged)
        z.append(zmerged)
        sigmas.append(sigmamerged)
    return (
        np.array(ra),
        np.array(dec),
        np.asarray(z),
        np.array(zqso),
        np.array(ids),
        np.asarray(sigmas),
        np.asarray(deltas),
    )

get_id_list

get_id_list(namefile)

Group delta objects by their primary key (line-of-sight id).

Parameters:

Name	Type	Description	Default
namefile	list[str]	Delta file paths.	required

Returns:

Name	Type	Description
tuple		(Deltas, ids) — a dict mapping each id to its list of delta
		objects, and the list of unique ids.

Source code in lelantos/cosmology.py

def get_id_list(namefile):
    """Group delta objects by their primary key (line-of-sight id).

    Args:
        namefile (list[str]): Delta file paths.

    Returns:
        tuple: ``(Deltas, ids)`` — a dict mapping each id to its list of delta
        objects, and the list of unique ids.
    """
    ids = []
    for i in range(len(namefile)):
        delta_tomo = tomographic_objects.Delta(name=namefile[i], pk1d_type=True)
        delta_tomo.read()
        id = []
        for i in range(len(delta_tomo.delta_array)):
            id.append(tomographic_objects.Delta.primary_key(delta_tomo.delta_array[i]))
        ids.append(id)
    ids = np.concatenate(ids)
    ids = list(set(ids))
    Deltas = {}
    for i in range(len(ids)):
        Deltas[ids[i]] = []
        for j in range(len(namefile)):
            delta_tomo = tomographic_objects.Delta(name=namefile[j], pk1d_type=True)
            delta_tomo.read()
            for k in range(len(delta_tomo.delta_array)):
                if (
                    tomographic_objects.Delta.primary_key(delta_tomo.delta_array[k])
                    == ids[i]
                ):
                    Deltas[ids[i]].append(delta_tomo.delta_array[k])
    return (Deltas, ids)

delete_los_extrema

delete_los_extrema(listz, listsigmas, listdelta)

Trim repeated lines of sight to their common redshift range.

Parameters:

Name	Type	Description	Default
listz	list[array]	Per-exposure redshift arrays of one LOS.	required
listsigmas	list[array]	Matching sigma arrays.	required
listdelta	list[array]	Matching delta arrays.	required

Returns:

Name	Type	Description
tuple		(listz, listsigmas, listdelta) trimmed to the shared z range.

Source code in lelantos/cosmology.py

def delete_los_extrema(listz, listsigmas, listdelta):
    """Trim repeated lines of sight to their common redshift range.

    Args:
        listz (list[array]): Per-exposure redshift arrays of one LOS.
        listsigmas (list[array]): Matching sigma arrays.
        listdelta (list[array]): Matching delta arrays.

    Returns:
        tuple: ``(listz, listsigmas, listdelta)`` trimmed to the shared z range.
    """
    # Get the list of common elements in the pack to have minimum and maximum redshifts
    zmin, zmax = 0, 10**10
    zcommon = listz[0]
    for j in range(1, len(listz)):
        zcommon = list(set(zcommon).intersection(listz[j]))
    zmin = np.min(zcommon)
    zmax = np.max(zcommon)

    # Deleting pixels at the beginning and at the end of LOSs
    for j in range(len(listz)):
        lineToDeleteFirst = 0
        k = 0
        while listz[j][k] != zmin:
            lineToDeleteFirst = lineToDeleteFirst + 1
            k = k + 1
        lineToDeleteLast = 0
        k = -1
        while listz[j][k] != zmax:
            lineToDeleteLast = lineToDeleteLast + 1
            k = k - 1
        listz[j] = listz[j][lineToDeleteFirst : len(listz[j]) - lineToDeleteLast]
        listsigmas[j] = listsigmas[j][
            lineToDeleteFirst : len(listsigmas[j]) - lineToDeleteLast
        ]
        listdelta[j] = listdelta[j][
            lineToDeleteFirst : len(listdelta[j]) - lineToDeleteLast
        ]
    return (listz, listsigmas, listdelta)

delete_missing_pixels

delete_missing_pixels(listz, listsigmas, listdelta)

Align repeated lines of sight by dropping non-shared pixels.

Parameters:

Name	Type	Description	Default
listz	list[array]	Per-exposure redshift arrays of one LOS.	required
listsigmas	list[array]	Matching sigma arrays.	required
listdelta	list[array]	Matching delta arrays.	required

Returns:

Name	Type	Description
tuple		(listz, listsigmas, listdelta, lenlists) with all exposures
		sharing the same pixels, plus their common lengths.

Source code in lelantos/cosmology.py

def delete_missing_pixels(listz, listsigmas, listdelta):
    """Align repeated lines of sight by dropping non-shared pixels.

    Args:
        listz (list[array]): Per-exposure redshift arrays of one LOS.
        listsigmas (list[array]): Matching sigma arrays.
        listdelta (list[array]): Matching delta arrays.

    Returns:
        tuple: ``(listz, listsigmas, listdelta, lenlists)`` with all exposures
        sharing the same pixels, plus their common lengths.
    """
    # Ensuring that all LOS have the same lenght
    lenlists = []
    for j in range(len(listz)):
        lenlists.append(len(listz[j]))

    # Selection of the pixels to delete in case of a missing pixel along one LOS + Deletion
    while np.max(lenlists) != np.min(lenlists):
        eltTodelete = [[] for j in range(len(listz))]
        mi = 10**10
        mins = []
        for j in range(len(lenlists)):
            if lenlists[j] < mi:
                mi = lenlists[j]
                mins = [j]
            elif lenlists[j] == mi:
                mins.append(j)
        for j in range(len(mins)):
            for k in range(len(listz)):
                for m in range(len(listz[k])):
                    if (k != mins[j]) & (
                        np.isin([listz[k][m]], listz[mins[j]]) == False
                    ):
                        eltTodelete[k].append(m)
        newlistz, newlistsigma, newlistdelta = (
            [[] for n in range(len(listz))],
            [[] for n in range(len(listz))],
            [[] for n in range(len(listz))],
        )
        for j in range(len(eltTodelete)):
            eltTodeletej = list(set(eltTodelete[j]))
            for k in range(len(listz[j])):
                if np.isin([k], eltTodeletej) == False:
                    newlistz[j].append(listz[j][k])
                    newlistsigma[j].append(listsigmas[j][k])
                    newlistdelta[j].append(listdelta[j][k])
        listz = newlistz
        listsigmas = newlistsigma
        listdelta = newlistdelta
        lenlists = []
        for j in range(len(listz)):
            lenlists.append(len(listz[j]))
        return (listz, listsigmas, listdelta, lenlists)

compute_shape_size_parallel

compute_shape_size_parallel(extremum_coord, number_chunks, overlaping, shape_sub_map)

Compute the full-map pixel shape and physical size for a tiled run.

Accounts for the pixels shared between overlapping chunks so the merged map has a consistent resolution.

Parameters:

Name	Type	Description	Default
extremum_coord	sequence	(Xmin, Xmax, Ymin, Ymax, Zmin, Zmax).	required
number_chunks	tuple[int]	Chunk grid (nx, ny).	required
overlaping	float	Chunk overlap (Mpc.h^-1).	required
shape_sub_map	tuple[int]	Pixel shape of one sub-map.	required

Returns:

Name	Type	Description
tuple		(shape, size) of the full merged map.

Source code in lelantos/cosmology.py

def compute_shape_size_parallel(
    extremum_coord, number_chunks, overlaping, shape_sub_map
):
    """Compute the full-map pixel shape and physical size for a tiled run.

    Accounts for the pixels shared between overlapping chunks so the merged map
    has a consistent resolution.

    Args:
        extremum_coord (sequence): ``(Xmin, Xmax, Ymin, Ymax, Zmin, Zmax)``.
        number_chunks (tuple[int]): Chunk grid ``(nx, ny)``.
        overlaping (float): Chunk overlap (Mpc.h^-1).
        shape_sub_map (tuple[int]): Pixel shape of one sub-map.

    Returns:
        tuple: ``(shape, size)`` of the full merged map.
    """
    if overlaping is None:
        overlaping = 0.0
    minx, maxx = extremum_coord[0], extremum_coord[1]
    miny, maxy = extremum_coord[2], extremum_coord[3]
    minz, maxz = extremum_coord[4], extremum_coord[5]
    intervalx = maxx - minx
    intervaly = maxy - miny
    intervalz = maxz - minz
    subIntervalx = intervalx / number_chunks[0]
    subIntervaly = intervaly / number_chunks[1]
    shape_x = number_chunks[0] * shape_sub_map[0]
    shape_y = number_chunks[1] * shape_sub_map[1]
    remove_shape_x, remove_shape_y = 0, 0
    for i in range(number_chunks[0]):
        for j in range(number_chunks[1]):
            if (i == number_chunks[0] - 1) & (i == 0):
                intervalxChunk = [i * subIntervalx, (i + 1) * subIntervalx]
            elif i == 0:
                intervalxChunk = [
                    i * subIntervalx,
                    (i + 1) * subIntervalx + overlaping,
                ]
            elif i == number_chunks[0] - 1:
                intervalxChunk = [i * subIntervalx - overlaping, intervalx]
            else:
                intervalxChunk = [
                    i * subIntervalx - overlaping,
                    (i + 1) * subIntervalx + overlaping,
                ]
            if (j == number_chunks[1] - 1) & (j == 0):
                intervalyChunk = [j * subIntervaly, (j + 1) * subIntervaly]
            elif j == 0:
                intervalyChunk = [
                    j * subIntervaly,
                    (j + 1) * subIntervaly + overlaping,
                ]
            elif j == number_chunks[1] - 1:
                intervalyChunk = [j * subIntervaly - overlaping, intervaly]
            else:
                intervalyChunk = [
                    j * subIntervaly - overlaping,
                    (j + 1) * subIntervaly + overlaping,
                ]
            size = (
                intervalxChunk[1] - intervalxChunk[0],
                intervalyChunk[1] - intervalyChunk[0],
                intervalz,
            )
            pixel_to_remove = np.around(
                utils.pixel_per_mpc(size, shape_sub_map) * overlaping, 0
            ).astype(int)
            if number_chunks[0] != 1:
                if (i == 0) | (i == number_chunks[0] - 1):
                    remove_shape_x = remove_shape_x + pixel_to_remove[0]
                else:
                    remove_shape_x = remove_shape_x + 2 * pixel_to_remove[0]
            if number_chunks[1] != 1:
                if (j == 0) | (j == number_chunks[1] - 1):
                    remove_shape_y = remove_shape_y + pixel_to_remove[1]
                else:
                    remove_shape_y = remove_shape_y + 2 * pixel_to_remove[1]
    shape_x = shape_x - remove_shape_x // number_chunks[1]
    shape_y = shape_y - remove_shape_y // number_chunks[0]

    size = (maxx - minx, maxy - miny, maxz - minz)
    shape = (shape_x, shape_y, shape_sub_map[2])
    return (shape, size)

compute_shape_size_parallel_from_interface

compute_shape_size_parallel_from_interface(ramin, ramax, decmin, decmax, zmin, zmax, Omega_m, coordinate_transform, number_chunks, overlaping, shape_sub_map, N_coord_edge=100)

Estimate the map shape/size from a footprint and chunking config.

Samples the edges of the (RA, Dec, z) footprint, converts them to cartesian coordinates, and feeds the bounding box to :func:compute_shape_size_parallel.

Parameters:

Name	Type	Description	Default
ramin, ramax, decmin, decmax	float	Angular footprint (degrees).	required
zmin, zmax	float	Redshift window.	required
Omega_m	float	Fiducial matter density.	required
coordinate_transform	str	Sky<->cartesian transform mode.	required
number_chunks	tuple[int]	Chunk grid.	required
overlaping	float	Chunk overlap (Mpc.h^-1).	required
shape_sub_map	tuple[int]	Pixel shape of one sub-map.	required
N_coord_edge	int	Number of samples per footprint edge.	100

Returns:

Name	Type	Description
tuple		(shape, size) of the full merged map.

Source code in lelantos/cosmology.py

def compute_shape_size_parallel_from_interface(
    ramin,
    ramax,
    decmin,
    decmax,
    zmin,
    zmax,
    Omega_m,
    coordinate_transform,
    number_chunks,
    overlaping,
    shape_sub_map,
    N_coord_edge=100,
):
    """Estimate the map shape/size from a footprint and chunking config.

    Samples the edges of the (RA, Dec, z) footprint, converts them to cartesian
    coordinates, and feeds the bounding box to
    :func:`compute_shape_size_parallel`.

    Args:
        ramin, ramax, decmin, decmax (float): Angular footprint (degrees).
        zmin, zmax (float): Redshift window.
        Omega_m (float): Fiducial matter density.
        coordinate_transform (str): Sky<->cartesian transform mode.
        number_chunks (tuple[int]): Chunk grid.
        overlaping (float): Chunk overlap (Mpc.h^-1).
        shape_sub_map (tuple[int]): Pixel shape of one sub-map.
        N_coord_edge (int, optional): Number of samples per footprint edge.

    Returns:
        tuple: ``(shape, size)`` of the full merged map.
    """
    ra = np.linspace(ramin, ramax, N_coord_edge)
    dec = np.linspace(decmin, decmax, N_coord_edge)
    z = np.linspace(zmin, zmax, N_coord_edge)

    cube_edge_coord = np.concatenate(
        [
            np.array([[ra[i], decmin, zmin] for i in range(N_coord_edge)]),
            np.array([[ramin, dec[i], zmin] for i in range(N_coord_edge)]),
            np.array([[ramin, decmin, z[i]] for i in range(N_coord_edge)]),
            np.array([[ramax, dec[i], zmin] for i in range(N_coord_edge)]),
            np.array([[ramax, decmin, z[i]] for i in range(N_coord_edge)]),
            np.array([[ra[i], decmax, zmin] for i in range(N_coord_edge)]),
            np.array([[ramin, decmax, z[i]] for i in range(N_coord_edge)]),
            np.array([[ra[i], decmin, zmax] for i in range(N_coord_edge)]),
            np.array([[ramin, dec[i], zmax] for i in range(N_coord_edge)]),
            np.array([[ra[i], decmax, zmax] for i in range(N_coord_edge)]),
            np.array([[ramax, dec[i], zmax] for i in range(N_coord_edge)]),
            np.array([[ramax, decmax, z[i]] for i in range(N_coord_edge)]),
        ]
    )

    (rcomov, distang, _, _) = utils.get_cosmo_function(Omega_m)
    suplementary_parameters = utils.return_suplementary_parameters(
        coordinate_transform, zmin=zmin, zmax=zmax
    )

    cartesian_cube_edge_coord = np.zeros_like(cube_edge_coord)

    (
        cartesian_cube_edge_coord[:, 0],
        cartesian_cube_edge_coord[:, 1],
        cartesian_cube_edge_coord[:, 2],
    ) = utils.convert_sky_to_cartesian(
        np.radians(cube_edge_coord[:, 0]),
        np.radians(cube_edge_coord[:, 1]),
        cube_edge_coord[:, 2],
        coordinate_transform,
        rcomov=rcomov,
        distang=distang,
        suplementary_parameters=suplementary_parameters,
    )

    Xmin = np.min(cartesian_cube_edge_coord[:, 0])
    Xmax = np.max(cartesian_cube_edge_coord[:, 0])
    Ymin = np.min(cartesian_cube_edge_coord[:, 1])
    Ymax = np.max(cartesian_cube_edge_coord[:, 1])
    Zmin = np.min(cartesian_cube_edge_coord[:, 2])
    Zmax = np.max(cartesian_cube_edge_coord[:, 2])
    extremum_coord = [Xmin, Xmax, Ymin, Ymax, Zmin, Zmax]

    shape, size = compute_shape_size_parallel(
        extremum_coord,
        number_chunks,
        overlaping,
        shape_sub_map,
    )
    return shape, size