selfisys package
Subpackages
- selfisys.utils package
- Submodules
- selfisys.utils.examples_utils module
- selfisys.utils.logger module
- selfisys.utils.low_level module
- selfisys.utils.parser module
- selfisys.utils.path_utils module
- selfisys.utils.plot_examples module
- selfisys.utils.plot_params module
- selfisys.utils.plot_utils module
- selfisys.utils.timestepping module
- selfisys.utils.tools module
- selfisys.utils.workers module
- Module contents
Submodules
selfisys.global_parameters module
selfisys.grf module
- selfisys.grf.primordial_grf(L, N, seedphases, fname_powerspectrum, fname_outputinitialdensity, force_sim=False, return_g=False, verbose=0)[source]
Generate a Gaussian random field from a specified input power spectrum.
- Parameters:
L (float) – Side length of the simulation box in Mpc/h.
N (int) – Grid resolution (number of cells per dimension).
seedphases (int) – Seed for random phase generation (for reproducibility).
fname_powerspectrum (str) – File path to the input power spectrum.
fname_outputinitialdensity (str) – File path to store the generated initial density field.
force_sim (bool, optional) – If True, regenerate the GRF even if the output file exists. Default is False.
return_g (bool, optional) – If True, return the GRF as a numpy array. Default is False.
verbose (int, optional) – Verbosity level (0 = silent, 1 = progress, 2 = detailed). Default is 0.
- Raises:
OSError – If the power spectrum file cannot be read.
RuntimeError – If an unexpected error occurs during power spectrum reading.
- Returns:
The GRF data if return_g is True, otherwise None.
- Return type:
numpy.ndarray or None
selfisys.normalise_hb module
- selfisys.normalise_hb.define_normalisation(hidden_box: HiddenBox, Pbins: ndarray, cosmo: Dict[str, Any], N: int, min_k_norma: float = 0.04, npar: int = 1, force: bool = False) ndarray [source]
Define the normalisation constants for the HiddenBox instance.
- Parameters:
hidden_box (HiddenBox) – Instance of the HiddenBox class.
Pbins (ndarray) – Array of P bin values.
cosmo (dict) – Cosmological and infrastructure parameters.
N (int) – Number of realisations required.
min_k_norma (float, optional) – Minimum k value to compute the normalisation constants.
npar (int, optional) – Number of parallel processes to use. Default is 1.
force (bool, optional) – If True, force recomputation. Default is False.
- Returns:
norm_csts – Normalisation constants for the HiddenBox instance.
- Return type:
ndarray
- selfisys.normalise_hb.worker_normalisation(hidden_box: HiddenBox, params: Tuple[Dict[str, Any], list, list, bool]) ndarray [source]
Worker function to compute the normalisation constants, compatible with Python multiprocessing.
- Parameters:
hidden_box (HiddenBox) – Instance of the HiddenBox class.
params (tuple) – A tuple containing (cosmo, seedphase, seednoise, force).
- Returns:
phi – Computed summary statistics.
- Return type:
ndarray
- selfisys.normalise_hb.worker_normalisation_public(hidden_box, cosmo: Dict[str, Any], N: int, i: int)[source]
Run the i-th simulation required to compute the normalisation constants.
- Parameters:
hidden_box (HiddenBox) – Instance of the HiddenBox class.
cosmo (dict) – Cosmological and some infrastructure parameters.
N (int) – Total number of realisations required.
i (int) – Index of the simulation to be computed.
selfisys.prior module
- selfisys.prior.get_summary(x, bins, normalisation=None, kmax=1.4)[source]
Compute a power-spectrum summary for given cosmological parameters.
- Parameters:
x (array-like) – Cosmological parameters [h, Omega_b, Omega_m, n_s, sigma_8].
bins (array-like) – Wavenumber bins.
normalisation (float or None, optional) – Normalisation constant to scale the resulting spectrum.
kmax (float, optional) – Maximum wavenumber for get_Pk.
- Returns:
theta – The computed power-spectrum values, optionally normalised.
- Return type:
ndarray
- Raises:
RuntimeError – If the power-spectrum computation fails unexpectedly.
- selfisys.prior.logposterior_hyperparameters_parallel(selfi, theta_fiducial, Nbin_min, Nbin_max, theta_norm, k_corr, alpha_cv)[source]
Compute the log-posterior for the hyperparameters of the prior from [leclercq2019primordial], for use within the SelfiSys pipeline.
- Parameters:
selfi (object) – The selfi object.
theta_fiducial (ndarray) – Fiducial spectrum.
Nbin_min (int) – Minimum bin index for the wavenumber range.
Nbin_max (int) – Maximum bin index for the wavenumber range.
theta_norm (float) – Hyperparameter controlling the overall uncertainty.
k_corr (float) – Hyperparameter controlling correlation scale.
alpha_cv (float) – Cosmic variance strength.
- Returns:
The log-posterior value for the given hyperparameters.
- Return type:
float
- Raises:
RuntimeError – If the log-posterior computation fails unexpectedly.
- selfisys.prior.perform_prior_optimisation_and_plot(selfi, theta_fiducial, theta_norm_mean=0.1, theta_norm_std=0.3, k_corr_mean=0.02, k_corr_std=0.015, k_opt_min=0.0, k_opt_max=1.4, theta_norm_min=0.04, theta_norm_max=0.12, k_corr_min=0.012, k_corr_max=0.02, meshsize=30, Nbin_min=0, Nbin_max=100, theta_norm=0.05, k_corr=0.015, alpha_cv=0.00065, plot=True, savepath=None)[source]
Optimise the hyperparameters for the selfi2019 prior (from [leclercq2019primordial]).
- Parameters:
selfi (object) – The selfi object.
theta_fiducial (ndarray) – Fiducial spectrum.
theta_norm_mean (float, optional) – Mean of the Gaussian hyperprior on theta_norm. Default 0.1.
theta_norm_std (float, optional) – Standard deviation of the hyperprior on theta_norm. Default 0.3.
k_corr_mean (float, optional) – Mean of the Gaussian hyperprior on k_corr. Default 0.020.
k_corr_std (float, optional) – Standard deviation of the hyperprior on k_corr. Default 0.015.
k_opt_min (float, optional) – Minimum wavenumber for the prior optimisation. Default 0.0.
k_opt_max (float, optional) – Maximum wavenumber for the prior optimisation. Default 1.4.
theta_norm_min (float, optional) – Lower bound for theta_norm in the mesh. Default 0.04.
theta_norm_max (float, optional) – Upper bound for theta_norm in the mesh. Default 0.12.
k_corr_min (float, optional) – Lower bound for k_corr in the mesh. Default 0.012.
k_corr_max (float, optional) – Upper bound for k_corr in the mesh. Default 0.02.
meshsize (int, optional) – Number of points in each dimension of the plot mesh. Default 30.
Nbin_min (int, optional) – Minimum bin index for restricting the prior. Default 0.
Nbin_max (int, optional) – Maximum bin index for restricting the prior. Default 100.
theta_norm (float, optional) – Initial or default guess of theta_norm. Default 0.05.
k_corr (float, optional) – Initial or default guess of k_corr. Default 0.015.
alpha_cv (float, optional) – Cosmic variance term or similar. Default 0.00065.
plot (bool, optional) – If True, generate and show/save a 2D contour plot. Default True.
savepath (str, optional) – File path to save the plot. If None, the plot is displayed.
- Returns:
(theta_norm, k_corr) after optimisation.
- Return type:
tuple
- Raises:
OSError – If file operations fail during saving the prior or posterior.
RuntimeError – If the optimisation fails unexpectedly.
- class selfisys.prior.planck_prior(Omega_mean, Omega_cov, bins, normalisation, kmax, nsamples=10000, nthreads=-1, EPS_K=1e-07, EPS_residual=0.001, filename=None)[source]
Bases:
object
Custom prior for the SelfiSys pipeline. This is the prior used in [hoellinger2024diagnosing], based on the Planck 2018 cosmological parameters.
This class provides methods to compute a power-spectrum prior from a prior distribution of cosmological parameters, using a Gaussian fit. See equation (7) in [hoellinger2024diagnosing].
- Parameters:
Omega_mean (array-like) – Mean of the prior distribution on cosmological parameters.
Omega_cov (array-like) – Covariance matrix of the prior distribution on cosmological parameters.
bins (array-like) – Wavenumbers where the power spectrum is evaluated.
normalisation (float or None) – If not None, divide the power spectra by the normalisation.
kmax (float) – Maximum wavenumber for computations.
nsamples (int, optional) – Number of samples drawn from the prior on the cosmological parameters. Default is 10,000.
nthreads (int, optional) – Number of CPU threads for parallel tasks. Default is -1, that is, auto-detect the number of available threads.
EPS_K (float, optional) – Regularisation parameter for covariance inversion. Default 1e-7.
EPS_residual (float, optional) – Additional cutoff for matrix inversion. Default 1e-3.
filename (str or None, optional) – Path to a .npy file to store or load precomputed power spectra.
- mean
Mean of the computed power spectra.
- Type:
ndarray
- covariance
Covariance matrix of the computed power spectra.
- Type:
ndarray
- inv_covariance
Inverse of the covariance matrix.
- Type:
ndarray
- Raises:
OSError – If file reading or writing fails.
RuntimeError – For unexpected HPC or multi-processing errors.
- property Nbin_max
Index of the maximal wavenumber given self.kmax.
- property Nbin_min
Index of the minimal wavenumber given k_min.
- compute()[source]
Compute the prior (mean, covariance, and inverse covariance).
If self.filename exists, tries to load the prior. Otherwise, samples from the prior distribution on cosmological parameters and evaluates the power spectra in parallel.
- Raises:
OSError – If self.filename is not writable/accessible.
RuntimeError – If multi-processing or power-spectra computations fail.
- classmethod load(fname)[source]
Load the prior from input file.
- Parameters:
fname (str) – Input HDF5 filename.
- Returns:
The prior object.
- Return type:
prior
- Raises:
OSError – If the file cannot be read or is invalid.
- logpdf(theta, theta_mean, theta_covariance, theta_icov)[source]
Return the log prior probability at a given point in parameter space.
- Parameters:
theta (ndarray) – Evaluation point in parameter space.
theta_mean (ndarray) – Prior mean vector.
theta_covariance (ndarray) – Prior covariance matrix.
theta_icov (ndarray) – Inverse of the prior covariance matrix.
- Returns:
Log prior probability value.
- Return type:
float
- selfisys.prior.worker_class(params)[source]
Worker function to compute power spectra with CLASS, compatible with Python multiprocessing.
- Parameters:
params (tuple) – (x, bins, normalisation, kmax) where x is an array-like of cosmological parameters, bins is the wavenumber array, normalisation is a float or None, and kmax is a float.
- Returns:
theta – Power-spectrum summary from get_summary.
- Return type:
ndarray
selfisys.sbmy_interface module
- selfisys.sbmy_interface.compute_Phi(G_ss_path, P_ss_path, g_obj, norm, AliasingCorr=True, verbosity=1)[source]
Compute the summary statistics from a field object, based on a provided summary-statistics Fourier grid and baseline spectrum.
- Parameters:
G_ss_path (str) – Path to the FourierGrid file used for summary-statistics.
P_ss_path (str) – Path to the baseline power spectrum file for normalisation.
g_obj (Field) – Input field object from which to compute summary statistics.
norm (ndarray) – Normalisation constants for the summary statistics.
AliasingCorr (bool, optional) – Whether to apply aliasing correction. Default is True.
verbosity (int, optional) – Verbosity level (0=quiet, 1=normal, 2=debug). Default 1.
- Returns:
Phi – Vector of summary statistics.
- Return type:
ndarray
- Raises:
OSError – If file reading fails at G_ss_path or P_ss_path.
RuntimeError – If unexpected issues occur during computation.
- selfisys.sbmy_interface.generate_white_noise_Field(L, size, seedphase, fname_whitenoise, seedname_whitenoise, force_phase=False)[source]
Generate a white noise realisation in physical space and write it to disk.
- Parameters:
L (float) – Size of the simulation box (in Mpc/h).
size (int) – Number of grid points along each axis.
seedphase (int or list of int) – User-provided seed to generate the initial white noise.
fname_whitenoise (str) – File path to write the white noise realisation.
seedname_whitenoise (str) – File path to write the seed state of the RNG.
force_phase (bool, optional) – If True, forces regeneration of the random phases. Default is False.
- Raises:
OSError – If file writing fails or directory paths are invalid.
RuntimeError – For unexpected issues.
- selfisys.sbmy_interface.get_power_spectrum_from_cosmo(L, size, cosmo, fname_power_spectrum, force=False)[source]
Compute a power spectrum from cosmological parameters and save it to disk.
- Parameters:
L (float) – Size of the simulation box (in Mpc/h).
size (int) – Number of grid points along each axis.
cosmo (dict) – Cosmological parameters (and infrastructure parameters).
fname_power_spectrum (str) – Name (including path) of the power spectrum file to read/write.
force (bool, optional) – If True, forces recomputation even if the file exists. Default is False.
- Raises:
OSError – If file writing fails or the directory path is invalid.
RuntimeError – For unexpected issues during power spectrum computation.
- selfisys.sbmy_interface.handle_time_stepping(aa: List[float], total_steps: int, modeldir: str, figuresdir: str, sim_params: str, force: bool = False) Tuple[List[int] | None, float | None] [source]
Create and merge individual time-stepping objects.
- Parameters:
aa (list of float) – List of scale factors in ascending order.
total_steps (int) – Total number of time steps to distribute among the provided scale factors.
modeldir (str) – Directory path to store generated time-stepping files.
figuresdir (str) – Directory path to store time-stepping plots.
sim_params (str) – Simulation parameter string (e.g., “custom”, “std”, “nograv”).
force (bool, optional) – Whether to force recompute the time-stepping files. Default is False.
- Returns:
merged_path (str) – Path to the merged time-stepping file.
indices_steps_cumul (list of int or None) – Cumulative indices for the distributed steps. Returns None if using splitLPT or if sim_params indicates an alternative strategy.
eff_redshifts (float or None) – Effective redshift derived from the final scale factor in ‘custom’ or ‘nograv’ mode. None otherwise.
- Raises:
NotImplementedError – If a unsupported time-stepping strategy is used.
OSError – If file or directory operations fail.
RuntimeError – If unexpected issues occur during time-stepping setup.
- selfisys.sbmy_interface.setup_sbmy_parfiles(d, cosmology, file_names, hiddenbox_params, force=False)[source]
Set up Simbelmynë parameter file (please refer to the Simbelmynë documentation for more details).
- Parameters:
d (int) – Index (from 1 to S) specifying a direction in parameter space, 0 for the expansion point, or -1 for mock data.
cosmology (array, double, dimension=5) – Cosmological parameters.
file_names (dict) – Dictionary containing the names of the input/output files for the simulation.
hiddenbox_params (dict) – See the HiddenBox class for more details.
force (bool, optional, default=False) – If True, forces recompute the simulation parameter files.
selfisys.selection_functions module
- class selfisys.selection_functions.LognormalSelection(L=None, selection_params=None, survey_mask_path=None, local_select_path=None, size=None)[source]
Bases:
object
Class to generate radial selection functions.
- init_selection(reset=False)[source]
Initialise the radial selection functions.
- Parameters:
reset (bool, optional) – Whether to reset the selection function.
- lognormals_z_to_x(xx, mask, params, spline)[source]
Convert log-normal distributions from redshift to distance.
- Parameters:
xx (array-like) – Comoving distances at which to evaluate the distributions.
mask (ndarray or None) – Survey mask C(n).
params (tuple of arrays) – Parameters for the distributions (ss, mm, rr).
spline (UnivariateSpline) – Linear interpolator for the distance-redshift relation.
- Returns:
Tuple containing redshifts and list of distributions.
- Return type:
tuple
- multiple_lognormal(x, mask, ss, ll, rr)[source]
Compute multiple log-normal distributions.
- Parameters:
x (ndarray) – Input array.
mask (ndarray or None) – Survey mask C(n).
ss (array_like) – Standard deviations for each distribution.
ll (array_like) – Means for each distribution.
rr (array_like) – Rescaling factors for each distribution.
- Returns:
List of log-normal distributions.
- Return type:
list of ndarray
- multiple_lognormal_z(x, mask, ss, mm, rr)[source]
Compute multiple rescaled lognormal distributions as functions of redshift.
- Parameters:
x (ndarray) – Input array (redshifts).
mask (ndarray or None) – Survey mask C(n).
ss (array_like) – Standard deviations of the lognormal distributions.
mm (array_like) – Means of the lognormal distributions.
rr (array_like) – Rescaling factors for each distribution.
- Returns:
List of log-normal distributions.
- Return type:
list of ndarray
- static one_lognormal(x, std, mean, rescale=None)[source]
Rescaled log-normal distribution.
- Parameters:
x (ndarray) – Input array.
std (float) – Standard deviation of the distribution.
mean (float) – Mean of the distribution.
rescale (float, optional) – Rescaling factor. If None, the distribution is normalised such that its maximum value is 1.
- Returns:
Log-normal distribution evaluated at x.
- Return type:
ndarray
- static one_lognormal_z(x, sig2, mu, rescale=None)[source]
Compute a log-normal distribution in redshift.
- Parameters:
x (ndarray) – Input array.
sig2 (float) – Variance of the distribution.
mu (float) – Mean of the distribution.
rescale (float, optional) – Rescaling factor.
- Returns:
Log-normal distribution evaluated at x.
- Return type:
ndarray
selfisys.selfi_interface module
- selfisys.selfi_interface.PrintMessage(required_verbosity: int, message: str, verbosity: int) None [source]
Print a message to standard output using pyselfi.utils.PrintMessage.
- Parameters:
required_verbosity (int) – The verbosity level required to display the message.
message (str) – The actual message to display.
verbosity (int) – The current verbosity level (0=quiet, 1=normal, 2=debug).
selfisys.setup_model module
- class selfisys.setup_model.ModelSetup(size, L, P, S, G_sim_path, G_ss_path, Pbins_bnd, Pbins, k_s, P_ss_obj_path, P_0, planck_Pk)[source]
Bases:
NamedTuple
- G_sim_path: str
Alias for field number 4
- G_ss_path: str
Alias for field number 5
- L: float
Alias for field number 1
- P: int
Alias for field number 2
- P_0: ndarray
Alias for field number 10
- P_ss_obj_path: str
Alias for field number 9
- Pbins: ndarray
Alias for field number 7
- Pbins_bnd: ndarray
Alias for field number 6
- S: int
Alias for field number 3
- k_s: ndarray
Alias for field number 8
- planck_Pk: ndarray
Alias for field number 11
- size: int
Alias for field number 0
- selfisys.setup_model.compute_alpha_cv(workdir: str, k_s: ndarray, size: int, L: float, window_fct_path: str | None = None, force: bool = False) None [source]
Compute the cosmic variance parameter alpha_cv.
- Parameters:
workdir (str) – Directory where the results will be stored.
k_s (np.ndarray) – Support wavenumbers.
size (int) – Number of elements in each direction of the box.
L (float) – Comoving length of the box in Mpc/h.
window_fct_path (str, optional) – Path to the window function file.
force (bool) – If True, forces recomputation of the inputs.
- selfisys.setup_model.setup_model(workdir: str, params_planck: dict, params_P0: dict, size: int = 256, L: float = 3600.0, S: int = 100, N_exact: int = 8, Pinit: int = 50, trim_threshold: int = 100, minval: float | None = None, maxval: float | None = None, force: bool = False) ModelSetup [source]
Set up the model by computing or loading necessary grids and parameters.
- Parameters:
workdir (str) – Directory where the results will be stored.
params_planck (dict) – Parameters for the Planck 2018 cosmology.
params_P0 (dict) – Parameters for the normalisation power spectrum.
size (int) – Number of elements in each direction of the box.
L (float) – Comoving length of the box in Mpc/h.
S (int) – Number of support wavenumbers for the input power spectra.
N_exact (int) – Number of support wavenumbers matching the Fourier grid.
Pinit (int) – Maximum number of bins for the summaries.
trim_threshold (int) – Minimum number of modes required per bin.
minval (float, optional) – Minimum k value for the summaries.
maxval (float, optional) – Maximum k value for the summaries.
force (bool) – If True, forces recomputation of the inputs.
- Returns:
A named tuple containing: - size (int): Number of elements in each direction of the box. - L (float): Comoving length of the box in Mpc/h. - P (int): Number of bins for the summaries. - S (int): Number of support wavenumbers for input powerspectra. - G_sim_path (str): Path to the full Fourier grid file. - G_ss_path (str): Path to the Fourier grid for summaries file. - Pbins_bnd (np.ndarray): Boundaries of summary bins. - Pbins (np.ndarray): Centres of the bins for the summaries. - k_s (np.ndarray): Support wavenumbers for input power spectra. - P_ss_obj_path (str): Path to the summary power spectrum file. - P_0 (np.ndarray): Normalisation power spectrum values. - planck_Pk (np.ndarray): Planck 2018 power spectrum values.
- Return type: