CorrHOD_cutsky#

Warning

This class is still in development and is not yet ready for use.

While the analysis code is stable, the creation of cutsky simulation is still in development.

What is a cutsky ?#

In theory#

While the cubic box is a snapshot of the Universe at a given time (i.e. a redshift), a cutsky aims to be a simulation of an acual survey. It is a simulation of the sky, with a given geometry.

Several effects are included in a cutsky, such as Redshift Space Distortions (RSD), a survey geometry (footprint), some observational contraints (like fiber assignment, incompleteness, …), etc.

As the further we look in space the further we look in time, a cutsky is not a snapshot of the Universe at a given time, as the redshift evolves with the distance. this causes the number density of galaxies to change with the redshift.

In theory, a cutsky and cubic mock of the same cosmology should have the same clustering properties.

In practice (computationaly)#

A cutsky is usually created from a cubic box, and is therefore a subset of the cubic box.

As it is not a periodic box, the computation method changes to become a bit more complicated. We will need random points matching the footprint of the survey, to be able to compute the analyis. The more points, the better the analysis will be, but the longer it will take to compute.

Some methods to handle the data and randoms catalogs can be found in the CorrHOD.catalogs module. (see Catalogs handling)

How to use the class#

from CorrHOD import CorrHOD_cutsky

Warning

The whole “Creating the simulation and turn it in a cutsky” patrt of the class is not yet implemented. This part of the documentation will be updated when the class will be ready.

However, the analysis part of the class is ready and can be used.

Tip

You can import the cubic_dict catalog in the class like this

Object = CorrHOD_cutsky(path2config) # Useless here, but we need it to create the object
Object.cutsky_dict = data_cutsky
Object.randoms_dict = randoms_cutsky

As the first steps of this class are the same as the CorrHOD_cubic class, please refer to the documentation of that class for the first steps of the analysis.

Object.get_tracer_positions()
Object.get_tracer_weights()

Note

There are actually three positionnal arrays created in the class : data_positions (resp. randoms_positions) are the positions in sky coordinates (Ra, Dec, Comoving distance), data_sky (resp. randoms_sky) are the positions in sky coordinates (Ra, Dec, Redshift), and data_cartesian (resp. randoms_cartesian) are the positions in cartesian coordinates (x, y, z).

Be careful when calling the positions in this class !

Tip

Also, the get_tracer_positions() method has a return_cartesian parameter, that allows to return the cartesian positions instead of the sky positions. This can be useful.

Weights and number density#

To normalize the power spectrum, some weights (called ‘FKP weights’ [1] ) need to be applied in the computation of the correlation functions.

They can be computed by the get_tracer_weights() method of the class, and are all set to 1 if not computed before calling the correlation functions.

Note

The densitysplit has not method to compute the weights yet, as it si not clear yet what those mean in this context.

However, the correlation functions will return wierd results if the weights of the catalogs are not evolving in the same way with the redshift. Therefore, setting all the weights to 1 for the quantiles and using the FKP weights for the data and randoms is not a good idea.

In the get_tracer_weights() method, the weights are computed for the data and randoms and splitted in the quantiles to their respectives galaxies (the same way as the densitysplit).

This is a temporary fix and doesn’t impact a lot the results, compared with setting all the weights to 1.

The number density functions can also be computed within the class, which can allow for some verifications. For example, for the densitysplit to agree with the cubic box, the number density of the quantiles should have roughly the same shape as the data and randoms catalogs.

The number density can be computed with the get_nz() method of the class.

Computing the analysis#

Computing the densitysplit and correlation functions is the same as for the CorrHOD_cubic class. (Under the hood, the process is a bit different but the use of the methods are the same).

Tip

After computing the densitysplit, if you want to downsample the data to reduce the computing time and memory usage, it is recommended to use the npoints or frac parameters of the downsample_data() method, rather than the new_n parameter, as the number density parameter is computed as the mean number density in the cutsky, and is less accurate than the one computed in the cubic box.

API#

class CorrHOD.cutsky.CorrHOD_cutsky(path2config: str, HOD_params: Optional[dict] = None, los: str = 'z', boxsize: float = 2000, cosmo: int = 0)[source][source]#

