CLASS MANUAL
|
Functions | |
int | lensing_cl_at_l (struct lensing *ple, int l, double *cl_lensed) |
int | lensing_init (struct precision *ppr, struct perturbations *ppt, struct harmonic *phr, struct fourier *pfo, struct lensing *ple) |
int | lensing_free (struct lensing *ple) |
int | lensing_indices (struct precision *ppr, struct harmonic *phr, struct lensing *ple) |
int | lensing_lensed_cl_tt (double *ksi, double **d00, double *w8, int nmu, struct lensing *ple) |
int | lensing_addback_cl_tt (struct lensing *ple, double *cl_tt) |
int | lensing_lensed_cl_te (double *ksiX, double **d20, double *w8, int nmu, struct lensing *ple) |
int | lensing_addback_cl_te (struct lensing *ple, double *cl_te) |
int | lensing_lensed_cl_ee_bb (double *ksip, double *ksim, double **d22, double **d2m2, double *w8, int nmu, struct lensing *ple) |
int | lensing_addback_cl_ee_bb (struct lensing *ple, double *cl_ee, double *cl_bb) |
int | lensing_d00 (double *mu, int num_mu, int lmax, double **d00) |
int | lensing_d11 (double *mu, int num_mu, int lmax, double **d11) |
int | lensing_d1m1 (double *mu, int num_mu, int lmax, double **d1m1) |
int | lensing_d2m2 (double *mu, int num_mu, int lmax, double **d2m2) |
int | lensing_d22 (double *mu, int num_mu, int lmax, double **d22) |
int | lensing_d20 (double *mu, int num_mu, int lmax, double **d20) |
int | lensing_d31 (double *mu, int num_mu, int lmax, double **d31) |
int | lensing_d3m1 (double *mu, int num_mu, int lmax, double **d3m1) |
int | lensing_d3m3 (double *mu, int num_mu, int lmax, double **d3m3) |
int | lensing_d40 (double *mu, int num_mu, int lmax, double **d40) |
int | lensing_d4m2 (double *mu, int num_mu, int lmax, double **d4m2) |
int | lensing_d4m4 (double *mu, int num_mu, int lmax, double **d4m4) |
Documented lensing module
Simon Prunet and Julien Lesgourgues, 6.12.2010
This module computes the lensed temperature and polarization anisotropy power spectra 's given the unlensed temperature, polarization and lensing potential spectra.
Follows Challinor and Lewis full-sky method, astro-ph/0502425
The following functions can be called from other modules:
int lensing_cl_at_l | ( | struct lensing * | ple, |
int | l, | ||
double * | cl_lensed | ||
) |
Anisotropy power spectra 's for all types, modes and initial conditions. SO FAR: ONLY SCALAR
This routine evaluates all the lensed 's at a given value of l by picking it in the pre-computed table. When relevant, it also sums over all initial conditions for each mode, and over all modes.
This function can be called from whatever module at whatever time, provided that lensing_init() has been called before, and lensing_free() has not been called yet.
ple | Input: pointer to lensing structure |
l | Input: multipole number |
cl_lensed | Output: lensed 's for all types (TT, TE, EE, etc..) |
int lensing_init | ( | struct precision * | ppr, |
struct perturbations * | ppt, | ||
struct harmonic * | phr, | ||
struct fourier * | pfo, | ||
struct lensing * | ple | ||
) |
This routine initializes the lensing structure (in particular, computes table of lensed anisotropy spectra )
ppr | Input: pointer to precision structure |
ppt | Input: pointer to perturbation structure (just in case, not used in current version...) |
phr | Input: pointer to harmonic structure |
pfo | Input: pointer to fourier structure |
ple | Output: pointer to initialized lensing structure |
Summary:
int lensing_free | ( | struct lensing * | ple | ) |
This routine frees all the memory space allocated by lensing_init().
To be called at the end of each run, only when no further calls to lensing_cl_at_l() are needed.
ple | Input: pointer to lensing structure (which fields must be freed) |
This routine defines indices and allocates tables in the lensing structure
ppr | Input: pointer to precision structure |
phr | Input: pointer to harmonic structure |
ple | Input/output: pointer to lensing structure |
int lensing_lensed_cl_tt | ( | double * | ksi, |
double ** | d00, | ||
double * | w8, | ||
int | nmu, | ||
struct lensing * | ple | ||
) |
This routine computes the lensed power spectra by Gaussian quadrature
ksi | Input: Lensed correlation function (ksi[index_mu]) |
d00 | Input: Legendre polynomials ( [l][index_mu]) |
w8 | Input: Legendre quadrature weights (w8[index_mu]) |
nmu | Input: Number of quadrature points (0<=index_mu<=nmu) |
ple | Input/output: Pointer to the lensing structure |
Integration by Gauss-Legendre quadrature.
int lensing_addback_cl_tt | ( | struct lensing * | ple, |
double * | cl_tt | ||
) |
This routine adds back the unlensed power spectrum Used in case of fast (and BB inaccurate) integration of correlation functions.
ple | Input/output: Pointer to the lensing structure |
cl_tt | Input: Array of unlensed power spectrum |
int lensing_lensed_cl_te | ( | double * | ksiX, |
double ** | d20, | ||
double * | w8, | ||
int | nmu, | ||
struct lensing * | ple | ||
) |
This routine computes the lensed power spectra by Gaussian quadrature
ksiX | Input: Lensed correlation function (ksiX[index_mu]) |
d20 | Input: Wigner d-function ( [l][index_mu]) |
w8 | Input: Legendre quadrature weights (w8[index_mu]) |
nmu | Input: Number of quadrature points (0<=index_mu<=nmu) |
ple | Input/output: Pointer to the lensing structure |
Integration by Gauss-Legendre quadrature.
int lensing_addback_cl_te | ( | struct lensing * | ple, |
double * | cl_te | ||
) |
This routine adds back the unlensed power spectrum Used in case of fast (and BB inaccurate) integration of correlation functions.
ple | Input/output: Pointer to the lensing structure |
cl_te | Input: Array of unlensed power spectrum |
int lensing_lensed_cl_ee_bb | ( | double * | ksip, |
double * | ksim, | ||
double ** | d22, | ||
double ** | d2m2, | ||
double * | w8, | ||
int | nmu, | ||
struct lensing * | ple | ||
) |
This routine computes the lensed power spectra by Gaussian quadrature
ksip | Input: Lensed correlation function (ksi+[index_mu]) |
ksim | Input: Lensed correlation function (ksi-[index_mu]) |
d22 | Input: Wigner d-function ( [l][index_mu]) |
d2m2 | Input: Wigner d-function ( [l][index_mu]) |
w8 | Input: Legendre quadrature weights (w8[index_mu]) |
nmu | Input: Number of quadrature points (0<=index_mu<=nmu) |
ple | Input/output: Pointer to the lensing structure |
Integration by Gauss-Legendre quadrature.
int lensing_addback_cl_ee_bb | ( | struct lensing * | ple, |
double * | cl_ee, | ||
double * | cl_bb | ||
) |
This routine adds back the unlensed , power spectra Used in case of fast (and BB inaccurate) integration of correlation functions.
ple | Input/output: Pointer to the lensing structure |
cl_ee | Input: Array of unlensed power spectrum |
cl_bb | Input: Array of unlensed power spectrum |
int lensing_d00 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d00 | ||
) |
This routine computes the d00 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d00 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d11 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d11 | ||
) |
This routine computes the d11 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d11 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d1m1 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d1m1 | ||
) |
This routine computes the d1m1 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d1m1 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d2m2 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d2m2 | ||
) |
This routine computes the d2m2 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d2m2 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d22 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d22 | ||
) |
This routine computes the d22 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d22 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d20 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d20 | ||
) |
This routine computes the d20 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d20 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d31 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d31 | ||
) |
This routine computes the d31 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d31 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d3m1 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d3m1 | ||
) |
This routine computes the d3m1 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d3m1 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d3m3 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d3m3 | ||
) |
This routine computes the d3m3 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d3m3 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d40 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d40 | ||
) |
This routine computes the d40 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d40 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d4m2 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d4m2 | ||
) |
This routine computes the d4m2 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d4m2 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003
int lensing_d4m4 | ( | double * | mu, |
int | num_mu, | ||
int | lmax, | ||
double ** | d4m4 | ||
) |
This routine computes the d4m4 term
mu | Input: Vector of cos(beta) values |
num_mu | Input: Number of cos(beta) values |
lmax | Input: maximum multipole |
d4m4 | Input/output: Result is stored here |
Wigner d-functions, computed by recurrence actual recurrence on for stability Formulae from Kostelec & Rockmore 2003