1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 // Library class for particle pt and y distributions used for
19 // muon spectrometer simulations.
20 // To be used with AliGenParam.
21 // The following particle typed can be simulated:
22 // pi, K, phi, omega, eta, J/Psi, Upsilon, charm and beauty mesons.
24 // andreas.morsch@cern.ch
29 #include "TDatabasePDG.h"
31 #include "AliGenMUONlib.h"
33 ClassImp(AliGenMUONlib)
36 Double_t AliGenMUONlib::PtPion(const Double_t *px, const Double_t* /*dummy*/)
39 // PT-PARAMETERIZATION CDF, PRL 61(88) 1819
40 // POWER LAW FOR PT > 500 MEV
41 // MT SCALING BELOW (T=160 MEV)
43 const Double_t kp0 = 1.3;
44 const Double_t kxn = 8.28;
45 const Double_t kxlim=0.5;
46 const Double_t kt=0.160;
47 const Double_t kxmpi=0.139;
49 Double_t y, y1, xmpi2, ynorm, a;
52 y1=TMath::Power(kp0/(kp0+kxlim),kxn);
54 ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
57 y=a*TMath::Power(kp0/(kp0+x),kxn);
59 y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
65 Double_t AliGenMUONlib::YPion( const Double_t *py, const Double_t */*dummy*/)
68 Double_t y=TMath::Abs(*py);
70 const Double_t ka = 7000.;
71 const Double_t kdy = 4.;
72 Double_t ex = y*y/(2*kdy*kdy);
73 return ka*TMath::Exp(-ex);
75 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
78 // particle composition
80 Int_t AliGenMUONlib::IpPion(TRandom *ran)
83 if (ran->Rndm() < 0.5) {
90 //____________________________________________________________
94 Double_t AliGenMUONlib::PtScal(Double_t pt, Int_t np)
96 // SCALING EN MASSE PAR RAPPORT A PTPI
97 // MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
98 const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
99 // VALUE MESON/PI AT 5 GEV
100 const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
102 Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
103 Double_t fmax2=f5/kfmax[np];
105 Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
106 Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
107 (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
108 return fmtscal*ptpion;
114 //____________________________________________________________
115 Double_t AliGenMUONlib::PtKaon( const Double_t *px, const Double_t */*dummy*/)
118 return PtScal(*px,2);
122 //____________________________________________________________
123 Double_t AliGenMUONlib::YKaon( const Double_t *py, const Double_t */*dummy*/)
126 Double_t y=TMath::Abs(*py);
128 const Double_t ka = 1000.;
129 const Double_t kdy = 4.;
131 Double_t ex = y*y/(2*kdy*kdy);
132 return ka*TMath::Exp(-ex);
135 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
138 // particle composition
140 Int_t AliGenMUONlib::IpKaon(TRandom *ran)
143 if (ran->Rndm() < 0.5) {
154 //____________________________________________________________
155 Double_t AliGenMUONlib::PtJpsiPPdummy(Double_t x, Double_t energy)
159 // from the fit of RHIC, CDF & LHC data, see arXiv:1103.2394
161 const Double_t kpt0 = 1.04*TMath::Power(energy,0.101);
162 const Double_t kxn = 3.9;
164 Double_t pass1 = 1.+0.363*(x/kpt0)*(x/kpt0);
165 return x/TMath::Power(pass1,kxn);
168 Double_t AliGenMUONlib::PtJpsiPP7000(const Double_t *px, const Double_t */*dummy*/)
173 return PtJpsiPPdummy(*px,7000);
176 Double_t AliGenMUONlib::PtJpsiPP2760(const Double_t *px, const Double_t */*dummy*/)
181 return PtJpsiPPdummy(*px,2760);
184 Double_t AliGenMUONlib::PtJpsiPP8800(const Double_t *px, const Double_t */*dummy*/)
189 return PtJpsiPPdummy(*px,8800);
192 Double_t AliGenMUONlib::PtJpsiPbPb2760ShFdummy(Double_t x, Int_t n)
194 // J/Psi shadowing factor vs pT for PbPb min. bias and 11 centr. bins (in 2.5<y<4)
196 // PbPb 2.76 TeV, for EKS98, minimum bias shadowing factor = 0.66 in 4pi
197 // S.Grigoryan, details presented at the PWG3-Muon meeting (05.10.2011)
198 // https://indico.cern.ch/conferenceDisplay.py?confId=157367
200 const Double_t f1[12] = {1, 1.128, 1.097, 1.037, 0.937, 0.821, 0.693, 0.558,
201 0.428, 0.317, 0.231, 0.156};
202 const Double_t f2[12] = {1, 1.313, 1.202, 1.039, 0.814, 0.593, 0.391, 0.224,
203 0.106, 0.041, 0.013, 0.002};
204 const Double_t c1[7] = {1.6077e+00, 7.6300e-02,-7.1880e-03, 3.4067e-04,
205 -9.2776e-06,1.5138e-07, 1.4652e-09};
206 const Double_t c2[7] = {6.2047e-01, 5.7653e-02,-4.1414e-03, 1.0301e-04,
207 9.6205e-07,-7.4098e-08, 5.0946e-09};
210 y1 = c1[j = 6]; y2 = c2[6];
211 while (j > 0) {y1 = y1 * x + c1[--j]; y2 = y2 * x + c2[j];}
213 y1 /= 1.+c1[6]*TMath::Power(x,6);
214 y2 /= 1.+c2[6]*TMath::Power(x,6);
216 y1 = 1 + (y1-2)*f1[n] + (y2+1-y1)*f2[n];
221 Double_t AliGenMUONlib::PtJpsiPbPb2760(const Double_t *px, const Double_t *dummy)
224 // PbPb 2.76 TeV, minimum bias 0-100 %
226 return PtJpsiPbPb2760ShFdummy(*px, 0) * PtJpsiPP2760(px, dummy);
229 Double_t AliGenMUONlib::PtJpsiPbPb2760c1(const Double_t *px, const Double_t *dummy)
232 // PbPb 2.76 TeV, 1st centrality bin 0-5 %
234 return PtJpsiPbPb2760ShFdummy(*px, 1) * PtJpsiPP2760(px, dummy);
237 Double_t AliGenMUONlib::PtJpsiPbPb2760c2(const Double_t *px, const Double_t *dummy)
240 // PbPb 2.76 TeV, 2nd centrality bin 5-10 %
242 return PtJpsiPbPb2760ShFdummy(*px, 2) * PtJpsiPP2760(px, dummy);
245 Double_t AliGenMUONlib::PtJpsiPbPb2760c3(const Double_t *px, const Double_t *dummy)
248 // PbPb 2.76 TeV, 3rd centrality bin 10-20 %
250 return PtJpsiPbPb2760ShFdummy(*px, 3) * PtJpsiPP2760(px, dummy);
253 Double_t AliGenMUONlib::PtJpsiPbPb2760c4(const Double_t *px, const Double_t *dummy)
256 // PbPb 2.76 TeV, 4th centrality bin 20-30 %
258 return PtJpsiPbPb2760ShFdummy(*px, 4) * PtJpsiPP2760(px, dummy);
261 Double_t AliGenMUONlib::PtJpsiPbPb2760c5(const Double_t *px, const Double_t *dummy)
264 // PbPb 2.76 TeV, 5th centrality bin 30-40 %
266 return PtJpsiPbPb2760ShFdummy(*px, 5) * PtJpsiPP2760(px, dummy);
269 Double_t AliGenMUONlib::PtJpsiPbPb2760c6(const Double_t *px, const Double_t *dummy)
272 // PbPb 2.76 TeV, 6th centrality bin 40-50 %
274 return PtJpsiPbPb2760ShFdummy(*px, 6) * PtJpsiPP2760(px, dummy);
277 Double_t AliGenMUONlib::PtJpsiPbPb2760c7(const Double_t *px, const Double_t *dummy)
280 // PbPb 2.76 TeV, 7th centrality bin 50-60 %
282 return PtJpsiPbPb2760ShFdummy(*px, 7) * PtJpsiPP2760(px, dummy);
285 Double_t AliGenMUONlib::PtJpsiPbPb2760c8(const Double_t *px, const Double_t *dummy)
288 // PbPb 2.76 TeV, 8th centrality bin 60-70 %
290 return PtJpsiPbPb2760ShFdummy(*px, 8) * PtJpsiPP2760(px, dummy);
293 Double_t AliGenMUONlib::PtJpsiPbPb2760c9(const Double_t *px, const Double_t *dummy)
296 // PbPb 2.76 TeV, 9th centrality bin 70-80 %
298 return PtJpsiPbPb2760ShFdummy(*px, 9) * PtJpsiPP2760(px, dummy);
301 Double_t AliGenMUONlib::PtJpsiPbPb2760c10(const Double_t *px, const Double_t *dummy)
304 // PbPb 2.76 TeV, 10th centrality bin 80-90 %
306 return PtJpsiPbPb2760ShFdummy(*px, 10) * PtJpsiPP2760(px, dummy);
309 Double_t AliGenMUONlib::PtJpsiPbPb2760c11(const Double_t *px, const Double_t *dummy)
312 // PbPb 2.76 TeV, 11th centrality bin 90-100 %
314 return PtJpsiPbPb2760ShFdummy(*px, 11) * PtJpsiPP2760(px, dummy);
317 Double_t AliGenMUONlib::PtJpsiPPb8800ShFdummy(Double_t x, Int_t n)
319 // J/Psi shadowing factor vs pT for pPb min. bias and 4 centr. bins (in 2.5<y<4)
321 // pPb 8.8 TeV, for EKS98, minimum bias shadowing factor = 0.78 in 4pi
323 const Double_t f[5] = {1, 1.33, 1.05, 0.67, 0.23};
324 const Double_t c[7] = {6.4922e-01, 6.4857e-03, 4.7268e-03,-9.5075e-04,
325 8.4075e-05,-4.2006e-06, 4.9970e-07};
329 while (j > 0) y = y * x + c[--j];
330 y /= 1 + c[6]*TMath::Power(x,6);
332 return 1 + (y-1)*f[n];
335 Double_t AliGenMUONlib::PtJpsiPPb8800(const Double_t *px, const Double_t *dummy)
338 // pPb 8.8 TeV, minimum bias 0-100 %
340 return PtJpsiPPb8800ShFdummy(*px, 0) * PtJpsiPP8800(px, dummy);
343 Double_t AliGenMUONlib::PtJpsiPPb8800c1(const Double_t *px, const Double_t *dummy)
346 // pPb 8.8 TeV, 1st centrality bin 0-20 %
348 return PtJpsiPPb8800ShFdummy(*px, 1) * PtJpsiPP8800(px, dummy);
351 Double_t AliGenMUONlib::PtJpsiPPb8800c2(const Double_t *px, const Double_t *dummy)
354 // pPb 8.8 TeV, 2nd centrality bin 20-40 %
356 return PtJpsiPPb8800ShFdummy(*px, 2) * PtJpsiPP8800(px, dummy);
359 Double_t AliGenMUONlib::PtJpsiPPb8800c3(const Double_t *px, const Double_t *dummy)
362 // pPb 8.8 TeV, 3rd centrality bin 40-60 %
364 return PtJpsiPPb8800ShFdummy(*px, 3) * PtJpsiPP8800(px, dummy);
367 Double_t AliGenMUONlib::PtJpsiPPb8800c4(const Double_t *px, const Double_t *dummy)
370 // pPb 8.8 TeV, 4th centrality bin 60-100 %
372 return PtJpsiPPb8800ShFdummy(*px, 4) * PtJpsiPP8800(px, dummy);
375 Double_t AliGenMUONlib::PtJpsiPbP8800ShFdummy(Double_t x, Int_t n)
377 // J/Psi shadowing factor vs pT for Pbp min. bias and 4 centr. bins (in 2.5<y<4)
379 // Pbp 8.8 TeV, for EKS98, minimum bias shadowing factor = 0.78 in 4pi
381 const Double_t f[5] = {1, 1.33, 1.05, 0.67, 0.23};
382 const Double_t c[7] = {8.7562e-01, 2.1944e-02, 7.8509e-03,-1.3979e-03,
383 3.8513e-05, 4.2008e-06, 1.7088e-06};
387 while (j > 0) y = y * x + c[--j];
388 y /= 1 + c[6]*TMath::Power(x,6);
390 return 1 + (y-1)*f[n];
393 Double_t AliGenMUONlib::PtJpsiPbP8800(const Double_t *px, const Double_t *dummy)
396 // Pbp 8.8 TeV, minimum bias 0-100 %
398 return PtJpsiPbP8800ShFdummy(*px, 0) * PtJpsiPP8800(px, dummy);
401 Double_t AliGenMUONlib::PtJpsiPbP8800c1(const Double_t *px, const Double_t *dummy)
404 // Pbp 8.8 TeV, 1st centrality bin 0-20 %
406 return PtJpsiPbP8800ShFdummy(*px, 1) * PtJpsiPP8800(px, dummy);
409 Double_t AliGenMUONlib::PtJpsiPbP8800c2(const Double_t *px, const Double_t *dummy)
412 // Pbp 8.8 TeV, 2nd centrality bin 20-40 %
414 return PtJpsiPbP8800ShFdummy(*px, 2) * PtJpsiPP8800(px, dummy);
417 Double_t AliGenMUONlib::PtJpsiPbP8800c3(const Double_t *px, const Double_t *dummy)
420 // Pbp 8.8 TeV, 3rd centrality bin 40-60 %
422 return PtJpsiPbP8800ShFdummy(*px, 3) * PtJpsiPP8800(px, dummy);
425 Double_t AliGenMUONlib::PtJpsiPbP8800c4(const Double_t *px, const Double_t *dummy)
428 // Pbp 8.8 TeV, 4th centrality bin 60-100 %
430 return PtJpsiPbP8800ShFdummy(*px, 4) * PtJpsiPP8800(px, dummy);
433 Double_t AliGenMUONlib::PtJpsi( const Double_t *px, const Double_t */*dummy*/)
436 const Double_t kpt0 = 4.;
437 const Double_t kxn = 3.6;
440 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
441 return x/TMath::Power(pass1,kxn);
444 Double_t AliGenMUONlib::PtJpsiCDFscaled( const Double_t *px, const Double_t */*dummy*/)
449 // scaled from CDF data at 2 TeV
450 // see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
452 const Double_t kpt0 = 5.100;
453 const Double_t kxn = 4.102;
456 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
457 return x/TMath::Power(pass1,kxn);
460 Double_t AliGenMUONlib::PtJpsiCDFscaledPP( const Double_t *px, const Double_t */*dummy*/)
465 // scaled from CDF data at 2 TeV
467 const Double_t kpt0 = 5.630;
468 const Double_t kxn = 4.071;
471 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
472 return x/TMath::Power(pass1,kxn);
475 Double_t AliGenMUONlib::PtJpsiCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
480 // scaled from CDF data at 2 TeV
482 const Double_t kpt0 = 5.334;
483 const Double_t kxn = 4.071;
486 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
487 return x/TMath::Power(pass1,kxn);
490 Double_t AliGenMUONlib::PtJpsiCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
495 // scaled from CDF data at 2 TeV
497 const Double_t kpt0 = 5.245;
498 const Double_t kxn = 4.071;
501 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
502 return x/TMath::Power(pass1,kxn);
505 Double_t AliGenMUONlib::PtJpsiCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
510 // scaled from CDF data at 2 TeV
512 const Double_t kpt0 = 5.072;
513 const Double_t kxn = 4.071;
516 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
517 return x/TMath::Power(pass1,kxn);
520 Double_t AliGenMUONlib::PtJpsiCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
525 // scaled from CDF data at 2 TeV
527 const Double_t kpt0 = 4.647;
528 const Double_t kxn = 4.071;
531 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
532 return x/TMath::Power(pass1,kxn);
535 Double_t AliGenMUONlib::PtJpsiCDFscaledPP3( const Double_t *px, const Double_t */*dummy*/)
540 // scaled from CDF data at 1.9 TeV
542 const Double_t kpt0 = 4.435;
543 const Double_t kxn = 4.071;
546 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
547 return x/TMath::Power(pass1,kxn);
550 Double_t AliGenMUONlib::PtJpsiCDFscaledPP2( const Double_t *px, const Double_t */*dummy*/)
555 // fit of the CDF data at 1.96 TeV
557 const Double_t kpt0 = 4.233;
558 const Double_t kxn = 4.071;
561 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
562 return x/TMath::Power(pass1,kxn);
565 Double_t AliGenMUONlib::PtJpsiCDFscaledPPb9( const Double_t *px, const Double_t *dummy)
569 // pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.79
571 Double_t c[5] = {6.42774e-01, 1.86168e-02, -6.77296e-04, 8.93512e-06, 1.31586e-07};
576 while (j > 0) y = y * x + c[--j];
578 Double_t d = 1.+c[4]*TMath::Power(x,4);
579 return y/d * AliGenMUONlib::PtJpsiCDFscaledPP9(px,dummy);
582 Double_t AliGenMUONlib::PtJpsiCDFscaledPbP9( const Double_t *px, const Double_t *dummy)
586 // Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.79
588 Double_t c[5] = {8.58557e-01, 5.39791e-02, -4.75180e-03, 2.49463e-04, 5.52396e-05};
593 while (j > 0) y = y * x + c[--j];
595 Double_t d = 1.+c[4]*TMath::Power(x,4);
596 return y/d * AliGenMUONlib::PtJpsiCDFscaledPP9(px,dummy);
599 Double_t AliGenMUONlib::PtJpsiCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
603 // PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.66
605 Double_t c[5] = {6.01022e-01, 4.70988e-02, -2.27917e-03, 3.09885e-05, 1.31955e-06};
610 while (j > 0) y = y * x + c[--j];
612 Double_t d = 1.+c[4]*TMath::Power(x,4);
613 return y/d * AliGenMUONlib::PtJpsiCDFscaledPP4(px,dummy);
616 Double_t AliGenMUONlib::PtJpsiFlat( const Double_t */*px*/, const Double_t */*dummy*/ )
621 Double_t AliGenMUONlib::PtJpsiPbPb( const Double_t *px, const Double_t */*dummy*/)
628 // mc = 1.4 GeV, pt-kick 1 GeV
632 -2.13098e+00, 9.46552e+00, -5.06799e+00, 1.27260e+00,
633 -1.83806e-01, 1.55853e-02, -7.23241e-04, 1.42105e-05
640 while (j > 0) y = y * x +c[--j];
641 y = x * TMath::Exp(y);
648 Double_t AliGenMUONlib::PtJpsiBPbPb( const Double_t *px, const Double_t */*dummy*/)
652 Double_t x0 = 4.0384;
656 Double_t y = x / TMath::Power((1. + (x/x0)*(x/x0)), n);
662 Double_t AliGenMUONlib::PtJpsiPP( const Double_t *px, const Double_t */*dummy*/)
669 // mc = 1.4 GeV, pt-kick 1 GeV
672 Float_t c[4] = {8.47471e+00, -1.93567e+00, 1.50271e-01, -5.51212e-03};
678 while (j > 0) y = y * x +c[--j];
679 y = x * TMath::Exp(y);
688 //____________________________________________________________
689 Double_t AliGenMUONlib::YJpsiPPdummy(Double_t x, Double_t energy)
693 // from the fit of RHIC + LHC data, see arXiv:1103.2394
695 x = x/TMath::Log(energy/3.097);
697 Double_t y = TMath::Exp(-x/0.4/0.4/2);
702 Double_t AliGenMUONlib::YJpsiPPpoly(Double_t x, Double_t energy)
706 // from the fit of RHIC + LHC data, see arXiv:1103.2394
708 x = x/TMath::Log(energy/3.097);
710 Double_t y = 1 - 6.9*x*x;
715 Double_t AliGenMUONlib::YJpsiPP7000(const Double_t *px, const Double_t */*dummy*/)
720 return YJpsiPPdummy(*px, 7000);
723 Double_t AliGenMUONlib::YJpsiPP2760(const Double_t *px, const Double_t */*dummy*/)
728 return YJpsiPPdummy(*px, 2760);
731 Double_t AliGenMUONlib::YJpsiPP8800(const Double_t *px, const Double_t */*dummy*/)
736 return YJpsiPPdummy(*px, 8800);
739 Double_t AliGenMUONlib::YJpsiPPpoly7000(const Double_t *px, const Double_t */*dummy*/)
744 return YJpsiPPpoly(*px, 7000);
747 Double_t AliGenMUONlib::YJpsiPPpoly2760(const Double_t *px, const Double_t */*dummy*/)
752 return YJpsiPPpoly(*px, 2760);
755 Double_t AliGenMUONlib::YJpsiPbPb2760ShFdummy(Double_t x, Int_t n)
757 // J/Psi shadowing factor vs y for PbPb min. bias and 11 centr. bins
759 // PbPb 2.76 TeV, for EKS98, minimum bias shadowing factor = 0.66 in 4pi
761 const Double_t f1[12] = {1, 1.128, 1.097, 1.037, 0.937, 0.821, 0.693, 0.558,
762 0.428, 0.317, 0.231, 0.156};
763 const Double_t f2[12] = {1, 1.313, 1.202, 1.039, 0.814, 0.593, 0.391, 0.224,
764 0.106, 0.041, 0.013, 0.002};
765 const Double_t c1[5] = {1.5591e+00, 7.5827e-03, 2.0676e-03,-1.1717e-04, 1.5237e-06};
766 const Double_t c2[5] = {6.0861e-01, 4.8854e-03, 1.3685e-03,-7.9182e-05, 1.0475e-06};
771 y1 = c1[j = 4]; y2 = c2[4];
772 while (j > 0) {y1 = y1 * x + c1[--j]; y2 = y2 * x + c2[j];}
774 y1 = 1 + (y1-2)*f1[n] + (y2+1-y1)*f2[n];
779 Double_t AliGenMUONlib::YJpsiPbPb2760(const Double_t *px, const Double_t *dummy)
782 // PbPb 2.76 TeV, minimum bias 0-100 %
784 return YJpsiPbPb2760ShFdummy(*px, 0) * YJpsiPP2760(px, dummy);
787 Double_t AliGenMUONlib::YJpsiPbPb2760c1(const Double_t *px, const Double_t *dummy)
790 // PbPb 2.76 TeV, 1st centrality bin 0-5 %
792 return YJpsiPbPb2760ShFdummy(*px, 1) * YJpsiPP2760(px, dummy);
795 Double_t AliGenMUONlib::YJpsiPbPb2760c2(const Double_t *px, const Double_t *dummy)
798 // PbPb 2.76 TeV, 2nd centrality bin 5-10 %
800 return YJpsiPbPb2760ShFdummy(*px, 2) * YJpsiPP2760(px, dummy);
803 Double_t AliGenMUONlib::YJpsiPbPb2760c3(const Double_t *px, const Double_t *dummy)
806 // PbPb 2.76 TeV, 3rd centrality bin 10-20 %
808 return YJpsiPbPb2760ShFdummy(*px, 3) * YJpsiPP2760(px, dummy);
811 Double_t AliGenMUONlib::YJpsiPbPb2760c4(const Double_t *px, const Double_t *dummy)
814 // PbPb 2.76 TeV, 4th centrality bin 20-30 %
816 return YJpsiPbPb2760ShFdummy(*px, 4) * YJpsiPP2760(px, dummy);
819 Double_t AliGenMUONlib::YJpsiPbPb2760c5(const Double_t *px, const Double_t *dummy)
822 // PbPb 2.76 TeV, 5th centrality bin 30-40 %
824 return YJpsiPbPb2760ShFdummy(*px, 5) * YJpsiPP2760(px, dummy);
827 Double_t AliGenMUONlib::YJpsiPbPb2760c6(const Double_t *px, const Double_t *dummy)
830 // PbPb 2.76 TeV, 6th centrality bin 40-50 %
832 return YJpsiPbPb2760ShFdummy(*px, 6) * YJpsiPP2760(px, dummy);
835 Double_t AliGenMUONlib::YJpsiPbPb2760c7(const Double_t *px, const Double_t *dummy)
838 // PbPb 2.76 TeV, 7th centrality bin 50-60 %
840 return YJpsiPbPb2760ShFdummy(*px, 7) * YJpsiPP2760(px, dummy);
843 Double_t AliGenMUONlib::YJpsiPbPb2760c8(const Double_t *px, const Double_t *dummy)
846 // PbPb 2.76 TeV, 8th centrality bin 60-70 %
848 return YJpsiPbPb2760ShFdummy(*px, 8) * YJpsiPP2760(px, dummy);
851 Double_t AliGenMUONlib::YJpsiPbPb2760c9(const Double_t *px, const Double_t *dummy)
854 // PbPb 2.76 TeV, 9th centrality bin 70-80 %
856 return YJpsiPbPb2760ShFdummy(*px, 9) * YJpsiPP2760(px, dummy);
859 Double_t AliGenMUONlib::YJpsiPbPb2760c10(const Double_t *px, const Double_t *dummy)
862 // PbPb 2.76 TeV, 10th centrality bin 80-90 %
864 return YJpsiPbPb2760ShFdummy(*px, 10) * YJpsiPP2760(px, dummy);
867 Double_t AliGenMUONlib::YJpsiPbPb2760c11(const Double_t *px, const Double_t *dummy)
870 // PbPb 2.76 TeV, 11th centrality bin 90-100 %
872 return YJpsiPbPb2760ShFdummy(*px, 11) * YJpsiPP2760(px, dummy);
875 Double_t AliGenMUONlib::YJpsiPP8800dummy(Double_t px)
877 return AliGenMUONlib::YJpsiPP8800(&px, (Double_t*) 0);
880 Double_t AliGenMUONlib::YJpsiPPb8800ShFdummy(Double_t x, Int_t n)
882 // J/Psi shadowing factor vs y for pPb min. bias and 4 centr. bins
884 // pPb 8.8 TeV, for EKS98, minimum bias shadowing factor = 0.78 in 4pi
886 const Double_t f[5] = {1, 1.33, 1.05, 0.67, 0.23};
887 const Double_t c[7] = {7.4372e-01, 2.3299e-02, 2.8678e-03, 1.9595e-03,
888 3.2849e-04,-4.0547e-05,-7.9732e-06};
892 while (j > 0) y = y * x + c[--j];
895 return 1 +(y-1)*f[n];
898 Double_t AliGenMUONlib::YJpsiPPb8800(const Double_t *px, const Double_t */*dummy*/)
901 // pPb 8.8 TeV, minimum bias 0-100 %
903 Double_t x = px[0] + 0.47; // rapidity shift
904 return YJpsiPPb8800ShFdummy(x, 0) * YJpsiPP8800dummy(x);
907 Double_t AliGenMUONlib::YJpsiPPb8800c1(const Double_t *px, const Double_t */*dummy*/)
910 // pPb 8.8 TeV, 1st centrality bin 0-20 %
912 Double_t x = px[0] + 0.47; // rapidity shift
913 return YJpsiPPb8800ShFdummy(x, 1) * YJpsiPP8800dummy(x);
916 Double_t AliGenMUONlib::YJpsiPPb8800c2(const Double_t *px, const Double_t */*dummy*/)
919 // pPb 8.8 TeV, 2nd centrality bin 20-40 %
921 Double_t x = px[0] + 0.47; // rapidity shift
922 return YJpsiPPb8800ShFdummy(x, 2) * YJpsiPP8800dummy(x);
925 Double_t AliGenMUONlib::YJpsiPPb8800c3(const Double_t *px, const Double_t */*dummy*/)
928 // pPb 8.8 TeV, 3rd centrality bin 40-60 %
930 Double_t x = px[0] + 0.47; // rapidity shift
931 return YJpsiPPb8800ShFdummy(x, 3) * YJpsiPP8800dummy(x);
934 Double_t AliGenMUONlib::YJpsiPPb8800c4(const Double_t *px, const Double_t */*dummy*/)
937 // pPb 8.8 TeV, 4th centrality bin 60-100 %
939 Double_t x = px[0] + 0.47; // rapidity shift
940 return YJpsiPPb8800ShFdummy(x, 4) * YJpsiPP8800dummy(x);
943 Double_t AliGenMUONlib::YJpsiPbP8800(const Double_t *px, const Double_t */*dummy*/)
946 // Pbp 8.8 TeV, minimum bias 0-100 %
948 Double_t x = -px[0] + 0.47; // rapidity shift
949 return YJpsiPPb8800ShFdummy(x, 0) * YJpsiPP8800dummy(x);
952 Double_t AliGenMUONlib::YJpsiPbP8800c1(const Double_t *px, const Double_t */*dummy*/)
955 // Pbp 8.8 TeV, 1st centrality bin 0-20 %
957 Double_t x = -px[0] + 0.47; // rapidity shift
958 return YJpsiPPb8800ShFdummy(x, 1) * YJpsiPP8800dummy(x);
961 Double_t AliGenMUONlib::YJpsiPbP8800c2(const Double_t *px, const Double_t */*dummy*/)
964 // Pbp 8.8 TeV, 2nd centrality bin 20-40 %
966 Double_t x = -px[0] + 0.47; // rapidity shift
967 return YJpsiPPb8800ShFdummy(x, 2) * YJpsiPP8800dummy(x);
970 Double_t AliGenMUONlib::YJpsiPbP8800c3(const Double_t *px, const Double_t */*dummy*/)
973 // Pbp 8.8 TeV, 3rd centrality bin 40-60 %
975 Double_t x = -px[0] + 0.47; // rapidity shift
976 return YJpsiPPb8800ShFdummy(x, 3) * YJpsiPP8800dummy(x);
979 Double_t AliGenMUONlib::YJpsiPbP8800c4(const Double_t *px, const Double_t */*dummy*/)
982 // Pbp 8.8 TeV, 4th centrality bin 60-100 %
984 Double_t x = -px[0] + 0.47; // rapidity shift
985 return YJpsiPPb8800ShFdummy(x, 4) * YJpsiPP8800dummy(x);
988 Double_t AliGenMUONlib::YJpsi(const Double_t *py, const Double_t */*dummy*/)
991 const Double_t ky0 = 4.;
992 const Double_t kb=1.;
994 Double_t y=TMath::Abs(*py);
999 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
1003 Double_t AliGenMUONlib::YJpsiFlat( const Double_t */*py*/, const Double_t */*dummy*/ )
1009 Double_t AliGenMUONlib::YJpsiPbPb( const Double_t *px, const Double_t */*dummy*/)
1019 // mc = 1.4 GeV, pt-kick 1 GeV
1021 Double_t c[5] = {-6.03425e+02, 4.98257e+02, -1.38794e+02, 1.62209e+01, -6.85955e-01};
1022 Double_t x = TMath::Abs(px[0]);
1030 while (j > 0) y = y * x + c[--j];
1038 Double_t AliGenMUONlib::YJpsiCDFscaled( const Double_t *px, const Double_t* dummy)
1041 return AliGenMUONlib::YJpsiPbPb(px, dummy);
1044 Double_t AliGenMUONlib::YJpsiCDFscaledPP( const Double_t *px, const Double_t* dummy)
1047 return AliGenMUONlib::YJpsiPP(px, dummy);
1050 Double_t AliGenMUONlib::YJpsiCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
1055 // scaled from YJpsiPP(14 TeV) using 10 TeV / 14 TeV ratio of y-spectra in LO pQCD.
1056 // see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
1059 Double_t c[5] = {2.46681e+01, 8.91486e+01, -3.21227e+01, 3.63075e+00, -1.32047e-01};
1061 Double_t x = TMath::Abs(px[0]);
1065 y = 98.523 - 1.3664 * x * x;
1066 } else if (x < 7.5) {
1069 while (j > 0) y = y * x + c[--j];
1079 Double_t AliGenMUONlib::YJpsiCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
1084 // rescaling of YJpsiPP(14 TeV) using 8.8 TeV / 14 TeV ratio of y-spectra in LO QCD
1086 Double_t c[5] = {3.33882e+02, -1.30980e+02, 2.59082e+01, -3.08935e+00, 1.56375e-01};
1087 Double_t x = TMath::Abs(px[0]);
1091 y = 99.236 - 1.5498 * x * x;
1092 } else if (x < 7.4) {
1095 while (j > 0) y = y * x + c[--j];
1105 Double_t AliGenMUONlib::YJpsiCDFscaledPP9dummy(Double_t px)
1107 return AliGenMUONlib::YJpsiCDFscaledPP9(&px, (Double_t*) 0);
1110 Double_t AliGenMUONlib::YJpsiCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
1115 // scaled from YJpsiPP(14 TeV) using 7 TeV / 14 TeV ratio of y-spectra in LO pQCD.
1118 Double_t c[5] = {6.71181e+02, -3.69240e+02, 8.89644e+01, -1.04937e+01, 4.80959e-01};
1120 Double_t x = TMath::Abs(px[0]);
1124 y = 100.78 - 1.8353 * x * x;
1125 } else if (x < 7.3) {
1128 while (j > 0) y = y * x + c[--j];
1138 Double_t AliGenMUONlib::YJpsiCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
1143 // rescaling of YJpsiPP(14 TeV) using 3.94 TeV / 14 TeV ratio of y-spectra in LO QCD
1145 Double_t c[5] = {4.00785e+02, -1.41159e+01, -3.28599e+01, 5.53048e+00, -2.45151e-01};
1146 Double_t x = TMath::Abs(px[0]);
1150 y = 107.389 - 2.7454 * x * x;
1151 } else if (x < 7.0) {
1154 while (j > 0) y = y * x + c[--j];
1164 Double_t AliGenMUONlib::YJpsiCDFscaledPP3( const Double_t *px, const Double_t *dummy)
1167 return AliGenMUONlib::YJpsiPP2760(px, dummy);
1170 Double_t AliGenMUONlib::YJpsiCDFscaledPP2( const Double_t *px, const Double_t */*dummy*/)
1175 return YJpsiPPdummy(*px, 1960);
1178 Double_t AliGenMUONlib::YJpsiPP( const Double_t *px, const Double_t */*dummy*/)
1188 // mc = 1.4 GeV, pt-kick 1 GeV
1191 Double_t c[5] = {1.38532e+00, 1.00596e+02, -3.46378e+01, 3.94172e+00, -1.48319e-01};
1192 Double_t x = TMath::Abs(px[0]);
1196 y = 96.455 - 0.8483 * x * x;
1197 } else if (x < 7.9) {
1200 while (j > 0) y = y * x + c[--j];
1208 Double_t AliGenMUONlib::YJpsiCDFscaledPPb9( const Double_t *px, const Double_t */*dummy*/)
1212 // pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.79
1214 Double_t c[7] = {7.52296e-01, 2.49917e-02, 3.36500e-03, 1.91187e-03, 2.92154e-04,
1215 -4.16509e-05,-7.62709e-06};
1217 Double_t x = px[0] + 0.47; // rapidity shift
1220 while (j > 0) y = y * x + c[--j];
1223 return y * AliGenMUONlib::YJpsiCDFscaledPP9dummy(x);
1226 Double_t AliGenMUONlib::YJpsiCDFscaledPbP9( const Double_t *px, const Double_t */*dummy*/)
1230 // Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.79
1232 Double_t c[7] = {7.52296e-01, 2.49917e-02, 3.36500e-03, 1.91187e-03, 2.92154e-04,
1233 -4.16509e-05,-7.62709e-06};
1235 Double_t x = -px[0] + 0.47; // rapidity shift
1238 while (j > 0) y = y * x + c[--j];
1241 return y * AliGenMUONlib::YJpsiCDFscaledPP9dummy(x);
1244 Double_t AliGenMUONlib::YJpsiCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
1248 // PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.66
1250 Double_t c[4] = {5.95228e-01, 9.45069e-03, 2.44710e-04, -1.32894e-05};
1251 Double_t x = px[0]*px[0];
1255 while (j > 0) y = y * x + c[--j];
1258 return y * AliGenMUONlib::YJpsiCDFscaledPP4(px,dummy);
1261 Double_t AliGenMUONlib::YJpsiBPbPb( const Double_t *px, const Double_t */*dummy*/)
1265 // J/Psi from B->J/Psi X
1270 Double_t c[7] = {7.37025e-02, 0., -2.94487e-03, 0., 6.07953e-06, 0., 5.39219e-07};
1272 Double_t x = TMath::Abs(px[0]);
1280 while (j > 0) y = y * x + c[--j];
1288 // particle composition
1290 Int_t AliGenMUONlib::IpJpsi(TRandom *)
1292 // J/Psi composition
1295 Int_t AliGenMUONlib::IpPsiP(TRandom *)
1297 // Psi prime composition
1300 Int_t AliGenMUONlib::IpJpsiFamily(TRandom *)
1302 // J/Psi composition
1304 Float_t r = gRandom->Rndm();
1319 //____________________________________________________________
1320 Double_t AliGenMUONlib::PtUpsilonPPdummy(Double_t x, Double_t energy)
1324 // from the fit of CDF & LHC data, like for J/Psi in arXiv:1103.2394
1326 const Double_t kpt0 = 1.96*TMath::Power(energy,0.095);
1327 const Double_t kxn = 3.4;
1329 Double_t pass1 = 1.+0.471*(x/kpt0)*(x/kpt0);
1330 return x/TMath::Power(pass1,kxn);
1333 Double_t AliGenMUONlib::PtUpsilonPP7000(const Double_t *px, const Double_t */*dummy*/)
1338 return PtUpsilonPPdummy(*px,7000);
1341 Double_t AliGenMUONlib::PtUpsilonPP2760(const Double_t *px, const Double_t */*dummy*/)
1346 return PtUpsilonPPdummy(*px,2760);
1349 Double_t AliGenMUONlib::PtUpsilonPP8800(const Double_t *px, const Double_t */*dummy*/)
1354 return PtUpsilonPPdummy(*px,8800);
1357 Double_t AliGenMUONlib::PtUpsilonPbPb2760ShFdummy(Double_t x, Int_t n)
1359 // Usilon shadowing factor vs pT for PbPb min. bias and 11 centr. bins (in 2.5<y<4)
1361 // PbPb 2.76 TeV, for EKS98, minimum bias shadowing factor = 0.87 in 4pi
1363 const Double_t f1[12] = {1, 1.128, 1.097, 1.037, 0.937, 0.821, 0.693, 0.558,
1364 0.428, 0.317, 0.231, 0.156};
1365 const Double_t f2[12] = {1, 1.313, 1.202, 1.039, 0.814, 0.593, 0.391, 0.224,
1366 0.106, 0.041, 0.013, 0.002};
1367 const Double_t c1[7] = {1.9089e+00, 1.2969e-03, 8.9786e-05,-5.3062e-06,
1368 -1.0046e-06,6.1446e-08, 1.0885e-09};
1369 const Double_t c2[7] = {8.8423e-01,-8.7488e-05, 5.9857e-04,-5.7959e-05,
1370 2.0059e-06,-2.7343e-08, 6.6053e-10};
1373 y1 = c1[j = 6]; y2 = c2[6];
1374 while (j > 0) {y1 = y1 * x + c1[--j]; y2 = y2 * x + c2[j];}
1376 y1 /= 1.+c1[6]*TMath::Power(x,6);
1377 y2 /= 1.+c2[6]*TMath::Power(x,6);
1379 y1 = 1 + (y1-2)*f1[n] + (y2+1-y1)*f2[n];
1384 Double_t AliGenMUONlib::PtUpsilonPbPb2760(const Double_t *px, const Double_t *dummy)
1387 // PbPb 2.76 TeV, minimum bias 0-100 %
1389 return PtUpsilonPbPb2760ShFdummy(*px, 0) * PtUpsilonPP2760(px, dummy);
1392 Double_t AliGenMUONlib::PtUpsilonPbPb2760c1(const Double_t *px, const Double_t *dummy)
1395 // PbPb 2.76 TeV, 1st centrality bin 0-5 %
1397 return PtUpsilonPbPb2760ShFdummy(*px, 1) * PtUpsilonPP2760(px, dummy);
1400 Double_t AliGenMUONlib::PtUpsilonPbPb2760c2(const Double_t *px, const Double_t *dummy)
1403 // PbPb 2.76 TeV, 2nd centrality bin 5-10 %
1405 return PtUpsilonPbPb2760ShFdummy(*px, 2) * PtUpsilonPP2760(px, dummy);
1408 Double_t AliGenMUONlib::PtUpsilonPbPb2760c3(const Double_t *px, const Double_t *dummy)
1411 // PbPb 2.76 TeV, 3rd centrality bin 10-20 %
1413 return PtUpsilonPbPb2760ShFdummy(*px, 3) * PtUpsilonPP2760(px, dummy);
1416 Double_t AliGenMUONlib::PtUpsilonPbPb2760c4(const Double_t *px, const Double_t *dummy)
1419 // PbPb 2.76 TeV, 4th centrality bin 20-30 %
1421 return PtUpsilonPbPb2760ShFdummy(*px, 4) * PtUpsilonPP2760(px, dummy);
1424 Double_t AliGenMUONlib::PtUpsilonPbPb2760c5(const Double_t *px, const Double_t *dummy)
1427 // PbPb 2.76 TeV, 5th centrality bin 30-40 %
1429 return PtUpsilonPbPb2760ShFdummy(*px, 5) * PtUpsilonPP2760(px, dummy);
1432 Double_t AliGenMUONlib::PtUpsilonPbPb2760c6(const Double_t *px, const Double_t *dummy)
1435 // PbPb 2.76 TeV, 6th centrality bin 40-50 %
1437 return PtUpsilonPbPb2760ShFdummy(*px, 6) * PtUpsilonPP2760(px, dummy);
1440 Double_t AliGenMUONlib::PtUpsilonPbPb2760c7(const Double_t *px, const Double_t *dummy)
1443 // PbPb 2.76 TeV, 7th centrality bin 50-60 %
1445 return PtUpsilonPbPb2760ShFdummy(*px, 7) * PtUpsilonPP2760(px, dummy);
1448 Double_t AliGenMUONlib::PtUpsilonPbPb2760c8(const Double_t *px, const Double_t *dummy)
1451 // PbPb 2.76 TeV, 8th centrality bin 60-70 %
1453 return PtUpsilonPbPb2760ShFdummy(*px, 8) * PtUpsilonPP2760(px, dummy);
1456 Double_t AliGenMUONlib::PtUpsilonPbPb2760c9(const Double_t *px, const Double_t *dummy)
1459 // PbPb 2.76 TeV, 9th centrality bin 70-80 %
1461 return PtUpsilonPbPb2760ShFdummy(*px, 9) * PtUpsilonPP2760(px, dummy);
1464 Double_t AliGenMUONlib::PtUpsilonPbPb2760c10(const Double_t *px, const Double_t *dummy)
1467 // PbPb 2.76 TeV, 10th centrality bin 80-90 %
1469 return PtUpsilonPbPb2760ShFdummy(*px, 10) * PtUpsilonPP2760(px, dummy);
1472 Double_t AliGenMUONlib::PtUpsilonPbPb2760c11(const Double_t *px, const Double_t *dummy)
1475 // PbPb 2.76 TeV, 11th centrality bin 90-100 %
1477 return PtUpsilonPbPb2760ShFdummy(*px, 11) * PtUpsilonPP2760(px, dummy);
1480 Double_t AliGenMUONlib::PtUpsilonPPb8800ShFdummy(Double_t x, Int_t n)
1482 // Upsilon shadowing factor vs pT for pPb min. bias and 4 centr. bins (in 2.5<y<4)
1484 // pPb 8.8 TeV, for EKS98, minimum bias shadowing factor = 0.89 in 4pi
1486 const Double_t f[5] = {1, 1.33, 1.05, 0.67, 0.23};
1487 const Double_t c[5] = {7.6561e-01, 1.1360e-04, 4.9596e-04,-3.0287e-05, 3.7555e-06};
1491 while (j > 0) y = y * x + c[--j];
1492 y /= 1 + c[4]*TMath::Power(x,4);
1494 return 1 + (y-1)*f[n];
1497 Double_t AliGenMUONlib::PtUpsilonPPb8800(const Double_t *px, const Double_t *dummy)
1500 // pPb 8.8 TeV, minimum bias 0-100 %
1502 return PtUpsilonPPb8800ShFdummy(*px, 0) * PtUpsilonPP8800(px, dummy);
1505 Double_t AliGenMUONlib::PtUpsilonPPb8800c1(const Double_t *px, const Double_t *dummy)
1508 // pPb 8.8 TeV, 1st centrality bin 0-20 %
1510 return PtUpsilonPPb8800ShFdummy(*px, 1) * PtUpsilonPP8800(px, dummy);
1513 Double_t AliGenMUONlib::PtUpsilonPPb8800c2(const Double_t *px, const Double_t *dummy)
1516 // pPb 8.8 TeV, 2nd centrality bin 20-40 %
1518 return PtUpsilonPPb8800ShFdummy(*px, 2) * PtUpsilonPP8800(px, dummy);
1521 Double_t AliGenMUONlib::PtUpsilonPPb8800c3(const Double_t *px, const Double_t *dummy)
1524 // pPb 8.8 TeV, 3rd centrality bin 40-60 %
1526 return PtUpsilonPPb8800ShFdummy(*px, 3) * PtUpsilonPP8800(px, dummy);
1529 Double_t AliGenMUONlib::PtUpsilonPPb8800c4(const Double_t *px, const Double_t *dummy)
1532 // pPb 8.8 TeV, 4th centrality bin 60-100 %
1534 return PtUpsilonPPb8800ShFdummy(*px, 4) * PtUpsilonPP8800(px, dummy);
1537 Double_t AliGenMUONlib::PtUpsilonPbP8800ShFdummy(Double_t x, Int_t n)
1539 // Upsilon shadowing factor vs pT for Pbp min. bias and 4 centr. bins (in 2.5<y<4)
1541 // Pbp 8.8 TeV, for EKS98, minimum bias shadowing factor = 0.89 in 4pi
1543 const Double_t f[5] = {1, 1.33, 1.05, 0.67, 0.23};
1544 const Double_t c[5] = {1.0975, 3.1905e-03,-2.0477e-04, 8.5270e-06, 2.5343e-06};
1548 while (j > 0) y = y * x + c[--j];
1549 y /= 1 + c[4]*TMath::Power(x,4);
1551 return 1 + (y-1)*f[n];
1554 Double_t AliGenMUONlib::PtUpsilonPbP8800(const Double_t *px, const Double_t *dummy)
1557 // Pbp 8.8 TeV, minimum bias 0-100 %
1559 return PtUpsilonPbP8800ShFdummy(*px, 0) * PtUpsilonPP8800(px, dummy);
1562 Double_t AliGenMUONlib::PtUpsilonPbP8800c1(const Double_t *px, const Double_t *dummy)
1565 // Pbp 8.8 TeV, 1st centrality bin 0-20 %
1567 return PtUpsilonPbP8800ShFdummy(*px, 1) * PtUpsilonPP8800(px, dummy);
1570 Double_t AliGenMUONlib::PtUpsilonPbP8800c2(const Double_t *px, const Double_t *dummy)
1573 // Pbp 8.8 TeV, 2nd centrality bin 20-40 %
1575 return PtUpsilonPbP8800ShFdummy(*px, 2) * PtUpsilonPP8800(px, dummy);
1578 Double_t AliGenMUONlib::PtUpsilonPbP8800c3(const Double_t *px, const Double_t *dummy)
1581 // Pbp 8.8 TeV, 3rd centrality bin 40-60 %
1583 return PtUpsilonPbP8800ShFdummy(*px, 3) * PtUpsilonPP8800(px, dummy);
1586 Double_t AliGenMUONlib::PtUpsilonPbP8800c4(const Double_t *px, const Double_t *dummy)
1589 // Pbp 8.8 TeV, 4th centrality bin 60-100 %
1591 return PtUpsilonPbP8800ShFdummy(*px, 4) * PtUpsilonPP8800(px, dummy);
1594 Double_t AliGenMUONlib::PtUpsilon( const Double_t *px, const Double_t */*dummy*/ )
1597 const Double_t kpt0 = 5.3;
1598 const Double_t kxn = 2.5;
1601 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1602 return x/TMath::Power(pass1,kxn);
1605 Double_t AliGenMUONlib::PtUpsilonCDFscaled( const Double_t *px, const Double_t */*dummy*/ )
1608 const Double_t kpt0 = 7.753;
1609 const Double_t kxn = 3.042;
1612 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1613 return x/TMath::Power(pass1,kxn);
1616 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP( const Double_t *px, const Double_t */*dummy*/ )
1622 // scaled from CDF data at 2 TeV
1624 const Double_t kpt0 = 8.610;
1625 const Double_t kxn = 3.051;
1628 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1629 return x/TMath::Power(pass1,kxn);
1632 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
1638 // scaled from CDF data at 2 TeV
1640 const Double_t kpt0 = 8.235;
1641 const Double_t kxn = 3.051;
1644 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1645 return x/TMath::Power(pass1,kxn);
1648 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
1653 // scaled from CDF data at 2 TeV
1655 const Double_t kpt0 = 8.048;
1656 const Double_t kxn = 3.051;
1659 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1660 return x/TMath::Power(pass1,kxn);
1663 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
1669 // scaled from CDF data at 2 TeV
1671 const Double_t kpt0 = 7.817;
1672 const Double_t kxn = 3.051;
1675 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1676 return x/TMath::Power(pass1,kxn);
1679 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
1684 // scaled from CDF data at 2 TeV
1686 const Double_t kpt0 = 7.189;
1687 const Double_t kxn = 3.051;
1690 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1691 return x/TMath::Power(pass1,kxn);
1694 Double_t AliGenMUONlib::PtUpsilonCDFscaledPPb9( const Double_t *px, const Double_t *dummy)
1698 // pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.89
1700 Double_t c[5] = {7.64952e-01, 1.12501e-04, 4.96038e-04, -3.03198e-05, 3.74035e-06};
1705 while (j > 0) y = y * x + c[--j];
1707 Double_t d = 1.+c[4]*TMath::Power(x,4);
1708 return y/d * AliGenMUONlib::PtUpsilonCDFscaledPP9(px,dummy);
1711 Double_t AliGenMUONlib::PtUpsilonCDFscaledPbP9( const Double_t *px, const Double_t *dummy)
1715 // Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.89
1717 Double_t c[5] = {1.09881e+00, 3.08329e-03, -2.00356e-04, 8.28991e-06, 2.52576e-06};
1722 while (j > 0) y = y * x + c[--j];
1724 Double_t d = 1.+c[4]*TMath::Power(x,4);
1725 return y/d * AliGenMUONlib::PtUpsilonCDFscaledPP9(px,dummy);
1728 Double_t AliGenMUONlib::PtUpsilonCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
1732 // PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.85
1734 Double_t c[5] = {8.65872e-01, 2.05465e-03, 2.56063e-04, -1.65598e-05, 2.29209e-06};
1739 while (j > 0) y = y * x + c[--j];
1741 Double_t d = 1.+c[4]*TMath::Power(x,4);
1742 return y/d * AliGenMUONlib::PtUpsilonCDFscaledPP4(px,dummy);
1745 Double_t AliGenMUONlib::PtUpsilonFlat( const Double_t */*px*/, const Double_t */*dummy*/ )
1750 Double_t AliGenMUONlib::PtUpsilonPbPb( const Double_t *px, const Double_t */*dummy*/)
1759 // mc = 1.4 GeV, pt-kick 1 GeV
1763 -1.03488e+01, 1.28065e+01, -6.60500e+00, 1.66140e+00,
1764 -2.34293e-01, 1.86925e-02, -7.80708e-04, 1.30610e-05
1770 while (j > 0) y = y * x +c[--j];
1771 y = x * TMath::Exp(y);
1778 Double_t AliGenMUONlib::PtUpsilonPP( const Double_t *px, const Double_t */*dummy*/)
1787 // mc = 1.4 GeV, pt-kick 1 GeV
1790 Double_t c[8] = {-7.93955e+00, 1.06306e+01, -5.21392e+00, 1.19703e+00,
1791 -1.45718e-01, 8.95151e-03, -2.04806e-04, -1.13053e-06};
1797 while (j > 0) y = y * x +c[--j];
1798 y = x * TMath::Exp(y);
1808 //____________________________________________________________
1809 Double_t AliGenMUONlib::YUpsilonPPdummy(Double_t x, Double_t energy)
1813 // from the fit of CDF & LHC data, like for J/Psi in arXiv:1103.2394
1815 x = x/TMath::Log(energy/9.46);
1817 Double_t y = TMath::Exp(-x/0.4/0.4/2);
1822 Double_t AliGenMUONlib::YUpsilonPPpoly(Double_t x, Double_t energy)
1826 // from the fit of CDF & LHC data, like for J/Psi in arXiv:1103.2394
1828 x = x/TMath::Log(energy/9.46);
1830 Double_t y = 1 - 6.9*x*x;
1835 Double_t AliGenMUONlib::YUpsilonPP7000(const Double_t *px, const Double_t */*dummy*/)
1840 return YUpsilonPPdummy(*px, 7000);
1843 Double_t AliGenMUONlib::YUpsilonPP2760(const Double_t *px, const Double_t */*dummy*/)
1848 return YUpsilonPPdummy(*px, 2760);
1851 Double_t AliGenMUONlib::YUpsilonPP8800(const Double_t *px, const Double_t */*dummy*/)
1856 return YUpsilonPPdummy(*px, 8800);
1859 Double_t AliGenMUONlib::YUpsilonPPpoly7000(const Double_t *px, const Double_t */*dummy*/)
1864 return YUpsilonPPpoly(*px, 7000);
1867 Double_t AliGenMUONlib::YUpsilonPPpoly2760(const Double_t *px, const Double_t */*dummy*/)
1872 return YUpsilonPPpoly(*px, 2760);
1875 Double_t AliGenMUONlib::YUpsilonPbPb2760ShFdummy(Double_t x, Int_t n)
1877 // Upsilon shadowing factor vs y for PbPb min. bias and 11 centr. bins
1879 // PbPb 2.76 TeV, for EKS98, minimum bias shadowing factor = 0.87 in 4pi
1881 const Double_t f1[12] = {1, 1.128, 1.097, 1.037, 0.937, 0.821, 0.693, 0.558,
1882 0.428, 0.317, 0.231, 0.156};
1883 const Double_t f2[12] = {1, 1.313, 1.202, 1.039, 0.814, 0.593, 0.391, 0.224,
1884 0.106, 0.041, 0.013, 0.002};
1885 const Double_t c1[5] = {1.8547e+00, 1.6717e-02,-2.1285e-04,-9.7662e-05, 2.5768e-06};
1886 const Double_t c2[5] = {8.6029e-01, 1.1742e-02,-2.7944e-04,-6.7973e-05, 1.8838e-06};
1891 y1 = c1[j = 4]; y2 = c2[4];
1892 while (j > 0) {y1 = y1 * x + c1[--j]; y2 = y2 * x + c2[j];}
1894 y1 = 1 + (y1-2)*f1[n] + (y2+1-y1)*f2[n];
1899 Double_t AliGenMUONlib::YUpsilonPbPb2760(const Double_t *px, const Double_t *dummy)
1902 // PbPb 2.76 TeV, minimum bias 0-100 %
1904 return YUpsilonPbPb2760ShFdummy(*px, 0) * YUpsilonPP2760(px, dummy);
1907 Double_t AliGenMUONlib::YUpsilonPbPb2760c1(const Double_t *px, const Double_t *dummy)
1910 // PbPb 2.76 TeV, 1st centrality bin 0-5 %
1912 return YUpsilonPbPb2760ShFdummy(*px, 1) * YUpsilonPP2760(px, dummy);
1915 Double_t AliGenMUONlib::YUpsilonPbPb2760c2(const Double_t *px, const Double_t *dummy)
1918 // PbPb 2.76 TeV, 2nd centrality bin 5-10 %
1920 return YUpsilonPbPb2760ShFdummy(*px, 2) * YUpsilonPP2760(px, dummy);
1923 Double_t AliGenMUONlib::YUpsilonPbPb2760c3(const Double_t *px, const Double_t *dummy)
1926 // PbPb 2.76 TeV, 3rd centrality bin 10-20 %
1928 return YUpsilonPbPb2760ShFdummy(*px, 3) * YUpsilonPP2760(px, dummy);
1931 Double_t AliGenMUONlib::YUpsilonPbPb2760c4(const Double_t *px, const Double_t *dummy)
1934 // PbPb 2.76 TeV, 4th centrality bin 20-30 %
1936 return YUpsilonPbPb2760ShFdummy(*px, 4) * YUpsilonPP2760(px, dummy);
1939 Double_t AliGenMUONlib::YUpsilonPbPb2760c5(const Double_t *px, const Double_t *dummy)
1942 // PbPb 2.76 TeV, 5th centrality bin 30-40 %
1944 return YUpsilonPbPb2760ShFdummy(*px, 5) * YUpsilonPP2760(px, dummy);
1947 Double_t AliGenMUONlib::YUpsilonPbPb2760c6(const Double_t *px, const Double_t *dummy)
1950 // PbPb 2.76 TeV, 6th centrality bin 40-50 %
1952 return YUpsilonPbPb2760ShFdummy(*px, 6) * YUpsilonPP2760(px, dummy);
1955 Double_t AliGenMUONlib::YUpsilonPbPb2760c7(const Double_t *px, const Double_t *dummy)
1958 // PbPb 2.76 TeV, 7th centrality bin 50-60 %
1960 return YUpsilonPbPb2760ShFdummy(*px, 7) * YUpsilonPP2760(px, dummy);
1963 Double_t AliGenMUONlib::YUpsilonPbPb2760c8(const Double_t *px, const Double_t *dummy)
1966 // PbPb 2.76 TeV, 8th centrality bin 60-70 %
1968 return YUpsilonPbPb2760ShFdummy(*px, 8) * YUpsilonPP2760(px, dummy);
1971 Double_t AliGenMUONlib::YUpsilonPbPb2760c9(const Double_t *px, const Double_t *dummy)
1974 // PbPb 2.76 TeV, 9th centrality bin 70-80 %
1976 return YUpsilonPbPb2760ShFdummy(*px, 9) * YUpsilonPP2760(px, dummy);
1979 Double_t AliGenMUONlib::YUpsilonPbPb2760c10(const Double_t *px, const Double_t *dummy)
1982 // PbPb 2.76 TeV, 10th centrality bin 80-90 %
1984 return YUpsilonPbPb2760ShFdummy(*px, 10) * YUpsilonPP2760(px, dummy);
1987 Double_t AliGenMUONlib::YUpsilonPbPb2760c11(const Double_t *px, const Double_t *dummy)
1990 // PbPb 2.76 TeV, 11th centrality bin 90-100 %
1992 return YUpsilonPbPb2760ShFdummy(*px, 11) * YUpsilonPP2760(px, dummy);
1995 Double_t AliGenMUONlib::YUpsilonPP8800dummy(Double_t px)
1997 return AliGenMUONlib::YUpsilonPP8800(&px, (Double_t*) 0);
2000 Double_t AliGenMUONlib::YUpsilonPPb8800ShFdummy(Double_t x, Int_t n)
2002 // Upsilon shadowing factor vs y for pPb min. bias and 4 centr. bins
2004 // pPb 8.8 TeV, for EKS98, minimum bias shadowing factor = 0.89 in 4pi
2006 const Double_t f[5] = {1, 1.33, 1.05, 0.67, 0.23};
2007 const Double_t c[7] = {8.6581e-01, 4.6111e-02, 7.6911e-03, 8.7313e-04,
2008 -1.4700e-04,-5.0975e-05,-3.5718e-06};
2012 while (j > 0) y = y * x + c[--j];
2015 return 1 +(y-1)*f[n];
2018 Double_t AliGenMUONlib::YUpsilonPPb8800(const Double_t *px, const Double_t */*dummy*/)
2021 // pPb 8.8 TeV, minimum bias 0-100 %
2023 Double_t x = px[0] + 0.47; // rapidity shift
2024 return YUpsilonPPb8800ShFdummy(x, 0) * YUpsilonPP8800dummy(x);
2027 Double_t AliGenMUONlib::YUpsilonPPb8800c1(const Double_t *px, const Double_t */*dummy*/)
2030 // pPb 8.8 TeV, 1st centrality bin 0-20 %
2032 Double_t x = px[0] + 0.47; // rapidity shift
2033 return YUpsilonPPb8800ShFdummy(x, 1) * YUpsilonPP8800dummy(x);
2036 Double_t AliGenMUONlib::YUpsilonPPb8800c2(const Double_t *px, const Double_t */*dummy*/)
2039 // pPb 8.8 TeV, 2nd centrality bin 20-40 %
2041 Double_t x = px[0] + 0.47; // rapidity shift
2042 return YUpsilonPPb8800ShFdummy(x, 2) * YUpsilonPP8800dummy(x);
2045 Double_t AliGenMUONlib::YUpsilonPPb8800c3(const Double_t *px, const Double_t */*dummy*/)
2048 // pPb 8.8 TeV, 3rd centrality bin 40-60 %
2050 Double_t x = px[0] + 0.47; // rapidity shift
2051 return YUpsilonPPb8800ShFdummy(x, 3) * YUpsilonPP8800dummy(x);
2054 Double_t AliGenMUONlib::YUpsilonPPb8800c4(const Double_t *px, const Double_t */*dummy*/)
2057 // pPb 8.8 TeV, 4th centrality bin 60-100 %
2059 Double_t x = px[0] + 0.47; // rapidity shift
2060 return YUpsilonPPb8800ShFdummy(x, 4) * YUpsilonPP8800dummy(x);
2063 Double_t AliGenMUONlib::YUpsilonPbP8800(const Double_t *px, const Double_t */*dummy*/)
2066 // Pbp 8.8 TeV, minimum bias 0-100 %
2068 Double_t x = -px[0] + 0.47; // rapidity shift
2069 return YUpsilonPPb8800ShFdummy(x, 0) * YUpsilonPP8800dummy(x);
2072 Double_t AliGenMUONlib::YUpsilonPbP8800c1(const Double_t *px, const Double_t */*dummy*/)
2075 // Pbp 8.8 TeV, 1st centrality bin 0-20 %
2077 Double_t x = -px[0] + 0.47; // rapidity shift
2078 return YUpsilonPPb8800ShFdummy(x, 1) * YUpsilonPP8800dummy(x);
2081 Double_t AliGenMUONlib::YUpsilonPbP8800c2(const Double_t *px, const Double_t */*dummy*/)
2084 // Pbp 8.8 TeV, 2nd centrality bin 20-40 %
2086 Double_t x = -px[0] + 0.47; // rapidity shift
2087 return YUpsilonPPb8800ShFdummy(x, 2) * YUpsilonPP8800dummy(x);
2090 Double_t AliGenMUONlib::YUpsilonPbP8800c3(const Double_t *px, const Double_t */*dummy*/)
2093 // Pbp 8.8 TeV, 3rd centrality bin 40-60 %
2095 Double_t x = -px[0] + 0.47; // rapidity shift
2096 return YUpsilonPPb8800ShFdummy(x, 3) * YUpsilonPP8800dummy(x);
2099 Double_t AliGenMUONlib::YUpsilonPbP8800c4(const Double_t *px, const Double_t */*dummy*/)
2102 // Pbp 8.8 TeV, 4th centrality bin 60-100 %
2104 Double_t x = -px[0] + 0.47; // rapidity shift
2105 return YUpsilonPPb8800ShFdummy(x, 4) * YUpsilonPP8800dummy(x);
2108 Double_t AliGenMUONlib::YUpsilon(const Double_t *py, const Double_t */*dummy*/)
2111 const Double_t ky0 = 3.;
2112 const Double_t kb=1.;
2114 Double_t y=TMath::Abs(*py);
2119 yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2124 Double_t AliGenMUONlib::YUpsilonPbPb( const Double_t *px, const Double_t */*dummy*/)
2134 // mc = 1.4 GeV, pt-kick 1 GeV
2137 Double_t c[7] = {3.40036e-01, -3.98882e-07, -4.48398e-03, 8.46411e-08, -6.10854e-04,
2138 -2.99753e-09, 1.28895e-05};
2139 Double_t x = TMath::Abs(px[0]);
2140 if (x > 5.55) return 0.;
2142 Double_t y = c[j = 6];
2143 while (j > 0) y = y * x +c[--j];
2147 Double_t AliGenMUONlib::YUpsilonCDFscaled( const Double_t *px, const Double_t *dummy)
2150 return AliGenMUONlib::YUpsilonPbPb(px, dummy);
2154 Double_t AliGenMUONlib::YUpsilonCDFscaledPP( const Double_t *px, const Double_t *dummy)
2157 return AliGenMUONlib::YUpsilonPP(px, dummy);
2161 Double_t AliGenMUONlib::YUpsilonFlat( const Double_t */*px*/, const Double_t */*dummy*/)
2168 Double_t AliGenMUONlib::YUpsilonCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
2173 // scaled from YUpsilonPP(14 TeV) using 10 TeV / 14 TeV ratio of y-spectra in LO pQCD.
2174 // see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
2176 Double_t c[4] = {1., -2.17877e-02, -6.52830e-04, 1.40578e-05};
2177 Double_t x = TMath::Abs(px[0]);
2178 if (x > 6.1) return 0.;
2180 Double_t y = c[j = 3];
2181 while (j > 0) y = y * x*x +c[--j];
2185 Double_t AliGenMUONlib::YUpsilonCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
2190 // rescaling of YUpsilonPP(14 TeV) using 8.8 TeV / 14 TeV ratio of y-spectra in LO QCD
2192 Double_t c[4] = {1., -2.37621e-02, -6.29610e-04, 1.47976e-05};
2193 Double_t x = TMath::Abs(px[0]);
2194 if (x > 6.1) return 0.;
2196 Double_t y = c[j = 3];
2197 while (j > 0) y = y * x*x +c[--j];
2201 Double_t AliGenMUONlib::YUpsilonCDFscaledPP9dummy(Double_t px)
2203 return AliGenMUONlib::YUpsilonCDFscaledPP9(&px, (Double_t*) 0);
2206 Double_t AliGenMUONlib::YUpsilonCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
2211 // scaled from YUpsilonPP(14 TeV) using 7 TeV / 14 TeV ratio of y-spectra in LO pQCD.
2213 Double_t c[4] = {1., -2.61009e-02, -6.83937e-04, 1.78451e-05};
2214 Double_t x = TMath::Abs(px[0]);
2215 if (x > 6.0) return 0.;
2217 Double_t y = c[j = 3];
2218 while (j > 0) y = y * x*x +c[--j];
2222 Double_t AliGenMUONlib::YUpsilonCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
2227 // rescaling of YUpsilonPP(14 TeV) using 3.94 TeV / 14 TeV ratio of y-spectra in LO QCD
2229 Double_t c[4] = {1., -3.91924e-02, -4.26184e-04, 2.10914e-05};
2230 Double_t x = TMath::Abs(px[0]);
2231 if (x > 5.7) return 0.;
2233 Double_t y = c[j = 3];
2234 while (j > 0) y = y * x*x +c[--j];
2239 Double_t AliGenMUONlib::YUpsilonPP( const Double_t *px, const Double_t */*dummy*/)
2249 // mc = 1.4 GeV, pt-kick 1 GeV
2251 Double_t c[7] = {8.91936e-01, -6.46645e-07, -1.52774e-02, 4.28677e-08, -7.01517e-04,
2252 -6.20539e-10, 1.29943e-05};
2253 Double_t x = TMath::Abs(px[0]);
2254 if (x > 6.2) return 0.;
2256 Double_t y = c[j = 6];
2257 while (j > 0) y = y * x +c[--j];
2261 Double_t AliGenMUONlib::YUpsilonCDFscaledPPb9( const Double_t *px, const Double_t */*dummy*/)
2265 // pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.89
2267 Double_t c[7] = {8.71829e-01, 4.77467e-02, 8.09671e-03, 6.45294e-04, -2.15730e-04,
2268 -4.67538e-05,-2.11683e-06};
2270 Double_t x = px[0] + 0.47; // rapidity shift
2273 while (j > 0) y = y * x + c[--j];
2276 return y * AliGenMUONlib::YUpsilonCDFscaledPP9dummy(x);
2279 Double_t AliGenMUONlib::YUpsilonCDFscaledPbP9( const Double_t *px, const Double_t */*dummy*/)
2283 // Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.89
2285 Double_t c[7] = {8.71829e-01, 4.77467e-02, 8.09671e-03, 6.45294e-04, -2.15730e-04,
2286 -4.67538e-05,-2.11683e-06};
2288 Double_t x = -px[0] + 0.47; // rapidity shift
2291 while (j > 0) y = y * x + c[--j];
2294 return y * AliGenMUONlib::YUpsilonCDFscaledPP9dummy(x);
2297 Double_t AliGenMUONlib::YUpsilonCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
2301 // PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.85
2303 Double_t c[4] = {8.27837e-01, 1.70115e-02, -1.26046e-03, 1.52091e-05};
2304 Double_t x = px[0]*px[0];
2308 while (j > 0) y = y * x + c[--j];
2311 return y * AliGenMUONlib::YUpsilonCDFscaledPP4(px,dummy);
2315 // particle composition
2317 Int_t AliGenMUONlib::IpUpsilon(TRandom *)
2322 Int_t AliGenMUONlib::IpUpsilonP(TRandom *)
2327 Int_t AliGenMUONlib::IpUpsilonPP(TRandom *)
2332 Int_t AliGenMUONlib::IpUpsilonFamily(TRandom *)
2336 Float_t r = gRandom->Rndm();
2340 } else if (r < 0.896) {
2353 // pt-distribution (by scaling of pion distribution)
2354 //____________________________________________________________
2355 Double_t AliGenMUONlib::PtPhi( const Double_t *px, const Double_t */*dummy*/)
2358 return PtScal(*px,7);
2361 Double_t AliGenMUONlib::YPhi( const Double_t *px, const Double_t */*dummy*/)
2365 return YJpsi(px,dum);
2367 // particle composition
2369 Int_t AliGenMUONlib::IpPhi(TRandom *)
2379 // pt-distribution (by scaling of pion distribution)
2380 //____________________________________________________________
2381 Double_t AliGenMUONlib::PtOmega( const Double_t *px, const Double_t */*dummy*/)
2384 return PtScal(*px,5);
2387 Double_t AliGenMUONlib::YOmega( const Double_t *px, const Double_t */*dummy*/)
2391 return YJpsi(px,dum);
2393 // particle composition
2395 Int_t AliGenMUONlib::IpOmega(TRandom *)
2397 // Omega composition
2406 // pt-distribution (by scaling of pion distribution)
2407 //____________________________________________________________
2408 Double_t AliGenMUONlib::PtEta( const Double_t *px, const Double_t */*dummy*/)
2411 return PtScal(*px,3);
2414 Double_t AliGenMUONlib::YEta( const Double_t *px, const Double_t */*dummy*/)
2418 return YJpsi(px,dum);
2420 // particle composition
2422 Int_t AliGenMUONlib::IpEta(TRandom *)
2433 //____________________________________________________________
2434 Double_t AliGenMUONlib::PtCharm( const Double_t *px, const Double_t */*dummy*/)
2437 const Double_t kpt0 = 2.25;
2438 const Double_t kxn = 3.17;
2441 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2442 return x/TMath::Power(pass1,kxn);
2445 Double_t AliGenMUONlib::PtCharmCentral( const Double_t *px, const Double_t */*dummy*/)
2448 const Double_t kpt0 = 2.12;
2449 const Double_t kxn = 2.78;
2452 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2453 return x/TMath::Power(pass1,kxn);
2455 Double_t AliGenMUONlib::PtCharmF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2457 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
2458 // PtCharmFiMjSkPP = PtCharmF0M0S0PP * (dN(i,j,k)/dpt / dN(0,0,0)/dpt)_MNR
2459 // i=0,1,2; j=0,1,2; k=0,1,...,6
2460 // dN(i,j,k)/dpt - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
2461 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
2462 // calculations for the following inputs:
2463 // Peterson fragmentation function (F) with \epsilon_c = 0.02, 0.002 & 0.11
2464 // for i=0,1 & 2 respectively; quark mass (M) of 1.5, 1.3 & 1.7 GeV
2465 // for j=0,1 & 2 respectively;
2466 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
2467 // with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2
2468 // for k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
2469 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
2470 // June 2008, Smbat.Grigoryan@cern.ch
2473 // Pythia6.214 (kCharmppMNRwmi, PDF = CTEQ5L, quark mass = 1.2 GeV, PtHard > 2.76 GeV/c)
2474 // for pp collisions at 14 TeV with one c-cbar pair per event.
2475 // Corresponding NLO total cross section is 5.68 mb
2478 const Double_t kpt0 = 2.2930;
2479 const Double_t kxn = 3.1196;
2480 Double_t c[3]={-5.2180e-01,1.8753e-01,2.8669e-02};
2483 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2484 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2486 Double_t AliGenMUONlib::PtCharmF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2489 // Corresponding NLO total cross section is 6.06 mb
2490 const Double_t kpt0 = 2.8669;
2491 const Double_t kxn = 3.1044;
2492 Double_t c[3]={-4.6714e-01,1.5005e-01,4.5003e-02};
2495 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2496 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2498 Double_t AliGenMUONlib::PtCharmF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2501 // Corresponding NLO total cross section is 6.06 mb
2502 const Double_t kpt0 = 1.8361;
2503 const Double_t kxn = 3.2966;
2504 Double_t c[3]={-6.1550e-01,2.6498e-01,1.0728e-02};
2507 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2508 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2510 Double_t AliGenMUONlib::PtCharmF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
2513 // Corresponding NLO total cross section is 7.69 mb
2514 const Double_t kpt0 = 2.1280;
2515 const Double_t kxn = 3.1397;
2516 Double_t c[3]={-5.4021e-01,2.0944e-01,2.5211e-02};
2519 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2520 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2522 Double_t AliGenMUONlib::PtCharmF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
2525 // Corresponding NLO total cross section is 4.81 mb
2526 const Double_t kpt0 = 2.4579;
2527 const Double_t kxn = 3.1095;
2528 Double_t c[3]={-5.1497e-01,1.7532e-01,3.2429e-02};
2531 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2532 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2534 Double_t AliGenMUONlib::PtCharmF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
2537 // Corresponding NLO total cross section is 14.09 mb
2538 const Double_t kpt0 = 2.1272;
2539 const Double_t kxn = 3.1904;
2540 Double_t c[3]={-4.6088e-01,2.1918e-01,2.3055e-02};
2543 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2544 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2546 Double_t AliGenMUONlib::PtCharmF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
2549 // Corresponding NLO total cross section is 1.52 mb
2550 const Double_t kpt0 = 2.8159;
2551 const Double_t kxn = 3.0857;
2552 Double_t c[3]={-6.4691e-01,2.0289e-01,2.4922e-02};
2555 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2556 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2558 Double_t AliGenMUONlib::PtCharmF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
2561 // Corresponding NLO total cross section is 3.67 mb
2562 const Double_t kpt0 = 2.7297;
2563 const Double_t kxn = 3.3019;
2564 Double_t c[3]={-6.2216e-01,1.9031e-01,1.5341e-02};
2567 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2568 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2570 Double_t AliGenMUONlib::PtCharmF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
2573 // Corresponding NLO total cross section is 3.38 mb
2574 const Double_t kpt0 = 2.3894;
2575 const Double_t kxn = 3.1075;
2576 Double_t c[3]={-4.9742e-01,1.7032e-01,2.5994e-02};
2579 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2580 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2582 Double_t AliGenMUONlib::PtCharmF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
2585 // Corresponding NLO total cross section is 10.37 mb
2586 const Double_t kpt0 = 2.0187;
2587 const Double_t kxn = 3.3011;
2588 Double_t c[3]={-3.9869e-01,2.9248e-01,1.1763e-02};
2591 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2592 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2594 Double_t AliGenMUONlib::PtCharmF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
2597 // Corresponding NLO total cross section is 7.22 mb
2598 const Double_t kpt0 = 2.1089;
2599 const Double_t kxn = 3.1848;
2600 Double_t c[3]={-4.6275e-01,1.8114e-01,2.1363e-02};
2603 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2604 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
2608 Double_t AliGenMUONlib::YCharm( const Double_t *px, const Double_t */*dummy*/)
2610 // Charm y :: Carrer & Dainese : ALICE-INT-2003-019 v.3 (hep-ph/0311225)
2611 // Pythia tuned to reproduce the distribution given by the HVQMNR program based on NLO calculations (pQCD)
2612 // shadowing + kt broadening
2615 Double_t c[2]={-2.42985e-03,-2.31001e-04};
2616 Double_t y=1+(c[0]*TMath::Power(x,2))+(c[1]*TMath::Power(x,4));
2619 if (TMath::Abs(x)>8) {
2623 ycharm=TMath::Power(y,3);
2628 Double_t AliGenMUONlib::YCharmF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2630 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
2631 // YCharmFiMjSkPP = YCharmF0M0S0PP * (dN(i,j,k)/dy / dN(0,0,0)/dy)_MNR
2632 // i=0,1,2; j=0,1,2; k=0,1,...,6
2633 // dN(i,j,k)/dy - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
2634 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
2635 // calculations for the following inputs:
2636 // Peterson fragmentation function (F) with \epsilon_c = 0.02, 0.002 & 0.11
2637 // for i=0,1 & 2 respectively; quark mass (M) of 1.5, 1.3 & 1.7 GeV
2638 // for j=0,1 & 2 respectively;
2639 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
2640 // with a/b = 1/1,1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2 for
2641 // k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
2642 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
2643 // June 2008, Smbat.Grigoryan@cern.ch
2646 // Pythia6.214 (kCharmppMNRwmi, PDF = CTEQ5L, quark mass = 1.2 GeV, PtHard > 2.76 GeV/c)
2647 // for pp collisions at 14 TeV with one c-cbar pair per event.
2648 // Corresponding NLO total cross section is 5.68 mb
2651 Double_t c[2]={7.0909e-03,6.1967e-05};
2652 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2655 if (TMath::Abs(x)>9) {
2659 ycharm=TMath::Power(y,3);
2664 Double_t AliGenMUONlib::YCharmF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2667 // Corresponding NLO total cross section is 6.06 mb
2669 Double_t c[2]={6.9707e-03,6.0971e-05};
2670 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2673 if (TMath::Abs(x)>9) {
2677 ycharm=TMath::Power(y,3);
2682 Double_t AliGenMUONlib::YCharmF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2685 // Corresponding NLO total cross section is 6.06 mb
2687 Double_t c[2]={7.1687e-03,6.5303e-05};
2688 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2691 if (TMath::Abs(x)>9) {
2695 ycharm=TMath::Power(y,3);
2700 Double_t AliGenMUONlib::YCharmF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
2703 // Corresponding NLO total cross section is 7.69 mb
2705 Double_t c[2]={5.9090e-03,7.1854e-05};
2706 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2709 if (TMath::Abs(x)>9) {
2713 ycharm=TMath::Power(y,3);
2718 Double_t AliGenMUONlib::YCharmF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
2721 // Corresponding NLO total cross section is 4.81 mb
2723 Double_t c[2]={8.0882e-03,5.5872e-05};
2724 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2727 if (TMath::Abs(x)>9) {
2731 ycharm=TMath::Power(y,3);
2736 Double_t AliGenMUONlib::YCharmF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
2739 // Corresponding NLO total cross section is 14.09 mb
2741 Double_t c[2]={7.2520e-03,6.2691e-05};
2742 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2745 if (TMath::Abs(x)>9) {
2749 ycharm=TMath::Power(y,3);
2754 Double_t AliGenMUONlib::YCharmF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
2757 // Corresponding NLO total cross section is 1.52 mb
2759 Double_t c[2]={1.1040e-04,1.4498e-04};
2760 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2763 if (TMath::Abs(x)>9) {
2767 ycharm=TMath::Power(y,3);
2772 Double_t AliGenMUONlib::YCharmF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
2775 // Corresponding NLO total cross section is 3.67 mb
2777 Double_t c[2]={-3.1328e-03,1.8270e-04};
2778 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2781 if (TMath::Abs(x)>9) {
2785 ycharm=TMath::Power(y,3);
2790 Double_t AliGenMUONlib::YCharmF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
2793 // Corresponding NLO total cross section is 3.38 mb
2795 Double_t c[2]={7.0865e-03,6.2532e-05};
2796 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2799 if (TMath::Abs(x)>9) {
2803 ycharm=TMath::Power(y,3);
2808 Double_t AliGenMUONlib::YCharmF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
2811 // Corresponding NLO total cross section is 10.37 mb
2813 Double_t c[2]={7.7070e-03,5.3533e-05};
2814 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2817 if (TMath::Abs(x)>9) {
2821 ycharm=TMath::Power(y,3);
2826 Double_t AliGenMUONlib::YCharmF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
2829 // Corresponding NLO total cross section is 7.22 mb
2831 Double_t c[2]={7.9195e-03,5.3823e-05};
2832 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
2835 if (TMath::Abs(x)>9) {
2839 ycharm=TMath::Power(y,3);
2846 Int_t AliGenMUONlib::IpCharm(TRandom *ran)
2848 // Charm composition
2852 random = ran->Rndm();
2853 // Taux de production Carrer & Dainese : ALICE-INT-2003-019 v.3
2854 // >>>>> cf. tab 4 p 11
2856 if (random < 0.30) {
2858 } else if (random < 0.60) {
2860 } else if (random < 0.70) {
2862 } else if (random < 0.80) {
2864 } else if (random < 0.86) {
2866 } else if (random < 0.92) {
2868 } else if (random < 0.96) {
2882 //____________________________________________________________
2883 Double_t AliGenMUONlib::PtBeauty( const Double_t *px, const Double_t */*dummy*/)
2886 const Double_t kpt0 = 6.53;
2887 const Double_t kxn = 3.59;
2890 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2891 return x/TMath::Power(pass1,kxn);
2894 Double_t AliGenMUONlib::PtBeautyCentral( const Double_t *px, const Double_t */*dummy*/)
2897 const Double_t kpt0 = 6.14;
2898 const Double_t kxn = 2.93;
2901 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2902 return x/TMath::Power(pass1,kxn);
2904 Double_t AliGenMUONlib::PtBeautyF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2906 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
2907 // PtBeautyFiMjSkPP = PtBeautyF0M0S0PP * (dN(i,j,k)/dpt / dN(0,0,0)/dpt)_MNR
2908 // i=0,1,2; j=0,1,2; k=0,1,...,6
2909 // dN(i,j,k)/dpt - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
2910 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
2911 // calculations for the following inputs:
2912 // Peterson fragmentation function (F) with \epsilon_b = 0.001, 0.0002 & 0.004
2913 // for i=0,1 & 2 respectively; quark mass (M) of 4.75, 4.5 & 5.0 GeV
2914 // for j=0,1 & 2 respectively;
2915 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
2916 // with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2 for
2917 // k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
2918 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
2919 // June 2008, Smbat.Grigoryan@cern.ch
2922 // Pythia6.214 (kBeautyppMNRwmi, PDF = CTEQ5L, quark mass = 4.75 GeV, PtHard > 2.76 GeV/c)
2923 // for pp collisions at 14 TeV with one b-bbar pair per event.
2924 // Corresponding NLO total cross section is 0.494 mb
2926 const Double_t kpt0 = 8.0575;
2927 const Double_t kxn = 3.1921;
2930 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2931 return x/TMath::Power(pass1,kxn);
2933 Double_t AliGenMUONlib::PtBeautyF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2936 // Corresponding NLO total cross section is 0.445 mb
2937 const Double_t kpt0 = 8.6239;
2938 const Double_t kxn = 3.2911;
2941 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2942 return x/TMath::Power(pass1,kxn);
2944 Double_t AliGenMUONlib::PtBeautyF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
2947 // Corresponding NLO total cross section is 0.445 mb
2948 const Double_t kpt0 = 7.3367;
2949 const Double_t kxn = 3.0692;
2952 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2953 return x/TMath::Power(pass1,kxn);
2955 Double_t AliGenMUONlib::PtBeautyF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
2958 // Corresponding NLO total cross section is 0.518 mb
2959 const Double_t kpt0 = 7.6409;
2960 const Double_t kxn = 3.1364;
2963 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2964 return x/TMath::Power(pass1,kxn);
2966 Double_t AliGenMUONlib::PtBeautyF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
2969 // Corresponding NLO total cross section is 0.384 mb
2970 const Double_t kpt0 = 8.4948;
2971 const Double_t kxn = 3.2546;
2974 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2975 return x/TMath::Power(pass1,kxn);
2977 Double_t AliGenMUONlib::PtBeautyF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
2980 // Corresponding NLO total cross section is 0.648 mb
2981 const Double_t kpt0 = 7.6631;
2982 const Double_t kxn = 3.1621;
2985 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2986 return x/TMath::Power(pass1,kxn);
2988 Double_t AliGenMUONlib::PtBeautyF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
2991 // Corresponding NLO total cross section is 0.294 mb
2992 const Double_t kpt0 = 8.7245;
2993 const Double_t kxn = 3.2213;
2996 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2997 return x/TMath::Power(pass1,kxn);
2999 Double_t AliGenMUONlib::PtBeautyF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
3002 // Corresponding NLO total cross section is 0.475 mb
3003 const Double_t kpt0 = 8.5296;
3004 const Double_t kxn = 3.2187;
3007 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3008 return x/TMath::Power(pass1,kxn);
3010 Double_t AliGenMUONlib::PtBeautyF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
3013 // Corresponding NLO total cross section is 0.324 mb
3014 const Double_t kpt0 = 7.9440;
3015 const Double_t kxn = 3.1614;
3018 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3019 return x/TMath::Power(pass1,kxn);
3021 Double_t AliGenMUONlib::PtBeautyF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
3024 // Corresponding NLO total cross section is 0.536 mb
3025 const Double_t kpt0 = 8.2408;
3026 const Double_t kxn = 3.3029;
3029 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3030 return x/TMath::Power(pass1,kxn);
3032 Double_t AliGenMUONlib::PtBeautyF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
3035 // Corresponding NLO total cross section is 0.420 mb
3036 const Double_t kpt0 = 7.8041;
3037 const Double_t kxn = 3.2094;
3040 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3041 return x/TMath::Power(pass1,kxn);
3045 Double_t AliGenMUONlib::YBeauty( const Double_t *px, const Double_t */*dummy*/)
3047 // Beauty y :: Carrer & Dainese : ALICE-INT-2003-019 v.3 (hep-ph/0311225)
3048 // Pythia tuned to reproduce the distribution given by the HVQMNR program based on NLO calculations (pQCD)
3049 // shadowing + kt broadening
3052 Double_t c[2]={-1.27590e-02,-2.42731e-04};
3053 Double_t y=1+c[0]*TMath::Power(x,2)+c[1]*TMath::Power(x,4);
3056 if (TMath::Abs(x)>6) {
3060 ybeauty=TMath::Power(y,3);
3065 Double_t AliGenMUONlib::YBeautyF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
3067 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
3068 // YBeautyFiMjSkPP = YBeautyF0M0S0PP * (dN(i,j,k)/dy / dN(0,0,0)/dy)_MNR
3069 // i=0,1,2; j=0,1,2; k=0,1,...,6
3070 // dN(i,j,k)/dy - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
3071 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
3072 // calculations for the following inputs:
3073 // Peterson fragmentation function (F) with \epsilon_b = 0.001, 0.0002 & 0.004
3074 // for i=0,1 & 2 respectively; quark mass (M) of 4.75, 4.5 & 5.0 GeV
3075 // for j=0,1 & 2 respectively;
3076 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
3077 // with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2
3078 // for k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
3079 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
3080 // June 2008, Smbat.Grigoryan@cern.ch
3083 // Pythia6.214 (kBeautyppMNRwmi, PDF = CTEQ5L, quark mass = 4.75 GeV, PtHard > 2.76 GeV/c)
3084 // for pp collisions at 14 TeV with one b-bbar pair per event.
3085 // Corresponding NLO total cross section is 0.494 mb
3089 Double_t c[2]={1.2350e-02,9.2667e-05};
3090 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3093 if (TMath::Abs(x)>7.6) {
3097 ybeauty=TMath::Power(y,3);
3102 Double_t AliGenMUONlib::YBeautyF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
3105 // Corresponding NLO total cross section is 0.445 mb
3107 Double_t c[2]={1.2292e-02,9.1847e-05};
3108 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3111 if (TMath::Abs(x)>7.6) {
3115 ybeauty=TMath::Power(y,3);
3120 Double_t AliGenMUONlib::YBeautyF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
3123 // Corresponding NLO total cross section is 0.445 mb
3125 Double_t c[2]={1.2436e-02,9.3709e-05};
3126 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3129 if (TMath::Abs(x)>7.6) {
3133 ybeauty=TMath::Power(y,3);
3138 Double_t AliGenMUONlib::YBeautyF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
3141 // Corresponding NLO total cross section is 0.518 mb
3143 Double_t c[2]={1.1714e-02,1.0068e-04};
3144 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3147 if (TMath::Abs(x)>7.6) {
3151 ybeauty=TMath::Power(y,3);
3156 Double_t AliGenMUONlib::YBeautyF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
3159 // Corresponding NLO total cross section is 0.384 mb
3161 Double_t c[2]={1.2944e-02,8.5500e-05};
3162 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3165 if (TMath::Abs(x)>7.6) {
3169 ybeauty=TMath::Power(y,3);
3174 Double_t AliGenMUONlib::YBeautyF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
3177 // Corresponding NLO total cross section is 0.648 mb
3179 Double_t c[2]={1.2455e-02,9.2713e-05};
3180 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3183 if (TMath::Abs(x)>7.6) {
3187 ybeauty=TMath::Power(y,3);
3192 Double_t AliGenMUONlib::YBeautyF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
3195 // Corresponding NLO total cross section is 0.294 mb
3197 Double_t c[2]={1.0897e-02,1.1878e-04};
3198 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3201 if (TMath::Abs(x)>7.6) {
3205 ybeauty=TMath::Power(y,3);
3210 Double_t AliGenMUONlib::YBeautyF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
3213 // Corresponding NLO total cross section is 0.475 mb
3215 Double_t c[2]={1.0912e-02,1.1858e-04};
3216 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3219 if (TMath::Abs(x)>7.6) {
3223 ybeauty=TMath::Power(y,3);
3228 Double_t AliGenMUONlib::YBeautyF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
3231 // Corresponding NLO total cross section is 0.324 mb
3233 Double_t c[2]={1.2378e-02,9.2490e-05};
3234 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3237 if (TMath::Abs(x)>7.6) {
3241 ybeauty=TMath::Power(y,3);
3246 Double_t AliGenMUONlib::YBeautyF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
3249 // Corresponding NLO total cross section is 0.536 mb
3251 Double_t c[2]={1.2886e-02,8.2912e-05};
3252 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3255 if (TMath::Abs(x)>7.6) {
3259 ybeauty=TMath::Power(y,3);
3264 Double_t AliGenMUONlib::YBeautyF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
3267 // Corresponding NLO total cross section is 0.420 mb
3269 Double_t c[2]={1.3106e-02,8.0115e-05};
3270 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
3273 if (TMath::Abs(x)>7.6) {
3277 ybeauty=TMath::Power(y,3);
3283 Int_t AliGenMUONlib::IpBeauty(TRandom *ran)
3285 // Beauty Composition
3288 random = ran->Rndm();
3290 // Taux de production Carrer & Dainese : ALICE-INT-2003-019 v.3
3291 // >>>>> cf. tab 4 p 11
3293 if (random < 0.20) {
3295 } else if (random < 0.40) {
3297 } else if (random < 0.605) {
3299 } else if (random < 0.81) {
3301 } else if (random < 0.87) {
3303 } else if (random < 0.93) {
3305 } else if (random < 0.965) {
3315 typedef Double_t (*GenFunc) (const Double_t*, const Double_t*);
3316 GenFunc AliGenMUONlib::GetPt(Int_t param, const char* tname) const
3318 // Return pointer to pT parameterisation
3319 TString sname = TString(tname);
3337 if (sname == "Vogt" || sname == "Vogt PbPb") {
3339 } else if (sname == "Vogt pp") {
3341 } else if (sname == "pp 7") {
3343 } else if (sname == "pp 2.76") {
3345 } else if (sname == "PbPb 2.76") {
3346 func=PtJpsiPbPb2760;
3347 } else if (sname == "PbPb 2.76c1") {
3348 func=PtJpsiPbPb2760c1;
3349 } else if (sname == "PbPb 2.76c2") {
3350 func=PtJpsiPbPb2760c2;
3351 } else if (sname == "PbPb 2.76c3") {
3352 func=PtJpsiPbPb2760c3;
3353 } else if (sname == "PbPb 2.76c4") {
3354 func=PtJpsiPbPb2760c4;
3355 } else if (sname == "PbPb 2.76c5") {
3356 func=PtJpsiPbPb2760c5;
3357 } else if (sname == "PbPb 2.76c6") {
3358 func=PtJpsiPbPb2760c6;
3359 } else if (sname == "PbPb 2.76c7") {
3360 func=PtJpsiPbPb2760c7;
3361 } else if (sname == "PbPb 2.76c8") {
3362 func=PtJpsiPbPb2760c8;
3363 } else if (sname == "PbPb 2.76c9") {
3364 func=PtJpsiPbPb2760c9;
3365 } else if (sname == "PbPb 2.76c10") {
3366 func=PtJpsiPbPb2760c10;
3367 } else if (sname == "PbPb 2.76c11") {
3368 func=PtJpsiPbPb2760c11;
3369 } else if (sname == "pp 7 poly") {
3371 } else if (sname == "pp 2.76 poly") {
3373 } else if (sname == "pp 8.8") {
3375 } else if (sname == "pPb 8.8") {
3377 } else if (sname == "pPb 8.8c1") {
3378 func=PtJpsiPPb8800c1;
3379 } else if (sname == "pPb 8.8c2") {
3380 func=PtJpsiPPb8800c2;
3381 } else if (sname == "pPb 8.8c3") {
3382 func=PtJpsiPPb8800c3;
3383 } else if (sname == "pPb 8.8c4") {
3384 func=PtJpsiPPb8800c4;
3385 } else if (sname == "Pbp 8.8") {
3387 } else if (sname == "Pbp 8.8c1") {
3388 func=PtJpsiPbP8800c1;
3389 } else if (sname == "Pbp 8.8c2") {
3390 func=PtJpsiPbP8800c2;
3391 } else if (sname == "Pbp 8.8c3") {
3392 func=PtJpsiPbP8800c3;
3393 } else if (sname == "Pbp 8.8c4") {
3394 func=PtJpsiPbP8800c4;
3395 } else if (sname == "CDF scaled") {
3396 func=PtJpsiCDFscaled;
3397 } else if (sname == "CDF pp") {
3398 func=PtJpsiCDFscaledPP;
3399 } else if (sname == "CDF pp 10") {
3400 func=PtJpsiCDFscaledPP10;
3401 } else if (sname == "CDF pp 8.8") {
3402 func=PtJpsiCDFscaledPP9;
3403 } else if (sname == "CDF pp 7" || sname == "CDF pp 7 flat y") {
3404 func=PtJpsiCDFscaledPP7;
3405 } else if (sname == "CDF pp 3.94") {
3406 func=PtJpsiCDFscaledPP4;
3407 } else if (sname == "CDF pp 2.76") {
3408 func=PtJpsiCDFscaledPP3;
3409 } else if (sname == "CDF pp 1.9") {
3410 func=PtJpsiCDFscaledPP2;
3411 } else if (sname == "CDF pPb 8.8") {
3412 func=PtJpsiCDFscaledPPb9;
3413 } else if (sname == "CDF Pbp 8.8") {
3414 func=PtJpsiCDFscaledPbP9;
3415 } else if (sname == "CDF PbPb 3.94") {
3416 func=PtJpsiCDFscaledPbPb4;
3417 } else if (sname == "Flat" || sname == "CDF pp 7 flat pt") {
3426 case kUpsilonFamily:
3430 if (sname == "Vogt" || sname == "Vogt PbPb") {
3432 } else if (sname == "Vogt pp") {
3434 } else if (sname == "pp 7") {
3435 func=PtUpsilonPP7000;
3436 } else if (sname == "pp 2.76") {
3437 func=PtUpsilonPP2760;
3438 } else if (sname == "PbPb 2.76") {
3439 func=PtUpsilonPbPb2760;
3440 } else if (sname == "PbPb 2.76c1") {
3441 func=PtUpsilonPbPb2760c1;
3442 } else if (sname == "PbPb 2.76c2") {
3443 func=PtUpsilonPbPb2760c2;
3444 } else if (sname == "PbPb 2.76c3") {
3445 func=PtUpsilonPbPb2760c3;
3446 } else if (sname == "PbPb 2.76c4") {
3447 func=PtUpsilonPbPb2760c4;
3448 } else if (sname == "PbPb 2.76c5") {
3449 func=PtUpsilonPbPb2760c5;
3450 } else if (sname == "PbPb 2.76c6") {
3451 func=PtUpsilonPbPb2760c6;
3452 } else if (sname == "PbPb 2.76c7") {
3453 func=PtUpsilonPbPb2760c7;
3454 } else if (sname == "PbPb 2.76c8") {
3455 func=PtUpsilonPbPb2760c8;
3456 } else if (sname == "PbPb 2.76c9") {
3457 func=PtUpsilonPbPb2760c9;
3458 } else if (sname == "PbPb 2.76c10") {
3459 func=PtUpsilonPbPb2760c10;
3460 } else if (sname == "PbPb 2.76c11") {
3461 func=PtUpsilonPbPb2760c11;
3462 } else if (sname == "pp 7 poly") {
3463 func=PtUpsilonPP7000;
3464 } else if (sname == "pp 2.76 poly") {
3465 func=PtUpsilonPP2760;
3466 } else if (sname == "pp 8.8") {
3467 func=PtUpsilonPP8800;
3468 } else if (sname == "pPb 8.8") {
3469 func=PtUpsilonPPb8800;
3470 } else if (sname == "pPb 8.8c1") {
3471 func=PtUpsilonPPb8800c1;
3472 } else if (sname == "pPb 8.8c2") {
3473 func=PtUpsilonPPb8800c2;
3474 } else if (sname == "pPb 8.8c3") {
3475 func=PtUpsilonPPb8800c3;
3476 } else if (sname == "pPb 8.8c4") {
3477 func=PtUpsilonPPb8800c4;
3478 } else if (sname == "Pbp 8.8") {
3479 func=PtUpsilonPbP8800;
3480 } else if (sname == "Pbp 8.8c1") {
3481 func=PtUpsilonPbP8800c1;
3482 } else if (sname == "Pbp 8.8c2") {
3483 func=PtUpsilonPbP8800c2;
3484 } else if (sname == "Pbp 8.8c3") {
3485 func=PtUpsilonPbP8800c3;
3486 } else if (sname == "Pbp 8.8c4") {
3487 func=PtUpsilonPbP8800c4;
3488 } else if (sname == "CDF scaled") {
3489 func=PtUpsilonCDFscaled;
3490 } else if (sname == "CDF pp") {
3491 func=PtUpsilonCDFscaledPP;
3492 } else if (sname == "CDF pp 10") {
3493 func=PtUpsilonCDFscaledPP10;
3494 } else if (sname == "CDF pp 8.8") {
3495 func=PtUpsilonCDFscaledPP9;
3496 } else if (sname == "CDF pp 7") {
3497 func=PtUpsilonCDFscaledPP7;
3498 } else if (sname == "CDF pp 3.94") {
3499 func=PtUpsilonCDFscaledPP4;
3500 } else if (sname == "CDF pPb 8.8") {
3501 func=PtUpsilonCDFscaledPPb9;
3502 } else if (sname == "CDF Pbp 8.8") {
3503 func=PtUpsilonCDFscaledPbP9;
3504 } else if (sname == "CDF PbPb 3.94") {
3505 func=PtUpsilonCDFscaledPbPb4;
3506 } else if (sname == "Flat") {
3513 if (sname == "F0M0S0 pp") {
3514 func=PtCharmF0M0S0PP;
3515 } else if (sname == "F1M0S0 pp") {
3516 func=PtCharmF1M0S0PP;
3517 } else if (sname == "F2M0S0 pp") {
3518 func=PtCharmF2M0S0PP;
3519 } else if (sname == "F0M1S0 pp") {
3520 func=PtCharmF0M1S0PP;
3521 } else if (sname == "F0M2S0 pp") {
3522 func=PtCharmF0M2S0PP;
3523 } else if (sname == "F0M0S1 pp") {
3524 func=PtCharmF0M0S1PP;
3525 } else if (sname == "F0M0S2 pp") {
3526 func=PtCharmF0M0S2PP;
3527 } else if (sname == "F0M0S3 pp") {
3528 func=PtCharmF0M0S3PP;
3529 } else if (sname == "F0M0S4 pp") {
3530 func=PtCharmF0M0S4PP;
3531 } else if (sname == "F0M0S5 pp") {
3532 func=PtCharmF0M0S5PP;
3533 } else if (sname == "F0M0S6 pp") {
3534 func=PtCharmF0M0S6PP;
3535 } else if (sname == "central") {
3536 func=PtCharmCentral;
3542 if (sname == "F0M0S0 pp") {
3543 func=PtBeautyF0M0S0PP;
3544 } else if (sname == "F1M0S0 pp") {
3545 func=PtBeautyF1M0S0PP;
3546 } else if (sname == "F2M0S0 pp") {
3547 func=PtBeautyF2M0S0PP;
3548 } else if (sname == "F0M1S0 pp") {
3549 func=PtBeautyF0M1S0PP;
3550 } else if (sname == "F0M2S0 pp") {
3551 func=PtBeautyF0M2S0PP;
3552 } else if (sname == "F0M0S1 pp") {
3553 func=PtBeautyF0M0S1PP;
3554 } else if (sname == "F0M0S2 pp") {
3555 func=PtBeautyF0M0S2PP;
3556 } else if (sname == "F0M0S3 pp") {
3557 func=PtBeautyF0M0S3PP;
3558 } else if (sname == "F0M0S4 pp") {
3559 func=PtBeautyF0M0S4PP;
3560 } else if (sname == "F0M0S5 pp") {
3561 func=PtBeautyF0M0S5PP;
3562 } else if (sname == "F0M0S6 pp") {
3563 func=PtBeautyF0M0S6PP;
3564 } else if (sname == "central") {
3565 func=PtBeautyCentral;
3571 if (sname == "2010 Pos PP") {
3572 func=PtPionPos2010PP;
3573 } else if (sname == "2010 Neg PP") {
3574 func=PtPionNeg2010PP;
3580 if (sname == "2010 Pos PP") {
3581 func=PtKaonPos2010PP;
3582 } else if (sname == "2010 Neg PP") {
3583 func=PtKaonNeg2010PP;
3596 printf("<AliGenMUONlib::GetPt> unknown parametrisation\n");
3601 GenFunc AliGenMUONlib::GetY(Int_t param, const char* tname) const
3604 // Return pointer to y- parameterisation
3606 TString sname = TString(tname);
3624 if (sname == "Vogt" || sname == "Vogt PbPb") {
3626 } else if (sname == "Vogt pp"){
3628 } else if (sname == "pp 7") {
3630 } else if (sname == "pp 2.76") {
3632 } else if (sname == "PbPb 2.76") {
3634 } else if (sname == "PbPb 2.76c1") {
3635 func=YJpsiPbPb2760c1;
3636 } else if (sname == "PbPb 2.76c2") {
3637 func=YJpsiPbPb2760c2;
3638 } else if (sname == "PbPb 2.76c3") {
3639 func=YJpsiPbPb2760c3;
3640 } else if (sname == "PbPb 2.76c4") {
3641 func=YJpsiPbPb2760c4;
3642 } else if (sname == "PbPb 2.76c5") {
3643 func=YJpsiPbPb2760c5;
3644 } else if (sname == "PbPb 2.76c6") {
3645 func=YJpsiPbPb2760c6;
3646 } else if (sname == "PbPb 2.76c7") {
3647 func=YJpsiPbPb2760c7;
3648 } else if (sname == "PbPb 2.76c8") {
3649 func=YJpsiPbPb2760c8;
3650 } else if (sname == "PbPb 2.76c9") {
3651 func=YJpsiPbPb2760c9;
3652 } else if (sname == "PbPb 2.76c10") {
3653 func=YJpsiPbPb2760c10;
3654 } else if (sname == "PbPb 2.76c11") {
3655 func=YJpsiPbPb2760c11;
3656 } else if (sname == "pp 7 poly") {
3657 func=YJpsiPPpoly7000;
3658 } else if (sname == "pp 2.76 poly") {
3659 func=YJpsiPPpoly2760;
3660 } else if (sname == "pp 8.8") {
3662 } else if (sname == "pPb 8.8") {
3664 } else if (sname == "pPb 8.8c1") {
3665 func=YJpsiPPb8800c1;
3666 } else if (sname == "pPb 8.8c2") {
3667 func=YJpsiPPb8800c2;
3668 } else if (sname == "pPb 8.8c3") {
3669 func=YJpsiPPb8800c3;
3670 } else if (sname == "pPb 8.8c4") {
3671 func=YJpsiPPb8800c4;
3672 } else if (sname == "Pbp 8.8") {
3674 } else if (sname == "Pbp 8.8c1") {
3675 func=YJpsiPbP8800c1;
3676 } else if (sname == "Pbp 8.8c2") {
3677 func=YJpsiPbP8800c2;
3678 } else if (sname == "Pbp 8.8c3") {
3679 func=YJpsiPbP8800c3;
3680 } else if (sname == "Pbp 8.8c4") {
3681 func=YJpsiPbP8800c4;
3682 } else if (sname == "CDF scaled") {
3683 func=YJpsiCDFscaled;
3684 } else if (sname == "CDF pp") {
3685 func=YJpsiCDFscaledPP;
3686 } else if (sname == "CDF pp 10") {
3687 func=YJpsiCDFscaledPP10;
3688 } else if (sname == "CDF pp 8.8") {
3689 func=YJpsiCDFscaledPP9;
3690 } else if (sname == "CDF pp 7" || sname == "CDF pp 7 flat pt") {
3691 func=YJpsiCDFscaledPP7;
3692 } else if (sname == "CDF pp 3.94") {
3693 func=YJpsiCDFscaledPP4;
3694 } else if (sname == "CDF pp 2.76") {
3695 func=YJpsiCDFscaledPP3;
3696 } else if (sname == "CDF pp 1.9") {
3697 func=YJpsiCDFscaledPP2;
3698 } else if (sname == "CDF pPb 8.8") {
3699 func=YJpsiCDFscaledPPb9;
3700 } else if (sname == "CDF Pbp 8.8") {
3701 func=YJpsiCDFscaledPbP9;
3702 } else if (sname == "CDF PbPb 3.94") {
3703 func=YJpsiCDFscaledPbPb4;
3704 } else if (sname == "Flat" || sname == "CDF pp 7 flat y") {
3713 case kUpsilonFamily:
3717 if (sname == "Vogt" || sname == "Vogt PbPb") {
3719 } else if (sname == "Vogt pp") {
3721 } else if (sname == "pp 7") {
3722 func=YUpsilonPP7000;
3723 } else if (sname == "pp 2.76") {
3724 func=YUpsilonPP2760;
3725 } else if (sname == "PbPb 2.76") {
3726 func=YUpsilonPbPb2760;
3727 } else if (sname == "PbPb 2.76c1") {
3728 func=YUpsilonPbPb2760c1;
3729 } else if (sname == "PbPb 2.76c2") {
3730 func=YUpsilonPbPb2760c2;
3731 } else if (sname == "PbPb 2.76c3") {
3732 func=YUpsilonPbPb2760c3;
3733 } else if (sname == "PbPb 2.76c4") {
3734 func=YUpsilonPbPb2760c4;
3735 } else if (sname == "PbPb 2.76c5") {
3736 func=YUpsilonPbPb2760c5;
3737 } else if (sname == "PbPb 2.76c6") {
3738 func=YUpsilonPbPb2760c6;
3739 } else if (sname == "PbPb 2.76c7") {
3740 func=YUpsilonPbPb2760c7;
3741 } else if (sname == "PbPb 2.76c8") {
3742 func=YUpsilonPbPb2760c8;
3743 } else if (sname == "PbPb 2.76c9") {
3744 func=YUpsilonPbPb2760c9;
3745 } else if (sname == "PbPb 2.76c10") {
3746 func=YUpsilonPbPb2760c10;
3747 } else if (sname == "PbPb 2.76c11") {
3748 func=YUpsilonPbPb2760c11;
3749 } else if (sname == "pp 7 poly") {
3750 func=YUpsilonPPpoly7000;
3751 } else if (sname == "pp 2.76 poly") {
3752 func=YUpsilonPPpoly2760;
3753 } else if (sname == "pp 8.8") {
3754 func=YUpsilonPP8800;
3755 } else if (sname == "pPb 8.8") {
3756 func=YUpsilonPPb8800;
3757 } else if (sname == "pPb 8.8c1") {
3758 func=YUpsilonPPb8800c1;
3759 } else if (sname == "pPb 8.8c2") {
3760 func=YUpsilonPPb8800c2;
3761 } else if (sname == "pPb 8.8c3") {
3762 func=YUpsilonPPb8800c3;
3763 } else if (sname == "pPb 8.8c4") {
3764 func=YUpsilonPPb8800c4;
3765 } else if (sname == "Pbp 8.8") {
3766 func=YUpsilonPbP8800;
3767 } else if (sname == "Pbp 8.8c1") {
3768 func=YUpsilonPbP8800c1;
3769 } else if (sname == "Pbp 8.8c2") {
3770 func=YUpsilonPbP8800c2;
3771 } else if (sname == "Pbp 8.8c3") {
3772 func=YUpsilonPbP8800c3;
3773 } else if (sname == "Pbp 8.8c4") {
3774 func=YUpsilonPbP8800c4;
3775 } else if (sname == "CDF scaled") {
3776 func=YUpsilonCDFscaled;
3777 } else if (sname == "CDF pp") {
3778 func=YUpsilonCDFscaledPP;
3779 } else if (sname == "CDF pp 10") {
3780 func=YUpsilonCDFscaledPP10;
3781 } else if (sname == "CDF pp 8.8") {
3782 func=YUpsilonCDFscaledPP9;
3783 } else if (sname == "CDF pp 7") {
3784 func=YUpsilonCDFscaledPP7;
3785 } else if (sname == "CDF pp 3.94") {
3786 func=YUpsilonCDFscaledPP4;
3787 } else if (sname == "CDF pPb 8.8") {
3788 func=YUpsilonCDFscaledPPb9;
3789 } else if (sname == "CDF Pbp 8.8") {
3790 func=YUpsilonCDFscaledPbP9;
3791 } else if (sname == "CDF PbPb 3.94") {
3792 func=YUpsilonCDFscaledPbPb4;
3793 } else if (sname == "Flat") {
3800 if (sname == "F0M0S0 pp") {
3801 func=YCharmF0M0S0PP;
3802 } else if (sname == "F1M0S0 pp") {
3803 func=YCharmF1M0S0PP;
3804 } else if (sname == "F2M0S0 pp") {
3805 func=YCharmF2M0S0PP;
3806 } else if (sname == "F0M1S0 pp") {
3807 func=YCharmF0M1S0PP;
3808 } else if (sname == "F0M2S0 pp") {
3809 func=YCharmF0M2S0PP;
3810 } else if (sname == "F0M0S1 pp") {
3811 func=YCharmF0M0S1PP;
3812 } else if (sname == "F0M0S2 pp") {
3813 func=YCharmF0M0S2PP;
3814 } else if (sname == "F0M0S3 pp") {
3815 func=YCharmF0M0S3PP;
3816 } else if (sname == "F0M0S4 pp") {
3817 func=YCharmF0M0S4PP;
3818 } else if (sname == "F0M0S5 pp") {
3819 func=YCharmF0M0S5PP;
3820 } else if (sname == "F0M0S6 pp") {
3821 func=YCharmF0M0S6PP;
3827 if (sname == "F0M0S0 pp") {
3828 func=YBeautyF0M0S0PP;
3829 } else if (sname == "F1M0S0 pp") {
3830 func=YBeautyF1M0S0PP;
3831 } else if (sname == "F2M0S0 pp") {
3832 func=YBeautyF2M0S0PP;
3833 } else if (sname == "F0M1S0 pp") {
3834 func=YBeautyF0M1S0PP;
3835 } else if (sname == "F0M2S0 pp") {
3836 func=YBeautyF0M2S0PP;
3837 } else if (sname == "F0M0S1 pp") {
3838 func=YBeautyF0M0S1PP;
3839 } else if (sname == "F0M0S2 pp") {
3840 func=YBeautyF0M0S2PP;
3841 } else if (sname == "F0M0S3 pp") {
3842 func=YBeautyF0M0S3PP;
3843 } else if (sname == "F0M0S4 pp") {
3844 func=YBeautyF0M0S4PP;
3845 } else if (sname == "F0M0S5 pp") {
3846 func=YBeautyF0M0S5PP;
3847 } else if (sname == "F0M0S6 pp") {
3848 func=YBeautyF0M0S6PP;
3854 if (sname == "2010 Pos PP") {
3855 func=YKaonPion2010PP;
3856 } else if (sname == "2010 Neg PP") {
3857 func=YKaonPion2010PP;
3863 if (sname == "2010 Pos PP") {
3864 func=YKaonPion2010PP;
3865 } else if (sname == "2010 Neg PP") {
3866 func=YKaonPion2010PP;
3879 printf("<AliGenMUONlib::GetY> unknown parametrisation\n");
3889 //____________________________________________________________
3890 Double_t AliGenMUONlib::PtChic0( const Double_t *px, const Double_t */*dummy*/)
3893 const Double_t kpt0 = 4.;
3894 const Double_t kxn = 3.6;
3897 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3898 return x/TMath::Power(pass1,kxn);
3900 Double_t AliGenMUONlib::PtChic1( const Double_t *px, const Double_t */*dummy*/)
3903 const Double_t kpt0 = 4.;
3904 const Double_t kxn = 3.6;
3907 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3908 return x/TMath::Power(pass1,kxn);
3910 Double_t AliGenMUONlib::PtChic2( const Double_t *px, const Double_t */*dummy*/)
3913 const Double_t kpt0 = 4.;
3914 const Double_t kxn = 3.6;
3917 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3918 return x/TMath::Power(pass1,kxn);
3920 Double_t AliGenMUONlib::PtChic( const Double_t *px, const Double_t */*dummy*/)
3923 const Double_t kpt0 = 4.;
3924 const Double_t kxn = 3.6;
3927 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
3928 return x/TMath::Power(pass1,kxn);
3933 //____________________________________________________________
3934 Double_t AliGenMUONlib::YChic0(const Double_t *py, const Double_t */*dummy*/)
3937 const Double_t ky0 = 4.;
3938 const Double_t kb=1.;
3940 Double_t y=TMath::Abs(*py);
3945 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
3949 Double_t AliGenMUONlib::YChic1(const Double_t *py, const Double_t */*dummy*/)
3952 const Double_t ky0 = 4.;
3953 const Double_t kb=1.;
3955 Double_t y=TMath::Abs(*py);
3960 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
3964 Double_t AliGenMUONlib::YChic2(const Double_t *py, const Double_t */*dummy*/)
3967 const Double_t ky0 = 4.;
3968 const Double_t kb=1.;
3970 Double_t y=TMath::Abs(*py);
3975 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
3979 Double_t AliGenMUONlib::YChic(const Double_t *py, const Double_t */*dummy*/)
3982 const Double_t ky0 = 4.;
3983 const Double_t kb=1.;
3985 Double_t y=TMath::Abs(*py);
3990 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
3994 // particle composition
3996 Int_t AliGenMUONlib::IpChic0(TRandom *)
4002 Int_t AliGenMUONlib::IpChic1(TRandom *)
4007 Int_t AliGenMUONlib::IpChic2(TRandom *)
4009 // Chi_c2 prime composition
4012 Int_t AliGenMUONlib::IpChic(TRandom *)
4016 Float_t r = gRandom->Rndm();
4019 } else if( r < 0.377 ) {
4028 //_____________________________________________________________
4030 typedef Int_t (*GenFuncIp) (TRandom *);
4031 GenFuncIp AliGenMUONlib::GetIp(Int_t param, const char* tname) const
4033 // Return pointer to particle type parameterisation
4034 TString sname = TString(tname);
4060 case kUpsilonFamily:
4061 func=IpUpsilonFamily;
4076 if (sname == "2010 Pos PP") {
4078 } else if (sname == "2010 Neg PP") {
4085 if (sname == "2010 Pos PP") {
4087 } else if (sname == "2010 Neg PP") {
4107 printf("<AliGenMUONlib::GetIp> unknown parametrisation\n");
4114 Float_t AliGenMUONlib::Interpolate(Float_t x, Float_t* y, Float_t x0,
4119 // Neville's alorithm for interpolation
4125 // n: number of data points
4126 // no: order of polynom
4128 Float_t* c = new Float_t[n];
4129 Float_t* d = new Float_t[n];
4131 for (i = 0; i < n; i++) {
4136 Int_t ns = int((x - x0)/dx);
4140 for (m = 0; m < no; m++) {
4141 for (i = 0; i < n-m; i++) {
4142 Float_t ho = x0 + Float_t(i) * dx - x;
4143 Float_t hp = x0 + Float_t(i+m+1) * dx - x;
4144 Float_t w = c[i+1] - d[i];
4145 Float_t den = ho-hp;
4152 if (2*ns < (n-m-1)) {
4164 //=============================================================================
4165 Double_t AliGenMUONlib::PtPionPos2010PP(const Double_t *px, const Double_t* /*dummy*/)
4168 const Double_t par[3] = {2.27501, 0.116141, 5.59591};
4169 Double_t pt = px[0];
4170 Double_t m0 = TDatabasePDG::Instance()->GetParticle(211)->Mass();
4171 Double_t mt = TMath::Sqrt(m0*m0 + pt*pt);
4172 Double_t nc = par[1]*par[2];
4173 Double_t t1 = (par[2]-1.)/nc/(nc/(par[2]-2.)+m0);
4174 Double_t t2 = TMath::Power(1.+(mt-m0)/nc, -1.*par[2]);
4175 Double_t fn = par[0] * pt * t1 * t2;
4179 //=============================================================================
4180 Double_t AliGenMUONlib::PtPionNeg2010PP(const Double_t *px, const Double_t* /*dummy*/)
4183 const Double_t par[3] = {2.25188, 0.12176, 5.91166};
4184 Double_t pt = px[0];
4185 Double_t m0 = TDatabasePDG::Instance()->GetParticle(211)->Mass();
4186 Double_t mt = TMath::Sqrt(m0*m0 + pt*pt);
4187 Double_t nc = par[1]*par[2];
4188 Double_t t1 = (par[2]-1.)/nc/(nc/(par[2]-2.)+m0);
4189 Double_t t2 = TMath::Power(1.+(mt-m0)/nc, -1.*par[2]);
4190 Double_t fn = par[0] * pt * t1 * t2;
4194 //=============================================================================
4195 Double_t AliGenMUONlib::PtKaonPos2010PP(const Double_t *px, const Double_t* /*dummy*/)
4198 const Double_t par[3] = {0.279386, 0.195466, 6.59587};
4199 Double_t pt = px[0];
4200 Double_t m0 = TDatabasePDG::Instance()->GetParticle(321)->Mass();
4201 Double_t mt = TMath::Sqrt(m0*m0 + pt*pt);
4202 Double_t nc = par[1]*par[2];
4203 Double_t t1 = (par[2]-1.)/nc/(nc/(par[2]-2.)+m0);
4204 Double_t t2 = TMath::Power(1.+(mt-m0)/nc, -1.*par[2]);
4205 Double_t fn = par[0] * pt * t1 * t2;
4209 //=============================================================================
4210 Double_t AliGenMUONlib::PtKaonNeg2010PP(const Double_t *px, const Double_t* /*dummy*/)
4213 const Double_t par[3] = {0.278927, 0.189049, 6.43006};
4214 Double_t pt = px[0];
4215 Double_t m0 = TDatabasePDG::Instance()->GetParticle(321)->Mass();
4216 Double_t mt = TMath::Sqrt(m0*m0 + pt*pt);
4217 Double_t nc = par[1]*par[2];
4218 Double_t t1 = (par[2]-1.)/nc/(nc/(par[2]-2.)+m0);
4219 Double_t t2 = TMath::Power(1.+(mt-m0)/nc, -1.*par[2]);
4220 Double_t fn = par[0] * pt * t1 * t2;
4224 //=============================================================================
4225 Double_t AliGenMUONlib::YKaonPion2010PP(const Double_t *px, const Double_t* /*dummy*/)
4229 Double_t sigma = 2.35;
4230 Double_t kernal = y/2./sigma;
4231 Double_t fxn = TMath::Exp(-1.*kernal*kernal);
4235 //=============================================================================
4236 Int_t AliGenMUONlib::IpPionPos(TRandom *)
4242 //=============================================================================
4243 Int_t AliGenMUONlib::IpPionNeg(TRandom *)
4249 //=============================================================================
4250 Int_t AliGenMUONlib::IpKaonPos(TRandom *)
4256 //=============================================================================
4257 Int_t AliGenMUONlib::IpKaonNeg(TRandom *)