Bases: CorrHOD_cubic

This class is used to compute the 2PCF and the autocorrelation and cross-correlation of the quantiles of the DensitySplit. It takes HOD parameters and a cosmology as input and uses AbacusHOD to generate a DESI cutsky mock. It then uses DensitySplit to compute the density field and the quantiles. Finally, it uses pycorr to compute the 2PCF and the autocorrelation and cross-correlation of the quantiles.

If a cutsky file is to be provided, the cutsky must be in the form of a dictionary in sky coordinates. Every computation in this class is done in sky coordinates, except for the densitysplit. The conversion is done using the densitysplit package. Unless the name of the variable indicates it, the positions are in (ra, dec, comoving distance) coordinates. If the name of the variable contains ‘sky’, the positions are in (ra, dec, z) coordinates.

Initialises the CorrHOD object. Note that the ‘tracer’ parameter is fixed to ‘LRG’ by default. It can be changed later by manually reassigning the parameter. This tracer is chosen to study LRG and BGS. It can be changed but the code is not guaranteed to work.

Parameters:

  • path2config (str) – Path to the config file of AbacusHOD. It is used to load the simulation. See the documentation for more details.

  • HOD_params (dict, optional) – Dictionary containing the HOD parameters. The keys must be the same as the ones used in AbacusHOD. If not provided, the default parameters of AbacusHOD are used. Defaults to None.

  • los (str, optional) – Line-of-sight direction in which to apply the RSD. Defaults to ‘z’.

  • boxsize (float, optional) – Size of the simulation box in Mpc/h. Defaults to 2000.

  • cosmo (int, optional) – Index of the cosmology to use. This index is used to load the cosmology from the AbacusSummit object. Defaults to 0. See the documentation of AbacusSummit for more details.

Raises:

ValueError – If the simulation is not precomputed. See the documentation of AbacusHOD for more details.

get_tracer_positions(return_cartesian: bool = False)[source][source]#

Get the positions of the tracers (data and randoms) in cartesian coordinates. Note : For now, we assume that the RSD are already applied to the data

Returns:

  • data_positions (array_like) – The positions of the galaxies in cartesian coordinates.

  • randoms_positions (array_like) – The positions of the randoms in cartesian coordinates.

get_tracer_weights(edges: Optional[list] = None, area: Optional[float] = None, fsky: Optional[float] = None, P0: float = 7000)[source][source]#

Compute the weights for the tracer and the randoms. Note : If neither area nor fsky are set, the sky fraction of the data will be computed, using the current data sky positions.

Parameters:

  • edges (list, optional) – The edges of the bins used to compute the number density. If set to None, the edges are computed using Scott’s rule. Defaults to None.

  • area (float, optional) – The area of the survey in square degrees. Defaults to None.

  • fsky (float, optional) – The fraction of the sky covered by the survey. Defaults to None.

  • P0 (float, optional) – The power spectrum normalization (TODO : Check this definition). Defaults to 7000 for the BGS.

Returns:

  • data_weights (array_like) – The weights of the galaxies.

  • randoms_weights (array_like) – The weights of the randoms.

compute_DensitySplit(smooth_radius: float = 10, cellsize: float = 5, nquantiles: int = 10, sampling: str = 'randoms', randoms=None, filter_shape: str = 'Gaussian', return_density: bool = True, nthread=16)[source][source]#

Compute the density field and the quantiles using DensitySplit. (see epaillas/densitysplit for more details)

Parameters:

  • smooth_radius (float, optional) – The radius of the Gaussian smoothing in Mpc/h. Defaults to 10.

  • cellsize (float, optional) – The size of the cells in the mesh. Defaults to 10.

  • nquantiles (int, optional) – The number of quantiles to compute. Defaults to 10.

  • sampling (str, optional) – The sampling to use. Can be ‘randoms’ or ‘data’. Defaults to ‘randoms’.

  • randoms (np.ndarray, optional) – The positions of the randoms. It is a (N,3) array of points in the range [0, boxsize] for each coordinate. If set to none, a the randoms positions of the data will be used. Defaults to None.

  • filter_shape (str, optional) – The shape of the filter to use. Can be ‘Gaussian’ or ‘TopHat’. Defaults to ‘Gaussian’.

  • return_density (bool, optional) – If True, the density field is returned. Defaults to True.

  • nthread (int, optional) – The number of threads to use. Defaults to 16.

