#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_integration.h>
#include "ccl.h"

Include dependency graph for ccl_cls.c:

Classes
struct	IntLensPar

struct	IntMagPar

struct	IntClPar

Macros
#define	CCL_FRAC_RELEVANT 5E-4

Functions
static void	get_support_interval (int n, double x, double y, double frac, double xmin_out, double xmax_out)

static double	speval_bis (double x, void *params)

void	ccl_cl_workspace_free (CCL_ClWorkspace *w)

CCL_ClWorkspace *	ccl_cl_workspace_new (int lmax, int l_limber, double l_logstep, int l_linstep, int *status)

CCL_ClWorkspace *	ccl_cl_workspace_new_limber (int lmax, double l_logstep, int l_linstep, int *status)

static double	integrand_wl (double chip, void *params)

static int	window_lensing (double chi, ccl_cosmology cosmo, SplPar spl_pz, double chi_max, double *win)

static double	integrand_mag (double chip, void *params)

static int	window_magnification (double chi, ccl_cosmology cosmo, SplPar spl_pz, SplPar spl_sz, double chi_max, double win)

static void	clt_init_nz (CCL_ClTracer clt, ccl_cosmology cosmo, int nz_n, double z_n, double n, int *status)

static void	clt_init_bz (CCL_ClTracer clt, ccl_cosmology cosmo, int nz_b, double z_b, double b, int *status)

static void	clt_init_wM (CCL_ClTracer clt, ccl_cosmology cosmo, int nz_s, double z_s, double s, int *status)

static void	clt_nc_init (CCL_ClTracer clt, ccl_cosmology cosmo, int has_rsd, int has_magnification, int nz_n, double z_n, double n, int nz_b, double z_b, double b, int nz_s, double z_s, double s, int *status)

static void	clt_init_wL (CCL_ClTracer clt, ccl_cosmology cosmo, int *status)

static void	clt_init_rf (CCL_ClTracer clt, ccl_cosmology cosmo, int nz_rf, double z_rf, double rf, int *status)

static void	clt_init_ba (CCL_ClTracer clt, ccl_cosmology cosmo, int nz_ba, double z_ba, double ba, int *status)

static void	clt_wl_init (CCL_ClTracer clt, ccl_cosmology cosmo, int has_intrinsic_alignment, int nz_n, double z_n, double n, int nz_ba, double z_ba, double ba, int nz_rf, double z_rf, double rf, int *status)

static CCL_ClTracer *	cl_tracer (ccl_cosmology cosmo, int tracer_type, int has_rsd, int has_magnification, int has_intrinsic_alignment, int nz_n, double z_n, double n, int nz_b, double z_b, double b, int nz_s, double z_s, double s, int nz_ba, double z_ba, double ba, int nz_rf, double z_rf, double rf, double z_source, int status)

CCL_ClTracer *	ccl_cl_tracer (ccl_cosmology cosmo, int tracer_type, int has_rsd, int has_magnification, int has_intrinsic_alignment, int nz_n, double z_n, double n, int nz_b, double z_b, double b, int nz_s, double z_s, double s, int nz_ba, double z_ba, double ba, int nz_rf, double z_rf, double rf, double z_source, int status)

void	ccl_cl_tracer_free (CCL_ClTracer *clt)

CCL_ClTracer *	ccl_cl_tracer_cmblens (ccl_cosmology cosmo, double z_source, int status)

CCL_ClTracer *	ccl_cl_tracer_number_counts (ccl_cosmology cosmo, int has_rsd, int has_magnification, int nz_n, double z_n, double n, int nz_b, double z_b, double b, int nz_s, double z_s, double s, int status)

CCL_ClTracer *	ccl_cl_tracer_number_counts_simple (ccl_cosmology cosmo, int nz_n, double z_n, double n, int nz_b, double z_b, double b, int status)

CCL_ClTracer *	ccl_cl_tracer_lensing (ccl_cosmology cosmo, int has_alignment, int nz_n, double z_n, double n, int nz_ba, double z_ba, double ba, int nz_rf, double z_rf, double rf, int status)

CCL_ClTracer *	ccl_cl_tracer_lensing_simple (ccl_cosmology cosmo, int nz_n, double z_n, double n, int status)

static double	f_dens (double a, ccl_cosmology cosmo, CCL_ClTracer clt, int *status)

static double	f_rsd (double a, ccl_cosmology cosmo, CCL_ClTracer clt, int *status)

static double	f_mag (double a, double chi, ccl_cosmology cosmo, CCL_ClTracer clt, int *status)

static double	transfer_nc (int l, double k, ccl_cosmology cosmo, CCL_ClWorkspace w, CCL_ClTracer clt, int status)

static double	f_lensing (double a, double chi, ccl_cosmology cosmo, CCL_ClTracer clt, int *status)

static double	f_IA_NLA (double a, double chi, ccl_cosmology cosmo, CCL_ClTracer clt, int *status)

static double	transfer_wl (int l, double k, ccl_cosmology cosmo, CCL_ClWorkspace w, CCL_ClTracer clt, int status)

static double	transfer_cmblens (int l, double k, ccl_cosmology cosmo, CCL_ClTracer clt, int *status)

static double	transfer_wrap (int il, double lk, ccl_cosmology cosmo, CCL_ClWorkspace w, CCL_ClTracer clt, int status)

static double	cl_integrand (double lk, void *params)

static void	get_k_interval (ccl_cosmology cosmo, CCL_ClWorkspace w, CCL_ClTracer clt1, CCL_ClTracer clt2, int l, double lkmin, double lkmax)

static double	ccl_angular_cl_native (ccl_cosmology cosmo, CCL_ClWorkspace cw, int il, CCL_ClTracer clt1, CCL_ClTracer clt2, int *status)

void	ccl_angular_cls (ccl_cosmology cosmo, CCL_ClWorkspace w, CCL_ClTracer clt1, CCL_ClTracer clt2, int nl_out, int l_out, double cl_out, int *status)

static int	check_clt_fa_inconsistency (CCL_ClTracer *clt, int func_code)

double	ccl_get_tracer_fa (ccl_cosmology cosmo, CCL_ClTracer clt, double a, int func_code, int *status)

int	ccl_get_tracer_fas (ccl_cosmology cosmo, CCL_ClTracer clt, int na, double a, double fa, int func_code, int *status)

Macro Definition Documentation

#define CCL_FRAC_RELEVANT 5E-4

Definition at line 15 of file ccl_cls.c.

Referenced by clt_init_nz().

Function Documentation

static double ccl_angular_cl_native	(	ccl_cosmology *	cosmo,
		CCL_ClWorkspace *	cw,
		int	il,
		CCL_ClTracer *	clt1,
		CCL_ClTracer *	clt2,
		int *	status
	)

static

Definition at line 889 of file ccl_cls.c.

References ccl_check_status(), ccl_cosmology_set_status_message(), CCL_ERROR_INTEG, ccl_gsl, ccl_raise_gsl_warning(), cl_integrand(), IntClPar::clt1, IntClPar::clt2, cl_cmbl_bm::cosmo, IntClPar::cosmo, get_k_interval(), IntClPar::il, ccl_gsl_params::INTEGRATION_LIMBER_EPSREL, ccl_gsl_params::INTEGRATION_LIMBER_GAUSS_KRONROD_POINTS, CCL_ClWorkspace::l_arr, ccl_gsl_params::N_ITERATION, IntClPar::status, w, and IntClPar::w.

