Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / EVGEN / AliGenEMlib.cxx
CommitLineData
e40b9538 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id: AliGenEMlib.cxx 30052 2008-11-25 14:54:18Z morsch $ */
17
18/////////////////////////////////////////////////////////////////////////////
19// //
20// Implementation of AliGenEMlib for electron, di-electron, and photon //
21// cocktail calculations. //
22// It is based on AliGenGSIlib. //
23// //
24// Responsible: R.Averbeck@gsi.de //
25// //
26/////////////////////////////////////////////////////////////////////////////
27
28
4ae1c9f0 29#include <Riostream.h>
e40b9538 30#include "TMath.h"
31#include "TRandom.h"
32#include "TString.h"
33#include "AliGenEMlib.h"
34
35
36ClassImp(AliGenEMlib)
37
4ae1c9f0 38//Initializers for static members
39Int_t AliGenEMlib::fgSelectedPtParam=AliGenEMlib::kPizero7TeVpp;
40Int_t AliGenEMlib::fgSelectedCentrality=AliGenEMlib::kpp;
41Int_t AliGenEMlib::fgSelectedV2Systematic=AliGenEMlib::kNoV2Sys;
6078e216 42
43Double_t AliGenEMlib::CrossOverLc(const double a, const double b, const double x){
4ae1c9f0 44 if(x<b-a/2) return 1.0;
6078e216 45 else if(x>b+a/2) return 0.0;
46 else return cos(((x-b)/a+0.5)*TMath::Pi())/2+0.5;
47}
48Double_t AliGenEMlib::CrossOverRc(const double a, const double b, const double x){
49 return 1-CrossOverLc(a,b,x);
50}
51
4ae1c9f0 52const Double_t AliGenEMlib::fgkV2param[16][15] = {
53 // charged pion cent, based on
54 { 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // pp no V2
55 ,{ 6.792545e-02, 1.324227e+00, 2.071866e-02, 3.385157e+00, 5.183718e+00, 2.783082e-02, 5.000000e+00, 4.078812e+00, 5.444857e-02, 1.555569e-03, 1.099253e+00, 0, 1, 8.823293e-04, 1.310387e+00 } // 0-5, Francesco
56 ,{ 1.219655e-01, 1.237698e+00, 2.595132e-02, 3.157085e+00, 3.697252e+00, 4.576262e-02, 1.103851e+00, 5.102592e+00, 1.072299e-01, 2.338602e-03, 8.475128e-01, 0, 1, 1.631554e-03, 9.432292e-01 } // 5-10, Francesco
57 ,{ 1.758981e-01, 1.253238e+00, 2.823197e-02, 3.455476e+00, 4.283204e+00, 6.053191e-02, 6.737159e-01, 4.420205e+00, 1.733567e-01, 3.292614e-03, 6.389341e-01, 0, 1, 3.234254e-03, 5.304084e-01 } // 10-20,Francesco
58 ,{ 2.231108e-01, 1.316781e+00, 3.537289e-02, 3.228243e+00, 3.697226e+00, 7.131631e-02, 5.902192e-01, 4.938186e+00, 2.123481e-01, 4.775233e-03, 5.707533e-01, 0, 1, 3.906720e-03, 5.261247e-01 } // 20-30,Francesco
59 ,{ 2.045341e-01, 1.324054e+00, 1.125395e-01, 2.179936e+00, 5.614603e+00, 8.946014e-02, 5.601938e-01, 2.959228e+00, 3.669442e-01, 5.318697e-03, 6.032741e-01, 0, 1, 4.231316e-03, 5.563485e-01 } // 30-40,Francesco
60 ,{ 2.074778e-01, 1.894322e+00, 6.019224e-02, 2.414500e+00, 3.822721e+00, 9.891310e-02, 4.976126e-01, 1.256500e+00, 5.210237e-01, 5.687290e-03, 6.785355e-01, 0, 1, 3.954051e-03, 6.913428e-01 } // 40-50,Francesco
61 ,{ 1.961099e-01, 2.122719e+00, 7.633636e-02, 1.731213e+00, 2.427536e+00, 1.026013e-01, 7.932261e-01, 4.082226e+00, 2.495494e-01, 5.677386e-03, 8.405024e-01, 0, 1, 2.961753e-03, 1.008892e+00 } // 50-60,Francesco
62 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 0-10
63 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 20-40
64 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 40-60
65 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 60-80
66 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 0-20
67 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 0-40
68 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 20-80
69 ,{ 0.0000000000, 0.0000000000, 0.0000000000,-1.0000000000, 1.0000000000, 0.0000000000, 0.0000000000, 0.0000000000, 2.0000000000, 0.0000000000, 1.0000000000, 0, 1, 0.0000000000, 1.0000000000 } // 40-80
70};
71
72// MASS 0=>PIZERO, 1=>ETA, 2=>RHO, 3=>OMEGA, 4=>ETAPRIME, 5=>PHI, 6=>JPSI
73const Double_t AliGenEMlib::fgkHM[7] = {0.13498, 0.54751, 0.7755, 0.78265, 0.95778, 1.01946, 3.0969};
74
75const Double_t AliGenEMlib::fgkMtFactor[2][7] = {
76 // {1.0, 0.5, 1.0, 0.9, 0.4, 0.23, 0.054}, // factor for pp from arXiv:1110.3929
77 // {1.0, 0.55, 1.0, 0.9, 0.4, 0.25, 0.004} // factor for PbPb from arXiv:1110.3929
78 //{1., 0.48, 1.0, 0.9, 0.25, 0.4}, (old values)
79 //{1., 0.48, 1.0, 0.9, 0.4, 0.25}, (nlo values)
80 //{1., 0.48, 1.0, 0.8, 0.4, 0.2, 0.06} (combination of nlo and LHC measurements)
81 //https://aliceinfo.cern.ch/Figure/node/2634
82 //https://aliceinfo.cern.ch/Figure/node/2788
83 //https://aliceinfo.cern.ch/Figure/node/4403
84 //J/Psi PbPb from Comparison with Julian Books J/Psi -> e+e-, might be contradicting with https://aliceinfo.cern.ch/Figure/node/3457
85 //best guess:
86 {1., 0.4, 1.0, 0.8, 0.4, 0.2, 0.0054}, //pp
87 {1., 0.4, 1.0, 0.8, 0.4, 0.25, 0.01} //PbPb
88};
6078e216 89
e40b9538 90//==========================================================================
91//
92// Definition of Particle Distributions
93//
94//==========================================================================
95
96//--------------------------------------------------------------------------
97//
4ae1c9f0 98// General functions
99//
100//--------------------------------------------------------------------------
101Double_t AliGenEMlib::PtModifiedHagedornThermal(const Double_t pt,
102 const Double_t c,
103 const Double_t p0,
104 const Double_t p1,
105 const Double_t n,
106 const Double_t cT,
107 const Double_t T)
108{
109 // Modified Hagedorn Thermal fit to Picharged for PbPb:
110 Double_t invYield;
111 invYield = c/TMath::Power(p0+pt/p1,n) + cT*exp(-1.0*pt/T);
112
113 return invYield*(2*TMath::Pi()*pt);
114}
115
116
117
118Double_t AliGenEMlib::PtModifiedHagedornExp(const Double_t pt,
119 const Double_t c,
120 const Double_t p1,
121 const Double_t p2,
122 const Double_t p0,
123 const Double_t n)
124{
125 // Modified Hagedorn exponentiel fit to Pizero for PbPb:
126 Double_t invYield;
127 invYield = c*TMath::Power(exp(-1*(p1*pt-p2*pt*pt))+pt/p0,-n);
128
129 return invYield*(2*TMath::Pi()*pt);
130}
131
132
133Double_t AliGenEMlib::PtModifiedHagedornExp2(const Double_t pt,
134 const Double_t c,
135 const Double_t a,
136 const Double_t b,
137 const Double_t p0,
138 const Double_t p1,
139 const Double_t d,
140 const Double_t n)
141{
142 // Modified Hagedorn exponential fit to charged pions for pPb:
143 Double_t invYield;
144 invYield = c*TMath::Power(exp(-a*pt-b*pt*pt)+pt/p0+TMath::Power(pt/p1,d),-n);
145
146 return invYield*(2*TMath::Pi()*pt);
147}
148
149Double_t AliGenEMlib::PtTsallis(const Double_t pt,
150 const Double_t m,
151 const Double_t c,
152 const Double_t T,
153 const Double_t n)
154{
155 // Tsallis fit to Pizero for pp:
156 Double_t mt;
157 Double_t invYield;
158
159 mt = sqrt(m*m + pt*pt);
160 invYield = c*((n-1.)*(n-2.))/(n*T*(n*T+m*(n-2.)))*pow(1.+(mt-m)/(n*T),-n);
161
162 return invYield*(2*TMath::Pi()*pt);
163}
164
165
166
167//--------------------------------------------------------------------------
168//
e40b9538 169// Pizero
170//
171//--------------------------------------------------------------------------
172Int_t AliGenEMlib::IpPizero(TRandom *)
173{
4ae1c9f0 174 // Return pizero pdg code
e40b9538 175 return 111;
176}
177
178Double_t AliGenEMlib::PtPizero( const Double_t *px, const Double_t */*dummy*/ )
179{
4ae1c9f0 180 // fit functions and corresponding parameter of Pizero pT for pp @ 2.76 TeV and @ 7 TeV and for PbPb @ 2.76 TeV
181
182 Double_t km=0.;
183 Double_t kc=0.;
184 Double_t kn=0.;
185 Double_t kcT=0.;
186 Double_t kT=0.;
187 Double_t kp0=0.;
188 Double_t kp1=0.;
189 Double_t kp2=0.;
190 Double_t ka=0.;
191 Double_t kb=0.;
192 Double_t kd=0.;
193
194 switch(fgSelectedPtParam|fgSelectedCentrality) {
195 // fit to pi charged v1
196 // charged pion from ToF, unidentified hadrons scaled with pion from TPC
197 // for Pb-Pb @ 2.76 TeV
198 case kPichargedPbPb|k0005:
199 kc=1347.5; kp0=0.9393; kp1=2.254; kn=11.294; kcT=0.002537; kT=2.414;
200 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
201 break;
202 case kPichargedPbPb|k0510:
203 kc=1256.1; kp0=0.9545; kp1=2.248; kn=11.291; kcT=0.002662; kT=2.326;
204 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
205 break;
206 case kPichargedPbPb|k2030:
207 kc=7421.6; kp0=1.2059; kp1=1.520; kn=10.220; kcT=0.002150; kT=2.196;
208 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
209 break;
210 case kPichargedPbPb|k3040:
211 kc=1183.2; kp0=1.0478; kp1=1.623; kn=9.8073; kcT=0.00198333; kT=2.073;
212 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
213 break;
214 // the following is what went into the Pb-Pb preliminary approval (0-10%)
215 case kPichargedPbPb|k0010:
216 kc=1296.0; kp0=0.968; kp1=2.567; kn=12.27; kcT=0.004219; kT=2.207;
217 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
218 break;
219 case kPichargedPbPb|k1020:
220 kc=986.0; kp0=0.9752; kp1=2.376; kn=11.62; kcT=0.003116; kT=2.213;
221 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
222 break;
223 case kPichargedPbPb|k2040:
224 kc=17337.0; kp0=1.337; kp1=1.507; kn=10.629; kcT=0.00184; kT=2.234;
225 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
226 break;
227 case kPichargedPbPb|k4050:
228 kc=6220.0; kp0=1.322; kp1=1.224; kn=9.378; kcT=0.000595; kT=2.383;
229 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
230 break;
231 case kPichargedPbPb|k5060:
232 kc=2319.0; kp0=1.267; kp1=1.188; kn=9.044; kcT=0.000437; kT=2.276;
233 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
234 break;
235 case kPichargedPbPb|k4060:
236 kc=4724.0; kp0=1.319; kp1=1.195; kn=9.255; kcT=0.000511; kT=2.344;
237 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
238 break;
239 case kPichargedPbPb|k6080:
240 kc=2842.0; kp0=1.465; kp1=0.8324; kn=8.167; kcT=0.0001049; kT=2.29;
241 return PtModifiedHagedornThermal(*px,kc,kp0,kp1,kn,kcT,kT);
242 break;
e40b9538 243
6078e216 244
4ae1c9f0 245 // fit to pizero from conversion analysis
246 // for PbPb @ 2.76 TeV
247 // Pi0 spectra --> not final results
248 case kPizeroPbPb|k0005:
249 kc=1952.832; kp1=0.264; kp2=0.069; kp0=1.206; kn=9.732;
250 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
251 break;
252 case kPizeroPbPb|k0010:
253 kc=1810.029; kp1=0.291; kp2=0.059; kp0=1.170; kn=9.447;
254 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
255 break;
256 case kPizeroPbPb|k0020:
257 kc=856.241; kp1=-0.409; kp2=-0.127; kp0=1.219; kn=9.030;
258 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
259 break;
260 case kPizeroPbPb|k1020:
261 kc=509.781; kp1=-0.784; kp2=-0.120; kp0=0.931; kn=7.299;
262 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
263 break;
264 case kPizeroPbPb|k2040:
265 kc=541.049; kp1=0.542; kp2=-0.069; kp0=0.972; kn=7.866;
266 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
267 break;
268 case kPizeroPbPb|k2080:
269 kc=222.577; kp1=0.634; kp2=0.009; kp0=0.915; kn=7.431;
270 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
271 break;
272 case kPizeroPbPb|k4080:
273 kc=120.153; kp1=0.7; kp2=-0.14; kp0=0.835; kn=6.980;
274 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
275 break;
276 case kPizeroPbPb|k0040:
277 kc=560.532; kp1=0.548; kp2=-0.048; kp0=1.055; kn=8.132;
278 return PtModifiedHagedornExp(*px,kc,kp1,kp2,kp0,kn);
279 break;
280
281
282 // fit to charged pions for p-Pb @ 5.02TeV
283 case kPichargedPPb:
284 kc=235.5; ka=0.6903; kb=0.06864; kp0=2.289; kp1=0.5872; kd=0.6474; kn=7.842;
285 return PtModifiedHagedornExp2(*px,kc,ka,kb,kp0,kp1,kd,kn);
286 break;
287
288
289 // Tsallis fit to final pizero (PHOS+PCM) -> used for publication
290 // for pp @ 7 TeV
291 case kPizero7TeVpp:
292 case kPizeroEta7TeVpp:
293 km=0.13498; kc=28.01; kT=0.139; kn=6.875;
294 return PtTsallis(*px,km,kc,kT,kn);
295 break;
296 case kPizero7TeVpplow:
297 case kPizeroEta7TeVpplow:
298 km=0.13498; kc=23.84; kT=0.147; kn=7.025;
299 return PtTsallis(*px,km,kc,kT,kn);
300 break;
301 case kPizero7TeVpphigh:
302 case kPizeroEta7TeVpphigh:
303 km=0.13498; kc=32.47; kT=0.132; kn=6.749;
304 return PtTsallis(*px,km,kc,kT,kn);
305 break;
306 // Tsallis fit to pizero: preliminary result from PCM and PHOS (QM'11)
307 // for pp @ 2.76 TeV
308 case kPizero2760GeVpp:
309 case kPizeroEta2760GeVpp:
310 km = 0.13498; kc = 19.75; kT = 0.130; kn = 7.051;
311 return PtTsallis(*px,km,kc,kT,kn);
312 break;
313 case kPizero2760GeVpplow:
314 case kPizeroEta2760GeVpplow:
315 km = 0.13498; kc = 16.12; kT = 0.142; kn = 7.327;
316 return PtTsallis(*px,km,kc,kT,kn);
317 break;
318 case kPizero2760GeVpphigh:
319 case kPizeroEta2760GeVpphigh:
320 km = 0.13498; kc = 25.18; kT = 0.118; kn = 6.782;
321 return PtTsallis(*px,km,kc,kT,kn);
322 break;
323
324 default:
325 return NULL;
326 }
327
e40b9538 328}
329
330Double_t AliGenEMlib::YPizero( const Double_t *py, const Double_t */*dummy*/ )
331{
332 return YFlat(*py);
333
334}
335
6078e216 336Double_t AliGenEMlib::V2Pizero( const Double_t *px, const Double_t */*dummy*/ )
337{
4ae1c9f0 338 int selectedCentrality;
339 double n1,n2,v1,v2;
340 switch(fgSelectedCentrality) {
341 case k0010:
342 selectedCentrality=fgSelectedCentrality;
343 fgSelectedCentrality=k0005;
344 n1=PtPizero(px,(Double_t*) 0);
345 v1=V2Param(px,fgkV2param[fgSelectedCentrality]);
346 fgSelectedCentrality=k0510;
347 n2=PtPizero(px,(Double_t*) 0);
348 v2=V2Param(px,fgkV2param[fgSelectedCentrality]);
349 fgSelectedCentrality=selectedCentrality;
350 return (n1*v1+n2*v2)/(n1+n2);
351 break;
352 case k2040:
353 selectedCentrality=fgSelectedCentrality;
354 fgSelectedCentrality=k2030;
355 n1=PtPizero(px,(Double_t*) 0);
356 v1=V2Param(px,fgkV2param[fgSelectedCentrality]);
357 fgSelectedCentrality=k3040;
358 n2=PtPizero(px,(Double_t*) 0);
359 v2=V2Param(px,fgkV2param[fgSelectedCentrality]);
360 fgSelectedCentrality=selectedCentrality;
361 return (n1*v1+n2*v2)/(n1+n2);
362 break;
363
364 default:
365 return V2Param(px,fgkV2param[fgSelectedCentrality]);
366 }
6078e216 367}
368
4ae1c9f0 369
370
e40b9538 371//--------------------------------------------------------------------------
372//
373// Eta
374//
375//--------------------------------------------------------------------------
376Int_t AliGenEMlib::IpEta(TRandom *)
377{
4ae1c9f0 378 // Return eta pdg code
e40b9538 379 return 221;
380}
381
382Double_t AliGenEMlib::PtEta( const Double_t *px, const Double_t */*dummy*/ )
383{
4ae1c9f0 384
385 // fit functions and corresponding parameter of Eta pT for pp @ 2.76 TeV and @ 7 TeV
386 // and mtscaled pT
387
388 Double_t km = 0.;
389 Double_t kc = 0.;
390 Double_t kT = 0.;
391 Double_t kn = 0.;
392
393 switch(fgSelectedPtParam){
394 // Tsallis fit to final eta (PHOS+PCM) -> used for final publication
395 // for pp @ 7 TeV
396 case kPizeroEta7TeVpp:
397 km = 0.547853; kc = 2.496; kT = 0.229; kn = 6.985;
398 return PtTsallis(*px,km,kc,kT,kn);
399 break;
400 case kPizeroEta7TeVpplow:
401 km = 0.547853; kc = 1.970; kT = 0.253; kn = 7.591;
402 return PtTsallis(*px,km,kc,kT,kn);
403 break;
404 case kPizeroEta7TeVpphigh:
405 km = 0.547853; kc = 3.060; kT = 0.212; kn = 6.578;
406 return PtTsallis(*px,km,kc,kT,kn);
407 break;
408 // Tsallis fit to preliminary eta (QM'11)
409 // for pp @ 2.76 TeV
410 case kPizeroEta2760GeVpp:
411 km = 0.547853; kc = 1.971; kT = 0.188; kn = 6.308;
412 return PtTsallis(*px,km,kc,kT,kn);
413 case kPizeroEta2760GeVpplow:
414 km = 0.547853; kc = 1.228; kT = 0.220; kn = 7.030;
415 return PtTsallis(*px,km,kc,kT,kn);
416 break;
417 case kPizeroEta2760GeVpphigh:
418 km = 0.547853; kc = 2.802; kT = 0.164; kn = 5.815;
419 return PtTsallis(*px,km,kc,kT,kn);
420 break;
421
422 default:
0db2f441 423 return MtScal(*px,1);
4ae1c9f0 424 break;
425
426 }
427
e40b9538 428}
429
430Double_t AliGenEMlib::YEta( const Double_t *py, const Double_t */*dummy*/ )
431{
432 return YFlat(*py);
6078e216 433}
e40b9538 434
6078e216 435Double_t AliGenEMlib::V2Eta( const Double_t *px, const Double_t */*dummy*/ )
436{
4ae1c9f0 437 return KEtScal(*px,1); //V2Param(px,fgkV2param[1][fgSelectedV2Param]);
e40b9538 438}
439
440//--------------------------------------------------------------------------
441//
442// Rho
443//
444//--------------------------------------------------------------------------
445Int_t AliGenEMlib::IpRho(TRandom *)
446{
4ae1c9f0 447 // Return rho pdg code
e40b9538 448 return 113;
449}
450
451Double_t AliGenEMlib::PtRho( const Double_t *px, const Double_t */*dummy*/ )
452{
4ae1c9f0 453 // Rho pT
0db2f441 454 return MtScal(*px,2);
e40b9538 455}
456
457Double_t AliGenEMlib::YRho( const Double_t *py, const Double_t */*dummy*/ )
458{
459 return YFlat(*py);
4ae1c9f0 460}
e40b9538 461
4ae1c9f0 462Double_t AliGenEMlib::V2Rho( const Double_t *px, const Double_t */*dummy*/ )
463{
464 return KEtScal(*px,2);
e40b9538 465}
466
467//--------------------------------------------------------------------------
468//
469// Omega
470//
471//--------------------------------------------------------------------------
472Int_t AliGenEMlib::IpOmega(TRandom *)
473{
4ae1c9f0 474 // Return omega pdg code
e40b9538 475 return 223;
476}
477
478Double_t AliGenEMlib::PtOmega( const Double_t *px, const Double_t */*dummy*/ )
479{
4ae1c9f0 480 // Omega pT
0db2f441 481 return MtScal(*px,3);
e40b9538 482}
483
484Double_t AliGenEMlib::YOmega( const Double_t *py, const Double_t */*dummy*/ )
485{
486 return YFlat(*py);
4ae1c9f0 487}
e40b9538 488
4ae1c9f0 489Double_t AliGenEMlib::V2Omega( const Double_t *px, const Double_t */*dummy*/ )
490{
491 return KEtScal(*px,3);
e40b9538 492}
493
4ae1c9f0 494
e40b9538 495//--------------------------------------------------------------------------
496//
497// Etaprime
498//
499//--------------------------------------------------------------------------
500Int_t AliGenEMlib::IpEtaprime(TRandom *)
501{
4ae1c9f0 502 // Return etaprime pdg code
e40b9538 503 return 331;
504}
505
506Double_t AliGenEMlib::PtEtaprime( const Double_t *px, const Double_t */*dummy*/ )
507{
4ae1c9f0 508 // Eta pT
0db2f441 509 return MtScal(*px,4);
e40b9538 510}
511
512Double_t AliGenEMlib::YEtaprime( const Double_t *py, const Double_t */*dummy*/ )
513{
514 return YFlat(*py);
515
516}
517
4ae1c9f0 518Double_t AliGenEMlib::V2Etaprime( const Double_t *px, const Double_t */*dummy*/ )
519{
520 return KEtScal(*px,4);
521}
522
e40b9538 523//--------------------------------------------------------------------------
524//
525// Phi
526//
527//--------------------------------------------------------------------------
528Int_t AliGenEMlib::IpPhi(TRandom *)
529{
4ae1c9f0 530 // Return phi pdg code
e40b9538 531 return 333;
532}
533
534Double_t AliGenEMlib::PtPhi( const Double_t *px, const Double_t */*dummy*/ )
535{
4ae1c9f0 536 // Phi pT
0db2f441 537 return MtScal(*px,5);
e40b9538 538}
539
540Double_t AliGenEMlib::YPhi( const Double_t *py, const Double_t */*dummy*/ )
541{
542 return YFlat(*py);
e40b9538 543}
544
4ae1c9f0 545Double_t AliGenEMlib::V2Phi( const Double_t *px, const Double_t */*dummy*/ )
e40b9538 546{
4ae1c9f0 547 return KEtScal(*px,5);
548}
549
e40b9538 550//--------------------------------------------------------------------------
551//
4ae1c9f0 552// Jpsi
e40b9538 553//
554//--------------------------------------------------------------------------
4ae1c9f0 555Int_t AliGenEMlib::IpJpsi(TRandom *)
556{
557 // Return phi pdg code
558 return 443;
559}
560
561Double_t AliGenEMlib::PtJpsi( const Double_t *px, const Double_t */*dummy*/ )
562{
563 // Jpsi pT
564 return MtScal(*px,6);
565}
566
567Double_t AliGenEMlib::YJpsi( const Double_t *py, const Double_t */*dummy*/ )
568{
569 return YFlat(*py);
570}
571
572Double_t AliGenEMlib::V2Jpsi( const Double_t *px, const Double_t */*dummy*/ )
573{
574 const int oldSys=fgSelectedV2Systematic;
575 fgSelectedV2Systematic=kNoV2Sys;
576 double ret=0;
577
578 switch(oldSys){
579 case kLoV2Sys: ret=0; break;
580 case kNoV2Sys: ret=KEtScal(*px,6)/2; break;
581 case kUpV2Sys: ret=KEtScal(*px,6); break;
582 }
583
584 fgSelectedV2Systematic=oldSys;
585 return ret;
586}
587
588Double_t AliGenEMlib::YFlat(Double_t /*y*/)
589{
590 //--------------------------------------------------------------------------
591 //
592 // flat rapidity distribution
593 //
594 //--------------------------------------------------------------------------
e40b9538 595
596 Double_t dNdy = 1.;
597
598 return dNdy;
599
600}
601
602//=============================================================
603//
604// Mt-scaling
605//
606//=============================================================
607//
4ae1c9f0 608Double_t AliGenEMlib::MtScal(Double_t pt, Int_t np)
e40b9538 609{
4ae1c9f0 610 // Function for the calculation of the Pt distribution for a
611 // given particle np, from the pizero Pt distribution using
612 // mt scaling.
e40b9538 613
4ae1c9f0 614 Double_t scaledPt = sqrt(pt*pt + fgkHM[np]*fgkHM[np] - fgkHM[0]*fgkHM[0]);
e40b9538 615 Double_t scaledYield = PtPizero(&scaledPt, (Double_t*) 0);
616
4ae1c9f0 617 // VALUE MESON/PI AT 5 GeV/c
618 Double_t NormPt = 5.;
619 Double_t scaledNormPt = sqrt(NormPt*NormPt + fgkHM[np]*fgkHM[np] - fgkHM[0]*fgkHM[0]);
e40b9538 620
4ae1c9f0 621 Double_t norm = fgkMtFactor[int(bool(fgSelectedCentrality))][np] * (PtPizero(&NormPt, (Double_t*) 0) / PtPizero(&scaledNormPt, (Double_t*) 0));
e40b9538 622
623 return norm*(pt/scaledPt)*scaledYield;
624}
625
4ae1c9f0 626Double_t AliGenEMlib::KEtScal(const Double_t pt, Int_t np)
627{
628 const int nq=2; //number of quarks for particle np, here always 2
629 Double_t scaledPt = sqrt(pow(2/nq*(sqrt(pt*pt+fgkHM[np]*fgkHM[np])-fgkHM[np])+fgkHM[0],2)-fgkHM[0]*fgkHM[0]);
630 double val=V2Pizero(&scaledPt, (Double_t*) 0);
631 static const double syserr[12]={0., 0.09, 0.07, 0.06, 0.04, 0.04, 0.04, 0.05, 0., 0., 0., 0.}; //based on pi vs kaon
632 double sys=fgSelectedV2Systematic*min(fgkV2param[fgSelectedCentrality][0],fgkV2param[fgSelectedCentrality][8])*syserr[fgSelectedCentrality];
633 return TMath::Max(val+sys,0.0);
634}
635
6078e216 636Double_t AliGenEMlib::V2Param(const Double_t *px, const Double_t *par)
637{
638 // Very general parametrization of the v2
639
4ae1c9f0 640 const double &pt=px[0];
641 double val=CrossOverLc(par[4],par[3],pt)*(2*par[0]/(1+TMath::Exp(par[1]*(par[2]-pt)))-par[0])+CrossOverRc(par[4],par[3],pt)*((par[8]-par[5])/(1+TMath::Exp(par[6]*(pt-par[7])))+par[5]);
642 double sys=fgSelectedV2Systematic*par[11+fgSelectedV2Systematic*2]*pow(pt,par[12+fgSelectedV2Systematic*2]);
643 return TMath::Max(val+sys,0.0);
6078e216 644}
645
646Double_t AliGenEMlib::V2Flat(const Double_t */*px*/, const Double_t */*param*/)
647{
648 // Flat v2
649
650 return 1.0;
651}
652
4ae1c9f0 653
654
e40b9538 655//==========================================================================
656//
657// Set Getters
658//
659//==========================================================================
660
661typedef Double_t (*GenFunc) (const Double_t*, const Double_t*);
662
663typedef Int_t (*GenFuncIp) (TRandom *);
664
665GenFunc AliGenEMlib::GetPt(Int_t param, const char * tname) const
666{
4ae1c9f0 667 // Return pointer to pT parameterisation
e40b9538 668 GenFunc func=0;
669 TString sname(tname);
670
671 switch (param)
672 {
673 case kPizero:
674 func=PtPizero;
675 break;
676 case kEta:
677 func=PtEta;
678 break;
679 case kRho:
680 func=PtRho;
681 break;
682 case kOmega:
683 func=PtOmega;
684 break;
685 case kEtaprime:
686 func=PtEtaprime;
687 break;
688 case kPhi:
689 func=PtPhi;
690 break;
4ae1c9f0 691 case kJpsi:
692 func=PtJpsi;
693 break;
e40b9538 694
695 default:
696 func=0;
697 printf("<AliGenEMlib::GetPt> unknown parametrisation\n");
698 }
699 return func;
700}
701
702GenFunc AliGenEMlib::GetY(Int_t param, const char * tname) const
703{
4ae1c9f0 704 // Return pointer to y- parameterisation
e40b9538 705 GenFunc func=0;
706 TString sname(tname);
707
708 switch (param)
709 {
710 case kPizero:
711 func=YPizero;
712 break;
713 case kEta:
714 func=YEta;
715 break;
716 case kRho:
717 func=YRho;
718 break;
719 case kOmega:
720 func=YOmega;
721 break;
722 case kEtaprime:
723 func=YEtaprime;
724 break;
725 case kPhi:
726 func=YPhi;
727 break;
4ae1c9f0 728 case kJpsi:
729 func=YJpsi;
730 break;
e40b9538 731
732 default:
733 func=0;
734 printf("<AliGenEMlib::GetY> unknown parametrisation\n");
735 }
736 return func;
737}
738
739GenFuncIp AliGenEMlib::GetIp(Int_t param, const char * tname) const
740{
4ae1c9f0 741 // Return pointer to particle type parameterisation
e40b9538 742 GenFuncIp func=0;
743 TString sname(tname);
744
745 switch (param)
746 {
747 case kPizero:
748 func=IpPizero;
749 break;
750 case kEta:
751 func=IpEta;
752 break;
753 case kRho:
754 func=IpRho;
755 break;
756 case kOmega:
757 func=IpOmega;
758 break;
759 case kEtaprime:
760 func=IpEtaprime;
761 break;
762 case kPhi:
763 func=IpPhi;
764 break;
4ae1c9f0 765 case kJpsi:
766 func=IpJpsi;
767 break;
e40b9538 768
769 default:
770 func=0;
771 printf("<AliGenEMlib::GetIp> unknown parametrisation\n");
772 }
773 return func;
774}
775
6078e216 776GenFunc AliGenEMlib::GetV2(Int_t param, const char * tname) const
777{
778 // Return pointer to v2-parameterisation
779 GenFunc func=0;
780 TString sname(tname);
e40b9538 781
6078e216 782 switch (param)
783 {
784 case kPizero:
785 func=V2Pizero;
786 break;
787 case kEta:
788 func=V2Eta;
789 break;
790 case kRho:
791 func=V2Pizero;
792 break;
793 case kOmega:
794 func=V2Pizero;
795 break;
796 case kEtaprime:
797 func=V2Pizero;
798 break;
799 case kPhi:
800 func=V2Pizero;
801 break;
4ae1c9f0 802 case kJpsi:
803 func=V2Jpsi;
804 break;
e40b9538 805
6078e216 806 default:
807 func=0;
808 printf("<AliGenEMlib::GetV2> unknown parametrisation\n");
809 }
810 return func;
811}
4ae1c9f0 812