Returns:

  • quantiles (np.ndarray) – The quantiles. It is a (nquantiles,N,3) array.

  • density (np.ndarray, optional) – The density field. It is a (N,3) array.

get_nz(edges: Optional[list] = None, area: Optional[float] = None, fsky: Optional[float] = None)[source][source]#

Computes the n(z) function of the data and the quantiles (if they exist). Note : If neither area nor fsky are set, the sky fraction of the data will be computed, using the current data sky positions.

Parameters:

  • edges (list, optional) – The edges of the bins used to compute the number density. If set to None, the edges are computed using Scott’s rule. Defaults to None.

  • area (float, optional) – The area of the survey in square degrees. Defaults to None.

  • fsky (float, optional) – The fraction of the sky covered by the survey. Defaults to None.

Returns:

  • nz_data (InterpolatedUnivariateSpline) – The n(z) function of the data.

  • nz_randoms (InterpolatedUnivariateSpline) – The n(z) function of the randoms.

  • nz_functions (list, optional) – The n(z) function for each quantile. It is a list of InterpolatedUnivariateSpline objects.

downsample_data(frac: Optional[float] = None, new_n: Optional[float] = None, npoints: Optional[int] = None)[source][source]#

Downsamples the data to a given number density. The number of points in the randoms and the quantiles will be adjusted accordingly. (i.e. if the number density of the data is 10 times smaller, the number of points in the randoms and the quantiles will be 10 times smaller, while keeping the proportionality to the data)

Since the number density varies with the redshift, the downsampling to the new number density is done by taking the mean n(z) of the data as reference. This should have n(z) remaining approx. the same. (TODO : Test this ?)

Parameters:

new_n (float) – The wanted number density of the data after downsampling. The number density of the randoms and the quantiles will be adjusted accordingly.

Returns:

  • data_positions (np.ndarray) – The positions of the data after downsampling.

  • randoms_positions (np.ndarray) – The positions of the randoms after downsampling.

  • quantiles (np.ndarray, optional) – The quantiles after downsampling.

compute_auto_corr(quantile: int, mode: str = 'smu', edges: list = [np.linspace(0.1, 200, 50), np.linspace(-1, 1, 60)], mpicomm=None, mpiroot=None, nthread: int = 16)[source][source]#

Compute the auto-correlation of a quantile. The result will also be saved in the class, as self.CF[‘Auto’][f’DS{quantile}’]

Parameters:

  • quantile (int) – Index of the quantile to compute. It must be between 0 and nquantiles-1.

  • mode (str, optional) – The mode to use for the 2PCF. Defaults to ‘smu’. See pycorr for more details .

  • edges (list, optional) – The edges to use for the 2PCF. Defaults to [np.linspace(0.1, 200, 50), np.linspace(-1, 1, 60)]. See pycorr for more details.

  • mpicomm (_type_, optional) – The MPI communicator. Defaults to None.

  • mpiroot (_type_, optional) – The MPI root. Defaults to None.

  • nthread (int, optional) – The number of threads to use. Defaults to 16.

Returns:

xi_quantile – The auto-correlation of the quantile.

Return type:

pycorr.TwoPointCorrelationFunction

compute_cross_corr(quantile: int, mode: str = 'smu', edges: list = [np.linspace(0.1, 200, 50), np.linspace(-1, 1, 60)], mpicomm=None, mpiroot=None, nthread: int = 16)[source][source]#

Compute the cross-correlation of a quantile with the galaxies. The result will also be saved in the class, as self.CF[‘Cross’][f’DS{quantile}’]

Parameters:

  • quantile (int) – Index of the quantile to compute. It must be between 0 and nquantiles-1.

  • mode (str, optional) – The mode to use for the 2PCF. Defaults to ‘smu’. See pycorr for more details .

  • edges (list, optional) – The edges to use for the 2PCF. Defaults to [np.linspace(0.1, 200, 50), np.linspace(-1, 1, 60)]. See pycorr for more details.

  • mpicomm (_type_, optional) – The MPI communicator. Defaults to None.

  • mpiroot (_type_, optional) – The MPI root. Defaults to None.

  • nthread (int, optional) – The number of threads to use. Defaults to 16.