Referenced by ccl_angular_cls().

 {
   int clastatus=0, gslstatus;
   IntClPar ipar;
   double result=0,eresult;
   double lkmin,lkmax;
   gsl_function F;
   gsl_integration_workspace *w=gsl_integration_workspace_alloc(ccl_gsl->N_ITERATION);
 
   ipar.il=il;
   ipar.cosmo=cosmo;
   ipar.w=cw;
   ipar.clt1=clt1;
   ipar.clt2=clt2;
   ipar.status = &clastatus;
   F.function=&cl_integrand;
   F.params=&ipar;
   get_k_interval(cosmo,cw,clt1,clt2,cw->l_arr[il],&lkmin,&lkmax);
   // This computes the angular power spectra in the Limber approximation between two quantities a and b:
   //  C_ell^ab = 2/(2*ell+1) * Integral[ Delta^a_ell(k) Delta^b_ell(k) * P(k) , k_min < k < k_max ]
   // Note that we use log10(k) as an integration variable, and the ell-dependent prefactor is included
   // at the end of this function.
   gslstatus=gsl_integration_qag(&F, lkmin, lkmax, 0,
                                 ccl_gsl->INTEGRATION_LIMBER_EPSREL, ccl_gsl->N_ITERATION,
                                 ccl_gsl->INTEGRATION_LIMBER_GAUSS_KRONROD_POINTS,
                                 w, &result, &eresult);
   gsl_integration_workspace_free(w);
 
   // Test if a round-off error occured in the evaluation of the integral
   // If so, try another integration function, more robust but potentially slower
   if(gslstatus == GSL_EROUND) {
     ccl_raise_gsl_warning(gslstatus, "ccl_cls.c: ccl_angular_cl_native(): Default GSL integration failure, attempting backup method.");
     gsl_integration_cquad_workspace *w_cquad= gsl_integration_cquad_workspace_alloc (ccl_gsl->N_ITERATION);
     size_t nevals=0;
     gslstatus=gsl_integration_cquad(&F, lkmin, lkmax, 0,
                     ccl_gsl->INTEGRATION_LIMBER_EPSREL,
                     w_cquad, &result, &eresult, &nevals);
     gsl_integration_cquad_workspace_free(w_cquad);
   }
   if(gslstatus!=GSL_SUCCESS || *ipar.status) {
     ccl_raise_gsl_warning(gslstatus, "ccl_cls.c: ccl_angular_cl_native():");
     // If an error status was already set, don't overwrite it.
     if(*status == 0){
         *status=CCL_ERROR_INTEG;
         ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_angular_cl_native(): error integrating over k\n");
     }
     return -1;
   }
   ccl_check_status(cosmo,status);
 
   return result*M_LN10/(cw->l_arr[il]+0.5);
 }

void ccl_angular_cls	(	ccl_cosmology *	cosmo,
		CCL_ClWorkspace *	w,
		CCL_ClTracer *	clt1,
		CCL_ClTracer *	clt2,
		int	nl_out,
		int *	l,
		double *	cl,
		int *	status
	)

Computes limber or non-limber power spectrum for two different tracers

Parameters

cosmo	Cosmological parameters
w	a ClWorkspace
clt1	a Cltracer
clt2	a Cltracer
nl_out	the maximum to ell to compute C_ell
l	an array of ell values
cl	the C_ell output array
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: void

Definition at line 943 of file ccl_cls.c.

References ccl_angular_cl_native(), ccl_check_status(), ccl_cosmology_set_status_message(), CCL_ERROR_MEMORY, CCL_ERROR_SPLINE_EV, ccl_spline_eval(), ccl_spline_free(), ccl_spline_init(), ccl_weak_lensing_tracer, CCL_ClTracer::has_magnification, CCL_ClWorkspace::l_arr, CCL_ClWorkspace::l_limber, CCL_ClWorkspace::lmax, CCL_ClWorkspace::n_ls, and CCL_ClTracer::tracer_type.

Referenced by compare_cls(), compare_corr(), main(), and test_angpow_precision().

 {
   int ii,do_angpow;
   double *l_nodes,*cl_nodes;
   SplPar *spcl_nodes;
   
   //First check if ell range is within workspace
   for(ii=0;ii<nl_out;ii++) {
     if(l_out[ii]>w->lmax) {
       *status=CCL_ERROR_SPLINE_EV;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_angular_cls(); "
          "requested l beyond range allowed by workspace\n");
       return;
     }
   }
 
   if(*status==0) {
     //Allocate array for power spectrum at interpolation nodes
     l_nodes=(double *)malloc(w->n_ls*sizeof(double));
     if(l_nodes==NULL) {
       *status=CCL_ERROR_MEMORY;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_angular_cls(); memory allocation\n");
     }
   }
 
   if(*status==0) {
     cl_nodes=(double *)malloc(w->n_ls*sizeof(double));
     if(cl_nodes==NULL) {
       *status=CCL_ERROR_MEMORY;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_cl_angular_cls(); memory allocation\n");
     }
   }
 
   if(*status==0) {
     for(ii=0;ii<w->n_ls;ii++)
       l_nodes[ii]=(double)(w->l_arr[ii]);
 
     do_angpow=0;
     //Now check if angpow is needed at all
     if(w->l_limber>0) {
       for(ii=0;ii<w->n_ls;ii++) {
     if(w->l_arr[ii]<=w->l_limber)
       do_angpow=1;
       }
     }
 #ifndef HAVE_ANGPOW
     do_angpow=0;
 #endif //HAVE_ANGPOW
   
     //Resort to Limber if we have lensing (this will hopefully only be temporary)
     if(clt1->tracer_type==ccl_weak_lensing_tracer || clt2->tracer_type==ccl_weak_lensing_tracer ||
        clt1->has_magnification || clt2->has_magnification) {
       do_angpow=0;
     }
 
     //Use angpow if non-limber is needed
     if(do_angpow)
       ccl_angular_cls_angpow(cosmo,w,clt1,clt2,cl_nodes,status);
     ccl_check_status(cosmo,status);
   }
 
   if(*status==0) {
     //Compute limber nodes
     for(ii=0;ii<w->n_ls;ii++) {
       if((!do_angpow) || (w->l_arr[ii]>w->l_limber))
     cl_nodes[ii]=ccl_angular_cl_native(cosmo,w,ii,clt1,clt2,status);
     }
 
     //Interpolate into ells requested by user
     spcl_nodes=ccl_spline_init(w->n_ls,l_nodes,cl_nodes,0,0);
     if(spcl_nodes==NULL) {
       *status=CCL_ERROR_MEMORY;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_cl_angular_cls(); memory allocation\n");
     }
   }
   
   if(*status==0) {
     for(ii=0;ii<nl_out;ii++)
       cl_out[ii]=ccl_spline_eval((double)(l_out[ii]),spcl_nodes);
   }
   
   //Cleanup
   ccl_spline_free(spcl_nodes);
   free(cl_nodes);
   free(l_nodes);
 }

CCL_ClTracer* ccl_cl_tracer	(	ccl_cosmology *	cosmo,
		int	tracer_type,
		int	has_rsd,
		int	has_magnification,
		int	has_intrinsic_alignment,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_b,
		double *	z_b,
		double *	b,
		int	nz_s,
		double *	z_s,
		double *	s,
		int	nz_ba,
		double *	z_ba,
		double *	ba,
		int	nz_rf,
		double *	z_rf,
		double *	rf,
		double	z_source,
		int *	status
	)

Constructor for a ClTracer.

Parameters

Tracer_type	pass ccl_number_counts_tracer (number counts), ccl_weak_lensing_tracer (weak lensing) or ccl_cmb_lensing_tracer (CMB lensing)
has_rsd	Set to 1 if you want to compute the RSD contribution to number counts (0 otherwise)
has_magnification	Set to 1 if you want to compute the magnification contribution to number counts (0 otherwise)
has_intrinsic_alignment	Set to 1 if you want to compute the IA contribution to shear
nz_n	Number of bins in z_n and n
z_n	Redshifts for each redshift interval of n
n	Number count of objects per redshift interval (Note: arbitrary normalization - renormalized inside)
nz_b	Number of bins in z_b and b
z_b	Redshifts for each redshift interval of b
b	Clustering bias in each redshift bin
nz_s	Number of bins in z_s and s
z_s	Redshifts for each redshift interval of s
s	Magnification bias in each redshift bin
nz_ba	Number of bins in z_ba and ba
z_ba	Redshifts for each redshift interval of ba
ba	Alignment bias in each redshift bin
nz_rf	Number of bins in z_f and f
z_rf	Redshifts for each redshift interval of rf
rf	Aligned red fraction in each redshift bin
z_source	Redshift of source plane for CMB lensing (z~1100 for CMB lensing).
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: CCL_ClTracer object

Definition at line 575 of file ccl_cls.c.

References ccl_check_status(), and cl_tracer().

Referenced by ccl_cl_tracer_cmblens(), ccl_cl_tracer_lensing(), ccl_cl_tracer_lensing_simple(), ccl_cl_tracer_number_counts(), and ccl_cl_tracer_number_counts_simple().

 {
   CCL_ClTracer *clt=cl_tracer(cosmo,tracer_type,has_rsd,has_magnification,has_intrinsic_alignment,
                   nz_n,z_n,n,nz_b,z_b,b,nz_s,z_s,s,
                   nz_ba,z_ba,ba,nz_rf,z_rf,rf,z_source,status);
   ccl_check_status(cosmo,status);
   return clt;
 }

CCL_ClTracer* ccl_cl_tracer_cmblens	(	ccl_cosmology *	cosmo,
		double	z_source,
		int *	status
	)

Simplified constructor for a CMB lensing ClTracer.

Parameters

z_source	Redshift of source plane (z~1100 for CMB lensing).
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: CCL_ClTracer object

Definition at line 614 of file ccl_cls.c.

References ccl_cl_tracer(), and ccl_cmb_lensing_tracer.

