camb

Python CAMB interface (http://camb.info)

camb.get_results(params)[source]

Calculate results for specified parameters and return CAMBdata instance for getting results.

Parameters:paramsmodel.CAMBparams instance
Returns:CAMBdata instance
camb.get_transfer_functions(params)[source]

Calculate transfer functions for specified parameters and return CAMBdata instance for getting results and subsequently calculating power spectra.

Parameters:paramsmodel.CAMBparams instance
Returns:CAMBdata instance
camb.get_matter_power_interpolator(params, zmin=0, zmax=10, nz_step=100, zs=None, kmax=10, nonlinear=True, var1=None, var2=None, hubble_units=True, k_hunit=True, return_z_k=False, k_per_logint=None, log_interp=True, extrap_kmax=None)[source]

Return a 2D spline interpolation object to evaluate matter power spectrum as function of z and k/h e.g:

PK = get_matter_power_evaluator(params);
print 'Power spectrum at z=0.5, k/h=0.1 is %s (Mpc/h)^3 '%(PK.P(0.5, 0.1))
Parameters:
  • paramsmodel.CAMBparams instance
  • zmin – minimum z (use 0 or smaller than you want for good interpolation)
  • zmax – maximum z (use larger than you want for good interpolation)
  • nz_step – number of steps to sample in z (default max allowed is 100)
  • zs – instead of zmin,zmax, nz_step, can specific explicit array of z values to spline from
  • kmax – maximum k
  • nonlinear – include non-linear correction from halo model
  • var1 – variable i (index, or name of variable; default delta_tot)
  • var2 – variable j (index, or name of variable; default delta_tot)
  • hubble_units – if true, output power spectrum in (Mpc/h)^{3} units, otherwise Mpc^{3}
  • k_hunit – if true, matter power is a function of k/h, if false, just k (both Mpc^{-1} units)
  • return_z_k – if true, return interpolator, z, k where z, k are the grid used
  • k_per_logint – specific uniform sampling over log k (if not set, uses optimized irregular sampling)
  • log_interp – if true, interpolate log of power spectrum (unless any values are negative in which case ignored)
  • extrap_kmax – if set, use power law extrapolation beyond kmax to extrap_kmax (useful for tails of integrals)
Returns:

RectBivariateSpline object PK, that can be called with PK(z,log(kh)) to get log matter power values. if return_z_k=True, instead return interpolator, z, k where z, k are the grid used

camb.get_age(params)[source]

Get age of universe for given set of parameters

Parameters:paramsmodel.CAMBparams instance
Returns:age of universe in gigayears
camb.get_zre_from_tau(params, tau)[source]

Get reionization redshift given optical depth tau

Parameters:
Returns:

reionization redshift (or negative number if error)

camb.set_z_outputs(z_outputs)[source]

Set the redshifts for calculating BAO parameters at

Parameters:z_outputs – array of redshifts
camb.set_feedback_level(level=1)[source]

Set the feedback level for internal CAMB calls :param level: zero for nothing, >1 for more

camb.set_params(cp=None, verbose=False, **params)[source]

Set all CAMB parameters at once, including parameters which are part of the CAMBparams structure, as well as global parameters.

E.g.

cp = camb.set_params(ns=1, omch2=0.1, ALens=1.2, lmax=2000)

This is equivalent to:

cp = model.CAMBparams() cp.set_cosmology(omch2=0.1) cp.set_for_lmax(lmax=2000) cp.InitPower.set_params(ns=1) lensing.ALens.value = 1.2

Parameters:
  • **params

    the values of the parameters

  • cp – use this CAMBparams instead of creating a new one
  • verbose – print out the equivalent set of commands
camb.set_halofit_version(version='takahashi')[source]

Set the halofit model for non-linear corrections.

Parameters:version

One of

camb.set_custom_scalar_sources(custom_sources, source_names=None, source_ell_scales=None, frame='CDM', code_path=None)[source]

Set custom sources for angular power spectrum using camb.symbolic sympy expressions.

Parameters:
  • custom_sources – list of sympy expressions for the angular power spectrum sources
  • source_names – optional list of string naes for the sources
  • source_ell_scales – list or dictionary of scalings for each source name, where for integer entry n, the source for multipole ell is scalled by sqrt((ell+n)!/(ell-n)!), i.e. n=2 for a new polarization-like source.
  • frame – if the source is not gauge invariant, frame in which to interpret result
  • code_path – optional path for output of source code for CAMB f90 source function
class camb.CAMBdata[source]

An object for storing transfer function data and parameters for CAMB. Not that it only stores transfer functions. If you want to get power spectra or background functions, you must have called one of the calculation functions for the parameters of interest more recently than any other call to these functions. You can can make multiple instances of CAMBdata and then later call power_spectra_from_transfer to calculate other quantities.

To quickly make a fully calculated CAMBdata instance for a set of parameters you can call get_results().

Variables:Params – the model.CAMBparams parameters being used
angular_diameter_distance(z)[source]

Get (non-comoving) angular diameter distance to redshift z.

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Parameters:z – redshift or array of redshifts
Returns:angular diameter distances, matching rank of z
angular_diameter_distance2(z1, z2)[source]

Get angular diameter distance between two redshifts

[r/(1+z2)]*sin_K([chi(z_1) - chi(z_1)]/r) where r is curvature radius and chi is the comoving radial distance

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Parameters:
  • z1 – redshift 1
  • z2 – redshift 2
Returns:

result

calc_background(params)[source]

Calculate the background evolution and thermal history. e.g. call this if you want to get derived parameters and call background functions :param params: model.CAMBparams instance to use

calc_background_no_thermo(params)[source]

Calculate the background evolution without calculating thermal history. e.g. call this if you want to just use angular_diameter_distance and similar background functions

Parameters:paramsmodel.CAMBparams instance to use
calc_power_spectra(params=None)[source]

Calculates transfer functions and power spectra.

Parameters:params – optional model.CAMBparams instance with parameters to use
calc_transfers(params, only_transfers=True)[source]

Calculate the transfer functions (for CMB and matter power, as determined by params.WantCls, params.WantTransfer)

Parameters:
  • paramsmodel.CAMBparams instance with parameters to use
  • only_transfers – only calculate transfer functions, no power spectra
Returns:

non-zero if error, zero if OK

comoving_radial_distance(z, tol=0.0001)[source]

Get comoving radial distance from us to redshift z in Mpc. This is efficient for arrays.

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Parameters:z – redshift
Returns:comoving radial distance (Mpc)
conformal_time(z)[source]

Conformal time from hot big bang to redshift z in Megaparsec. Use comoving_radial_distance for faster result for arrays.

Parameters:z – redshift or array of redshifts
Returns:eta(z)/Mpc
conformal_time_a1_a2(a1, a2)[source]

Get conformal time between two scale factors (=comoving radial distance travelled by light on light cone)

Parameters:
  • a1 – scale factor 1
  • a2 – scale factor 2
Returns:

eta(a2)-eta(a1) = chi(a1)-chi(a2) in Megaparsec

cosmomc_theta()[source]

Get theta_MC, an approximation of the radio of the sound horizon to the angular diameter distance at recombination.

Returns:theta_MC
get_BAO(redshifts, params)[source]

Get BAO parameters at given redshifts, using parameters in params

Parameters:
  • redshifts – list of redshifts
  • params – optional model.CAMBparams instance to use
Returns:

array of rs/DV, H, DA, F_AP for each redshift as 2D array

get_background_outputs()[source]

Get BAO values for redshifts (set redshifts using ‘set_z_outputs’)

Returns:rs/DV, H, DA, F_AP for each requested redshift (as 2D array)
get_background_redshift_evolution(z, vars=['x_e', 'opacity', 'visibility', 'cs2b'], format='dict')[source]

Get the evolution of background variables a function of redshift. For the moment a and H are rather perversely only available via get_time_evolution()

Parameters:
  • z – array of requested redshifts to output
  • vars – list of variable names to output
  • format – ‘dict’ or ‘array’, for either dict of 1D arrays indexed by name, or 2D array
Returns:

n_eta x len(vars) 2D numpy array of outputs or dict of 1D arrays

get_background_time_evolution(eta, vars=['x_e', 'opacity', 'visibility', 'cs2b'], format='dict')[source]

Get the evolution of background variables a function of conformal time. For the moment a and H are rather perversely only available via get_time_evolution()

Parameters:
  • eta – array of requested conformal times to output
  • vars – list of variable names to output
  • format – ‘dict’ or ‘array’, for either dict of 1D arrays indexed by name, or 2D array
Returns:

n_eta x len(vars) 2D numpy array of outputs or dict of 1D arrays

get_cmb_correlation_functions(params=None, lmax=None, spectrum='lensed_scalar', xvals=None, sampling_factor=1)[source]

Get the CMB correlation functions from the power spectra. By default evaluated at points cos(theta) = xvals that are roots of Legendre polynomials, for accurate back integration with correlations.corr2cl(). If xvals is explicitly given, instead calculates correlations at provided cos(theta) values.

Parameters:
  • params – optional model.CAMBparams instance with parameters to use. If None, must have previously set parameters and called calc_power_spectra (e.g. if you got this instance using camb.get_results),
  • lmax – optional maximum L to use from the cls arrays
  • spectrum – type of CMB power spectrum to get; default ‘lensed_scalar’, one of [‘total’, ‘unlensed_scalar’, ‘unlensed_total’, ‘lensed_scalar’, ‘tensor’]
  • xvals – optional array of cos(theta) values at which to calculate correlation function.
  • sampling_factor – multiple of lmax for the Gauss-Legendre order if xvals not given (default 1)
Returns:

if xvals not given: corrs, xvals, weights; if xvals specified, just corrs. corrs is 2D array corrs[i, ix], where ix=0,1,2,3 are T, Q+U, Q-U and cross, and i indexes xvals

get_cmb_power_spectra(params=None, lmax=None, spectra=['total', 'unlensed_scalar', 'unlensed_total', 'lensed_scalar', 'tensor', 'lens_potential'], CMB_unit=None, raw_cl=False)[source]

Get CMB power spectra, as requested by the ‘spectra’ argument. All power spectra are l(l+1)C_l/2pi self-owned numpy arrays (0..lmax, 0..3), where 0..3 index are TT, EE, BB TT, unless raw_cl is True in which case return just C_l. For the lens_potential the power spectrum returned is that of the deflection.

Parameters:
  • params – optional model.CAMBparams instance with parameters to use. If None, must have previously set parameters and called calc_power_spectra (e.g. if you got this instance using camb.get_results),
  • lmax – maximum l
  • spectra – list of names of spectra to get
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

dictionary of power spectrum arrays, indexed by names of requested spectra

get_cmb_transfer_data(tp='scalar')[source]

Get C_l transfer functions

Returns:class:.ClTransferData instance holding output arrays (copies, not pointers)
get_cmb_unlensed_scalar_array_dict(params=None, lmax=None, CMB_unit=None, raw_cl=False)[source]

Get all unlensed auto and cross power spectra, including any custom source functions set using set_custom_scalar_sources().

Parameters:
  • params – optional model.CAMBparams instance with parameters to use. If None, must have previously set parameters and called calc_power_spectra (e.g. if you got this instance using camb.get_results),
  • lmax – maximum l
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

dictionary of power spectrum arrays, index as TxT, TxE, W1xW2, custom_name_1xT... etc.

get_derived_params()[source]
Returns:dictionary of derived parameter values, indexed by name (‘kd’, ‘age’, etc..)
get_fsigma8()[source]

Get f*sigma_8 growth values (must previously have calculated power spectra). For general models f*sigma_8 is defined as in the Planck 2015 parameter paper in terms of the velocity-density correlation: sigma^2_{vd}/sigma_{dd} for 8 h^{-1} Mpc spheres.

Returns:array of f*sigma_8 values, in order of increasing time (decreasing redshift)
get_lens_potential_cls(lmax, CMB_unit=None, raw_cl=False)[source]

Get lensing deflection angle potential power spectrum, and cross-correlation with T and E. Must have already calculated power spectra. Power spectra are [l(l+1)]^2C_l^{phi phi}/2/pi and corresponding deflection cross-correlations.

Parameters:
  • lmax – lmax to output to
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return lensing potential C_L rather than [L(L+1)]^2C_L/2pi
Returns:

numpy array CL[0:lmax+1,0:3], where 0..2 indexes PP, PT, PE.

get_lensed_scalar_cls(lmax, CMB_unit=None, raw_cl=False)[source]

Get lensed scalar CMB power spectra. Must have already calculated power spectra.

Parameters:
  • lmax – lmax to output to
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

numpy array CL[0:lmax+1,0:4], where 0..3 indexes TT, EE, BB, TE.

get_linear_matter_power_spectrum(var1=None, var2=None, hubble_units=True, have_power_spectra=False, params=None, nonlinear=False)[source]

Calculates P_{xy}(k/h), where x, y are one of model.Transfer_cdm, model.Transfer_xx etc. The output k values are not regularly spaced, and not interpolated.

Parameters:
  • var1 – variable i (index, or name of variable; default delta_tot)
  • var2 – variable j (index, or name of variable; default delta_tot)
  • hubble_units – if true, output power spectrum in (Mpc/h) units, otherwise Mpc
  • have_power_spectra – set to True if already computed power spectra
  • params – if have_power_spectra=False, optional model.CAMBparams instance to specify new parameters
  • nonlinear – include non-linear correction from halo model
Returns:

kh, z, PK, where kz an z are arrays of k/h and z respectively, and PK[i,j] is value at z[i], k/h[j]

get_matter_power_spectrum(minkh=0.0001, maxkh=1.0, npoints=100, var1=None, var2=None, have_power_spectra=False, params=None)[source]

Calculates P_{xy}(k/h), where x, y are one of Transfer_cdm, Transfer_xx etc defined in ModelParams. The output k values are regularly log spaced and interpolated. If NonLinear is set, the result is non-linear.

Parameters:
  • minkh – minimum value of k/h for output grid (very low values < 1e-4 may not be calculated)
  • maxkh – maximum value of k/h (check consistent with input params.Transfer.kmax)
  • npoints – number of points equally spaced in log k
  • var1 – variable i (index, or name of variable; default delta_tot)
  • var2 – variable j (index, or name of variable; default delta_tot)
  • have_power_spectra – set to True if already computed power spectra
  • params – if have_power_spectra=False and want to specify new parameters, a model.CAMBparams instance
Returns:

kh, z, PK, where kz an z are arrays of k/h and z respectively, and PK[i,j] is value at z[i], k/h[j]

get_matter_transfer_data()[source]

Get matter transfer function data and sigma8 for calculated results.

Returns:MatterTransferData instance holding output arrays (copies, not pointers)
get_nonlinear_matter_power_spectrum(**kwargs)[source]

Calculates P_{xy}(k/h), where x, y are one of model.Transfer_cdm, model.Transfer_xx etc. The output k values are not regularly spaced, and not interpolated.

Parameters:
  • var1 – variable i (index, or name of variable; default delta_tot)
  • var2 – variable j (index, or name of variable; default delta_tot)
  • hubble_units – if true, output power spectrum in (Mpc/h)^{3} units, otherwise Mpc^{3}
  • have_power_spectra – set to True if already computed power spectra
  • params – if have_power_spectra=False, optional model.CAMBparams instance to specify new parameters
Returns:

kh, z, PK, where kz an z are arrays of k/h and z respectively, and PK[i,j] is value at z[i], k/h[j]

get_params()[source]

Get the parameters currently set. Returned object references stored data, so elements can be modified without calling set_params again.

Returns:model.CAMBparams instance pointing to the underlying parameters used by CAMB.
get_redshift_evolution(q, z, vars=['k/h', 'delta_cdm', 'delta_baryon', 'delta_photon', 'delta_neutrino', 'delta_nu', 'delta_tot', 'delta_nonu', 'delta_tot_de', 'Weyl', 'v_newtonian_cdm', 'v_newtonian_baryon', 'v_baryon_cdm', 'a', 'etak', 'H', 'growth', 'v_photon', 'pi_photon', 'E_2', 'v_neutrino', 'T_source', 'E_source', 'lens_potential_source'], lAccuracyBoost=4)[source]

Get the mode evolution as a function of redshift for some k values.

Parameters:
  • q – wavenumber values to calculate (or array of k values)
  • z – array of redshifts to output
  • vars – list of variable names or camb.symbolic sympy expressions to output
  • lAccuracyBoost – boost factor for ell accuracy (e.g. to get nice smooth curves for plotting)
Returns:

nd array, A_{qti}, size(q) x size(times) x len(vars), or 2d array if q is scalar

get_sigma8()[source]

Get sigma_8 values (must previously have calculated power spectra)

Returns:array of sigma_8 values, in order of increasing time (decreasing redshift)
get_tensor_cls(lmax, CMB_unit=None, raw_cl=False)[source]

Get tensor CMB power spectra. Must have already calculated power spectra.

Parameters:
  • lmax – lmax to output to
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

numpy array CL[0:lmax+1,0:4], where 0..3 indexes TT, EE, BB, TE

get_time_evolution(q, eta, vars=['k/h', 'delta_cdm', 'delta_baryon', 'delta_photon', 'delta_neutrino', 'delta_nu', 'delta_tot', 'delta_nonu', 'delta_tot_de', 'Weyl', 'v_newtonian_cdm', 'v_newtonian_baryon', 'v_baryon_cdm', 'a', 'etak', 'H', 'growth', 'v_photon', 'pi_photon', 'E_2', 'v_neutrino', 'T_source', 'E_source', 'lens_potential_source'], lAccuracyBoost=4, frame='CDM')[source]

Get the mode evolution as a function of conformal time for some k values.

Parameters:
  • q – wavenumber values to calculate (or array of k values)
  • eta – array of requested conformal times to output
  • vars – list of variable names or sympy symbolic expressions to output (using camb.symbolic)
  • lAccuracyBoost – factor by which to increase l_max in hierarchies compared to default - often needed to get nice smooth curves of acoustic oscillations for plotting.
  • frame – for symbolic expressions, can specify frame name if the variable is not gauge invariant. e.g. specifying Delta_g and frame=’Newtonian’ would give the Newtonian gauge photon density perturbation.
Returns:

nd array, A_{qti}, size(q) x size(times) x len(vars), or 2d array if q is scalar

get_total_cls(lmax, CMB_unit=None, raw_cl=False)[source]

Get lensed-scalar + tensor CMB power spectra. Must have already calculated power spectra.

Parameters:
  • lmax – lmax to output to
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

numpy array CL[0:lmax+1,0:4], where 0..3 indexes TT, EE, BB, TE

get_unlensed_scalar_array_cls(lmax)[source]

Get array of all cross power spectra. Must have already calculated power spectra. Results are dimensionless, and not scaled by custom_scaled_ell_fac.

Parameters:lmax – lmax to output to
Returns:numpy array CL[0:, 0:,0:lmax+1], where 0.. index T, E, deflection angle, source window functions
get_unlensed_scalar_cls(lmax, CMB_unit=None, raw_cl=False)[source]

Get unlensed scalar CMB power spectra. Must have already calculated power spectra.

Parameters:
  • lmax – lmax to output to
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

numpy array CL[0:lmax+1,0:4], where 0..3 indexes TT, EE, BB, TE. CL[:,2] will be zero.

get_unlensed_total_cls(lmax, CMB_unit=None, raw_cl=False)[source]

Get unlensed CMB power spectra, including tensors if relevant. Must have already calculated power spectra.

Parameters:
  • lmax – lmax to output to
  • CMB_unit – scale results from dimensionless. Use ‘muK’ for muK^2 units for CMB CL and muK units for lensing cross.
  • raw_cl – return C_L rather than L(L+1)C_L/2pi
Returns:

numpy array CL[0:lmax+1,0:4], where 0..3 indexes TT, EE, BB, TE.

h_of_z(z)[source]

Get Hubble rate at redshift z, in Mpc^{-1} units, scalar or array

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Use hubble_parameter instead if you want in [km/s/Mpc] units.

Parameters:z – redshift
Returns:H(z)
hubble_parameter(z)[source]

Get Hubble rate at redshift z, in km/s/Mpc units. Scalar or array.

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Parameters:z – redshift
Returns:H(z)/[km/s/Mpc]
luminosity_distance(z)[source]

Get luminosity distance from to redshift z.

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Parameters:z – redshift or array of redshifts
Returns:luminosity distance (matches rank of z)
physical_time(z)[source]

Get physical time from hot big bang to redshift z in Gigayears.

Parameters:z – redshift
Returns:t(z)/Gigayear
physical_time_a1_a2(a1, a2)[source]

Get physical time between two scalar factors in Gigayears

Must have called calc_background, calc_background_no_thermo or calculated transfer functions or power spectra.

Parameters:
  • a1 – scale factor 1
  • a2 – scale factor 2
Returns:

(age(a2)-age(a1))/Gigayear

power_spectra_from_transfer(initial_power_params)[source]

Assuming calc_transfers or calc_power_spectra have already been used, re-calculate the power spectra using a new set of initial power spectrum parameters with otherwise the same cosmology. This is typically much faster that re-calculating everything, as the transfer functions can be re-used.

Parameters:initial_power_paramsinitialpower.InitialPowerParams instance with new primordial power spectrum parameters
redshift_at_comoving_radial_distance(chi, nz_step=150, zmax=10000)[source]

Convert comoving radial distance array to redshift array. This is not calculated directly, but fit from a spline to a forward calculation of chi from z This is a utility routine, not optimized to be fast (though it can work on a large vector efficiently)

Parameters:
  • chi – comoving radial distance (in Mpc), scalar or array
  • nz_step – number of redshifts calculated internally for solving grid
  • zmax – maximum redshift in internal solving grid
Returns:

redshift at chi, scalar or array

set_params(params)[source]

Set parameters from params

Parameters:params – a model.CAMBparams instance
class camb.MatterTransferData[source]

MatterTransferData is the base class for storing matter power transfer function data for various q values. In a flat universe q=k, in a closed universe q is quantised.

To get an instance of this data, call camb.CAMBdata.get_matter_transfer_data()

Variables:
  • nq – number of q modes calculated
  • q – array of q values calculated
  • sigma_8 – array of sigma8 values for each redshift for each power spectrum
  • sigma2_vdelta_8 – array of v-delta8 correlation, so sigma2_vdelta_8/sigma_8 can define growth
  • transfer_data

    numpy array T[entry, q_index, z_index] storing transfer functions for each redshift and q; entry+1 can be

    • Transfer_kh = 1 (k/h)
    • Transfer_cdm = 2 (cdm)
    • Transfer_b = 3 (baryons)
    • Transfer_g = 4 (photons)
    • Transfer_r = 5 (massless neutrinos)
    • Transfer_nu = 6 (massive neutrinos)
    • Transfer_tot = 7 (total matter)
    • Transfer_nonu = 8 (total matter excluding neutrinos)
    • Transfer_tot_de = 9 (total including dark energy perturbations)
    • Transfer_Weyl = 10 (Weyl potential)
    • Transfer_Newt_vel_cdm = 11 (Newtonian CDM velocity)
    • Transfer_Newt_vel_baryon = 12 (Newtonian baryon velocity)
    • Transfer_vel_baryon_cdm = 13 (relative baryon-cdm velocity)
transfer_z(name, z_index=0)[source]

Get transfer function (function of q, for each q in self.q_trans) by name for given redshift index

Parameters:
  • name – parameter name
  • z_index – which redshift
Returns:

array of transfer function values for each calculated k