CosmoBolognaLib
Free Software C++/Python libraries for cosmological calculations
ModelFunction_TwoPointCorrelation.cpp
Go to the documentation of this file.
1 /*******************************************************************
2  * Copyright (C) 2016 by Federico Marulli and Alfonso Veropalumbo *
3  * federico.marulli3@unibo.it *
4  * *
5  * This program is free software; you can redistribute it and/or *
6  * modify it under the terms of the GNU General Public License as *
7  * published by the Free Software Foundation; either version 2 of *
8  * the License, or (at your option) any later version. *
9  * *
10  * This program is distributed in the hope that it will be useful,*
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
13  * GNU General Public License for more details. *
14  * *
15  * You should have received a copy of the GNU General Public *
16  * License along with this program; if not, write to the Free *
17  * Software Foundation, Inc., *
18  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
19  *******************************************************************/
20 
37 
38 using namespace std;
39 
40 using namespace cbl;
41 
42 
43 // ============================================================================================
44 
45 
46 std::vector<double> cbl::modelling::twopt::true_k_mu_AP (const double kk, const double mu, const double alpha_perp, const double alpha_par)
47 {
48  const double FF = alpha_par/alpha_perp;
49  const double fact = sqrt(1.+mu*mu*(1./(FF*FF)-1.));
50 
51  return {kk/alpha_perp*fact, mu/FF/fact};
52 }
53 
54 // ============================================================================================
55 
56 
57 double cbl::modelling::twopt::Pkmu_DeWiggled (const double kk, const double mu, const double sigmaNL_perp, const double sigmaNL_par, const double linear_growth_rate, const double bias, const double SigmaS, const std::shared_ptr<cbl::glob::FuncGrid> Pk, const std::shared_ptr<cbl::glob::FuncGrid> Pk_NW)
58 {
59  const double beta = linear_growth_rate/bias;
60 
61  const double KaiserBoost = pow(1.+mu*mu*beta, 2);
62  const double sigmaNL2 = 0.5*((1.-mu*mu)*sigmaNL_perp*sigmaNL_perp+mu*mu*sigmaNL_par*sigmaNL_par);
63 
64  const double Fstreaming = pow(1+kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*SigmaS*SigmaS, -2);
65 
66  const double sigmaNL = sqrt(sigmaNL_perp*sigmaNL_perp+sigmaNL_par*sigmaNL_par);
67  const double _Pk = (sigmaNL<1.e-5) ? Pk->operator()(kk) : (Pk->operator()(kk)-Pk_NW->operator()(kk))*exp(-kk*kk*sigmaNL2)+Pk_NW->operator()(kk);
68 
69  return bias*bias*KaiserBoost*Fstreaming*_Pk;
70 }
71 
72 
73 // ============================================================================================
74 
75 
76 double cbl::modelling::twopt::Pkmu_ModeCoupling (const double kk, const double mu, const double linear_growth_rate, const double bias, const double SigmaV, const double AMC, const std::shared_ptr<cbl::glob::FuncGrid> Pk, const std::shared_ptr<cbl::glob::FuncGrid> Pk_1loop)
77 {
78  const double beta = linear_growth_rate/bias;
79 
80  const double KaiserBoost = pow(1.+mu*mu*beta, 2);
81  const double Fstreaming = pow(1+kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*SigmaV*SigmaV, -2);
82  double Pk_NL = bias*bias*(Pk->operator()(kk)*exp(-kk*kk*SigmaV*SigmaV));
83 
84  if (kk<5)
85  Pk_NL += bias*bias*AMC*Pk_1loop->operator()(kk)*pow(2.*par::pi, -3);
86 
87  return KaiserBoost*Fstreaming*Pk_NL;
88 }
89 
90 
91 // ============================================================================================
92 
93 
94 double cbl::modelling::twopt::Pkmu_dispersion (const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const std::shared_ptr<cbl::glob::FuncGrid> Pklin)
95 {
96  double DispFactor;
97 
98  if (DFoG=="Gaussian")
99  DispFactor = exp(-kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav);
100 
101  else if (DFoG=="Lorentzian")
102  DispFactor = 1./(1.+(kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav));
103 
104  else
105  ErrorCBL("the chosen DFoG ("+DFoG+") is not currently implemented!", "Pkmu_dispersion", "ModelFunction_TwoPointCorrelation.cpp");
106 
107  return DispFactor*bias*bias*pow(1+linear_growth_rate/bias*mu*mu, 2)*Pklin->operator()(kk);
108 }
109 
110 // ============================================================================================
111 
112 
113 double cbl::modelling::twopt::Pkmu_Scoccimarro (const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const std::shared_ptr<cbl::glob::FuncGrid> Pk_DeltaDelta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_DeltaTheta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_ThetaTheta)
114 {
115  const double Pk_dd = Pk_DeltaDelta->operator()(kk);
116  const double Pk_dt = Pk_DeltaTheta->operator()(kk);
117  const double Pk_tt = Pk_ThetaTheta->operator()(kk);
118 
119  double DispFactor;
120 
121  if (DFoG=="Gaussian")
122  DispFactor = exp(-kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav);
123 
124  else if (DFoG=="Lorentzian")
125  DispFactor = 1./(1.+(kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav));
126 
127  else
128  ErrorCBL("the chosen DFoG ("+DFoG+") is not currently implemented!", "Pkmu_Scoccimarro", "ModelFunction_TwoPointCorrelation.cpp");
129 
130  const double beta = linear_growth_rate/bias;
131 
132  return DispFactor*bias*bias*(Pk_dd+2.*beta*pow(mu, 2)*Pk_dt+pow(beta, 2)*pow(mu, 4)*Pk_tt);
133 }
134 
135 
136 // ============================================================================================
137 
138 
139 double cbl::modelling::twopt::Pkmu_Scoccimarro_fitPezzotta (const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const double kd, const double kt, const std::shared_ptr<cbl::glob::FuncGrid> Pklin, const std::shared_ptr<cbl::glob::FuncGrid> Pknonlin)
140 {
141  const double beta = linear_growth_rate/bias;
142  double DispFactor;
143  double Pk_dd = Pknonlin->operator()(kk);
144  double Pk_dt = sqrt(Pknonlin->operator()(kk)*Pklin->operator()(kk)*exp(-kk/kd));
145  double Pk_tt = Pklin->operator()(kk)*exp(-kk/kt);
146 
147  if (DFoG=="Gaussian")
148  DispFactor = exp(-kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav);
149 
150  else if (DFoG=="Lorentzian")
151  DispFactor = 1./(1.+(kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav));
152 
153  else
154  ErrorCBL("the chosen DFoG ("+DFoG+") is not currently implemented!", "Pkmu_Scoccimarro_fitPezzotta", "ModelFunction_TwoPointCorrelation.cpp");
155 
156  return DispFactor*bias*bias*(Pk_dd+2.*beta*pow(mu, 2)*Pk_dt+pow(beta, 2)*pow(mu, 4)*Pk_tt);
157 }
158 
159 
160 // ============================================================================================
161 
162 
163 double cbl::modelling::twopt::Pkmu_Scoccimarro_fitBel (const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const double kd, const double bb, const double a1, const double a2, const double a3, const std::shared_ptr<cbl::glob::FuncGrid> Pklin, const std::shared_ptr<cbl::glob::FuncGrid> Pknonlin)
164 {
165  const double beta = linear_growth_rate/bias;
166  double DispFactor;
167  double Pk_dd = Pknonlin->operator()(kk);
168  double Pk_dt = sqrt(Pknonlin->operator()(kk)*Pklin->operator()(kk))*exp(-kk/kd-bb*pow(kk, 6.0));
169  double Pk_tt = Pklin->operator()(kk)*exp(-kk*(a1 + a2*kk + a3*kk*kk));
170 
171  if (DFoG=="Gaussian")
172  DispFactor = exp(-kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav);
173  else if (DFoG=="Lorentzian")
174  DispFactor = 1./(1.+(kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav));
175  else
176  ErrorCBL("the chosen DFoG ("+DFoG+") is not currently implemented!", "Pkmu_Scoccimarro_fitBel", "ModelFunction_TwoPointCorrelation.cpp");
177 
178  return DispFactor*bias*bias*(Pk_dd+2.*beta*pow(mu, 2)*Pk_dt+pow(beta, 2)*pow(mu, 4)*Pk_tt);
179 }
180 
181 
182 // ============================================================================================
183 
184 
185 double cbl::modelling::twopt::Pkmu_TNS (const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const std::shared_ptr<cbl::glob::FuncGrid> Pk_DeltaDelta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_DeltaTheta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_ThetaTheta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A11, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A12, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A22, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A23, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A33, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B12, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B13, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B14, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B22, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B23, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B24, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B33, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B34, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B44)
186 {
187  const double beta = linear_growth_rate/bias;
188 
189  const double Pk_dd = Pk_DeltaDelta->operator()(kk);
190  const double Pk_dt = Pk_DeltaTheta->operator()(kk);
191  const double Pk_tt = Pk_ThetaTheta->operator()(kk);
192 
193  const double pkA11 = Pk_A11->operator()(kk);
194  const double pkA12 = Pk_A12->operator()(kk);
195  const double pkA22 = Pk_A22->operator()(kk);
196  const double pkA23 = Pk_A23->operator()(kk);
197  const double pkA33 = Pk_A33->operator()(kk);
198  const double pkB12 = Pk_B12->operator()(kk);
199  const double pkB13 = Pk_B13->operator()(kk);
200  const double pkB14 = Pk_B14->operator()(kk);
201  const double pkB22 = Pk_B22->operator()(kk);
202  const double pkB23 = Pk_B23->operator()(kk);
203  const double pkB24 = Pk_B24->operator()(kk);
204  const double pkB33 = Pk_B33->operator()(kk);
205  const double pkB34 = Pk_B34->operator()(kk);
206  const double pkB44 = Pk_B44->operator()(kk);
207 
208  double DispFactor;
209  if (DFoG=="Gaussian")
210  DispFactor = exp(-kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav);
211  else if (DFoG=="Lorentzian")
212  DispFactor = 1./(1.+(kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav));
213  else
214  ErrorCBL("the chosen DFoG ("+DFoG+") is not currently implemented!", "Pkmu_TNS", "ModelFunction_TwoPointCorrelation.cpp");
215 
216  const double Pk_Scoccimarro_large_scales = bias*bias*(Pk_dd+2*beta*pow(mu, 2)*Pk_dt+pow(beta, 2)*pow(mu, 4)*Pk_tt);
217 
218  const double A2 = pow(mu, 2)*(beta*pkA11+pow(beta, 2)*pkA12);
219  const double A4 = pow(mu, 4)*pow(beta, 2)*(pkA22+beta*pkA23);
220  const double A6 = pow(mu, 6)*pow(beta, 3)*pkA33;
221  const double B2 = pow(mu, 2)*(pow(beta, 2)*pkB12+pow(beta, 3)*pkB13+pow(beta, 4)*pkB14);
222  const double B4 = pow(mu, 4)*(pow(beta, 2)*pkB22+pow(beta, 3)*pkB23+pow(beta, 4)*pkB24);
223  const double B6 = pow(mu, 6)*(pow(beta, 3)*pkB33+pow(beta, 4)*pkB34);
224  const double B8 = pow(mu, 8)*pow(beta, 4)*pkB44;
225 
226  const double Pk_A = pow(bias, 3)*(A2+A4+A6);
227  const double Pk_B = pow(bias, 4)*(B2+B4+B6+B8);
228 
229  return DispFactor*(Pk_Scoccimarro_large_scales+Pk_A+Pk_B);
230 }
231 
232 
233 // ============================================================================================
234 
235 
236 double cbl::modelling::twopt::Pkmu_eTNS (const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double bias2, const double sigmav, const double Ncorr, const std::shared_ptr<cbl::glob::FuncGrid> Pk_DeltaDelta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_DeltaTheta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_ThetaTheta, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A11, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A12, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A22, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A23, const std::shared_ptr<cbl::glob::FuncGrid> Pk_A33, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B12, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B13, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B14, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B22, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B23, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B24, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B33, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B34, const std::shared_ptr<cbl::glob::FuncGrid> Pk_B44, const std::shared_ptr<cbl::glob::FuncGrid> Pk_b2d, const std::shared_ptr<cbl::glob::FuncGrid> Pk_b2v, const std::shared_ptr<cbl::glob::FuncGrid> Pk_b22, const std::shared_ptr<cbl::glob::FuncGrid> Pk_bs2d, const std::shared_ptr<cbl::glob::FuncGrid> Pk_bs2v, const std::shared_ptr<cbl::glob::FuncGrid> Pk_b2s2, const std::shared_ptr<cbl::glob::FuncGrid> Pk_bs22, const std::shared_ptr<cbl::glob::FuncGrid> sigma32Pklin)
237 {
238  const double beta = linear_growth_rate/bias;
239 
240  const double pk_dd = Pk_DeltaDelta->operator()(kk);
241  const double pk_dt = Pk_DeltaTheta->operator()(kk);
242  const double pk_tt = Pk_ThetaTheta->operator()(kk);
243 
244  const double pkA11 = Pk_A11->operator()(kk);
245  const double pkA12 = Pk_A12->operator()(kk);
246  const double pkA22 = Pk_A22->operator()(kk);
247  const double pkA23 = Pk_A23->operator()(kk);
248  const double pkA33 = Pk_A33->operator()(kk);
249  const double pkB12 = Pk_B12->operator()(kk);
250  const double pkB13 = Pk_B13->operator()(kk);
251  const double pkB14 = Pk_B14->operator()(kk);
252  const double pkB22 = Pk_B22->operator()(kk);
253  const double pkB23 = Pk_B23->operator()(kk);
254  const double pkB24 = Pk_B24->operator()(kk);
255  const double pkB33 = Pk_B33->operator()(kk);
256  const double pkB34 = Pk_B34->operator()(kk);
257  const double pkB44 = Pk_B44->operator()(kk);
258 
259  const double pk_b2d = Pk_b2d->operator()(kk);
260  const double pk_b2v = Pk_b2v->operator()(kk);
261  const double pk_b22 = Pk_b22->operator()(kk);
262  const double pk_bs2d = Pk_bs2d->operator()(kk);
263  const double pk_bs2v = Pk_bs2v->operator()(kk);
264  const double pk_b2s2 = Pk_b2s2->operator()(kk);
265  const double pk_bs22 = Pk_bs22->operator()(kk);
266  const double sigma32pklin = sigma32Pklin->operator()(kk);
267 
268  double DispFactor;
269  if (DFoG=="Gaussian")
270  DispFactor = exp(-kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav);
271  else if (DFoG=="Lorentzian")
272  DispFactor = 1./(1.+(kk*kk*mu*mu*linear_growth_rate*linear_growth_rate*sigmav*sigmav));
273  else
274  ErrorCBL("the chosen DFoG ("+DFoG+") is not currently implemented!", "Pkmu_TNS", "ModelFunction_TwoPointCorrelation.cpp");
275 
276  const double A2 = pow(mu, 2)*(beta*pkA11+pow(beta, 2)*pkA12);
277  const double A4 = pow(mu, 4)*pow(beta, 2)*(pkA22+beta*pkA23);
278  const double A6 = pow(mu, 6)*pow(beta, 3)*pkA33;
279  const double B2 = pow(mu, 2)*(pow(beta, 2)*pkB12+pow(beta, 3)*pkB13+pow(beta, 4)*pkB14);
280  const double B4 = pow(mu, 4)*(pow(beta, 2)*pkB22+pow(beta, 3)*pkB23+pow(beta, 4)*pkB24);
281  const double B6 = pow(mu, 6)*(pow(beta, 3)*pkB33+pow(beta, 4)*pkB34);
282  const double B8 = pow(mu, 8)*pow(beta, 4)*pkB44;
283 
284  const double Pk_A = pow(bias, 3)*(A2+A4+A6);
285  const double Pk_B = pow(bias, 4)*(B2+B4+B6+B8);
286 
287  //----- NL bias terms -----
288  const double bs2 = (-4./7.)*(bias-1.);
289  const double b3nl = (32./315.)*(bias-1.);
290 
291  const double Pgdd = bias*bias*pk_dd+2*bias2*bias*pk_b2d+2*bs2*bias*pk_bs2d+2*b3nl*bias*sigma32pklin+bias2*bias2*pk_b22+2*bias2*bs2*pk_b2s2+bs2*bs2*pk_bs22+Ncorr;
292  const double Pgdv = bias*pk_dt+bias2*pk_b2v+bs2*pk_bs2v+b3nl*sigma32pklin;
293  const double Pk_Scoccimarro_large_scales = Pgdd+2.*linear_growth_rate*pow(mu, 2)*Pgdv+pow(linear_growth_rate, 2)*pow(mu, 4)*pk_tt;
294 
295  return DispFactor*(Pk_Scoccimarro_large_scales+Pk_A+Pk_B);
296 }
297 
298 
299 // ============================================================================================
300 
301 
302 double cbl::modelling::twopt::Pkmu (const double kk, const double mu, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double alpha_perp, const double alpha_par)
303 {
304  const double kp = true_k_mu_AP(kk, mu, alpha_perp, alpha_par)[0];
305  const double mup = true_k_mu_AP(kk, mu, alpha_perp, alpha_par)[1];
306 
307  double pkmu = 0;
308 
309  if (model=="dispersion_dewiggled") {
310  if (parameter.size()!=5)
311  ErrorCBL("the "+model+" model has 5 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
312  pkmu = Pkmu_DeWiggled(kp, mup, parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], pk_interp[0], pk_interp[1]);
313  }
314 
315  else if (model=="dispersion_modecoupling") {
316  if (parameter.size()!=4)
317  ErrorCBL("the "+model+" model has 4 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
318  pkmu = Pkmu_ModeCoupling(kp, mup, parameter[0], parameter[1], parameter[2], parameter[3], pk_interp[0], pk_interp[1]);
319  }
320 
321  else if (model=="dispersion_Gauss") {
322  if (parameter.size()!=3)
323  ErrorCBL("the "+model+" model has 3 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
324  pkmu = Pkmu_dispersion(kp, mup, "Gaussian", parameter[0], parameter[1], parameter[2], pk_interp[0]);
325  }
326 
327  else if (model=="dispersion_Lorentz") {
328  if (parameter.size()!=3)
329  ErrorCBL("the "+model+" model has 3 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
330  pkmu = Pkmu_dispersion(kp, mup, "Lorentzian", parameter[0], parameter[1], parameter[2], pk_interp[0]);
331  }
332 
333  else if (model=="Scoccimarro_Gauss") {
334  if (parameter.size()!=3)
335  ErrorCBL("the "+model+" model has 3 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
336  pkmu = Pkmu_Scoccimarro(kp, mup, "Gaussian", parameter[0], parameter[1], parameter[2], pk_interp[0], pk_interp[1], pk_interp[2]);
337  }
338 
339  else if (model=="Scoccimarro_Lorentz") {
340  if (parameter.size()!=3)
341  ErrorCBL("the "+model+" model has 3 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
342  pkmu = Pkmu_Scoccimarro(kp, mup, "Lorentzian", parameter[0], parameter[1], parameter[2], pk_interp[0], pk_interp[1], pk_interp[2]);
343  }
344 
345  else if (model=="Scoccimarro_Pezzotta_Gauss") {
346  if (parameter.size()!=5)
347  ErrorCBL("the "+model+" model has 5 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
348  pkmu = Pkmu_Scoccimarro_fitPezzotta(kp, mup, "Gaussian", parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], pk_interp[0], pk_interp[1]);
349  }
350 
351  else if (model=="Scoccimarro_Pezzotta_Lorentz") {
352  if (parameter.size()!=5)
353  ErrorCBL("the "+model+" model has 5 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
354  pkmu = Pkmu_Scoccimarro_fitPezzotta(kp, mup, "Lorentzian", parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], pk_interp[0], pk_interp[1]);
355  }
356 
357  else if (model=="Scoccimarro_Bel_Gauss") {
358  if (parameter.size()!=8)
359  ErrorCBL("the "+model+" model has 8 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
360  pkmu = Pkmu_Scoccimarro_fitBel(kp, mup, "Gaussian", parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], parameter[5], parameter[6], parameter[7], pk_interp[0], pk_interp[1]);
361  }
362 
363  else if (model=="Scoccimarro_Bel_Lorentz") {
364  if (parameter.size()!=8)
365  ErrorCBL("the "+model+" model has 8 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
366  pkmu = Pkmu_Scoccimarro_fitBel(kp, mup, "Lorentzian", parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], parameter[5], parameter[6], parameter[7], pk_interp[0], pk_interp[1]);
367  }
368 
369  else if (model=="TNS_Gauss") {
370  if (parameter.size()!=3)
371  ErrorCBL("the "+model+" model has 3 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
372  pkmu = Pkmu_TNS(kp, mup, "Gaussian", parameter[0], parameter[1], parameter[2], pk_interp[0], pk_interp[1], pk_interp[2], pk_interp[3], pk_interp[4], pk_interp[5], pk_interp[6], pk_interp[7], pk_interp[8], pk_interp[9], pk_interp[10], pk_interp[11], pk_interp[12], pk_interp[13], pk_interp[14], pk_interp[15], pk_interp[16]);
373  }
374 
375  else if (model=="TNS_Lorentz") {
376  if (parameter.size()!=3)
377  ErrorCBL("the "+model+" model has 3 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
378  pkmu = Pkmu_TNS(kp, mup, "Lorentzian", parameter[0], parameter[1], parameter[2], pk_interp[0], pk_interp[1], pk_interp[2], pk_interp[3], pk_interp[4], pk_interp[5], pk_interp[6], pk_interp[7], pk_interp[8], pk_interp[9], pk_interp[10], pk_interp[11], pk_interp[12], pk_interp[13], pk_interp[14], pk_interp[15], pk_interp[16]);
379  }
380 
381  else if (model=="eTNS_Gauss") {
382  if (parameter.size()!=5)
383  ErrorCBL("the "+model+" model has 5 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
384  pkmu = Pkmu_eTNS(kp, mup, "Gaussian", parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], pk_interp[0], pk_interp[1], pk_interp[2], pk_interp[3], pk_interp[4], pk_interp[5], pk_interp[6], pk_interp[7], pk_interp[8], pk_interp[9], pk_interp[10], pk_interp[11], pk_interp[12], pk_interp[13], pk_interp[14], pk_interp[15], pk_interp[16], pk_interp[17], pk_interp[18], pk_interp[19], pk_interp[20], pk_interp[21], pk_interp[22], pk_interp[23], pk_interp[24]);
385  }
386 
387  else if (model=="eTNS_Lorentz") {
388  if (parameter.size()!=5)
389  ErrorCBL("the "+model+" model has 5 parameters, while in the parameter vector in input has "+conv(parameter.size(), par::fINT)+" parameters!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
390  pkmu = Pkmu_eTNS(kp, mup, "Lorentzian", parameter[0], parameter[1], parameter[2], parameter[3], parameter[4], pk_interp[0], pk_interp[1], pk_interp[2], pk_interp[3], pk_interp[4], pk_interp[5], pk_interp[6], pk_interp[7], pk_interp[8], pk_interp[9], pk_interp[10], pk_interp[11], pk_interp[12], pk_interp[13], pk_interp[14], pk_interp[15], pk_interp[16], pk_interp[17], pk_interp[18], pk_interp[19], pk_interp[20], pk_interp[21], pk_interp[22], pk_interp[23], pk_interp[24]);
391  }
392 
393  else
394  ErrorCBL("the chosen model ("+model+") is not currently implemented!", "Pkmu", "ModelFunction_TwoPointCorrelation.cpp");
395 
396  return pkmu/(alpha_perp*alpha_perp*alpha_par);
397 }
398 
399 
400 // ============================================================================================
401 
402 
403 double cbl::modelling::twopt::Pk_l (const double kk, const int l, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par)
404 {
405  auto integrand = [&] (const double mu)
406  {
407  return modelling::twopt::Pkmu(kk, mu, model, parameter, pk_interp, alpha_perp, alpha_par)*legendre_polynomial(mu, l);
408  };
409 
410  return 0.5*(2*l+1)*wrapper::gsl::GSL_integrate_cquad(integrand, -1., 1., prec);
411 }
412 
413 // ============================================================================================
414 
415 
416 std::vector<double> cbl::modelling::twopt::Pk_l (const std::vector<double> kk, const int l, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par)
417 {
418  vector<double> Pkl(kk.size(), 0);
419 
420  for (size_t i=0; i<kk.size(); i++) {
421 
422  auto integrand = [&] (const double mu)
423  {
424  return modelling::twopt::Pkmu(kk[i], mu, model, parameter, pk_interp, alpha_perp, alpha_par)*legendre_polynomial(mu, l);
425  };
426 
427  Pkl[i] = 0.5*(2*l+1)*wrapper::gsl::GSL_integrate_cquad(integrand, -1., 1., prec);
428  }
429 
430  return Pkl;
431 }
432 
433 
434 // ============================================================================================
435 
436 
437 cbl::glob::FuncGrid cbl::modelling::twopt::Xil_interp (const std::vector<double> kk, const int l, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par)
438 {
439  vector<double> Pkl = Pk_l(kk, l, model, parameter, pk_interp, prec, alpha_perp, alpha_par);
440  vector<double> rr, Xil;
441 
442  cbl::wrapper::fftlog::transform_FFTlog(rr, Xil, 1, kk, Pkl, l, 0, 2.*par::pi, 1);
443 
444  cbl::glob::FuncGrid interp(rr, Xil, "Spline");
445 
446  return interp;
447 }
448 
449 
450 // ============================================================================================
451 
452 
453 std::vector<std::vector<double>> cbl::modelling::twopt::Xi_l (const std::vector<double> rr, const int nmultipoles, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par)
454 {
455  vector<vector<double>> Xil(3);
456 
457  for (int i=0; i<nmultipoles; i++) {
458  double sign = (i%2==0) ? 1 : -1;
459  Xil[i] = Xil_interp(pk_interp[0]->x(), 2*i, model, parameter, pk_interp, prec, alpha_perp, alpha_par).eval_func(rr);
460  for (size_t j=0; j<rr.size(); j++)
461  Xil[i][j] *= sign;
462  }
463 
464  return Xil;
465 }
466 
467 
468 // ============================================================================================
469 
470 
471 std::vector<double> cbl::modelling::twopt::Xi_l (const std::vector<double> rr, const std::vector<int> dataset_order, const std::vector<bool> use_pole, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par)
472 {
473  vector<cbl::glob::FuncGrid> interp_Xil(3);
474  vector<double> sign={1., -1., 1.};
475 
476  for (size_t i=0; i<3; i++)
477  if (use_pole[i])
478  interp_Xil[i] = Xil_interp(pk_interp[0]->x(), 2*i, model, parameter, pk_interp, prec, alpha_perp, alpha_par);
479 
480  vector<double> Xil(rr.size());
481 
482  for (size_t i=0; i<rr.size(); i++)
483  Xil[i] = sign[dataset_order[i]]*interp_Xil[dataset_order[i]](rr[i]);
484 
485  return Xil;
486 }
487 
488 
489 // ============================================================================================
490 
491 
492 double cbl::modelling::twopt::Xi_polar(const double rad_fid, const double mu_fid, const double alpha_perpendicular, const double alpha_parallel, const std::vector<std::shared_ptr<cbl::glob::FuncGrid>> xi_multipoles)
493 {
494  const double apar2 = alpha_parallel*alpha_parallel;
495  const double aperp2 = alpha_perpendicular*alpha_perpendicular;
496 
497  double mu_fid_sq = mu_fid*mu_fid;
498  const double factor = sqrt(apar2*mu_fid_sq+(1.-mu_fid_sq)*aperp2);
499  const double mu_true = mu_fid * alpha_parallel/factor;
500  const double mu_true_sq = mu_true*mu_true;
501  const double st = rad_fid*factor;
502 
503  return xi_multipoles[0]->operator()(st)+
504  xi_multipoles[1]->operator()(st)*0.5*(3*mu_true_sq-1)+
505  xi_multipoles[2]->operator()(st)*0.125*(35*mu_true_sq*mu_true_sq-30*mu_true_sq+3);
506 }
507 
508 
509 // ============================================================================================
510 
511 
512 std::vector<std::vector<double>> cbl::modelling::twopt::Xi_rppi (const std::vector<double> rp, const std::vector<double> pi, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par)
513 {
514  const int nr=200;
515  const size_t nmultipoles = pk_interp.size();
516  vector<double> sign = {1., -1., 1., 1, -1};
517  vector<double> rr = linear_bin_vector(nr, min(Min(rp), Min(pi))*0.999, sqrt(Max(rp)*Max(rp)+Max(pi)*Max(pi))*1.001);
518 
519  vector<glob::FuncGrid> interp_Xil(nmultipoles);
520 
521  for (size_t i=0; i<nmultipoles; i++)
522  interp_Xil[i] = Xil_interp(pk_interp[0]->x(), 2*i, model, parameter, pk_interp, prec, alpha_perp, alpha_par);
523 
524  vector<vector<double>> xi_rppi(rp.size(), vector<double>(pi.size(), 0));
525 
526  for (size_t i =0; i<rp.size(); i++)
527  for (size_t j =0; j<pi.size(); j++) {
528  double s = sqrt(rp[i]*rp[i]+pi[j]*pi[j]);
529  double mu = pi[j]/s;
530  for (size_t l=0; l<nmultipoles; l++)
531  xi_rppi[i][j] += sign[l]*interp_Xil[l](s)*legendre_polynomial (mu, l*2);
532  }
533 
534  return xi_rppi;
535 }
536 
537 
538 // ============================================================================================
539 
540 
541 std::vector<double> cbl::modelling::twopt::wp_from_Xi_rppi (const std::vector<double> rp, const double pimax, const std::string model, const std::vector<double> parameter, const std::vector<std::shared_ptr<glob::FuncGrid>> pk_interp, const double prec, const double alpha_perp, const double alpha_par) {
542 
543  vector<double> pi = linear_bin_vector(100, 1.e-4, pimax*1.001);
544 
545  vector<vector<double>> xi_rppi = modelling::twopt::Xi_rppi(rp, pi, model, parameter, pk_interp, prec, alpha_perp, alpha_par);
546  vector<double> wp(rp.size());
547 
548  for (size_t i=0; i<rp.size(); i++) {
549  glob::FuncGrid func(pi, xi_rppi[i], "Spline");
550  wp[i] = func.integrate_qag(0., pimax);
551  }
552 
553  return wp;
554 }
555 
556 
557 // ============================================================================================
558 
559 
560 std::vector<std::vector<double>> cbl::modelling::twopt::damped_Pk_terms (const std::vector<double> kk, const double linear_growth_rate, const double SigmaS, const std::shared_ptr<cbl::glob::FuncGrid> PkDM)
561 {
562  (void)PkDM;
563  double sqrt_pi = sqrt(par::pi);
564  vector<vector<double>> pk(3, vector<double>(kk.size(), 0));
565 
566  for (size_t i=0; i< kk.size(); i++) {
567  double kSigma = kk[i]*SigmaS;
568  double pk_dm = PkDM->operator()(kk[i]);
569  double erf_kSigma = erf(kSigma);
570 
571  pk[0][i] = pk_dm*sqrt_pi/(2.*kSigma)*erf_kSigma;
572  pk[1][i] = linear_growth_rate*pow(kSigma, -3.)*pk_dm*(0.5*sqrt_pi*erf_kSigma-kSigma*exp(-kSigma*kSigma));
573  pk[2][i] = linear_growth_rate*linear_growth_rate*pow(kSigma, -5.)*pk_dm*(3.*sqrt_pi/8.*erf_kSigma-0.25*kSigma*(2*kSigma*kSigma+3)*exp(-kSigma*kSigma));
574  }
575 
576 
577  return pk;
578 }
579 
580 // ============================================================================================
581 
582 
583 std::vector<double> cbl::modelling::twopt::damped_Xi (const std::vector<double> ss, const double bias, const double linear_growth_rate, const double SigmaS, const std::vector<double> kk, const std::shared_ptr<cbl::glob::FuncGrid> PkDM)
584 {
585  vector<vector<double>> pk_terms = modelling::twopt::damped_Pk_terms(kk, linear_growth_rate, SigmaS, PkDM);
586 
587  vector<double> xi(ss.size(), 0);
588 
589  for (size_t i=0; i<pk_terms.size(); i++) {
590  vector<double> xi_term = wrapper::fftlog::transform_FFTlog(ss, 1, kk, pk_terms[i], 0);
591  for (size_t j=0; j<ss.size(); j++)
592  xi[j] += pow(bias, 2-i)*xi_term[j];
593 
594  }
595 
596  return xi;
597 }
Global functions to model two-point correlation functions of any type.
The class FuncGrid.
Definition: FuncGrid.h:55
std::vector< double > eval_func(const std::vector< double > xx) const
evaluate the function at the xx points
Definition: FuncGrid.cpp:153
double integrate_qag(const double a, const double b, const double rel_err=1.e-2, const double abs_err=1.e-6, const int limit_size=1000, const int rule=6)
compute the definite integral with GSL qag method
Definition: FuncGrid.cpp:187
static const char fINT[]
conversion symbol for: int -> std::string
Definition: Constants.h:121
static const double pi
: the ratio of a circle's circumference to its diameter
Definition: Constants.h:199
std::vector< double > damped_Xi(const std::vector< double > ss, const double bias, const double linear_growth_rate, const double SigmaS, const std::vector< double > kk, const std::shared_ptr< cbl::glob::FuncGrid > PkDM)
the damped two-point correlation monopole; from Sereno et al. 2015
std::vector< double > wp_from_Xi_rppi(const std::vector< double > rp, const double pimax, const std::string model, const std::vector< double > parameter, const std::vector< std::shared_ptr< glob::FuncGrid >> pk_interp, const double prec=1.e-5, const double alpha_perp=1., const double alpha_par=1.)
the projected two-point correlation function
double Pkmu_Scoccimarro(const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const std::shared_ptr< cbl::glob::FuncGrid > Pk_DeltaDelta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_DeltaTheta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_ThetaTheta)
the redshift-space galaxy power spectrum, as a function of and , predicted by the Scoccimarro model
double Pkmu_DeWiggled(const double kk, const double mu, const double sigmaNL_perp, const double sigmaNL_par, const double linear_growth_rate, const double bias, const double SigmaS, const std::shared_ptr< cbl::glob::FuncGrid > Pk, const std::shared_ptr< cbl::glob::FuncGrid > Pk_NW)
the redshift-space galaxy power spectrum, as a function of and , predicted by the de-wiggled model
double Pkmu_dispersion(const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const std::shared_ptr< cbl::glob::FuncGrid > Pklin)
the redshift-space galaxy power spectrum, as a function of and , predicted by the dispersion model
std::vector< std::vector< double > > Xi_l(const std::vector< double > rr, const int nmultipoles, const std::string model, const std::vector< double > parameter, const std::vector< std::shared_ptr< glob::FuncGrid >> pk_interp, const double prec=1.e-5, const double alpha_perp=1., const double alpha_par=1.)
the multipole of order l of the two-point correlation function
double Pkmu_eTNS(const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double bias2, const double sigmav, const double Ncorr, const std::shared_ptr< cbl::glob::FuncGrid > Pk_DeltaDelta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_DeltaTheta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_ThetaTheta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A11, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A12, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A22, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A23, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A33, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B12, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B13, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B14, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B22, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B23, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B24, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B33, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B34, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B44, const std::shared_ptr< cbl::glob::FuncGrid > Pk_b2d, const std::shared_ptr< cbl::glob::FuncGrid > Pk_b2v, const std::shared_ptr< cbl::glob::FuncGrid > Pk_b22, const std::shared_ptr< cbl::glob::FuncGrid > Pk_bs2d, const std::shared_ptr< cbl::glob::FuncGrid > Pk_bs2v, const std::shared_ptr< cbl::glob::FuncGrid > Pk_b2s2, const std::shared_ptr< cbl::glob::FuncGrid > Pk_bs22, const std::shared_ptr< cbl::glob::FuncGrid > sigma32Pklin)
the redshift-space galaxy power spectrum, as a function of and , predicted by the extended TNS (Taru...
double Pkmu_Scoccimarro_fitBel(const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const double kd, const double bb, const double a1, const double a2, const double a3, const std::shared_ptr< cbl::glob::FuncGrid > Pklin, const std::shared_ptr< cbl::glob::FuncGrid > Pknonlin)
the redshift-space galaxy power spectrum, as a function of and , predicted by the Scoccimarro model
std::vector< std::vector< double > > Xi_rppi(const std::vector< double > rp, const std::vector< double > pi, const std::string model, const std::vector< double > parameter, const std::vector< std::shared_ptr< glob::FuncGrid >> pk_interp, const double prec=1.e-5, const double alpha_perp=1., const double alpha_par=1.)
the cartesian two-point correlation function
std::vector< std::vector< double > > damped_Pk_terms(const std::vector< double > kk, const double linear_growth_rate, const double SigmaS, const std::shared_ptr< cbl::glob::FuncGrid > PkDM)
the power spectrum terms obtained integrating the redshift space 2D power spectrum
cbl::glob::FuncGrid Xil_interp(const std::vector< double > kk, const int l, const std::string model, const std::vector< double > parameter, const std::vector< std::shared_ptr< glob::FuncGrid >> pk_interp, const double prec=1.e-5, const double alpha_perp=1., const double alpha_par=1.)
the interpolating function of multipole expansion of the two-point correlation function at a given or...
double Pkmu_ModeCoupling(const double kk, const double mu, const double linear_growth_rate, const double bias, const double sigmav, const double AMC, const std::shared_ptr< cbl::glob::FuncGrid > PkLin, const std::shared_ptr< cbl::glob::FuncGrid > PkMC)
the redshift-space galaxy power spectrum, as a function of and , predicted by the mode-coupling mode...
double Pkmu(const double kk, const double mu, const std::string model, const std::vector< double > parameter, const std::vector< std::shared_ptr< glob::FuncGrid >> pk_interp, const double alpha_perp=1., const double alpha_par=1.)
the power spectrum as a function of k and
double Xi_polar(const double rad_fid, const double mu_fid, const double alpha_perpendicular, const double alpha_parallel, const std::vector< std::shared_ptr< cbl::glob::FuncGrid >> xi_multipoles)
the polar two-point correlation function
double bias(const double Mmin, const double sigmalgM, const double M0, const double M1, const double alpha, const std::shared_ptr< void > inputs)
the mean galaxy bias
double Pkmu_Scoccimarro_fitPezzotta(const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const double kd, const double kt, const std::shared_ptr< cbl::glob::FuncGrid > Pklin, const std::shared_ptr< cbl::glob::FuncGrid > Pknonlin)
the redshift-space galaxy power spectrum, as a function of and , predicted by the Scoccimarro model
std::vector< double > true_k_mu_AP(const double kk, const double mu, const double alpha_perp, const double alpha_par)
true k and power spectrum coordinates as a function of observed ones
double Pk_l(const double kk, const int l, const std::string model, const std::vector< double > parameter, const std::vector< std::shared_ptr< glob::FuncGrid >> pk_interp, const double prec=1.e-5, const double alpha_perp=1., const double alpha_par=1.)
the multipole of order l of the power spectrum
double Pkmu_TNS(const double kk, const double mu, const std::string DFoG, const double linear_growth_rate, const double bias, const double sigmav, const std::shared_ptr< cbl::glob::FuncGrid > Pk_DeltaDelta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_DeltaTheta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_ThetaTheta, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A11, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A12, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A22, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A23, const std::shared_ptr< cbl::glob::FuncGrid > Pk_A33, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B12, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B13, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B14, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B22, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B23, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B24, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B33, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B34, const std::shared_ptr< cbl::glob::FuncGrid > Pk_B44)
the redshift-space galaxy power spectrum, as a function of and , predicted by the TNS (Taruya,...
std::vector< double > transform_FFTlog(const std::vector< double > yy, const int dir, const std::vector< double > xx, const std::vector< double > fx, const double mu=0, const double q=0, const double kr=1, const int kropt=0)
wrapper of the FFTlog to compute the discrete Fourier sine or cosine transform of a logarithmically...
Definition: FFTlog.cpp:43
double GSL_integrate_cquad(gsl_function func, const double a, const double b, const double rel_err=1.e-3, const double abs_err=0, const int nevals=100)
integral, using the gsl cquad method
Definition: GSLwrapper.cpp:159
The global namespace of the CosmoBolognaLib
Definition: CAMB.h:38
T Min(const std::vector< T > vect)
minimum element of a std::vector
Definition: Kernel.h:1324
std::string conv(const T val, const char *fact)
convert a number to a std::string
Definition: Kernel.h:903
std::vector< T > linear_bin_vector(const size_t nn, const T min, const T max)
fill a std::vector with linearly spaced values
Definition: Kernel.h:1604
int ErrorCBL(const std::string msg, const std::string functionCBL, const std::string fileCBL, const cbl::glob::ExitCode exitCode=cbl::glob::ExitCode::_error_)
throw an exception: it is used for handling exceptions inside the CosmoBolognaLib
Definition: Kernel.h:780
T Max(const std::vector< T > vect)
maximum element of a std::vector
Definition: Kernel.h:1336
double legendre_polynomial(const double mu, const int l)
the order l Legendre polynomial
Definition: Func.cpp:1786
double wp(const double rp, const std::vector< double > rr, const std::vector< double > xi, const double r_max=100.)
the projected two-point correlation function
Definition: FuncXi.cpp:192