Referenced by compare_cls(), and main().

 {
   return ccl_cl_tracer(cosmo,ccl_cmb_lensing_tracer,
                0,0,0,
                0,NULL,NULL,0,NULL,NULL,0,NULL,NULL,
                0,NULL,NULL,0,NULL,NULL,z_source,status);
 }

void ccl_cl_tracer_free ( CCL_ClTracer * clt )

Destructor for a Cltracer

Parameters

clt	a Cltracer

Returns: void

Definition at line 592 of file ccl_cls.c.

References ccl_number_counts_tracer, ccl_spline_free(), ccl_weak_lensing_tracer, CCL_ClTracer::has_intrinsic_alignment, CCL_ClTracer::has_magnification, CCL_ClTracer::spl_ba, CCL_ClTracer::spl_bz, CCL_ClTracer::spl_nz, CCL_ClTracer::spl_rf, CCL_ClTracer::spl_sz, CCL_ClTracer::spl_wL, CCL_ClTracer::spl_wM, and CCL_ClTracer::tracer_type.

Referenced by compare_cls(), main(), and test_angpow_precision().

 {
   if((clt->tracer_type==ccl_number_counts_tracer) || (clt->tracer_type==ccl_weak_lensing_tracer))
     ccl_spline_free(clt->spl_nz);
 
   if(clt->tracer_type==ccl_number_counts_tracer) {
     ccl_spline_free(clt->spl_bz);
     if(clt->has_magnification) {
       ccl_spline_free(clt->spl_sz);
       ccl_spline_free(clt->spl_wM);
     }
   }
   else if(clt->tracer_type==ccl_weak_lensing_tracer) {
     ccl_spline_free(clt->spl_wL);
     if(clt->has_intrinsic_alignment) {
       ccl_spline_free(clt->spl_ba);
       ccl_spline_free(clt->spl_rf);
     }
   }
   free(clt);
 }

CCL_ClTracer* ccl_cl_tracer_lensing	(	ccl_cosmology *	cosmo,
		int	has_alignment,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_ba,
		double *	z_ba,
		double *	ba,
		int	nz_rf,
		double *	z_rf,
		double *	rf,
		int *	status
	)

Simplified constructor for a lensing ClTracer.

Parameters

has_intrinsic_alignment	Set to 1 if you want to compute the IA contribution to shear
nz_n	Number of bins in z_n and n
z_n	Redshifts for each redshift interval of n
n	Number count of objects per redshift interval (Note: arbitrary normalization - renormalized inside)
nz_ba	Number of bins in z_ba and ba
z_ba	Redshifts for each redshift interval of ba
ba	Alignment bias in each redshift bin
nz_rf	Number of bins in z_f and f
z_rf	Redshifts for each redshift interval of rf
rf	Aligned red fraction in each redshift bin
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: CCL_ClTracer object

Definition at line 642 of file ccl_cls.c.

References ccl_cl_tracer(), and ccl_weak_lensing_tracer.

Referenced by compare_cls(), and compare_corr().

 {
   return ccl_cl_tracer(cosmo,ccl_weak_lensing_tracer,0,0,has_alignment,
                nz_n,z_n,n,-1,NULL,NULL,-1,NULL,NULL,
                nz_ba,z_ba,ba,nz_rf,z_rf,rf,0, status);
 }

CCL_ClTracer* ccl_cl_tracer_lensing_simple	(	ccl_cosmology *	cosmo,
		int	nz_n,
		double *	z_n,
		double *	n,
		int *	status
	)

Simplified constructor for a lensing ClTracer without intrinsic alignment.

Parameters

nz_n	Number of bins in z_n and n
z_n	Redshifts for each redshift interval of n
n	Number count of objects per redshift interval (Note: arbitrary normalization - renormalized inside)
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: CCL_ClTracer object

Definition at line 653 of file ccl_cls.c.

References ccl_cl_tracer(), and ccl_weak_lensing_tracer.

Referenced by compare_cls(), compare_corr(), and main().

 {
   return ccl_cl_tracer(cosmo,ccl_weak_lensing_tracer,0,0,0,
                nz_n,z_n,n,-1,NULL,NULL,-1,NULL,NULL,
                -1,NULL,NULL,-1,NULL,NULL,0, status);
 }

CCL_ClTracer* ccl_cl_tracer_number_counts	(	ccl_cosmology *	cosmo,
		int	has_rsd,
		int	has_magnification,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_b,
		double *	z_b,
		double *	b,
		int	nz_s,
		double *	z_s,
		double *	s,
		int *	status
	)

Simplified constructor for a clustering ClTracer.

Parameters

cosmo	Cosmological parameters
has_rsd	Set to 1 if you want to compute the RSD contribution to number counts (0 otherwise)
has_magnification	Set to 1 if you want to compute the magnification contribution to number counts (0 otherwise)
nz_n	Number of bins in z_n and n
z_n	Redshifts for each redshift interval of n
n	Number count of objects per redshift interval (Note: arbitrary normalization - renormalized inside)
nz_b	Number of bins in z_b and b
z_b	Redshifts for each redshift interval of b
b	Clustering bias in each redshift bin
nz_s	Number of bins in z_s and s
z_s	Redshifts for each redshift interval of s
s	Magnification bias in each redshift bin
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: CCL_ClTracer object

Definition at line 622 of file ccl_cls.c.

References ccl_cl_tracer(), and ccl_number_counts_tracer.

Referenced by main(), and test_angpow_precision().

 {
   return ccl_cl_tracer(cosmo,ccl_number_counts_tracer,has_rsd,has_magnification,0,
                nz_n,z_n,n,nz_b,z_b,b,nz_s,z_s,s,
                -1,NULL,NULL,-1,NULL,NULL,0, status);
 }

CCL_ClTracer* ccl_cl_tracer_number_counts_simple	(	ccl_cosmology *	cosmo,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_b,
		double *	z_b,
		double *	b,
		int *	status
	)

Simplified constructor for a ClTracer without magnification nor RSD.

Parameters

nz_n	Number of bins in z_n and n
z_n	Redshifts for each redshift interval of n
n	Number count of objects per redshift interval (Note: arbitrary normalization - renormalized inside)
nz_b	Number of bins in z_b and b
z_b	Redshifts for each redshift interval of b
b	Clustering bias in each redshift bin
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: CCL_ClTracer object

Definition at line 633 of file ccl_cls.c.

References ccl_cl_tracer(), and ccl_number_counts_tracer.

Referenced by compare_cls(), compare_corr(), and main().

 {
   return ccl_cl_tracer(cosmo,ccl_number_counts_tracer,0,0,0,
                nz_n,z_n,n,nz_b,z_b,b,-1,NULL,NULL,
                -1,NULL,NULL,-1,NULL,NULL,0, status);
 }

void ccl_cl_workspace_free ( CCL_ClWorkspace * w )

Definition at line 59 of file ccl_cls.c.

References CCL_ClWorkspace::l_arr.

Referenced by compare_cls(), compare_corr(), main(), and test_angpow_precision().

 {
   free(w->l_arr);
   free(w);
 }

CCL_ClWorkspace* ccl_cl_workspace_new	(	int	lmax,
		int	l_limber,
		double	l_logstep,
		int	l_linstep,
		int *	status
	)

Definition at line 65 of file ccl_cls.c.

References CCL_ERROR_MEMORY, CCL_MAX, CCL_ClWorkspace::l_arr, CCL_ClWorkspace::l_limber, CCL_ClWorkspace::l_linstep, CCL_ClWorkspace::l_logstep, CCL_ClWorkspace::lmax, CCL_ClWorkspace::n_ls, and w.

Referenced by ccl_cl_workspace_new_limber(), and test_angpow_precision().

 {
   int i_l,l0,increment;
   CCL_ClWorkspace *w=(CCL_ClWorkspace *)malloc(sizeof(CCL_ClWorkspace));
   if(w==NULL)
     *status=CCL_ERROR_MEMORY;
 
   if(*status==0) {
     //Set params
     w->lmax=lmax;
     w->l_limber=l_limber;
     w->l_logstep=l_logstep;
     w->l_linstep=l_linstep;
 
     //Compute number of multipoles
     i_l=0; l0=0;
     increment=CCL_MAX(((int)(l0*(w->l_logstep-1.))),1);
     while((l0 < w->lmax) && (increment < w->l_linstep)) {
       i_l++;
       l0+=increment;
       increment=CCL_MAX(((int)(l0*(w->l_logstep-1))),1);
     }
     increment=w->l_linstep;
     while(l0 < w->lmax) {
       i_l++;
       l0+=increment;
     }
 
     //Allocate array of multipoles
     w->n_ls=i_l+1;
     w->l_arr=(int *)malloc(w->n_ls*sizeof(int));
     if(w->l_arr==NULL)
       *status=CCL_ERROR_MEMORY;
   }
 
   if(*status==0) {
     //Redo the computation above and store values of ell
     i_l=0; l0=0;
     increment=CCL_MAX(((int)(l0*(w->l_logstep-1.))),1);
     while((l0 < w->lmax) && (increment < w->l_linstep)) {
       w->l_arr[i_l]=l0;
       i_l++;
       l0+=increment;
       increment=CCL_MAX(((int)(l0*(w->l_logstep-1))),1);
     }
     increment=w->l_linstep;
     while(l0 < w->lmax) {
       w->l_arr[i_l]=l0;
       i_l++;
       l0+=increment;
     }
     //Don't go further than lmaw
     w->l_arr[w->n_ls-1]=w->lmax;
   }
   
   return w;
 }

