abacusnbody.hod.zcv#
abacusnbody.hod.zcv.advect_fields module#
- abacusnbody.hod.zcv.advect_fields.main(path2config, want_rsd=False, alt_simname=None, save_3D_power=False, only_requested_fields=False)[source]#
Advect the initial conditions fields (1cb, delta, delta^2, s^2, nabla^2) to some desired redshift and saving the 3D Fourier fields and power spectra in ASDF files along the way.
- Parameters:
path2config (str) – name of the yaml containing parameter specifications.
want_rsd (bool, optional) – compute the advected fields and power spectra in redshift space? Default is False.
alt_simname (str, optional) – specify simulation name if different from yaml file.
save_3D_power (bool, optional) – save the 3D power spectra in individual ASDF files. Default is False.
only_requested_fields (bool, optional) – instead of all 5 fields, use only the fields specified in zcv_params. Default is False.
- Returns:
pk_ij_dict – dictionary containing the auto- and cross-power spectra of the 5 fields.
- Return type:
dict
abacusnbody.hod.zcv.ic_fields module#
- abacusnbody.hod.zcv.ic_fields.compress_asdf(asdf_fn, table, header)[source]#
Compress the dictionaries table and header using blsc into an ASDF file, asdf_fn.
- abacusnbody.hod.zcv.ic_fields.load_dens(ic_dir, sim_name, nmesh)[source]#
Load initial condition density field for the given AbacusSummit simulation.
- abacusnbody.hod.zcv.ic_fields.load_disp(ic_dir, sim_name, nmesh)[source]#
Load initial condition displacement fields for the given AbacusSummit simulation.
- abacusnbody.hod.zcv.ic_fields.gaussian_filter(field, nmesh, lbox, kcut)[source]#
Apply a fourier space gaussian filter to a field.
- Parameters:
field (array_like) – the field to filter.
nmesh (int) – size of the mesh.
lbox (float) – size of the box.
kcut (float) – the exponential cutoff to use in the gaussian filter
- Returns:
f_filt – Gaussian filtered version of field
- Return type:
array_like
- abacusnbody.hod.zcv.ic_fields.filter_field(delta_k, n1d, L, kcut, dtype=<class 'numpy.float32'>)[source]#
Compute nabla^2 delta in Fourier space.
- Parameters:
delta_k (array_like) – Fourier 3D field.
n1d (int) – size of the 3d array along x and y dimension.
L (float) – box size of the simulation.
kcut (float) – smoothing scale in Fourier space.
dtype (np.dtype) – float type (32 or 64) to use in calculations.
- Returns:
n2_fft – Fourier 3D field.
- Return type:
array_like
- abacusnbody.hod.zcv.ic_fields.get_n2_fft(delta_k, n1d, L, dtype=<class 'numpy.float32'>)[source]#
Compute nabla^2 delta in Fourier space.
- Parameters:
delta_k (array_like) – Fourier 3D field.
n1d (int) – size of the 3d array along x and y dimension.
L (float) – box size of the simulation.
dtype (np.dtype) – float type (32 or 64) to use in calculations.
- Returns:
n2_fft – Fourier 3D field.
- Return type:
array_like
- abacusnbody.hod.zcv.ic_fields.get_sij_fft(i_comp, j_comp, delta_k, n1d, L, dtype=<class 'numpy.float32'>)[source]#
Compute ijth component of the tidal tensor in Fourier space.
- Parameters:
i_comp (int) – ith component of the tensor.
j_comp (int) – jth component of the tensor.
delta_k (array_like) – Fourier 3D field.
n1d (int) – size of the 3d array along x and y dimension.
L (float) – box size of the simulation.
dtype (np.dtype) – float type (32 or 64) to use in calculations.
- Returns:
s_ij_fft – Fourier 3D field.
- Return type:
array_like
- abacusnbody.hod.zcv.ic_fields.add_ij(final_field, field_to_add, n1d, factor=1.0, dtype=<class 'numpy.float32'>)[source]#
Add field field_to_add to final_field with a constant factor.
- abacusnbody.hod.zcv.ic_fields.get_dk_to_s2(delta_k, nmesh, lbox)[source]#
Computes the square tidal field from the density FFT s^2 = s_ij s_ij, where s_ij = (k_i k_j / k^2 - delta_ij / 3 ) * delta_k.
- Parameters:
delta_k (array_like) – Fourier transformed density field.
nmesh (int) – size of the mesh.
lbox (float) – size of the box.
- Returns:
the tidal field (s^2).
- Return type:
tidesq
- abacusnbody.hod.zcv.ic_fields.get_dk_to_n2(delta_k, nmesh, lbox)[source]#
Computes the density curvature from the density field: nabla^2 delta = IFFT(-k^2 delta_k) :param delta_k: Fourier transformed density field. :type delta_k: array_like :param nmesh: size of the mesh. :type nmesh: int :param lbox: size of the box. :type lbox: float
- Returns:
real_gradsqdelta – the nabla^2 delta field
- Return type:
array_like
- abacusnbody.hod.zcv.ic_fields.get_fields(delta_lin, Lbox, nmesh)[source]#
Return the fields delta, delta^2, s^2, nabla^2 given the linear density field.
- abacusnbody.hod.zcv.ic_fields.main(path2config, alt_simname=None, verbose=False)[source]#
Save the initial conditions fields (1cb, delta, delta^2, s^2, nabla^2) as ASDF files.
Note: you can save the fields separately if they are using too much memory. TODO: the multiplications in Fourier space can be sped up with numba.
- Parameters:
path2config (str) – name of the yaml containing parameter specifications.
alt_simname (str, optional) – specify simulation name if different from yaml file.
verbose (bool, optional) – print some useful benchmark statements. Default is False.
abacusnbody.hod.zcv.linear_fields module#
Script for saving the linear density field power spectrum and cross correlations needed for Kaiser effects TODO: change to pyfftw
- abacusnbody.hod.zcv.linear_fields.main(path2config, alt_simname=None, save_3D_power=False)[source]#
Advect the initial conditions density field to some desired redshift and saving the 3D Fourier fields (delta, delta*mu^2) and power spectra in ASDF files along the way.
- Parameters:
path2config (str) – name of the yaml containing parameter specifications.
alt_simname (str, optional) – specify simulation name if different from yaml file.
save_3D_power (bool, optional) – save the 3D power spectra in individual ASDF files. Default is False.
- Returns:
pk_lin_dict – dictionary containing the auto- and cross-power spectra of the 2 fields.
- Return type:
dict
abacusnbody.hod.zcv.tools_cv module#
Tools for applying variance reduction (ZCV and LCV) (based on scripts from Joe DeRose).
- abacusnbody.hod.zcv.tools_cv.combine_spectra(k, spectra, bias_params, rsd=False, numerical_nabla=False)[source]#
ZCV: Given some bias parameters, compute the model power spectra.
- abacusnbody.hod.zcv.tools_cv.combine_cross_spectra(k, spectra, bias_params, rsd=False)[source]#
ZCV: Given some bias parameters, compute the model power spectra (no shotnoise in cross).
- abacusnbody.hod.zcv.tools_cv.combine_cross_kaiser_spectra(k, spectra_dict, D, bias, f_growth, rec_algo, R, rsd=False)[source]#
LCV: Convert measured templates into tracer-model power spectrum given bias parameters assuming the Kaiser approximation (Chen et al. 2019, 1907.00043).
- RecSym equation:
< D (b delta + f mu2 delta, tr > = D * (b < delta, tr> + f < mu2 delta, tr >)
- RecIso equation:
< D ((b+f mu^2)(1-S) + bS) delta, tr > = D * (b < delta, tr > + f (1-S) < delta, mu^2 >)
- abacusnbody.hod.zcv.tools_cv.combine_kaiser_spectra(k, spectra_dict, D, bias, f_growth, rec_algo, R, rsd=False)[source]#
LCV: Convert measured templates into model-model power spectrum given bias parameters assuming the Kaiser approximation (Chen et al. 2019, 1907.00043).
- RecSym equation:
< D (b delta + f mu2 delta, D (b delta + f mu2 delta) > = = D^2 (b^2 < delta, delta> + f^2 < mu2 delta, mu2 delta > + 2 b f < delta, mu2 delta >)
- RecIso Equation:
((b+f mu^2)(1-S) + bS) delta = (b + f (1-S) mu2) delta < D ((b+f mu^2)(1-S) + bS) delta, D ((b+f mu^2)(1-S) + bS) delta > = = D^2 (b^2 < delta, delta> + fefff^2 < mu2 delta, mu2 delta > + 2 b feff < delta, mu2 delta >)
- abacusnbody.hod.zcv.tools_cv.get_poles(k, pk, D, bias, f_growth, poles=[0, 2, 4])[source]#
Compute the len(poles) multipoles given the linear power spectrum, pk, the growth function, the growth factor and the bias.
- abacusnbody.hod.zcv.tools_cv.multipole_cov(pell, ell)[source]#
Factors appearing in the covariance matrix that couple the multipoles.
- abacusnbody.hod.zcv.tools_cv.measure_2pt_bias(k, pk_ij, pk_tt, kmax, keynames, kmin=0.0, rsd=False)[source]#
ZCV: Infer the bias based on the template power spectrum and tracer measurements. Note that the bias parameter ordering corresponds to rsd == False.
- abacusnbody.hod.zcv.tools_cv.combine_field_spectra_k3D_lcv(bias, f_growth, D, power_lin_fns, power_rsd_tr_fns, nmesh, Lbox, R, rec_algo)[source]#
LCV: Given bias parameters, compute the model-model auto and cross correlation for the 3D k-vector.
- abacusnbody.hod.zcv.tools_cv.combine_field_spectra_k3D(bias, power_ij_fns, keynames)[source]#
ZCV: Given bias parameters, compute the model-model cross correlation for the 3D k-vector.
- abacusnbody.hod.zcv.tools_cv.combine_field_cross_spectra_k3D(bias, power_tr_fns, keynames)[source]#
ZCV: Given bias parameters, compute the model-tracer cross correlation for the 3D k-vector.
- abacusnbody.hod.zcv.tools_cv.measure_2pt_bias_lcv(k, power_dict, power_rsd_tr_dict, D, f_growth, kmax, rsd, rec_algo, R, ellmax=2, kmin=0.0)[source]#
LCV: Function for getting the linear bias in the Kaiser approximation.
- abacusnbody.hod.zcv.tools_cv.read_power_dict(power_tr_dict, power_ij_dict, want_rsd, keynames, poles)[source]#
ZCV: Function for reading the power spectra and saving them in the same format as Zenbu.
- abacusnbody.hod.zcv.tools_cv.get_cfg(sim_name, z_this, nmesh)[source]#
ZCV: Configuration parameters.
- abacusnbody.hod.zcv.tools_cv.run_zcv(power_rsd_tr_dict, power_rsd_ij_dict, power_tr_dict, power_ij_dict, config)[source]#
Apply Zel’dovich control variates (ZCV) reduction to some measured power spectrum.
- abacusnbody.hod.zcv.tools_cv.run_zcv_field(power_rsd_tr_fns, power_rsd_ij_fns, power_tr_fns, power_ij_fns, config)[source]#
Apply Zel’dovich control variates (ZCV) reduction to some measured 3D power spectrum.
abacusnbody.hod.zcv.tracer_power module#
- abacusnbody.hod.zcv.tracer_power.get_tracer_power(tracer_pos, want_rsd, config, want_save=True, save_3D_power=False)[source]#
Compute the auto- and cross-correlation between a galaxy catalog (tracer_pos) and the advected Zel’dovich fields.
- Parameters:
tracer_pos (array_like) – galaxy positions with shape (N, 3)
want_rsd (bool) – compute the power spectra in redshift space?
config (str) – name of the yaml containing parameter specifications.
save_3D_power (bool, optional) – save the 3D power spectra in individual ASDF files. Default is False.
- Returns:
pk_tr_dict – dictionary containing the auto- and cross-power spectra of the tracer with the 5 fields.
- Return type:
dict
- abacusnbody.hod.zcv.tracer_power.get_recon_power(tracer_pos, random_pos, want_rsd, config, want_save=True, save_3D_power=False, want_load_tr_fft=False)[source]#
Compute the auto- and cross-correlation between a galaxy catalog (tracer_pos) and the advected initial conditions fields (delta, delta*mu^2).
- Parameters:
tracer_pos (array_like) – galaxy positions with shape (N, 3)
random_pos (array_like) – randoms positions with shape (M, 3)
want_rsd (bool) – compute the power spectra in redshift space?
config (str) – name of the yaml containing parameter specifications.
save_3D_power (bool, optional) – save the 3D power spectra in individual ASDF files. Default is False.
want_load_tr_fft (bool, optional) – want to load provided 3D Fourier tracer field? Default is False.
- Returns:
pk_tr_dict – dictionary containing the auto- and cross-power spectra of the tracer with the 2 fields.
- Return type:
dict
abacusnbody.hod.zcv.zenbu_window module#
Script for pre-saving zenbu for a given redshift and simulation.
- abacusnbody.hod.zcv.zenbu_window.periodic_window_function(nmesh, lbox, kout, kin, k2weight=True)[source]#
Compute matrix appropriate for convolving a finely evaluated theory prediction with the mode coupling matrix.
- Parameters:
nmesh (int) – size of the mesh used for power spectrum measurement.
lbox (float) – box size of the simulation.
kout (array_like) – k bins used for power spectrum measurement.
kin (array_like) – Defaults to None.
k2weight (bool) – Defaults to True.
- Returns:
window (array_like) – np.dot(window, pell_th) gives convovled theory
keff (array_like) – effective k value of each output k bin.
- abacusnbody.hod.zcv.zenbu_window.zenbu_spectra(k, z, cfg, kin, pin, pkclass=None, N=2700, jn=15, rsd=True, nmax=6, ngauss=6)[source]#
Compute the ZeNBu power spectra.
- abacusnbody.hod.zcv.zenbu_window.main(path2config, alt_simname=None, want_xi=False)[source]#
Save the mode-coupling window function and the ZeNBu power spectra as npz files given ZCV specs.
- Parameters:
path2config (str) – name of the yaml containing parameter specifications.
alt_simname (str, optional) – specify simulation name if different from yaml file.