Returns:

xi_quantile – The cross-correlation of the quantile.

Return type:

pycorr.TwoPointCorrelationFunction

compute_2pcf(mode: str = 'smu', edges: list = [np.linspace(0.1, 200, 50), np.linspace(-1, 1, 60)], mpicomm=None, mpiroot=None, nthread: int = 16)[source][source]#

Compute the 2PCF of the galaxies. The result will also be saved in the class, as self.CF[‘2PCF’]

Parameters:

  • mode (str, optional) – The mode to use for the 2PCF. Defaults to ‘smu’. See pycorr for more details .

  • edges (list, optional) – The edges to use for the 2PCF. Defaults to [np.linspace(0.1, 200, 50), np.linspace(-1, 1, 60)]. See pycorr for more details.

  • mpicomm (_type_, optional) – The MPI communicator. Defaults to None.

  • mpiroot (_type_, optional) – The MPI root. Defaults to None.

  • nthread (int, optional) – The number of threads to use. Defaults to 16.

Returns:

xi – The 2PCF of the galaxies.

Return type:

pycorr.TwoPointCorrelationFunction

save(hod_indice: int = 0, path: Optional[str] = None, save_HOD: bool = True, save_pos: bool = True, save_density: bool = True, save_quantiles: bool = True, save_CF: bool = True, save_xi: bool = True, los: str = 'average', save_all: bool = False)[source][source]#

Save the results of the CorrHOD object. Some of the results are saved as dictionaries in numpy files. They can be accessed using the np.load function, and the .item() method of the loaded object.

Parameters:

  • hod_indice (int, optional) – The indice of the set of HOD parameters provided to the CorrHOD Object. Defaults to 0. Be careful to set the right indice if you want to save the HOD parameters !

  • path (str, optional) – The path where to save the results. If None, the path is set to the output_dir of config file. Defaults to None.

  • save_HOD (bool, optional) – If True, the HOD parameters are saved. File saved as hod{hod_indice}_c{cosmo}_p{phase}.npy. Defaults to True.

  • save_pos (bool, optional) – If True, the cubic mock dictionary is saved. File saved as pos_hod{hod_indice}_c{cosmo}_p{phase}.npy. Defaults to True.

  • save_density (bool, optional) – If True, the density PDF is saved. File saved as density_hod{hod_indice}_c{cosmo}_p{phase}.npy. Defaults to True.

  • save_quantiles (bool, optional) – If True, the quantiles of the densitysplit are saved. File saved as quantiles_hod{hod_indice}_c{cosmo}_p{phase}.npy. Defaults to True.

  • save_CF (bool, optional) – If True, the 2PCF, the autocorrelation and cross-correlation of the quantiles are saved. The 2PCF is saved as tpcf_hod{hod_indice}_c{cosmo}_p{phase}.npy. It contains the separation bins in the s key and the poles in the 2PCF key. The Auto and Corr dictionaries are saved as ds_auto_hod{hod_indice}_c{cosmo}_p{phase}.npy and ds_cross_hod{hod_indice}_c{cosmo}_p{phase}.npy. Each dictionnary contains DS{quantile} keys with the CF of the quantile. The s key contains the separation bins. Defaults to True.

  • save_xi (bool, optional) – If True, the pycorr objects are saved in an unique dictionary. The 2PCF can be accessed using the xi[los][‘2PCF’] key. The auto and cross-correlations of the quantiles can be accessed using the xi[los][‘Auto’][‘DS{i}’] and xi[los][‘Cross’][‘DS{i}’] keys. Defaults to True.

  • los (str, optional) – The line of sight along which to save the 2PCF, the autocorrelation and cross-correlation of the quantiles. Can be ‘x’, ‘y’, ‘z’ or ‘average’. Defaults to ‘average’.

  • save_all (bool, optional) – If True, all the results are saved. This overrides the other options. Defaults to False.

Footnotes