CCL_ClWorkspace* ccl_cl_workspace_new_limber	(	int	lmax,
		double	l_logstep,
		int	l_linstep,
		int *	status
	)

Definition at line 124 of file ccl_cls.c.

References ccl_cl_workspace_new().

Referenced by compare_cls(), compare_corr(), and main().

 {
   return ccl_cl_workspace_new(lmax,-1,l_logstep,l_linstep,status);
 }

double ccl_get_tracer_fa	(	ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		double	a,
		int	func_code,
		int *	status
	)

Method to return certain redshift or distance-dependent internal quantities for a given tracer.

Parameters

cosmo	Cosmological parameters
clt	ClTracer object
a	scale factor at which the function is to be evaluated
func_code	integer defining which internal function to evaluate. Choose between: ccl_trf_nz (redshift distribution), ccl_trf_bz (clustering bias), ccl_trf_sz (magnification bias), ccl_trf_rf (aligned fraction), ccl_trf_ba (alignment bias), ccl_trf_wL (weak lensing window function), ccl_trf_wM (magnification window function)
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: interpolated value of the requested function

Definition at line 1048 of file ccl_cls.c.

References ccl_comoving_radial_distance(), ccl_cosmology_set_status_message(), CCL_ERROR_INCONSISTENT, ccl_spline_eval(), ccl_trf_ba, ccl_trf_bz, ccl_trf_nz, ccl_trf_rf, ccl_trf_sz, ccl_trf_wL, ccl_trf_wM, check_clt_fa_inconsistency(), CCL_ClTracer::spl_ba, CCL_ClTracer::spl_bz, CCL_ClTracer::spl_nz, CCL_ClTracer::spl_rf, CCL_ClTracer::spl_sz, CCL_ClTracer::spl_wL, CCL_ClTracer::spl_wM, and x.

 {
   SplPar *spl;
 
   if(check_clt_fa_inconsistency(clt,func_code)) {
     *status=CCL_ERROR_INCONSISTENT;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: inconsistent combination of tracer and internal function to be evaluated");
     return -1;
   }
 
   switch(func_code) {
   case ccl_trf_nz :
     spl=clt->spl_nz;
     break;
   case ccl_trf_bz :
     spl=clt->spl_bz;
     break;
   case ccl_trf_sz :
     spl=clt->spl_sz;
     break;
   case ccl_trf_rf :
     spl=clt->spl_rf;
     break;
   case ccl_trf_ba :
     spl=clt->spl_ba;
     break;
   case ccl_trf_wL :
     spl=clt->spl_wL;
     break;
   case ccl_trf_wM :
     spl=clt->spl_wM;
     break;
   }
 
   double x;
   if((func_code==ccl_trf_wL) || (func_code==ccl_trf_wM))
     x=ccl_comoving_radial_distance(cosmo,a,status); //x-variable is comoving distance for lensing kernels
   else
     x=1./a-1; //x-variable is redshift by default
   
   return ccl_spline_eval(x,spl);
 }

int ccl_get_tracer_fas	(	ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int	na,
		double *	a,
		double *	fa,
		int	func_code,
		int *	status
	)

Method to return certain redshift or distance-dependent internal quantities for a given tracer.

Parameters

cosmo	Cosmological parameters
clt	ClTracer object
na	number of points at which the function will be evaluated
a	na values of the scale factor at which the function is to be evaluated
fa	output array with na values that will store the interpolated function values
func_code	integer defining which internal function to evaluate. Choose between: ccl_trf_nz (redshift distribution), ccl_trf_bz (clustering bias), ccl_trf_sz (magnification bias), ccl_trf_rf (aligned fraction), ccl_trf_ba (alignment bias), ccl_trf_wL (weak lensing window function), ccl_trf_wM (magnification window function)
status	Status flag. 0 if there are no errors, nonzero otherwise. For specific cases see documentation for ccl_error.c

Returns: void

Definition at line 1091 of file ccl_cls.c.

References ccl_comoving_radial_distance(), ccl_cosmology_set_status_message(), CCL_ERROR_INCONSISTENT, ccl_spline_eval(), ccl_trf_ba, ccl_trf_bz, ccl_trf_nz, ccl_trf_rf, ccl_trf_sz, ccl_trf_wL, ccl_trf_wM, check_clt_fa_inconsistency(), CCL_ClTracer::spl_ba, CCL_ClTracer::spl_bz, CCL_ClTracer::spl_nz, CCL_ClTracer::spl_rf, CCL_ClTracer::spl_sz, CCL_ClTracer::spl_wL, CCL_ClTracer::spl_wM, and x.

Referenced by main().

 {
   SplPar *spl;
 
   if(check_clt_fa_inconsistency(clt,func_code)) {
     *status=CCL_ERROR_INCONSISTENT;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: inconsistent combination of tracer and internal function to be evaluated");
     return -1;
   }
   
   switch(func_code) {
   case ccl_trf_nz :
     spl=clt->spl_nz;
     break;
   case ccl_trf_bz :
     spl=clt->spl_bz;
     break;
   case ccl_trf_sz :
     spl=clt->spl_sz;
     break;
   case ccl_trf_rf :
     spl=clt->spl_rf;
     break;
   case ccl_trf_ba :
     spl=clt->spl_ba;
     break;
   case ccl_trf_wL :
     spl=clt->spl_wL;
     break;
   case ccl_trf_wM :
     spl=clt->spl_wM;
     break;
   }
   
   int compchi = (func_code==ccl_trf_wL) || (func_code==ccl_trf_wM);
 
   int ia;
   for(ia=0;ia<na;ia++) {
     double x;
     if(compchi) //x-variable is comoving distance for lensing kernels
       x=ccl_comoving_radial_distance(cosmo,a[ia],status);
     else //x-variable is redshift by default
       x=1./a[ia]-1;
     fa[ia]=ccl_spline_eval(x,spl);
   }
 
   return 0;
 }

static int check_clt_fa_inconsistency	(	CCL_ClTracer *	clt,
		int	func_code
	)

static

Definition at line 1032 of file ccl_cls.c.

References ccl_cmb_lensing_tracer, ccl_number_counts_tracer, ccl_trf_ba, ccl_trf_bz, ccl_trf_nz, ccl_trf_rf, ccl_trf_sz, ccl_trf_wL, ccl_trf_wM, ccl_weak_lensing_tracer, CCL_ClTracer::has_intrinsic_alignment, CCL_ClTracer::has_magnification, and CCL_ClTracer::tracer_type.

Referenced by ccl_get_tracer_fa(), and ccl_get_tracer_fas().

 {
   if(((func_code==ccl_trf_nz) && (clt->tracer_type==ccl_cmb_lensing_tracer)) || //lensing has no n(z)
      (((func_code==ccl_trf_bz) || (func_code==ccl_trf_sz) || (func_code==ccl_trf_wM)) &&
       (clt->tracer_type!=ccl_number_counts_tracer)) || //bias and magnification only for clustering
      (((func_code==ccl_trf_rf) || (func_code==ccl_trf_ba) || (func_code==ccl_trf_wL)) &&
       (clt->tracer_type!=ccl_weak_lensing_tracer))) //IAs only for weak lensing
     return 1;
   if((((func_code==ccl_trf_sz) || (func_code==ccl_trf_wM)) &&
       (clt->has_magnification==0)) || //correct combination, but no magnification
      (((func_code==ccl_trf_rf) || (func_code==ccl_trf_ba)) &&
       (clt->has_intrinsic_alignment==0))) //Correct combination, but no IAs
     return 1;
   return 0;
 }

static double cl_integrand	(	double	lk,
		void *	params
	)

static

Definition at line 823 of file ccl_cls.c.

References ccl_nonlin_matter_power(), ccl_scale_factor_of_chi(), IntClPar::clt1, IntClPar::clt2, IntClPar::cosmo, IntClPar::il, CCL_ClWorkspace::l_arr, p, pow(), IntClPar::status, transfer_wrap(), and IntClPar::w.

