pyccl.correlations module
- class pyccl.correlations.CorrelationMethods(value, names=None, *, module=None, qualname=None, type=None, start=1, boundary=None)[source]
Bases:
Enum
Choices of algorithms used to compute correlation functions:
‘Bessel’ is a direct integration using Bessel functions.
‘FFTLog’ is fast using a fast Fourier transform.
‘Legendre’ uses a sum over Legendre polynomials.
- FFTLOG = 'fftlog'
- BESSEL = 'bessel'
- LEGENDRE = 'legendre'
- class pyccl.correlations.CorrelationTypes(value, names=None, *, module=None, qualname=None, type=None, start=1, boundary=None)[source]
Bases:
Enum
Correlation function types.
- NN = 'NN'
- NG = 'NG'
- GG_PLUS = 'GG+'
- GG_MINUS = 'GG-'
- pyccl.correlations.correlation(cosmo, *, ell, C_ell, theta, type='NN', method='fftlog')[source]
Compute the angular correlation function.
\[\xi^{ab}_\pm(\theta) = \sum_\ell\frac{2\ell+1}{4\pi}\,(\pm1)^{s_b}\, C^{ab\pm}_\ell\,d^\ell_{s_a,\pm s_b}(\theta)\]where \(\theta\) is the angle between the two fields \(a\) and \(b\) with spins \(s_a\) and \(s_b\) after alignement of their tangential coordinate. \(d^\ell_{mm'}\) are the Wigner-d matrices and we have defined the power spectra
\[C^{ab\pm}_\ell \equiv (C^{a_Eb_E}_\ell \pm C^{a_Bb_B}_\ell)+i (C^{a_Bb_E}_\ell \mp C^{a_Eb_B}_\ell),\]which reduces to the \(EE\) power spectrum when all \(B\)-modes are 0.
The different spin combinations are:
\(s_a=s_b=0\) e.g. galaxy-galaxy, galaxy-\(\kappa\) and \(\kappa\)-\(\kappa\)
\(s_a=2\), \(s_b=0\) e.g. galaxy-shear, and \(\kappa\)-shear
\(s_a=s_b=2\) e.g. shear-shear.
Note
For scales smaller than \(\sim 0.1^{\circ}\), the input power spectrum should be sampled to sufficienly high \(\ell\) to ensure the Hankel transform is well-behaved. The following spline parameters, related to
FFTLog
-sampling may also be modified for accuracy:ccl.spline_params.ELL_MIN_CORR
ccl.spline_params.ELL_MAX_CORR
ccl.spline_params.N_ELL_CORR
.
- Parameters:
cosmo (
Cosmology
) – A Cosmology object.ell (array) – Multipoles corresponding to the input angular power spectrum.
C_ell (array) – Input angular power spectrum.
theta (
float
or array) – Angular separation(s) at which to calculate the angular correlation function (in degrees).type (
str
) – Type of correlation function. Choices:'NN'
(0x0),'NG'
(0x2),'GG+'
(2x2, \(\xi_+\)),'GG-'
(2x2, \(\xi_-\)), where numbers refer to the spins of the two quantities being cross-correlated (see Section 2.4.2 of the CCL paper). The naming system roughly follows the nomenclature used in TreeCorr.method (
str
) – Method to compute the correlation function. Choices:'Bessel'
(direct integration over Bessel function),'FFTLog'
(fast integration with FFTLog),'Legendre'
(brute-force sum over Legendre polynomials).
- Returns:
Value(s) of the correlation function at the input angular separations.
- Return type:
(
float
or array)
- pyccl.correlations.correlation_3d(cosmo, *, r, a, p_of_k_a='delta_matter:delta_matter')[source]
Compute the 3D correlation function:
\[\xi(r)\equiv\frac{1}{2\pi^2}\int dk\,k^2\,P(k)\,j_0(kr).\]- Parameters:
cosmo (
Cosmology
) – A Cosmology object.r (
float
or array) – distance(s) at which to calculate the 3D correlation function (in Mpc).a (
float
) – scale factor.p_of_k_a (
Pk2D
,str
orNone
) – 3D Power spectrum to integrate. If a string, it must correspond to one of the non-linear power spectra stored in cosmo (e.g. ‘delta_matter:delta_matter’).
- Returns:
Value(s) of the correlation function at the input distance(s).
- pyccl.correlations.correlation_multipole(cosmo, *, r, a, beta, ell, p_of_k_a='delta_matter:delta_matter')[source]
Compute the correlation function multipoles:
\[\xi_\ell(r)\equiv\frac{i^\ell}{2\pi^2}\int dk\,k^2\,P(k)\,j_\ell(kr).\]- Parameters:
cosmo (
Cosmology
) – A Cosmology object.r (
float
or array) – distance(s) at which to calculate the 3D correlation function (in Mpc).a (
float
) – scale factor.beta (
float
) – growth rate divided by galaxy bias.ell (
int
) – the desired multipolep_of_k_a (
Pk2D
,str
orNone
) – 3D Power spectrum to integrate. If a string, it must correspond to one of the non-linear power spectra stored in cosmo (e.g. ‘delta_matter:delta_matter’).
- Returns:
Value(s) of the correlation function at the input distance(s).
- pyccl.correlations.correlation_3dRsd(cosmo, *, r, a, mu, beta, p_of_k_a='delta_matter:delta_matter', use_spline=True)[source]
Compute the 3D correlation function with linear RSDs using multipoles:
\[\xi(r,\mu) = \sum_{\ell\in\{0,2,4\}}\xi_\ell(r)\,P_\ell(\mu)\]where \(P_\ell(\mu)\) are the Legendre polynomials, and \(\xi_\ell(r)\) are the correlation function multipoles.
- Parameters:
cosmo (
Cosmology
) – A Cosmology object.r (
float
or array) – distance(s) at which to calculate the 3D correlation function (in Mpc).a (
float
) – scale factor.mu (
float
) – cosine of the angle at which to calculate the 3D correlation function.beta (
float
) – growth rate divided by galaxy bias.p_of_k_a (
Pk2D
,str
orNone
) – 3D Power spectrum to integrate. If a string, it must correspond to one of the non-linear power spectra stored in cosmo (e.g. ‘delta_matter:delta_matter’).use_spline (
bool
) – switch that determines whether the RSD correlation function is calculated using global splines of multipoles.
- Returns:
Value(s) of the correlation function at the input distance(s) & angle.
- pyccl.correlations.correlation_3dRsd_avgmu(cosmo, *, r, a, beta, p_of_k_a='delta_matter:delta_matter')[source]
Compute the 3D correlation function averaged over angles with RSDs. Equivalent to calling
correlation_multipole()
withell=0
.- Parameters:
cosmo (
Cosmology
) – A Cosmology object.r (
float
or array) – distance(s) at which to calculate the 3D correlation function (in Mpc).a (
float
) – scale factor.beta (
float
) – growth rate divided by galaxy bias.p_of_k_a (
Pk2D
,str
orNone
) – 3D Power spectrum to integrate. If a string, it must correspond to one of the non-linear power spectra stored in cosmo (e.g. ‘delta_matter:delta_matter’).
- Returns:
Value(s) of the correlation function at the input distance(s) & angle.
- pyccl.correlations.correlation_pi_sigma(cosmo, *, pi, sigma, a, beta, use_spline=True, p_of_k_a='delta_matter:delta_matter')[source]
Compute the 3D correlation in \((\pi,\sigma)\) space. This is just
\[\xi(\pi,\sigma) = \xi(r=\sqrt{\pi^2+\sigma^2},\mu=\pi/r).\]- Parameters:
cosmo (
Cosmology
) – A Cosmology object.pi (
float
) – distance times cosine of the angle (in Mpc).sigma (
float
or array) – distance(s) times sine of the angle (in Mpc).a (
float
) – scale factor.beta (
float
) – growth rate divided by galaxy bias.p_of_k_a (
Pk2D
,str
orNone
) – 3D Power spectrum to integrate. If a string, it must correspond to one of the non-linear power spectra stored in cosmo (e.g. ‘delta_matter:delta_matter’).use_spline (
bool
) – switch that determines whether the RSD correlation function is calculated using global splines of multipoles.
- Returns:
Value(s) of the correlation function at the input pi and sigma.