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
30 #include "AliGenMUONlib.h"
32 ClassImp(AliGenMUONlib)
35 Double_t AliGenMUONlib::PtPion(const Double_t *px, const Double_t* /*dummy*/)
38 // PT-PARAMETERIZATION CDF, PRL 61(88) 1819
39 // POWER LAW FOR PT > 500 MEV
40 // MT SCALING BELOW (T=160 MEV)
42 const Double_t kp0 = 1.3;
43 const Double_t kxn = 8.28;
44 const Double_t kxlim=0.5;
45 const Double_t kt=0.160;
46 const Double_t kxmpi=0.139;
48 Double_t y, y1, xmpi2, ynorm, a;
51 y1=TMath::Power(kp0/(kp0+kxlim),kxn);
53 ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
56 y=a*TMath::Power(kp0/(kp0+x),kxn);
58 y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
64 Double_t AliGenMUONlib::YPion( const Double_t *py, const Double_t */*dummy*/)
67 Double_t y=TMath::Abs(*py);
69 const Double_t ka = 7000.;
70 const Double_t kdy = 4.;
71 Double_t ex = y*y/(2*kdy*kdy);
72 return ka*TMath::Exp(-ex);
74 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
77 // particle composition
79 Int_t AliGenMUONlib::IpPion(TRandom *ran)
82 if (ran->Rndm() < 0.5) {
89 //____________________________________________________________
93 Double_t AliGenMUONlib::PtScal(Double_t pt, Int_t np)
95 // SCALING EN MASSE PAR RAPPORT A PTPI
96 // MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
97 const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
98 // VALUE MESON/PI AT 5 GEV
99 const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
101 Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
102 Double_t fmax2=f5/kfmax[np];
104 Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
105 Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
106 (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
107 return fmtscal*ptpion;
113 //____________________________________________________________
114 Double_t AliGenMUONlib::PtKaon( const Double_t *px, const Double_t */*dummy*/)
117 return PtScal(*px,2);
121 //____________________________________________________________
122 Double_t AliGenMUONlib::YKaon( const Double_t *py, const Double_t */*dummy*/)
125 Double_t y=TMath::Abs(*py);
127 const Double_t ka = 1000.;
128 const Double_t kdy = 4.;
130 Double_t ex = y*y/(2*kdy*kdy);
131 return ka*TMath::Exp(-ex);
134 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
137 // particle composition
139 Int_t AliGenMUONlib::IpKaon(TRandom *ran)
142 if (ran->Rndm() < 0.5) {
153 //____________________________________________________________
154 Double_t AliGenMUONlib::PtJpsi( const Double_t *px, const Double_t */*dummy*/)
157 const Double_t kpt0 = 4.;
158 const Double_t kxn = 3.6;
161 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
162 return x/TMath::Power(pass1,kxn);
165 Double_t AliGenMUONlib::PtJpsiCDFscaled( const Double_t *px, const Double_t */*dummy*/)
170 // scaled from CDF data at 2 TeV
171 // see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
173 const Double_t kpt0 = 5.100;
174 const Double_t kxn = 4.102;
177 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
178 return x/TMath::Power(pass1,kxn);
181 Double_t AliGenMUONlib::PtJpsiCDFscaledPP( const Double_t *px, const Double_t */*dummy*/)
186 // scaled from CDF data at 2 TeV
188 const Double_t kpt0 = 5.630;
189 const Double_t kxn = 4.071;
192 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
193 return x/TMath::Power(pass1,kxn);
196 Double_t AliGenMUONlib::PtJpsiCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
201 // scaled from CDF data at 2 TeV
203 const Double_t kpt0 = 5.334;
204 const Double_t kxn = 4.071;
207 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
208 return x/TMath::Power(pass1,kxn);
211 Double_t AliGenMUONlib::PtJpsiFlat( const Double_t */*px*/, const Double_t */*dummy*/ )
216 Double_t AliGenMUONlib::PtJpsiPbPb( const Double_t *px, const Double_t */*dummy*/)
223 // mc = 1.4 GeV, pt-kick 1 GeV
227 -2.13098e+00, 9.46552e+00, -5.06799e+00, 1.27260e+00,
228 -1.83806e-01, 1.55853e-02, -7.23241e-04, 1.42105e-05
235 while (j > 0) y = y * x +c[--j];
236 y = x * TMath::Exp(y);
243 Double_t AliGenMUONlib::PtJpsiBPbPb( const Double_t *px, const Double_t */*dummy*/)
247 Double_t x0 = 4.0384;
251 Double_t y = x / TMath::Power((1. + (x/x0)*(x/x0)), n);
257 Double_t AliGenMUONlib::PtJpsiPP( const Double_t *px, const Double_t */*dummy*/)
264 // mc = 1.4 GeV, pt-kick 1 GeV
267 Float_t c[4] = {8.47471e+00, -1.93567e+00, 1.50271e-01, -5.51212e-03};
273 while (j > 0) y = y * x +c[--j];
274 y = x * TMath::Exp(y);
283 //____________________________________________________________
284 Double_t AliGenMUONlib::YJpsi(const Double_t *py, const Double_t */*dummy*/)
287 const Double_t ky0 = 4.;
288 const Double_t kb=1.;
290 Double_t y=TMath::Abs(*py);
295 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
299 Double_t AliGenMUONlib::YJpsiFlat( const Double_t */*py*/, const Double_t */*dummy*/ )
305 Double_t AliGenMUONlib::YJpsiPbPb( const Double_t *px, const Double_t */*dummy*/)
315 // mc = 1.4 GeV, pt-kick 1 GeV
317 Double_t c[5] = {-6.03425e+02, 4.98257e+02, -1.38794e+02, 1.62209e+01, -6.85955e-01};
318 Double_t x = TMath::Abs(px[0]);
326 while (j > 0) y = y * x + c[--j];
334 Double_t AliGenMUONlib::YJpsiCDFscaled( const Double_t *px, const Double_t* dummy)
337 return AliGenMUONlib::YJpsiPbPb(px, dummy);
340 Double_t AliGenMUONlib::YJpsiCDFscaledPP( const Double_t *px, const Double_t* dummy)
343 return AliGenMUONlib::YJpsiPP(px, dummy);
346 Double_t AliGenMUONlib::YJpsiCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
353 // scaled from YJpsiPP(14 TeV) using 10 TeV / 14 TeV ratio of y-spectra in LO pQCD.
354 // see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
357 Double_t c[5] = {2.46681e+01, 8.91486e+01, -3.21227e+01, 3.63075e+00, -1.32047e-01};
359 Double_t x = TMath::Abs(px[0]);
363 y = 98.523 - 1.3664 * x * x;
364 } else if (x < 7.5) {
367 while (j > 0) y = y * x + c[--j];
377 Double_t AliGenMUONlib::YJpsiPP( const Double_t *px, const Double_t */*dummy*/)
387 // mc = 1.4 GeV, pt-kick 1 GeV
390 Double_t c[5] = {1.38532e+00, 1.00596e+02, -3.46378e+01, 3.94172e+00, -1.48319e-01};
391 Double_t x = TMath::Abs(px[0]);
395 y = 96.455 - 0.8483 * x * x;
396 } else if (x < 7.9) {
399 while (j > 0) y = y * x + c[--j];
407 Double_t AliGenMUONlib::YJpsiBPbPb( const Double_t *px, const Double_t */*dummy*/)
411 // J/Psi from B->J/Psi X
416 Double_t c[7] = {7.37025e-02, 0., -2.94487e-03, 0., 6.07953e-06, 0., 5.39219e-07};
418 Double_t x = TMath::Abs(px[0]);
426 while (j > 0) y = y * x + c[--j];
434 // particle composition
436 Int_t AliGenMUONlib::IpJpsi(TRandom *)
441 Int_t AliGenMUONlib::IpPsiP(TRandom *)
443 // Psi prime composition
446 Int_t AliGenMUONlib::IpJpsiFamily(TRandom *)
450 Float_t r = gRandom->Rndm();
465 //____________________________________________________________
466 Double_t AliGenMUONlib::PtUpsilon( const Double_t *px, const Double_t */*dummy*/ )
469 const Double_t kpt0 = 5.3;
470 const Double_t kxn = 2.5;
473 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
474 return x/TMath::Power(pass1,kxn);
477 Double_t AliGenMUONlib::PtUpsilonCDFscaled( const Double_t *px, const Double_t */*dummy*/ )
480 const Double_t kpt0 = 7.753;
481 const Double_t kxn = 3.042;
484 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
485 return x/TMath::Power(pass1,kxn);
488 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP( const Double_t *px, const Double_t */*dummy*/ )
494 // scaled from CDF data at 2 TeV
496 const Double_t kpt0 = 8.610;
497 const Double_t kxn = 3.051;
500 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
501 return x/TMath::Power(pass1,kxn);
504 Double_t AliGenMUONlib::PtUpsilonCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
510 // scaled from CDF data at 2 TeV
512 const Double_t kpt0 = 8.235;
513 const Double_t kxn = 3.051;
516 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
517 return x/TMath::Power(pass1,kxn);
520 Double_t AliGenMUONlib::PtUpsilonFlat( const Double_t */*px*/, const Double_t */*dummy*/ )
525 Double_t AliGenMUONlib::PtUpsilonPbPb( const Double_t *px, const Double_t */*dummy*/)
534 // mc = 1.4 GeV, pt-kick 1 GeV
538 -1.03488e+01, 1.28065e+01, -6.60500e+00, 1.66140e+00,
539 -2.34293e-01, 1.86925e-02, -7.80708e-04, 1.30610e-05
545 while (j > 0) y = y * x +c[--j];
546 y = x * TMath::Exp(y);
553 Double_t AliGenMUONlib::PtUpsilonPP( const Double_t *px, const Double_t */*dummy*/)
562 // mc = 1.4 GeV, pt-kick 1 GeV
565 Double_t c[8] = {-7.93955e+00, 1.06306e+01, -5.21392e+00, 1.19703e+00,
566 -1.45718e-01, 8.95151e-03, -2.04806e-04, -1.13053e-06};
572 while (j > 0) y = y * x +c[--j];
573 y = x * TMath::Exp(y);
583 //____________________________________________________________
584 Double_t AliGenMUONlib::YUpsilon(const Double_t *py, const Double_t */*dummy*/)
587 const Double_t ky0 = 3.;
588 const Double_t kb=1.;
590 Double_t y=TMath::Abs(*py);
595 yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
600 Double_t AliGenMUONlib::YUpsilonPbPb( const Double_t *px, const Double_t */*dummy*/)
610 // mc = 1.4 GeV, pt-kick 1 GeV
613 Double_t c[7] = {3.40036e-01, -3.98882e-07, -4.48398e-03, 8.46411e-08, -6.10854e-04,
614 -2.99753e-09, 1.28895e-05};
615 Double_t x = TMath::Abs(px[0]);
616 if (x > 5.55) return 0.;
618 Double_t y = c[j = 6];
619 while (j > 0) y = y * x +c[--j];
623 Double_t AliGenMUONlib::YUpsilonCDFscaled( const Double_t *px, const Double_t *dummy)
626 return AliGenMUONlib::YUpsilonPbPb(px, dummy);
630 Double_t AliGenMUONlib::YUpsilonCDFscaledPP( const Double_t *px, const Double_t *dummy)
633 return AliGenMUONlib::YUpsilonPP(px, dummy);
637 Double_t AliGenMUONlib::YUpsilonFlat( const Double_t */*px*/, const Double_t */*dummy*/)
644 Double_t AliGenMUONlib::YUpsilonCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
651 // scaled from YUpsilonPP(14 TeV) using 10 TeV / 14 TeV ratio of y-spectra in LO pQCD.
652 // see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
654 Double_t c[4] = {1.12979e+00, -2.46155e-02, -7.37561e-04, 1.58824e-05};
655 Double_t x = TMath::Abs(px[0]);
656 if (x > 6.1) return 0.;
658 Double_t y = c[j = 3];
659 while (j > 0) y = y * x*x +c[--j];
663 Double_t AliGenMUONlib::YUpsilonPP( const Double_t *px, const Double_t */*dummy*/)
673 // mc = 1.4 GeV, pt-kick 1 GeV
675 Double_t c[7] = {8.91936e-01, -6.46645e-07, -1.52774e-02, 4.28677e-08, -7.01517e-04,
676 -6.20539e-10, 1.29943e-05};
677 Double_t x = TMath::Abs(px[0]);
678 if (x > 6.2) return 0.;
680 Double_t y = c[j = 6];
681 while (j > 0) y = y * x +c[--j];
685 // particle composition
687 Int_t AliGenMUONlib::IpUpsilon(TRandom *)
692 Int_t AliGenMUONlib::IpUpsilonP(TRandom *)
697 Int_t AliGenMUONlib::IpUpsilonPP(TRandom *)
702 Int_t AliGenMUONlib::IpUpsilonFamily(TRandom *)
706 Float_t r = gRandom->Rndm();
710 } else if (r < 0.896) {
723 // pt-distribution (by scaling of pion distribution)
724 //____________________________________________________________
725 Double_t AliGenMUONlib::PtPhi( const Double_t *px, const Double_t */*dummy*/)
728 return PtScal(*px,7);
731 Double_t AliGenMUONlib::YPhi( const Double_t *px, const Double_t */*dummy*/)
735 return YJpsi(px,dum);
737 // particle composition
739 Int_t AliGenMUONlib::IpPhi(TRandom *)
749 // pt-distribution (by scaling of pion distribution)
750 //____________________________________________________________
751 Double_t AliGenMUONlib::PtOmega( const Double_t *px, const Double_t */*dummy*/)
754 return PtScal(*px,5);
757 Double_t AliGenMUONlib::YOmega( const Double_t *px, const Double_t */*dummy*/)
761 return YJpsi(px,dum);
763 // particle composition
765 Int_t AliGenMUONlib::IpOmega(TRandom *)
776 // pt-distribution (by scaling of pion distribution)
777 //____________________________________________________________
778 Double_t AliGenMUONlib::PtEta( const Double_t *px, const Double_t */*dummy*/)
781 return PtScal(*px,3);
784 Double_t AliGenMUONlib::YEta( const Double_t *px, const Double_t */*dummy*/)
788 return YJpsi(px,dum);
790 // particle composition
792 Int_t AliGenMUONlib::IpEta(TRandom *)
803 //____________________________________________________________
804 Double_t AliGenMUONlib::PtCharm( const Double_t *px, const Double_t */*dummy*/)
807 const Double_t kpt0 = 2.25;
808 const Double_t kxn = 3.17;
811 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
812 return x/TMath::Power(pass1,kxn);
815 Double_t AliGenMUONlib::PtCharmCentral( const Double_t *px, const Double_t */*dummy*/)
818 const Double_t kpt0 = 2.12;
819 const Double_t kxn = 2.78;
822 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
823 return x/TMath::Power(pass1,kxn);
825 Double_t AliGenMUONlib::PtCharmF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
827 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
828 // PtCharmFiMjSkPP = PtCharmF0M0S0PP * (dN(i,j,k)/dpt / dN(0,0,0)/dpt)_MNR
829 // i=0,1,2; j=0,1,2; k=0,1,...,6
830 // dN(i,j,k)/dpt - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
831 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
832 // calculations for the following inputs:
833 // Peterson fragmentation function (F) with \epsilon_c = 0.02, 0.002 & 0.11
834 // for i=0,1 & 2 respectively; quark mass (M) of 1.5, 1.3 & 1.7 GeV
835 // for j=0,1 & 2 respectively;
836 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
837 // with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2
838 // for k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
839 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
840 // June 2008, Smbat.Grigoryan@cern.ch
843 // Pythia6.214 (kCharmppMNRwmi, PDF = CTEQ5L, quark mass = 1.2 GeV, PtHard > 2.76 GeV/c)
844 // for pp collisions at 14 TeV with one c-cbar pair per event.
845 // Corresponding NLO total cross section is 5.68 mb
848 const Double_t kpt0 = 2.2930;
849 const Double_t kxn = 3.1196;
850 Double_t c[3]={-5.2180e-01,1.8753e-01,2.8669e-02};
853 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
854 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
856 Double_t AliGenMUONlib::PtCharmF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
859 // Corresponding NLO total cross section is 6.06 mb
860 const Double_t kpt0 = 2.8669;
861 const Double_t kxn = 3.1044;
862 Double_t c[3]={-4.6714e-01,1.5005e-01,4.5003e-02};
865 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
866 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
868 Double_t AliGenMUONlib::PtCharmF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
871 // Corresponding NLO total cross section is 6.06 mb
872 const Double_t kpt0 = 1.8361;
873 const Double_t kxn = 3.2966;
874 Double_t c[3]={-6.1550e-01,2.6498e-01,1.0728e-02};
877 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
878 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
880 Double_t AliGenMUONlib::PtCharmF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
883 // Corresponding NLO total cross section is 7.69 mb
884 const Double_t kpt0 = 2.1280;
885 const Double_t kxn = 3.1397;
886 Double_t c[3]={-5.4021e-01,2.0944e-01,2.5211e-02};
889 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
890 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
892 Double_t AliGenMUONlib::PtCharmF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
895 // Corresponding NLO total cross section is 4.81 mb
896 const Double_t kpt0 = 2.4579;
897 const Double_t kxn = 3.1095;
898 Double_t c[3]={-5.1497e-01,1.7532e-01,3.2429e-02};
901 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
902 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
904 Double_t AliGenMUONlib::PtCharmF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
907 // Corresponding NLO total cross section is 14.09 mb
908 const Double_t kpt0 = 2.1272;
909 const Double_t kxn = 3.1904;
910 Double_t c[3]={-4.6088e-01,2.1918e-01,2.3055e-02};
913 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
914 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
916 Double_t AliGenMUONlib::PtCharmF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
919 // Corresponding NLO total cross section is 1.52 mb
920 const Double_t kpt0 = 2.8159;
921 const Double_t kxn = 3.0857;
922 Double_t c[3]={-6.4691e-01,2.0289e-01,2.4922e-02};
925 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
926 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
928 Double_t AliGenMUONlib::PtCharmF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
931 // Corresponding NLO total cross section is 3.67 mb
932 const Double_t kpt0 = 2.7297;
933 const Double_t kxn = 3.3019;
934 Double_t c[3]={-6.2216e-01,1.9031e-01,1.5341e-02};
937 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
938 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
940 Double_t AliGenMUONlib::PtCharmF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
943 // Corresponding NLO total cross section is 3.38 mb
944 const Double_t kpt0 = 2.3894;
945 const Double_t kxn = 3.1075;
946 Double_t c[3]={-4.9742e-01,1.7032e-01,2.5994e-02};
949 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
950 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
952 Double_t AliGenMUONlib::PtCharmF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
955 // Corresponding NLO total cross section is 10.37 mb
956 const Double_t kpt0 = 2.0187;
957 const Double_t kxn = 3.3011;
958 Double_t c[3]={-3.9869e-01,2.9248e-01,1.1763e-02};
961 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
962 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
964 Double_t AliGenMUONlib::PtCharmF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
967 // Corresponding NLO total cross section is 7.22 mb
968 const Double_t kpt0 = 2.1089;
969 const Double_t kxn = 3.1848;
970 Double_t c[3]={-4.6275e-01,1.8114e-01,2.1363e-02};
973 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
974 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
978 Double_t AliGenMUONlib::YCharm( const Double_t *px, const Double_t */*dummy*/)
980 // Charm y :: Carrer & Dainese : ALICE-INT-2003-019 v.3 (hep-ph/0311225)
981 // Pythia tuned to reproduce the distribution given by the HVQMNR program based on NLO calculations (pQCD)
982 // shadowing + kt broadening
985 Double_t c[2]={-2.42985e-03,-2.31001e-04};
986 Double_t y=1+(c[0]*TMath::Power(x,2))+(c[1]*TMath::Power(x,4));
989 if (TMath::Abs(x)>8) {
993 ycharm=TMath::Power(y,3);
998 Double_t AliGenMUONlib::YCharmF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1000 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1001 // YCharmFiMjSkPP = YCharmF0M0S0PP * (dN(i,j,k)/dy / dN(0,0,0)/dy)_MNR
1002 // i=0,1,2; j=0,1,2; k=0,1,...,6
1003 // dN(i,j,k)/dy - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1004 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1005 // calculations for the following inputs:
1006 // Peterson fragmentation function (F) with \epsilon_c = 0.02, 0.002 & 0.11
1007 // for i=0,1 & 2 respectively; quark mass (M) of 1.5, 1.3 & 1.7 GeV
1008 // for j=0,1 & 2 respectively;
1009 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1010 // with a/b = 1/1,1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2 for
1011 // k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1012 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1013 // June 2008, Smbat.Grigoryan@cern.ch
1016 // Pythia6.214 (kCharmppMNRwmi, PDF = CTEQ5L, quark mass = 1.2 GeV, PtHard > 2.76 GeV/c)
1017 // for pp collisions at 14 TeV with one c-cbar pair per event.
1018 // Corresponding NLO total cross section is 5.68 mb
1021 Double_t c[2]={7.0909e-03,6.1967e-05};
1022 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1025 if (TMath::Abs(x)>9) {
1029 ycharm=TMath::Power(y,3);
1034 Double_t AliGenMUONlib::YCharmF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1037 // Corresponding NLO total cross section is 6.06 mb
1039 Double_t c[2]={6.9707e-03,6.0971e-05};
1040 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1043 if (TMath::Abs(x)>9) {
1047 ycharm=TMath::Power(y,3);
1052 Double_t AliGenMUONlib::YCharmF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1055 // Corresponding NLO total cross section is 6.06 mb
1057 Double_t c[2]={7.1687e-03,6.5303e-05};
1058 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1061 if (TMath::Abs(x)>9) {
1065 ycharm=TMath::Power(y,3);
1070 Double_t AliGenMUONlib::YCharmF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
1073 // Corresponding NLO total cross section is 7.69 mb
1075 Double_t c[2]={5.9090e-03,7.1854e-05};
1076 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1079 if (TMath::Abs(x)>9) {
1083 ycharm=TMath::Power(y,3);
1088 Double_t AliGenMUONlib::YCharmF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
1091 // Corresponding NLO total cross section is 4.81 mb
1093 Double_t c[2]={8.0882e-03,5.5872e-05};
1094 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1097 if (TMath::Abs(x)>9) {
1101 ycharm=TMath::Power(y,3);
1106 Double_t AliGenMUONlib::YCharmF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
1109 // Corresponding NLO total cross section is 14.09 mb
1111 Double_t c[2]={7.2520e-03,6.2691e-05};
1112 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1115 if (TMath::Abs(x)>9) {
1119 ycharm=TMath::Power(y,3);
1124 Double_t AliGenMUONlib::YCharmF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
1127 // Corresponding NLO total cross section is 1.52 mb
1129 Double_t c[2]={1.1040e-04,1.4498e-04};
1130 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1133 if (TMath::Abs(x)>9) {
1137 ycharm=TMath::Power(y,3);
1142 Double_t AliGenMUONlib::YCharmF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
1145 // Corresponding NLO total cross section is 3.67 mb
1147 Double_t c[2]={-3.1328e-03,1.8270e-04};
1148 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1151 if (TMath::Abs(x)>9) {
1155 ycharm=TMath::Power(y,3);
1160 Double_t AliGenMUONlib::YCharmF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
1163 // Corresponding NLO total cross section is 3.38 mb
1165 Double_t c[2]={7.0865e-03,6.2532e-05};
1166 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1169 if (TMath::Abs(x)>9) {
1173 ycharm=TMath::Power(y,3);
1178 Double_t AliGenMUONlib::YCharmF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
1181 // Corresponding NLO total cross section is 10.37 mb
1183 Double_t c[2]={7.7070e-03,5.3533e-05};
1184 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1187 if (TMath::Abs(x)>9) {
1191 ycharm=TMath::Power(y,3);
1196 Double_t AliGenMUONlib::YCharmF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
1199 // Corresponding NLO total cross section is 7.22 mb
1201 Double_t c[2]={7.9195e-03,5.3823e-05};
1202 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1205 if (TMath::Abs(x)>9) {
1209 ycharm=TMath::Power(y,3);
1216 Int_t AliGenMUONlib::IpCharm(TRandom *ran)
1218 // Charm composition
1222 random = ran->Rndm();
1223 // Taux de production Carrer & Dainese : ALICE-INT-2003-019 v.3
1224 // >>>>> cf. tab 4 p 11
1226 if (random < 0.30) {
1228 } else if (random < 0.60) {
1230 } else if (random < 0.70) {
1232 } else if (random < 0.80) {
1234 } else if (random < 0.86) {
1236 } else if (random < 0.92) {
1238 } else if (random < 0.96) {
1252 //____________________________________________________________
1253 Double_t AliGenMUONlib::PtBeauty( const Double_t *px, const Double_t */*dummy*/)
1256 const Double_t kpt0 = 6.53;
1257 const Double_t kxn = 3.59;
1260 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1261 return x/TMath::Power(pass1,kxn);
1264 Double_t AliGenMUONlib::PtBeautyCentral( const Double_t *px, const Double_t */*dummy*/)
1267 const Double_t kpt0 = 6.14;
1268 const Double_t kxn = 2.93;
1271 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1272 return x/TMath::Power(pass1,kxn);
1274 Double_t AliGenMUONlib::PtBeautyF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1276 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1277 // PtBeautyFiMjSkPP = PtBeautyF0M0S0PP * (dN(i,j,k)/dpt / dN(0,0,0)/dpt)_MNR
1278 // i=0,1,2; j=0,1,2; k=0,1,...,6
1279 // dN(i,j,k)/dpt - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1280 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1281 // calculations for the following inputs:
1282 // Peterson fragmentation function (F) with \epsilon_b = 0.001, 0.0002 & 0.004
1283 // for i=0,1 & 2 respectively; quark mass (M) of 4.75, 4.5 & 5.0 GeV
1284 // for j=0,1 & 2 respectively;
1285 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1286 // with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2 for
1287 // k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1288 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1289 // June 2008, Smbat.Grigoryan@cern.ch
1292 // Pythia6.214 (kBeautyppMNRwmi, PDF = CTEQ5L, quark mass = 4.75 GeV, PtHard > 2.76 GeV/c)
1293 // for pp collisions at 14 TeV with one b-bbar pair per event.
1294 // Corresponding NLO total cross section is 0.494 mb
1296 const Double_t kpt0 = 8.0575;
1297 const Double_t kxn = 3.1921;
1300 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1301 return x/TMath::Power(pass1,kxn);
1303 Double_t AliGenMUONlib::PtBeautyF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1306 // Corresponding NLO total cross section is 0.445 mb
1307 const Double_t kpt0 = 8.6239;
1308 const Double_t kxn = 3.2911;
1311 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1312 return x/TMath::Power(pass1,kxn);
1314 Double_t AliGenMUONlib::PtBeautyF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1317 // Corresponding NLO total cross section is 0.445 mb
1318 const Double_t kpt0 = 7.3367;
1319 const Double_t kxn = 3.0692;
1322 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1323 return x/TMath::Power(pass1,kxn);
1325 Double_t AliGenMUONlib::PtBeautyF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
1328 // Corresponding NLO total cross section is 0.518 mb
1329 const Double_t kpt0 = 7.6409;
1330 const Double_t kxn = 3.1364;
1333 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1334 return x/TMath::Power(pass1,kxn);
1336 Double_t AliGenMUONlib::PtBeautyF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
1339 // Corresponding NLO total cross section is 0.384 mb
1340 const Double_t kpt0 = 8.4948;
1341 const Double_t kxn = 3.2546;
1344 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1345 return x/TMath::Power(pass1,kxn);
1347 Double_t AliGenMUONlib::PtBeautyF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
1350 // Corresponding NLO total cross section is 0.648 mb
1351 const Double_t kpt0 = 7.6631;
1352 const Double_t kxn = 3.1621;
1355 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1356 return x/TMath::Power(pass1,kxn);
1358 Double_t AliGenMUONlib::PtBeautyF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
1361 // Corresponding NLO total cross section is 0.294 mb
1362 const Double_t kpt0 = 8.7245;
1363 const Double_t kxn = 3.2213;
1366 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1367 return x/TMath::Power(pass1,kxn);
1369 Double_t AliGenMUONlib::PtBeautyF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
1372 // Corresponding NLO total cross section is 0.475 mb
1373 const Double_t kpt0 = 8.5296;
1374 const Double_t kxn = 3.2187;
1377 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1378 return x/TMath::Power(pass1,kxn);
1380 Double_t AliGenMUONlib::PtBeautyF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
1383 // Corresponding NLO total cross section is 0.324 mb
1384 const Double_t kpt0 = 7.9440;
1385 const Double_t kxn = 3.1614;
1388 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1389 return x/TMath::Power(pass1,kxn);
1391 Double_t AliGenMUONlib::PtBeautyF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
1394 // Corresponding NLO total cross section is 0.536 mb
1395 const Double_t kpt0 = 8.2408;
1396 const Double_t kxn = 3.3029;
1399 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1400 return x/TMath::Power(pass1,kxn);
1402 Double_t AliGenMUONlib::PtBeautyF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
1405 // Corresponding NLO total cross section is 0.420 mb
1406 const Double_t kpt0 = 7.8041;
1407 const Double_t kxn = 3.2094;
1410 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1411 return x/TMath::Power(pass1,kxn);
1415 Double_t AliGenMUONlib::YBeauty( const Double_t *px, const Double_t */*dummy*/)
1417 // Beauty y :: Carrer & Dainese : ALICE-INT-2003-019 v.3 (hep-ph/0311225)
1418 // Pythia tuned to reproduce the distribution given by the HVQMNR program based on NLO calculations (pQCD)
1419 // shadowing + kt broadening
1422 Double_t c[2]={-1.27590e-02,-2.42731e-04};
1423 Double_t y=1+c[0]*TMath::Power(x,2)+c[1]*TMath::Power(x,4);
1426 if (TMath::Abs(x)>6) {
1430 ybeauty=TMath::Power(y,3);
1435 Double_t AliGenMUONlib::YBeautyF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1437 // FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1438 // YBeautyFiMjSkPP = YBeautyF0M0S0PP * (dN(i,j,k)/dy / dN(0,0,0)/dy)_MNR
1439 // i=0,1,2; j=0,1,2; k=0,1,...,6
1440 // dN(i,j,k)/dy - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1441 // http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1442 // calculations for the following inputs:
1443 // Peterson fragmentation function (F) with \epsilon_b = 0.001, 0.0002 & 0.004
1444 // for i=0,1 & 2 respectively; quark mass (M) of 4.75, 4.5 & 5.0 GeV
1445 // for j=0,1 & 2 respectively;
1446 // factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1447 // with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2
1448 // for k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1449 // (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1450 // June 2008, Smbat.Grigoryan@cern.ch
1453 // Pythia6.214 (kBeautyppMNRwmi, PDF = CTEQ5L, quark mass = 4.75 GeV, PtHard > 2.76 GeV/c)
1454 // for pp collisions at 14 TeV with one b-bbar pair per event.
1455 // Corresponding NLO total cross section is 0.494 mb
1459 Double_t c[2]={1.2350e-02,9.2667e-05};
1460 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1463 if (TMath::Abs(x)>7.6) {
1467 ybeauty=TMath::Power(y,3);
1472 Double_t AliGenMUONlib::YBeautyF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1475 // Corresponding NLO total cross section is 0.445 mb
1477 Double_t c[2]={1.2292e-02,9.1847e-05};
1478 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1481 if (TMath::Abs(x)>7.6) {
1485 ybeauty=TMath::Power(y,3);
1490 Double_t AliGenMUONlib::YBeautyF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
1493 // Corresponding NLO total cross section is 0.445 mb
1495 Double_t c[2]={1.2436e-02,9.3709e-05};
1496 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1499 if (TMath::Abs(x)>7.6) {
1503 ybeauty=TMath::Power(y,3);
1508 Double_t AliGenMUONlib::YBeautyF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
1511 // Corresponding NLO total cross section is 0.518 mb
1513 Double_t c[2]={1.1714e-02,1.0068e-04};
1514 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1517 if (TMath::Abs(x)>7.6) {
1521 ybeauty=TMath::Power(y,3);
1526 Double_t AliGenMUONlib::YBeautyF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
1529 // Corresponding NLO total cross section is 0.384 mb
1531 Double_t c[2]={1.2944e-02,8.5500e-05};
1532 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1535 if (TMath::Abs(x)>7.6) {
1539 ybeauty=TMath::Power(y,3);
1544 Double_t AliGenMUONlib::YBeautyF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
1547 // Corresponding NLO total cross section is 0.648 mb
1549 Double_t c[2]={1.2455e-02,9.2713e-05};
1550 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1553 if (TMath::Abs(x)>7.6) {
1557 ybeauty=TMath::Power(y,3);
1562 Double_t AliGenMUONlib::YBeautyF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
1565 // Corresponding NLO total cross section is 0.294 mb
1567 Double_t c[2]={1.0897e-02,1.1878e-04};
1568 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1571 if (TMath::Abs(x)>7.6) {
1575 ybeauty=TMath::Power(y,3);
1580 Double_t AliGenMUONlib::YBeautyF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
1583 // Corresponding NLO total cross section is 0.475 mb
1585 Double_t c[2]={1.0912e-02,1.1858e-04};
1586 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1589 if (TMath::Abs(x)>7.6) {
1593 ybeauty=TMath::Power(y,3);
1598 Double_t AliGenMUONlib::YBeautyF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
1601 // Corresponding NLO total cross section is 0.324 mb
1603 Double_t c[2]={1.2378e-02,9.2490e-05};
1604 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1607 if (TMath::Abs(x)>7.6) {
1611 ybeauty=TMath::Power(y,3);
1616 Double_t AliGenMUONlib::YBeautyF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
1619 // Corresponding NLO total cross section is 0.536 mb
1621 Double_t c[2]={1.2886e-02,8.2912e-05};
1622 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1625 if (TMath::Abs(x)>7.6) {
1629 ybeauty=TMath::Power(y,3);
1634 Double_t AliGenMUONlib::YBeautyF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
1637 // Corresponding NLO total cross section is 0.420 mb
1639 Double_t c[2]={1.3106e-02,8.0115e-05};
1640 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1643 if (TMath::Abs(x)>7.6) {
1647 ybeauty=TMath::Power(y,3);
1653 Int_t AliGenMUONlib::IpBeauty(TRandom *ran)
1655 // Beauty Composition
1658 random = ran->Rndm();
1660 // Taux de production Carrer & Dainese : ALICE-INT-2003-019 v.3
1661 // >>>>> cf. tab 4 p 11
1663 if (random < 0.20) {
1665 } else if (random < 0.40) {
1667 } else if (random < 0.605) {
1669 } else if (random < 0.81) {
1671 } else if (random < 0.87) {
1673 } else if (random < 0.93) {
1675 } else if (random < 0.965) {
1685 typedef Double_t (*GenFunc) (const Double_t*, const Double_t*);
1686 GenFunc AliGenMUONlib::GetPt(Int_t param, const char* tname) const
1688 // Return pointer to pT parameterisation
1689 TString sname = TString(tname);
1705 if (sname == "Vogt" || sname == "Vogt PbPb") {
1707 } else if (sname == "Vogt pp") {
1709 } else if (sname == "CDF scaled") {
1710 func=PtJpsiCDFscaled;
1711 } else if (sname == "CDF pp") {
1712 func=PtJpsiCDFscaledPP;
1713 } else if (sname == "CDF pp 10") {
1714 func=PtJpsiCDFscaledPP10;
1715 } else if (sname == "Flat") {
1724 case kUpsilonFamily:
1728 if (sname == "Vogt" || sname == "Vogt PbPb") {
1730 } else if (sname == "Vogt pp") {
1732 } else if (sname == "CDF scaled") {
1733 func=PtUpsilonCDFscaled;
1734 } else if (sname == "CDF pp") {
1735 func=PtUpsilonCDFscaledPP;
1736 } else if (sname == "CDF pp 10") {
1737 func=PtUpsilonCDFscaledPP10;
1738 } else if (sname == "Flat") {
1745 if (sname == "F0M0S0 pp") {
1746 func=PtCharmF0M0S0PP;
1747 } else if (sname == "F1M0S0 pp") {
1748 func=PtCharmF1M0S0PP;
1749 } else if (sname == "F2M0S0 pp") {
1750 func=PtCharmF2M0S0PP;
1751 } else if (sname == "F0M1S0 pp") {
1752 func=PtCharmF0M1S0PP;
1753 } else if (sname == "F0M2S0 pp") {
1754 func=PtCharmF0M2S0PP;
1755 } else if (sname == "F0M0S1 pp") {
1756 func=PtCharmF0M0S1PP;
1757 } else if (sname == "F0M0S2 pp") {
1758 func=PtCharmF0M0S2PP;
1759 } else if (sname == "F0M0S3 pp") {
1760 func=PtCharmF0M0S3PP;
1761 } else if (sname == "F0M0S4 pp") {
1762 func=PtCharmF0M0S4PP;
1763 } else if (sname == "F0M0S5 pp") {
1764 func=PtCharmF0M0S5PP;
1765 } else if (sname == "F0M0S6 pp") {
1766 func=PtCharmF0M0S6PP;
1767 } else if (sname == "central") {
1768 func=PtCharmCentral;
1774 if (sname == "F0M0S0 pp") {
1775 func=PtBeautyF0M0S0PP;
1776 } else if (sname == "F1M0S0 pp") {
1777 func=PtBeautyF1M0S0PP;
1778 } else if (sname == "F2M0S0 pp") {
1779 func=PtBeautyF2M0S0PP;
1780 } else if (sname == "F0M1S0 pp") {
1781 func=PtBeautyF0M1S0PP;
1782 } else if (sname == "F0M2S0 pp") {
1783 func=PtBeautyF0M2S0PP;
1784 } else if (sname == "F0M0S1 pp") {
1785 func=PtBeautyF0M0S1PP;
1786 } else if (sname == "F0M0S2 pp") {
1787 func=PtBeautyF0M0S2PP;
1788 } else if (sname == "F0M0S3 pp") {
1789 func=PtBeautyF0M0S3PP;
1790 } else if (sname == "F0M0S4 pp") {
1791 func=PtBeautyF0M0S4PP;
1792 } else if (sname == "F0M0S5 pp") {
1793 func=PtBeautyF0M0S5PP;
1794 } else if (sname == "F0M0S6 pp") {
1795 func=PtBeautyF0M0S6PP;
1796 } else if (sname == "central") {
1797 func=PtBeautyCentral;
1822 printf("<AliGenMUONlib::GetPt> unknown parametrisation\n");
1827 GenFunc AliGenMUONlib::GetY(Int_t param, const char* tname) const
1830 // Return pointer to y- parameterisation
1832 TString sname = TString(tname);
1848 if (sname == "Vogt" || sname == "Vogt PbPb") {
1850 } else if (sname == "Vogt pp"){
1852 } else if (sname == "CDF scaled") {
1853 func=YJpsiCDFscaled;
1854 } else if (sname == "CDF pp") {
1855 func=YJpsiCDFscaledPP;
1856 } else if (sname == "CDF pp 10") {
1857 func=YJpsiCDFscaledPP10;
1858 } else if (sname == "Flat") {
1867 case kUpsilonFamily:
1871 if (sname == "Vogt" || sname == "Vogt PbPb") {
1873 } else if (sname == "Vogt pp") {
1875 } else if (sname == "CDF scaled") {
1876 func=YUpsilonCDFscaled;
1877 } else if (sname == "CDF pp") {
1878 func=YUpsilonCDFscaledPP;
1879 } else if (sname == "CDF pp 10") {
1880 func=YUpsilonCDFscaledPP10;
1881 } else if (sname == "Flat") {
1888 if (sname == "F0M0S0 pp") {
1889 func=YCharmF0M0S0PP;
1890 } else if (sname == "F1M0S0 pp") {
1891 func=YCharmF1M0S0PP;
1892 } else if (sname == "F2M0S0 pp") {
1893 func=YCharmF2M0S0PP;
1894 } else if (sname == "F0M1S0 pp") {
1895 func=YCharmF0M1S0PP;
1896 } else if (sname == "F0M2S0 pp") {
1897 func=YCharmF0M2S0PP;
1898 } else if (sname == "F0M0S1 pp") {
1899 func=YCharmF0M0S1PP;
1900 } else if (sname == "F0M0S2 pp") {
1901 func=YCharmF0M0S2PP;
1902 } else if (sname == "F0M0S3 pp") {
1903 func=YCharmF0M0S3PP;
1904 } else if (sname == "F0M0S4 pp") {
1905 func=YCharmF0M0S4PP;
1906 } else if (sname == "F0M0S5 pp") {
1907 func=YCharmF0M0S5PP;
1908 } else if (sname == "F0M0S6 pp") {
1909 func=YCharmF0M0S6PP;
1915 if (sname == "F0M0S0 pp") {
1916 func=YBeautyF0M0S0PP;
1917 } else if (sname == "F1M0S0 pp") {
1918 func=YBeautyF1M0S0PP;
1919 } else if (sname == "F2M0S0 pp") {
1920 func=YBeautyF2M0S0PP;
1921 } else if (sname == "F0M1S0 pp") {
1922 func=YBeautyF0M1S0PP;
1923 } else if (sname == "F0M2S0 pp") {
1924 func=YBeautyF0M2S0PP;
1925 } else if (sname == "F0M0S1 pp") {
1926 func=YBeautyF0M0S1PP;
1927 } else if (sname == "F0M0S2 pp") {
1928 func=YBeautyF0M0S2PP;
1929 } else if (sname == "F0M0S3 pp") {
1930 func=YBeautyF0M0S3PP;
1931 } else if (sname == "F0M0S4 pp") {
1932 func=YBeautyF0M0S4PP;
1933 } else if (sname == "F0M0S5 pp") {
1934 func=YBeautyF0M0S5PP;
1935 } else if (sname == "F0M0S6 pp") {
1936 func=YBeautyF0M0S6PP;
1961 printf("<AliGenMUONlib::GetY> unknown parametrisation\n");
1971 //____________________________________________________________
1972 Double_t AliGenMUONlib::PtChic0( const Double_t *px, const Double_t */*dummy*/)
1975 const Double_t kpt0 = 4.;
1976 const Double_t kxn = 3.6;
1979 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1980 return x/TMath::Power(pass1,kxn);
1982 Double_t AliGenMUONlib::PtChic1( const Double_t *px, const Double_t */*dummy*/)
1985 const Double_t kpt0 = 4.;
1986 const Double_t kxn = 3.6;
1989 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1990 return x/TMath::Power(pass1,kxn);
1992 Double_t AliGenMUONlib::PtChic2( const Double_t *px, const Double_t */*dummy*/)
1995 const Double_t kpt0 = 4.;
1996 const Double_t kxn = 3.6;
1999 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2000 return x/TMath::Power(pass1,kxn);
2002 Double_t AliGenMUONlib::PtChic( const Double_t *px, const Double_t */*dummy*/)
2005 const Double_t kpt0 = 4.;
2006 const Double_t kxn = 3.6;
2009 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2010 return x/TMath::Power(pass1,kxn);
2015 //____________________________________________________________
2016 Double_t AliGenMUONlib::YChic0(const Double_t *py, const Double_t */*dummy*/)
2019 const Double_t ky0 = 4.;
2020 const Double_t kb=1.;
2022 Double_t y=TMath::Abs(*py);
2027 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2031 Double_t AliGenMUONlib::YChic1(const Double_t *py, const Double_t */*dummy*/)
2034 const Double_t ky0 = 4.;
2035 const Double_t kb=1.;
2037 Double_t y=TMath::Abs(*py);
2042 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2046 Double_t AliGenMUONlib::YChic2(const Double_t *py, const Double_t */*dummy*/)
2049 const Double_t ky0 = 4.;
2050 const Double_t kb=1.;
2052 Double_t y=TMath::Abs(*py);
2057 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2061 Double_t AliGenMUONlib::YChic(const Double_t *py, const Double_t */*dummy*/)
2064 const Double_t ky0 = 4.;
2065 const Double_t kb=1.;
2067 Double_t y=TMath::Abs(*py);
2072 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2076 // particle composition
2078 Int_t AliGenMUONlib::IpChic0(TRandom *)
2084 Int_t AliGenMUONlib::IpChic1(TRandom *)
2089 Int_t AliGenMUONlib::IpChic2(TRandom *)
2091 // Chi_c2 prime composition
2094 Int_t AliGenMUONlib::IpChic(TRandom *)
2098 Float_t r = gRandom->Rndm();
2101 } else if( r < 0.377 ) {
2110 //_____________________________________________________________
2112 typedef Int_t (*GenFuncIp) (TRandom *);
2113 GenFuncIp AliGenMUONlib::GetIp(Int_t param, const char* /*tname*/) const
2115 // Return pointer to particle type parameterisation
2141 case kUpsilonFamily:
2142 func=IpUpsilonFamily;
2176 printf("<AliGenMUONlib::GetIp> unknown parametrisation\n");
2183 Float_t AliGenMUONlib::Interpolate(Float_t x, Float_t* y, Float_t x0,
2188 // Neville's alorithm for interpolation
2194 // n: number of data points
2195 // no: order of polynom
2197 Float_t* c = new Float_t[n];
2198 Float_t* d = new Float_t[n];
2200 for (i = 0; i < n; i++) {
2205 Int_t ns = int((x - x0)/dx);
2209 for (m = 0; m < no; m++) {
2210 for (i = 0; i < n-m; i++) {
2211 Float_t ho = x0 + Float_t(i) * dx - x;
2212 Float_t hp = x0 + Float_t(i+m+1) * dx - x;
2213 Float_t w = c[i+1] - d[i];
2214 Float_t den = ho-hp;
2221 if (2*ns < (n-m-1)) {