Referenced by ccl_angular_cl_native().

 {
   double d1,d2;
   IntClPar *p=(IntClPar *)params;
   d1=transfer_wrap(p->il,lk,p->cosmo,p->w,p->clt1,p->status);
   if(d1==0)
     return 0;
   d2=transfer_wrap(p->il,lk,p->cosmo,p->w,p->clt2,p->status);
   if(d2==0)
     return 0;
 
   double k=pow(10.,lk);
   double chi=(p->w->l_arr[p->il]+0.5)/k;
   double a=ccl_scale_factor_of_chi(p->cosmo,chi,p->status);
   double pk=ccl_nonlin_matter_power(p->cosmo,k,a,p->status);
   
   return k*pk*d1*d2;
 }

static CCL_ClTracer* cl_tracer	(	ccl_cosmology *	cosmo,
		int	tracer_type,
		int	has_rsd,
		int	has_magnification,
		int	has_intrinsic_alignment,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_b,
		double *	z_b,
		double *	b,
		int	nz_s,
		double *	z_s,
		double *	s,
		int	nz_ba,
		double *	z_ba,
		double *	ba,
		int	nz_rf,
		double *	z_rf,
		double *	rf,
		double	z_source,
		int *	status
	)

static

Definition at line 515 of file ccl_cls.c.

References ccl_cmb_lensing_tracer, ccl_comoving_radial_distance(), ccl_cosmology_set_status_message(), CCL_ERROR_INCONSISTENT, CCL_ERROR_MEMORY, ccl_h_over_h0(), ccl_number_counts_tracer, ccl_weak_lensing_tracer, CCL_ClTracer::chi_source, CCL_ClTracer::chimax, CCL_ClTracer::chimin, CLIGHT_HMPC, clt_nc_init(), clt_wl_init(), ccl_parameters::h, ccl_parameters::Omega_m, ccl_cosmology::params, CCL_ClTracer::prefac_lensing, and CCL_ClTracer::tracer_type.

Referenced by ccl_cl_tracer().

 {
   int clstatus=0;
   CCL_ClTracer *clt=(CCL_ClTracer *)malloc(sizeof(CCL_ClTracer));
   if(clt==NULL) {
     *status=CCL_ERROR_MEMORY;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_cl_tracer(): memory allocation\n");
   }
 
   if(*status==0) {
     clt->tracer_type=tracer_type;
     
     double hub=cosmo->params.h*ccl_h_over_h0(cosmo,1.,status)/CLIGHT_HMPC;
     clt->prefac_lensing=1.5*hub*hub*cosmo->params.Omega_m;
 
     if(tracer_type==ccl_number_counts_tracer)
       clt_nc_init(clt,cosmo,has_rsd,has_magnification,
           nz_n,z_n,n,nz_b,z_b,b,nz_s,z_s,s,status);
     else if(tracer_type==ccl_weak_lensing_tracer)
       clt_wl_init(clt,cosmo,has_intrinsic_alignment,
           nz_n,z_n,n,nz_ba,z_ba,ba,nz_rf,z_rf,rf,status);
     else if(tracer_type==ccl_cmb_lensing_tracer) {
       clt->chi_source=ccl_comoving_radial_distance(cosmo,1./(1+z_source),status);
       clt->chimax=clt->chi_source;
       clt->chimin=0;
     }
     else {
       free(clt);
       *status=CCL_ERROR_INCONSISTENT;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_cl_tracer(): unknown tracer type\n");
       return NULL;
     }
   }
 
   if(*status) {
     free(clt);
     clt=NULL;
   }
     
   return clt;
 }

static void clt_init_ba	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	nz_ba,
		double *	z_ba,
		double *	ba,
		int *	status
	)

static

Definition at line 477 of file ccl_cls.c.

References ccl_cosmology_set_status_message(), CCL_ERROR_SPLINE, ccl_spline_init(), and CCL_ClTracer::spl_ba.

Referenced by clt_wl_init().

 {
   //Initialize bias spline
   clt->spl_ba=ccl_spline_init(nz_ba,z_ba,ba,ba[0],ba[nz_ba-1]);
   if(clt->spl_ba==NULL) {
     *status=CCL_ERROR_SPLINE;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_ba(): error initializing spline for b(z)\n");
   }
 }

static void clt_init_bz	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	nz_b,
		double *	z_b,
		double *	b,
		int *	status
	)

static

Definition at line 320 of file ccl_cls.c.

References ccl_cosmology_set_status_message(), CCL_ERROR_SPLINE, ccl_spline_init(), and CCL_ClTracer::spl_bz.

Referenced by clt_nc_init().

 {
   //Initialize bias spline
   clt->spl_bz=ccl_spline_init(nz_b,z_b,b,b[0],b[nz_b-1]);
   if(clt->spl_bz==NULL) {
     *status=CCL_ERROR_SPLINE;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_bz(): error initializing spline for b(z)\n");
   }
 }

static void clt_init_nz	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	nz_n,
		double *	z_n,
		double *	n,
		int *	status
	)

static

Definition at line 261 of file ccl_cls.c.

References ccl_comoving_radial_distance(), ccl_cosmology_set_status_message(), CCL_ERROR_INTEG, CCL_ERROR_MEMORY, CCL_ERROR_SPLINE, CCL_FRAC_RELEVANT, ccl_gsl, ccl_raise_gsl_warning(), ccl_spline_free(), ccl_spline_init(), CCL_ClTracer::chimax, CCL_ClTracer::chimin, get_support_interval(), ccl_gsl_params::INTEGRATION_EPSREL, ccl_gsl_params::INTEGRATION_GAUSS_KRONROD_POINTS, ccl_gsl_params::N_ITERATION, speval_bis(), CCL_ClTracer::spl_nz, w, CCL_ClTracer::zmax, and CCL_ClTracer::zmin.

Referenced by clt_nc_init(), and clt_wl_init().

 {
   int gslstatus;
   gsl_function F;
   double nz_norm,nz_enorm;
   double *nz_normalized;
   
   //Find redshift range where the N(z) has support
   get_support_interval(nz_n,z_n,n,CCL_FRAC_RELEVANT,&(clt->zmin),&(clt->zmax));
   clt->chimax=ccl_comoving_radial_distance(cosmo,1./(1+clt->zmax),status);
   clt->chimin=ccl_comoving_radial_distance(cosmo,1./(1+clt->zmin),status);
   clt->spl_nz=ccl_spline_init(nz_n,z_n,n,0,0);
   if(clt->spl_nz==NULL) {
     *status=CCL_ERROR_SPLINE;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_nz(): error initializing spline for N(z)\n");
   }
 
   if(*status==0) {
     //Normalize n(z)
     nz_normalized=(double *)malloc(nz_n*sizeof(double));
     if(nz_normalized==NULL) {
       *status=CCL_ERROR_MEMORY;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_nz(): memory allocation\n");
       return;
     }
   }
   
   if(*status==0) {
     gsl_integration_workspace *w=gsl_integration_workspace_alloc(ccl_gsl->N_ITERATION);
     F.function=&speval_bis;
     F.params=clt->spl_nz;
     //Here we're just integrating the N(z) to normalize it to unit probability.
     gslstatus=gsl_integration_qag(&F, z_n[0], z_n[nz_n-1], 0,
                   ccl_gsl->INTEGRATION_EPSREL, ccl_gsl->N_ITERATION,
                   ccl_gsl->INTEGRATION_GAUSS_KRONROD_POINTS,
                   w, &nz_norm, &nz_enorm);
     gsl_integration_workspace_free(w);
     if(gslstatus!=GSL_SUCCESS) {
       ccl_raise_gsl_warning(gslstatus, "ccl_cls.c: clt_init_nz():");
       *status=CCL_ERROR_INTEG;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_nz(): integration error when normalizing N(z)\n");
     }
   }
   
   if(*status==0) {
     for(int ii=0;ii<nz_n;ii++)
       nz_normalized[ii]=n[ii]/nz_norm;
     ccl_spline_free(clt->spl_nz);
     clt->spl_nz=ccl_spline_init(nz_n,z_n,nz_normalized,0,0);
     if(clt->spl_nz==NULL) {
       *status=CCL_ERROR_SPLINE;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_nz(): error initializing normalized spline for N(z)\n");
     }
   }
   
   free(nz_normalized);
 }

static void clt_init_rf	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	nz_rf,
		double *	z_rf,
		double *	rf,
		int *	status
	)

static

Definition at line 466 of file ccl_cls.c.

References ccl_cosmology_set_status_message(), CCL_ERROR_SPLINE, ccl_spline_init(), and CCL_ClTracer::spl_rf.

Referenced by clt_wl_init().

 {
   //Initialize bias spline
   clt->spl_rf=ccl_spline_init(nz_rf,z_rf,rf,rf[0],rf[nz_rf-1]);
   if(clt->spl_rf==NULL) {
     *status=CCL_ERROR_SPLINE;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_rf(): error initializing spline for b(z)\n");
   }
 }

static void clt_init_wL	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int *	status
	)

static

Definition at line 414 of file ccl_cls.c.

References ccl_comoving_radial_distance(), ccl_cosmology_set_status_message(), CCL_ERROR_INTEG, CCL_ERROR_LINSPACE, CCL_ERROR_MEMORY, CCL_ERROR_SPLINE, ccl_linear_spacing(), ccl_spline_init(), CCL_ClTracer::chimin, CCL_ClTracer::spl_nz, CCL_ClTracer::spl_wL, window_lensing(), x, SplPar::xf, and CCL_ClTracer::zmin.

Referenced by clt_wl_init().

 {
   //Compute weak lensing kernel
   int nchi;
   double *x,*y;
   double dchi_here=5.;
   double zmax=clt->spl_nz->xf;
   double chimax=ccl_comoving_radial_distance(cosmo,1./(1+zmax),status);
   //TODO: The interval in chi (5. Mpc) should be made a macro
   
   //In this case we need to integrate all the way to z=0. Reset zmin and chimin
   clt->zmin=0;
   clt->chimin=0;
   nchi=(int)(chimax/dchi_here)+1;
   x=ccl_linear_spacing(0.,chimax,nchi);
   dchi_here=chimax/nchi;
   if(x==NULL || (fabs(x[0]-0)>1E-5) || (fabs(x[nchi-1]-chimax)>1e-5)) {
     *status=CCL_ERROR_LINSPACE;
     ccl_cosmology_set_status_message(cosmo,
                      "ccl_cls.c: clt_init_wL(): Error creating linear spacing in chi\n");
   }
   
   if(*status==0) {
     y=(double *)malloc(nchi*sizeof(double));
     if(y==NULL) {
       *status=CCL_ERROR_MEMORY;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_wL(): memory allocation\n");
     }
   }
 
   if(*status==0) {
     int clstatus=0;
     for(int j=0;j<nchi;j++)
       clstatus|=window_lensing(x[j],cosmo,clt->spl_nz,chimax,&(y[j]));
     if(clstatus) {
       *status=CCL_ERROR_INTEG;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_wL(): error computing lensing window\n");
     }
   }
   
   if(*status==0) {
     clt->spl_wL=ccl_spline_init(nchi,x,y,y[0],0);
     if(clt->spl_wL==NULL) {
       *status=CCL_ERROR_SPLINE;
       ccl_cosmology_set_status_message(cosmo,
                      "ccl_cls.c: clt_init_wL(): error initializing spline for lensing window\n");
     }
   }
   free(x); free(y);
 }

static void clt_init_wM	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	nz_s,
		double *	z_s,
		double *	s,
		int *	status
	)

static

Definition at line 331 of file ccl_cls.c.

References ccl_comoving_radial_distance(), ccl_cosmology_set_status_message(), CCL_ERROR_INTEG, CCL_ERROR_LINSPACE, CCL_ERROR_MEMORY, CCL_ERROR_SPLINE, ccl_linear_spacing(), ccl_spline_init(), CCL_ClTracer::chimin, CCL_ClTracer::spl_nz, CCL_ClTracer::spl_sz, CCL_ClTracer::spl_wM, window_magnification(), x, SplPar::xf, and CCL_ClTracer::zmin.

Referenced by clt_nc_init().

 {
   //Compute magnification kernel
   int nchi;
   double *x,*y;
   double dchi_here=5.;
   double zmax=clt->spl_nz->xf;
   double chimax=ccl_comoving_radial_distance(cosmo,1./(1+zmax),status);
   //TODO: The interval in chi (5. Mpc) should be made a macro
 
   //In this case we need to integrate all the way to z=0. Reset zmin and chimin
   clt->zmin=0;
   clt->chimin=0;
   clt->spl_sz=ccl_spline_init(nz_s,z_s,s,s[0],s[nz_s-1]);
   if(clt->spl_sz==NULL) {
     *status=CCL_ERROR_SPLINE;
     ccl_cosmology_set_status_message(cosmo,
                      "ccl_cls.c: clt_init_wM(): error initializing spline for s(z)\n");
   }
 
   if(*status==0) {
     nchi=(int)(chimax/dchi_here)+1;
     x=ccl_linear_spacing(0.,chimax,nchi);
     dchi_here=chimax/nchi;
     if(x==NULL || (fabs(x[0]-0)>1E-5) || (fabs(x[nchi-1]-chimax)>1e-5)) {
       *status=CCL_ERROR_LINSPACE;
       ccl_cosmology_set_status_message(cosmo,
                        "ccl_cls.c: clt_init_wM(): Error creating linear spacing in chi\n");
     }
   }
 
   if(*status==0) {
     y=(double *)malloc(nchi*sizeof(double));
     if(y==NULL) {
       *status=CCL_ERROR_MEMORY;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_wM(): memory allocation\n");
     }
   }
 
   if(*status==0) {
     int clstatus=0;
     for(int j=0;j<nchi;j++)
       clstatus|=window_magnification(x[j],cosmo,clt->spl_nz,clt->spl_sz,chimax,&(y[j]));
     if(clstatus) {
       *status=CCL_ERROR_INTEG;
       ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: clt_init_wM(): error computing lensing window\n");
     }
   }
 
   if(*status==0) {
     clt->spl_wM=ccl_spline_init(nchi,x,y,y[0],0);
     if(clt->spl_wM==NULL) {
       *status=CCL_ERROR_SPLINE;
       ccl_cosmology_set_status_message(cosmo,
                        "ccl_cls.c: clt_init_wM(): error initializing spline for lensing window\n");
     }
   }
   free(x); free(y);
 }

static void clt_nc_init	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	has_rsd,
		int	has_magnification,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_b,
		double *	z_b,
		double *	b,
		int	nz_s,
		double *	z_s,
		double *	s,
		int *	status
	)

static

Definition at line 393 of file ccl_cls.c.

References ccl_cosmology_set_status_message(), CCL_ERROR_NOT_IMPLEMENTED, clt_init_bz(), clt_init_nz(), clt_init_wM(), CCL_ClTracer::has_magnification, CCL_ClTracer::has_rsd, ccl_parameters::N_nu_mass, and ccl_cosmology::params.

Referenced by cl_tracer().

 {
   clt->has_rsd=has_rsd;
   clt->has_magnification=has_magnification;
 
   if ( ((cosmo->params.N_nu_mass)>0) && clt->has_rsd){
     *status=CCL_ERROR_NOT_IMPLEMENTED;
     ccl_cosmology_set_status_message(cosmo, "ccl_cls.c: ccl_cl_tracer_new(): Number counts tracers with RSD not yet implemented in cosmologies with massive neutrinos.");
     return;
   }
 
   clt_init_nz(clt,cosmo,nz_n,z_n,n,status);
   clt_init_bz(clt,cosmo,nz_b,z_b,b,status);
   if(clt->has_magnification)
     clt_init_wM(clt,cosmo,nz_s,z_s,s,status);
 }

static void clt_wl_init	(	CCL_ClTracer *	clt,
		ccl_cosmology *	cosmo,
		int	has_intrinsic_alignment,
		int	nz_n,
		double *	z_n,
		double *	n,
		int	nz_ba,
		double *	z_ba,
		double *	ba,
		int	nz_rf,
		double *	z_rf,
		double *	rf,
		int *	status
	)

static

Definition at line 488 of file ccl_cls.c.

References clt_init_ba(), clt_init_nz(), clt_init_rf(), clt_init_wL(), and CCL_ClTracer::has_intrinsic_alignment.

Referenced by cl_tracer().

 {
   clt->has_intrinsic_alignment=has_intrinsic_alignment;
 
   clt_init_nz(clt,cosmo,nz_n,z_n,n,status);
   clt_init_wL(clt,cosmo,status);
   if(clt->has_intrinsic_alignment) {
     clt_init_rf(clt,cosmo,nz_rf,z_rf,rf,status);
     clt_init_ba(clt,cosmo,nz_ba,z_ba,ba,status);
   }
 }

static double f_dens	(	double	a,
		ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 661 of file ccl_cls.c.

References ccl_h_over_h0(), ccl_spline_eval(), CLIGHT_HMPC, ccl_test_distances::h, ccl_parameters::h, ccl_cosmology::params, CCL_ClTracer::spl_bz, CCL_ClTracer::spl_nz, and z.

Referenced by transfer_nc().

 {
   double z=1./a-1;
   double pz=ccl_spline_eval(z,clt->spl_nz);
   double bz=ccl_spline_eval(z,clt->spl_bz);
   double h=cosmo->params.h*ccl_h_over_h0(cosmo,a,status)/CLIGHT_HMPC;
 
   return pz*bz*h;
 }

static double f_IA_NLA	(	double	a,
		double	chi,
		ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 736 of file ccl_cls.c.

References ccl_h_over_h0(), ccl_scale_factor_of_chi(), ccl_spline_eval(), CLIGHT_HMPC, ccl_test_distances::h, ccl_parameters::h, ccl_cosmology::params, CCL_ClTracer::spl_ba, CCL_ClTracer::spl_nz, CCL_ClTracer::spl_rf, and z.

Referenced by transfer_wl().

 {
   if(chi<=1E-10)
     return 0;
   else {
     double a=ccl_scale_factor_of_chi(cosmo,chi, status);
     double z=1./a-1;
     double pz=ccl_spline_eval(z,clt->spl_nz);
     double ba=ccl_spline_eval(z,clt->spl_ba);
     double rf=ccl_spline_eval(z,clt->spl_rf);
     double h=cosmo->params.h*ccl_h_over_h0(cosmo,a,status)/CLIGHT_HMPC;
 
     return pz*ba*rf*h/(chi*chi);
   }
 }

static double f_lensing	(	double	a,
		double	chi,
		ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 726 of file ccl_cls.c.

References ccl_spline_eval(), CCL_ClTracer::prefac_lensing, and CCL_ClTracer::spl_wL.

Referenced by transfer_wl().

 {
   double wL=ccl_spline_eval(chi,clt->spl_wL);
 
   if(wL<=0)
     return 0;
   else
     return clt->prefac_lensing*wL/(a*chi);
 }

static double f_mag	(	double	a,
		double	chi,
		ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 681 of file ccl_cls.c.

References ccl_spline_eval(), and CCL_ClTracer::spl_wM.

Referenced by transfer_nc().

 {
   double wM=ccl_spline_eval(chi,clt->spl_wM);
 
   if(wM<=0)
     return 0;
   else
     return wM/(a*chi);
 }

static double f_rsd	(	double	a,
		ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 671 of file ccl_cls.c.

References ccl_growth_rate(), ccl_h_over_h0(), ccl_spline_eval(), CLIGHT_HMPC, ccl_test_distances::h, ccl_parameters::h, ccl_cosmology::params, CCL_ClTracer::spl_nz, and z.

Referenced by transfer_nc().

 {
   double z=1./a-1;
   double pz=ccl_spline_eval(z,clt->spl_nz);
   double fg=ccl_growth_rate(cosmo,a,status);
   double h=cosmo->params.h*ccl_h_over_h0(cosmo,a,status)/CLIGHT_HMPC;
 
   return pz*fg*h;
 }

static void get_k_interval	(	ccl_cosmology *	cosmo,
		CCL_ClWorkspace *	w,
		CCL_ClTracer *	clt1,
		CCL_ClTracer *	clt2,
		int	l,
		double *	lkmin,
		double *	lkmax
	)

static

Definition at line 847 of file ccl_cls.c.

References ccl_number_counts_tracer, ccl_splines, CCL_ClTracer::chimax, CCL_ClTracer::chimin, CCL_ClTracer::has_magnification, ccl_spline_params::K_MAX, ccl_spline_params::K_MIN, and CCL_ClTracer::tracer_type.

Referenced by ccl_angular_cl_native().

 {
   double chimin,chimax;
   int cut_low_1=0,cut_low_2=0;
   
   //Define a minimum distance only if no lensing is needed
   if((clt1->tracer_type==ccl_number_counts_tracer) && (clt1->has_magnification==0)) cut_low_1=1;
   if((clt2->tracer_type==ccl_number_counts_tracer) && (clt2->has_magnification==0)) cut_low_2=1;
   
   if(cut_low_1) {
     if(cut_low_2) {
       chimin=fmax(clt1->chimin,clt2->chimin);
       chimax=fmin(clt1->chimax,clt2->chimax);
     }
     else {
       chimin=clt1->chimin;
       chimax=clt1->chimax;
     }
   }
   else if(cut_low_2) {
     chimin=clt2->chimin;
     chimax=clt2->chimax;
   }
   else {
     chimin=0.5*(l+0.5)/ccl_splines->K_MAX;
     chimax=2*(l+0.5)/ccl_splines->K_MIN;
   }
   
   if(chimin<=0)
     chimin=0.5*(l+0.5)/ccl_splines->K_MAX;
   
   *lkmax=log10(fmin( ccl_splines->K_MAX  ,2  *(l+0.5)/chimin));
   *lkmin=log10(fmax( ccl_splines->K_MIN  ,0.5*(l+0.5)/chimax));
 }

static void get_support_interval	(	int	n,
		double *	x,
		double *	y,
		double	frac,
		double *	xmin_out,
		double *	xmax_out
	)

static

Definition at line 19 of file ccl_cls.c.

Referenced by clt_init_nz().

 {
   int ix;
   double ythr=-1000;
 
   //Initialize as the original edges in case we don't find an interval
   *xmin_out=x[0];
   *xmax_out=x[n-1];
 
   //Find threshold
   for(ix=0;ix<n;ix++) {
     if(y[ix]>ythr) ythr=y[ix];
   }
   ythr*=frac;
 
   //Find minimum
   for(ix=0;ix<n;ix++) {
     if(y[ix]>=ythr) {
       *xmin_out=x[ix];
       break;
     }
   }
 
   //Find maximum
   for(ix=n-1;ix>=0;ix--) {
     if(y[ix]>=ythr) {
       *xmax_out=x[ix];
       break;
     }
   }
 }

static double integrand_mag	(	double	chip,
		void *	params
	)

static

Definition at line 202 of file ccl_cls.c.

References ccl_h_over_h0(), ccl_scale_factor_of_chi(), ccl_sinn(), ccl_spline_eval(), IntMagPar::chi, CLIGHT_HMPC, IntMagPar::cosmo, ccl_test_distances::h, ccl_parameters::h, p, ccl_cosmology::params, IntMagPar::spl_pz, IntMagPar::spl_sz, IntMagPar::status, mk_bins::sz, and z.

Referenced by window_magnification().

 {
   IntMagPar *p=(IntMagPar *)params;
   double chi=p->chi;
   double a=ccl_scale_factor_of_chi(p->cosmo,chip, p->status);
   double z=1./a-1;
   double pz=ccl_spline_eval(z,p->spl_pz);
   double sz=ccl_spline_eval(z,p->spl_sz);
   double h=p->cosmo->params.h*ccl_h_over_h0(p->cosmo,a, p->status)/CLIGHT_HMPC;
 
   if(chi==0)
     return h*pz*(1-2.5*sz);
   else
     return h*pz*(1-2.5*sz)*ccl_sinn(p->cosmo,chip-chi,p->status)/ccl_sinn(p->cosmo,chip,p->status);
 }

static double integrand_wl	(	double	chip,
		void *	params
	)

static

Definition at line 138 of file ccl_cls.c.

References ccl_h_over_h0(), ccl_scale_factor_of_chi(), ccl_sinn(), ccl_spline_eval(), IntLensPar::chi, CLIGHT_HMPC, IntLensPar::cosmo, ccl_test_distances::h, ccl_parameters::h, p, ccl_cosmology::params, IntLensPar::spl_pz, IntLensPar::status, and z.

Referenced by window_lensing().

 {
   IntLensPar *p=(IntLensPar *)params;
   double chi=p->chi;
   double a=ccl_scale_factor_of_chi(p->cosmo,chip, p->status);
   double z=1./a-1;
   double pz=ccl_spline_eval(z,p->spl_pz);
   double h=p->cosmo->params.h*ccl_h_over_h0(p->cosmo,a, p->status)/CLIGHT_HMPC;
 
   if(chi==0)
     return h*pz;
   else
     return h*pz*ccl_sinn(p->cosmo,chip-chi,p->status)/ccl_sinn(p->cosmo,chip,p->status);
 }

static double speval_bis	(	double	x,
		void *	params
	)

static

Definition at line 53 of file ccl_cls.c.

References ccl_spline_eval().

Referenced by clt_init_nz().

 {
   return ccl_spline_eval(x,(SplPar *)params);
 }

static double transfer_cmblens	(	int	l,
		double	k,
		ccl_cosmology *	cosmo,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 776 of file ccl_cls.c.

References ccl_scale_factor_of_chi(), CCL_ClTracer::chi_source, CCL_ClTracer::prefac_lensing, and w.

Referenced by transfer_wrap().

 {
   double chi=(l+0.5)/k;
   if(chi>=clt->chi_source)
     return 0;
 
   if(chi<=clt->chimax) {
     double a=ccl_scale_factor_of_chi(cosmo,chi,status);
     double w=1-chi/clt->chi_source;
     return clt->prefac_lensing*l*(l+1.)*w/(a*chi*k*k);
   }
   return 0;
 }

static double transfer_nc	(	int	l,
		double	k,
		ccl_cosmology *	cosmo,
		CCL_ClWorkspace *	w,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 697 of file ccl_cls.c.

References ccl_nonlin_matter_power(), ccl_scale_factor_of_chi(), f_dens(), f_mag(), f_rsd(), CCL_ClTracer::has_magnification, CCL_ClTracer::has_rsd, CCL_ClTracer::prefac_lensing, and sqrt().

Referenced by transfer_wrap().

 {
   double ret=0;
   double x0=(l+0.5);
   double chi0=x0/k;
   if(chi0<=clt->chimax) {
     double a0=ccl_scale_factor_of_chi(cosmo,chi0,status);
     double f_all=f_dens(a0,cosmo,clt,status);
     if(clt->has_rsd) {
       double x1=(l+1.5);
       double chi1=x1/k;
       if(chi1<=clt->chimax) {
     double a1=ccl_scale_factor_of_chi(cosmo,chi1,status);
     double pk0=ccl_nonlin_matter_power(cosmo,k,a0,status);
     double pk1=ccl_nonlin_matter_power(cosmo,k,a1,status);
     double fg0=f_rsd(a0,cosmo,clt,status);
     double fg1=f_rsd(a1,cosmo,clt,status);
     f_all+=fg0*(1.-l*(l-1.)/(x0*x0))-fg1*2.*sqrt((l+0.5)*pk1/((l+1.5)*pk0))/x1;
       }
     }
     if(clt->has_magnification)
       f_all+=-2*clt->prefac_lensing*l*(l+1)*f_mag(a0,chi0,cosmo,clt,status)/(k*k);
     ret=f_all;
   }
 
   return ret;
 }

static double transfer_wl	(	int	l,
		double	k,
		ccl_cosmology *	cosmo,
		CCL_ClWorkspace *	w,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 758 of file ccl_cls.c.

References ccl_scale_factor_of_chi(), f_IA_NLA(), f_lensing(), CCL_ClTracer::has_intrinsic_alignment, and sqrt().

Referenced by transfer_wrap().

 {
   double ret=0;
   double chi=(l+0.5)/k;
   if(chi<=clt->chimax) {
     double a=ccl_scale_factor_of_chi(cosmo,chi,status);
     double f_all=f_lensing(a,chi,cosmo,clt,status);
     if(clt->has_intrinsic_alignment)
       f_all+=f_IA_NLA(a,chi,cosmo,clt,status);
     
     ret=f_all;
   }
 
   return sqrt((l+2.)*(l+1.)*l*(l-1.))*ret/(k*k);
   //return (l+1.)*l*ret/(k*k);
 }

static double transfer_wrap	(	int	il,
		double	lk,
		ccl_cosmology *	cosmo,
		CCL_ClWorkspace *	w,
		CCL_ClTracer *	clt,
		int *	status
	)

static

Definition at line 795 of file ccl_cls.c.

References ccl_cmb_lensing_tracer, ccl_number_counts_tracer, ccl_weak_lensing_tracer, CCL_ClWorkspace::l_arr, pow(), CCL_ClTracer::tracer_type, transfer_cmblens(), transfer_nc(), and transfer_wl().

Referenced by cl_integrand().

 {
   double transfer_out=0;
   double k=pow(10.,lk);
 
   if(clt->tracer_type==ccl_number_counts_tracer)
     transfer_out=transfer_nc(w->l_arr[il],k,cosmo,w,clt,status);
   else if(clt->tracer_type==ccl_weak_lensing_tracer)
     transfer_out=transfer_wl(w->l_arr[il],k,cosmo,w,clt,status);
   else if(clt->tracer_type==ccl_cmb_lensing_tracer)
     transfer_out=transfer_cmblens(w->l_arr[il],k,cosmo,clt,status);
   else
     transfer_out=-1;
   return transfer_out;
 }

static int window_lensing	(	double	chi,
		ccl_cosmology *	cosmo,
		SplPar *	spl_pz,
		double	chi_max,
		double *	win
	)

static

Definition at line 159 of file ccl_cls.c.

References ccl_gsl, ccl_raise_gsl_warning(), IntLensPar::chi, cl_cmbl_bm::cosmo, IntLensPar::cosmo, integrand_wl(), ccl_gsl_params::INTEGRATION_EPSREL, ccl_gsl_params::INTEGRATION_GAUSS_KRONROD_POINTS, ccl_gsl_params::N_ITERATION, IntLensPar::spl_pz, IntLensPar::status, and w.

Referenced by clt_init_wL().

 {
   int gslstatus =0, status =0;
   double result,eresult;
   IntLensPar ip;
   gsl_function F;
   gsl_integration_workspace *w=gsl_integration_workspace_alloc(ccl_gsl->N_ITERATION);
 
   ip.chi=chi;
   ip.cosmo=cosmo;
   ip.spl_pz=spl_pz;
   ip.status = &status;
   F.function=&integrand_wl;
   F.params=&ip;
   // This conputes the lensing kernel:
   //   w_L(chi) = Integral[ dN/dchi(chi') * f(chi'-chi)/f(chi') , chi < chi' < chi_horizon ]
   // Where f(chi) is the comoving angular distance (which is just chi for zero curvature).
   gslstatus=gsl_integration_qag(&F, chi, chi_max, 0,
                                 ccl_gsl->INTEGRATION_EPSREL, ccl_gsl->N_ITERATION,
                                 ccl_gsl->INTEGRATION_GAUSS_KRONROD_POINTS,
                                 w, &result, &eresult);
   *win=result;
   gsl_integration_workspace_free(w);
   if(gslstatus!=GSL_SUCCESS || *ip.status) {
     ccl_raise_gsl_warning(gslstatus, "ccl_cls.c: window_lensing():");
     return 1;
   }
   //TODO: chi_max should be changed to chi_horizon
   //we should precompute this quantity and store it in cosmo by default
 
   return 0;
 }

static int window_magnification	(	double	chi,
		ccl_cosmology *	cosmo,
		SplPar *	spl_pz,
		SplPar *	spl_sz,
		double	chi_max,
		double *	win
	)

static

Definition at line 225 of file ccl_cls.c.

References ccl_gsl, ccl_raise_gsl_warning(), IntMagPar::chi, cl_cmbl_bm::cosmo, IntMagPar::cosmo, integrand_mag(), ccl_gsl_params::INTEGRATION_EPSREL, ccl_gsl_params::INTEGRATION_GAUSS_KRONROD_POINTS, ccl_gsl_params::N_ITERATION, IntMagPar::spl_pz, IntMagPar::spl_sz, IntMagPar::status, and w.

Referenced by clt_init_wM().

 {
   int gslstatus =0, status =0;
   double result,eresult;
   IntMagPar ip;
   gsl_function F;
   gsl_integration_workspace *w=gsl_integration_workspace_alloc(ccl_gsl->N_ITERATION);
 
   ip.chi=chi;
   ip.cosmo=cosmo;
   ip.spl_pz=spl_pz;
   ip.spl_sz=spl_sz;
   ip.status = &status;
   F.function=&integrand_mag;
   F.params=&ip;
   // This conputes the magnification lensing kernel:
   //   w_M(chi) = Integral[ dN/dchi(chi') * (1-5/2 * s(chi)) * f(chi'-chi)/f(chi') , chi < chi' < chi_horizon ]
   // Where f(chi) is the comoving angular distance (which is just chi for zero curvature)
   // and s(chi) is the magnification bias parameter.
   gslstatus=gsl_integration_qag(&F, chi, chi_max, 0,
                                 ccl_gsl->INTEGRATION_EPSREL, ccl_gsl->N_ITERATION,
                                 ccl_gsl->INTEGRATION_GAUSS_KRONROD_POINTS,
                                 w, &result, &eresult);
   *win=result;
   gsl_integration_workspace_free(w);
   if(gslstatus!=GSL_SUCCESS || *ip.status) {
     ccl_raise_gsl_warning(gslstatus, "ccl_cls.c: window_magnification():");
     return 1;
   }
   //TODO: chi_max should be changed to chi_horizon
   //we should precompute this quantity and store it in cosmo by default
 
   return 0;
 }

Classes

Macros

Functions

Macro Definition Documentation

Function Documentation