]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EVGEN/AliGenMUONlib.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fixes for coding violations
[u/mrichter/AliRoot.git] / EVGEN / AliGenMUONlib.cxx
CommitLineData
4c039060 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
88cb7938 16/* $Id$ */
4c039060 17
53904666 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.
23//
24// andreas.morsch@cern.ch
25//
26
65fb704d 27#include "TMath.h"
28#include "TRandom.h"
29
fe4da5cc 30#include "AliGenMUONlib.h"
5c3fd7ea 31
fe4da5cc 32ClassImp(AliGenMUONlib)
33//
34// Pions
75e0cc59 35Double_t AliGenMUONlib::PtPion(const Double_t *px, const Double_t* /*dummy*/)
fe4da5cc 36{
37//
38// PT-PARAMETERIZATION CDF, PRL 61(88) 1819
39// POWER LAW FOR PT > 500 MEV
40// MT SCALING BELOW (T=160 MEV)
41//
d90f80fd 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;
47 const Double_t kb=1.;
fe4da5cc 48 Double_t y, y1, xmpi2, ynorm, a;
49 Double_t x=*px;
50 //
d90f80fd 51 y1=TMath::Power(kp0/(kp0+kxlim),kxn);
52 xmpi2=kxmpi*kxmpi;
53 ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
fe4da5cc 54 a=ynorm/y1;
d90f80fd 55 if (x > kxlim)
56 y=a*TMath::Power(kp0/(kp0+x),kxn);
fe4da5cc 57 else
d90f80fd 58 y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
fe4da5cc 59 return y*x;
60}
753690b0 61//
62// y-distribution
63//
75e0cc59 64Double_t AliGenMUONlib::YPion( const Double_t *py, const Double_t */*dummy*/)
753690b0 65{
d90f80fd 66// Pion y
2280e6af 67 Double_t y=TMath::Abs(*py);
68/*
d90f80fd 69 const Double_t ka = 7000.;
70 const Double_t kdy = 4.;
d90f80fd 71 Double_t ex = y*y/(2*kdy*kdy);
72 return ka*TMath::Exp(-ex);
2280e6af 73*/
74 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
75
753690b0 76}
77// particle composition
78//
65fb704d 79Int_t AliGenMUONlib::IpPion(TRandom *ran)
753690b0 80{
d90f80fd 81// Pion composition
65fb704d 82 if (ran->Rndm() < 0.5) {
753690b0 83 return 211;
84 } else {
85 return -211;
86 }
87}
fe4da5cc 88
89//____________________________________________________________
90//
91// Mt-scaling
92
93Double_t AliGenMUONlib::PtScal(Double_t pt, Int_t np)
94{
95 // SCALING EN MASSE PAR RAPPORT A PTPI
96 // MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
d90f80fd 97 const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
fe4da5cc 98 // VALUE MESON/PI AT 5 GEV
d90f80fd 99 const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
d90f80fd 100 Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
101 Double_t fmax2=f5/kfmax[np];
fe4da5cc 102 // PIONS
103 Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
104 Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
d90f80fd 105 (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
fe4da5cc 106 return fmtscal*ptpion;
107}
108//
753690b0 109// kaon
110//
111// pt-distribution
112//____________________________________________________________
75e0cc59 113Double_t AliGenMUONlib::PtKaon( const Double_t *px, const Double_t */*dummy*/)
753690b0 114{
d90f80fd 115// Kaon pT
0db2f441 116 return PtScal(*px,1);
753690b0 117}
118
119// y-distribution
fe4da5cc 120//____________________________________________________________
75e0cc59 121Double_t AliGenMUONlib::YKaon( const Double_t *py, const Double_t */*dummy*/)
fe4da5cc 122{
d90f80fd 123// Kaon y
2280e6af 124 Double_t y=TMath::Abs(*py);
125/*
d90f80fd 126 const Double_t ka = 1000.;
127 const Double_t kdy = 4.;
fe4da5cc 128 //
d90f80fd 129 Double_t ex = y*y/(2*kdy*kdy);
130 return ka*TMath::Exp(-ex);
2280e6af 131*/
132
133 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
753690b0 134}
135
136// particle composition
137//
65fb704d 138Int_t AliGenMUONlib::IpKaon(TRandom *ran)
753690b0 139{
d90f80fd 140// Kaon composition
65fb704d 141 if (ran->Rndm() < 0.5) {
753690b0 142 return 321;
143 } else {
144 return -321;
145 }
fe4da5cc 146}
753690b0 147
fe4da5cc 148// J/Psi
149//
150//
151// pt-distribution
152//____________________________________________________________
75e0cc59 153Double_t AliGenMUONlib::PtJpsi( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 154{
d90f80fd 155// J/Psi pT
156 const Double_t kpt0 = 4.;
157 const Double_t kxn = 3.6;
fe4da5cc 158 Double_t x=*px;
159 //
d90f80fd 160 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
161 return x/TMath::Power(pass1,kxn);
fe4da5cc 162}
05932df6 163
75e0cc59 164Double_t AliGenMUONlib::PtJpsiCDFscaled( const Double_t *px, const Double_t */*dummy*/)
0e137c25 165{
166// J/Psi pT
9e9ae065 167//
168// PbPb 5.5 TeV
169// scaled from CDF data at 2 TeV
170// see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
171
172 const Double_t kpt0 = 5.100;
173 const Double_t kxn = 4.102;
174 Double_t x=*px;
175 //
176 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
177 return x/TMath::Power(pass1,kxn);
178}
179
75e0cc59 180Double_t AliGenMUONlib::PtJpsiCDFscaledPP( const Double_t *px, const Double_t */*dummy*/)
9ff13849 181{
182// J/Psi pT
183//
184// pp 14 TeV
9ff13849 185// scaled from CDF data at 2 TeV
186
9e9ae065 187 const Double_t kpt0 = 5.630;
188 const Double_t kxn = 4.071;
189 Double_t x=*px;
190 //
191 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
192 return x/TMath::Power(pass1,kxn);
193}
194
75e0cc59 195Double_t AliGenMUONlib::PtJpsiCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 196{
197// J/Psi pT
198//
199// pp 10 TeV
200// scaled from CDF data at 2 TeV
201
202 const Double_t kpt0 = 5.334;
203 const Double_t kxn = 4.071;
204 Double_t x=*px;
205 //
206 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
207 return x/TMath::Power(pass1,kxn);
208}
209
9fd56238 210Double_t AliGenMUONlib::PtJpsiCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
211{
212// J/Psi pT
213//
214// pp 8.8 TeV
215// scaled from CDF data at 2 TeV
216//
217 const Double_t kpt0 = 5.245;
218 const Double_t kxn = 4.071;
219 Double_t x=*px;
220 //
221 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
222 return x/TMath::Power(pass1,kxn);
223}
224
225Double_t AliGenMUONlib::PtJpsiCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
226{
227// J/Psi pT
228//
229// pp 7 TeV
230// scaled from CDF data at 2 TeV
231
232 const Double_t kpt0 = 5.072;
233 const Double_t kxn = 4.071;
234 Double_t x=*px;
235 //
236 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
237 return x/TMath::Power(pass1,kxn);
238}
239
240Double_t AliGenMUONlib::PtJpsiCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
241{
242// J/Psi pT
243//
244// pp 3.94 TeV
245// scaled from CDF data at 2 TeV
246//
247 const Double_t kpt0 = 4.647;
248 const Double_t kxn = 4.071;
249 Double_t x=*px;
250 //
251 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
252 return x/TMath::Power(pass1,kxn);
253}
254
255Double_t AliGenMUONlib::PtJpsiCDFscaledPPb9( const Double_t *px, const Double_t *dummy)
256{
257// J/Psi pT
258//
259// pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.80
260//
261 Double_t c[5] = {6.42774e-01, 1.86168e-02, -6.77296e-04, 8.93512e-06, 1.31586e-07};
262 Double_t x=*px;
263 Double_t y;
264 Int_t j;
265 y = c[j = 4];
266 while (j > 0) y = y * x + c[--j];
267 //
268 Double_t d = 1.+c[4]*TMath::Power(x,4);
269 return y/d * AliGenMUONlib::PtJpsiCDFscaledPP9(px,dummy);
270}
271
272Double_t AliGenMUONlib::PtJpsiCDFscaledPbP9( const Double_t *px, const Double_t *dummy)
273{
274// J/Psi pT
275//
276// Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.80
277//
278 Double_t c[5] = {8.58557e-01, 5.39791e-02, -4.75180e-03, 2.49463e-04, 5.52396e-05};
279 Double_t x=*px;
280 Double_t y;
281 Int_t j;
282 y = c[j = 4];
283 while (j > 0) y = y * x + c[--j];
284 //
285 Double_t d = 1.+c[4]*TMath::Power(x,4);
286 return y/d * AliGenMUONlib::PtJpsiCDFscaledPP9(px,dummy);
287}
288
289Double_t AliGenMUONlib::PtJpsiCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
290{
291// J/Psi pT
292//
293// PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.66
294//
295 Double_t c[5] = {6.01022e-01, 4.70988e-02, -2.27917e-03, 3.09885e-05, 1.31955e-06};
296 Double_t x=*px;
297 Double_t y;
298 Int_t j;
299 y = c[j = 4];
300 while (j > 0) y = y * x + c[--j];
301 //
302 Double_t d = 1.+c[4]*TMath::Power(x,4);
303 return y/d * AliGenMUONlib::PtJpsiCDFscaledPP4(px,dummy);
304}
305
75e0cc59 306Double_t AliGenMUONlib::PtJpsiFlat( const Double_t */*px*/, const Double_t */*dummy*/ )
4ca8d070 307{
308 return 1.;
309}
310
75e0cc59 311Double_t AliGenMUONlib::PtJpsiPbPb( const Double_t *px, const Double_t */*dummy*/)
05932df6 312{
1af7144e 313// J/Psi pT spectrum
05932df6 314//
315// R. Vogt 2002
316// PbPb 5.5 TeV
317// MRST HO
318// mc = 1.4 GeV, pt-kick 1 GeV
319//
1af7144e 320 Float_t x = px[0];
321 Float_t c[8] = {
322 -2.13098e+00, 9.46552e+00, -5.06799e+00, 1.27260e+00,
323 -1.83806e-01, 1.55853e-02, -7.23241e-04, 1.42105e-05
05932df6 324 };
1af7144e 325
3d905dd7 326 Double_t y;
36349df1 327 if (x < 10.) {
3d905dd7 328 Int_t j;
329 y = c[j = 7];
330 while (j > 0) y = y * x +c[--j];
331 y = x * TMath::Exp(y);
332 } else {
333 y = 0.;
334 }
1af7144e 335 return y;
05932df6 336}
17d28ba5 337
75e0cc59 338Double_t AliGenMUONlib::PtJpsiBPbPb( const Double_t *px, const Double_t */*dummy*/)
17d28ba5 339{
340// J/Psi pT spectrum
341// B -> J/Psi X
342 Double_t x0 = 4.0384;
343 Double_t n = 3.0288;
344
345 Double_t x = px[0];
346 Double_t y = x / TMath::Power((1. + (x/x0)*(x/x0)), n);
347
348 return y;
349}
350
351
75e0cc59 352Double_t AliGenMUONlib::PtJpsiPP( const Double_t *px, const Double_t */*dummy*/)
bb6e81ac 353{
354// J/Psi pT spectrum
355//
356// R. Vogt 2002
357// pp 14 TeV
358// MRST HO
359// mc = 1.4 GeV, pt-kick 1 GeV
360//
361 Float_t x = px[0];
362 Float_t c[4] = {8.47471e+00, -1.93567e+00, 1.50271e-01, -5.51212e-03};
363
364 Double_t y;
365 if (x < 10.) {
366 Int_t j;
367 y = c[j = 3];
368 while (j > 0) y = y * x +c[--j];
369 y = x * TMath::Exp(y);
370 } else {
371 y = 0.;
372 }
373 return y;
374}
375
fe4da5cc 376//
377// y-distribution
378//____________________________________________________________
75e0cc59 379Double_t AliGenMUONlib::YJpsi(const Double_t *py, const Double_t */*dummy*/)
fe4da5cc 380{
d90f80fd 381// J/psi y
382 const Double_t ky0 = 4.;
383 const Double_t kb=1.;
fe4da5cc 384 Double_t yj;
385 Double_t y=TMath::Abs(*py);
386 //
d90f80fd 387 if (y < ky0)
388 yj=kb;
fe4da5cc 389 else
d90f80fd 390 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
fe4da5cc 391 return yj;
392}
05932df6 393
75e0cc59 394Double_t AliGenMUONlib::YJpsiFlat( const Double_t */*py*/, const Double_t */*dummy*/ )
4ca8d070 395{
396 return 1.;
397}
398
05932df6 399
75e0cc59 400Double_t AliGenMUONlib::YJpsiPbPb( const Double_t *px, const Double_t */*dummy*/)
05932df6 401{
402
403//
404// J/Psi y
405//
406//
407// R. Vogt 2002
408// PbPb 5.5 TeV
409// MRST HO
410// mc = 1.4 GeV, pt-kick 1 GeV
411//
1af7144e 412 Double_t c[5] = {-6.03425e+02, 4.98257e+02, -1.38794e+02, 1.62209e+01, -6.85955e-01};
413 Double_t x = TMath::Abs(px[0]);
414 Double_t y;
415
416 if (x < 4.) {
417 y = 31.754;
418 } else if (x < 6) {
419 Int_t j;
420 y = c[j = 4];
421 while (j > 0) y = y * x + c[--j];
422 } else {
423 y =0.;
424 }
425
426 return y;
05932df6 427}
428
75e0cc59 429Double_t AliGenMUONlib::YJpsiCDFscaled( const Double_t *px, const Double_t* dummy)
0e137c25 430{
431 // J/Psi y
432 return AliGenMUONlib::YJpsiPbPb(px, dummy);
433}
434
75e0cc59 435Double_t AliGenMUONlib::YJpsiCDFscaledPP( const Double_t *px, const Double_t* dummy)
9ff13849 436{
437 // J/Psi y
438 return AliGenMUONlib::YJpsiPP(px, dummy);
439}
0e137c25 440
75e0cc59 441Double_t AliGenMUONlib::YJpsiCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 442{
9e9ae065 443// J/Psi y
444//
445// pp 10 TeV
446// scaled from YJpsiPP(14 TeV) using 10 TeV / 14 TeV ratio of y-spectra in LO pQCD.
447// see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
448//
449
450 Double_t c[5] = {2.46681e+01, 8.91486e+01, -3.21227e+01, 3.63075e+00, -1.32047e-01};
451
452 Double_t x = TMath::Abs(px[0]);
453 Double_t y;
454
455 if (x < 3.2) {
456 y = 98.523 - 1.3664 * x * x;
457 } else if (x < 7.5) {
458 Int_t j;
459 y = c[j = 4];
460 while (j > 0) y = y * x + c[--j];
461 } else {
462 y =0.;
463 }
464
465 if(y<0) y=0;
466
467 return y;
468}
469
9fd56238 470Double_t AliGenMUONlib::YJpsiCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
471{
472// J/Psi y
473//
474// pp 8.8 TeV
475// rescaling of YJpsiPP(14 TeV) using 8.8 TeV / 14 TeV ratio of y-spectra in LO QCD
476//
477 Double_t c[5] = {3.33882e+02, -1.30980e+02, 2.59082e+01, -3.08935e+00, 1.56375e-01};
478 Double_t x = TMath::Abs(px[0]);
479 Double_t y;
480
481 if (x < 3.7) {
482 y = 99.236 - 1.5498 * x * x;
483 } else if (x < 7.4) {
484 Int_t j;
485 y = c[j = 4];
486 while (j > 0) y = y * x + c[--j];
487 } else {
488 y =0.;
489 }
490
491 if(y<0) y=0;
492
493 return y;
494}
495
496Double_t AliGenMUONlib::YJpsiCDFscaledPP9dummy(Double_t px)
497{
498 return AliGenMUONlib::YJpsiCDFscaledPP9(&px, (Double_t*) 0);
499}
500
501Double_t AliGenMUONlib::YJpsiCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
502{
503// J/Psi y
504//
505// pp 7 TeV
506// scaled from YJpsiPP(14 TeV) using 7 TeV / 14 TeV ratio of y-spectra in LO pQCD.
507//
508
509 Double_t c[5] = {6.71181e+02, -3.69240e+02, 8.89644e+01, -1.04937e+01, 4.80959e-01};
510
511 Double_t x = TMath::Abs(px[0]);
512 Double_t y;
513
514 if (x < 4.0) {
515 y = 100.78 - 1.8353 * x * x;
516 } else if (x < 7.3) {
517 Int_t j;
518 y = c[j = 4];
519 while (j > 0) y = y * x + c[--j];
520 } else {
521 y =0.;
522 }
523
524 if(y<0) y=0;
525
526 return y;
527}
528
529Double_t AliGenMUONlib::YJpsiCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
530{
531// J/Psi y
532//
533// pp 3.94 TeV
534// rescaling of YJpsiPP(14 TeV) using 3.94 TeV / 14 TeV ratio of y-spectra in LO QCD
535//
536 Double_t c[5] = {4.00785e+02, -1.41159e+01, -3.28599e+01, 5.53048e+00, -2.45151e-01};
537 Double_t x = TMath::Abs(px[0]);
538 Double_t y;
539
540 if (x < 5.5) {
541 y = 107.389 - 2.7454 * x * x;
542 } else if (x < 7.0) {
543 Int_t j;
544 y = c[j = 4];
545 while (j > 0) y = y * x + c[--j];
546 } else {
547 y =0.;
548 }
549
550 if(y<0) y=0;
551
552 return y;
553}
554
75e0cc59 555Double_t AliGenMUONlib::YJpsiPP( const Double_t *px, const Double_t */*dummy*/)
bb6e81ac 556{
557
558//
559// J/Psi y
560//
561//
562// R. Vogt 2002
563// pp 14 TeV
564// MRST HO
565// mc = 1.4 GeV, pt-kick 1 GeV
566//
567
568 Double_t c[5] = {1.38532e+00, 1.00596e+02, -3.46378e+01, 3.94172e+00, -1.48319e-01};
569 Double_t x = TMath::Abs(px[0]);
570 Double_t y;
571
572 if (x < 2.5) {
573 y = 96.455 - 0.8483 * x * x;
574 } else if (x < 7.9) {
575 Int_t j;
576 y = c[j = 4];
577 while (j > 0) y = y * x + c[--j];
578 } else {
579 y =0.;
580 }
581
582 return y;
583}
584
9fd56238 585Double_t AliGenMUONlib::YJpsiCDFscaledPPb9( const Double_t *px, const Double_t */*dummy*/)
586{
587// J/Psi y
588//
589// pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.80
590//
591 Double_t c[7] = {7.52296e-01, 2.49917e-02, 3.36500e-03, 1.91187e-03, 2.92154e-04,
592 -4.16509e-05,-7.62709e-06};
593 Double_t y;
594 Double_t x = px[0] + 0.47; // rapidity shift
595 Int_t j;
596 y = c[j = 6];
597 while (j > 0) y = y * x + c[--j];
598 if(y<0) y=0;
599
600 return y * AliGenMUONlib::YJpsiCDFscaledPP9dummy(x);
601}
602
603Double_t AliGenMUONlib::YJpsiCDFscaledPbP9( const Double_t *px, const Double_t */*dummy*/)
604{
605// J/Psi y
606//
607// Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.80
608//
609 Double_t c[7] = {7.52296e-01, 2.49917e-02, 3.36500e-03, 1.91187e-03, 2.92154e-04,
610 -4.16509e-05,-7.62709e-06};
611 Double_t y;
612 Double_t x = -px[0] + 0.47; // rapidity shift
613 Int_t j;
614 y = c[j = 6];
615 while (j > 0) y = y * x + c[--j];
616 if(y<0) y=0;
617
618 return y * AliGenMUONlib::YJpsiCDFscaledPP9dummy(x);
619}
620
621Double_t AliGenMUONlib::YJpsiCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
622{
623// J/Psi y
624//
625// PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.66
626//
627 Double_t c[4] = {5.95228e-01, 9.45069e-03, 2.44710e-04, -1.32894e-05};
628 Double_t x = px[0]*px[0];
629 Double_t y;
630 Int_t j;
631 y = c[j = 3];
632 while (j > 0) y = y * x + c[--j];
633 if(y<0) y=0;
634
635 return y * AliGenMUONlib::YJpsiCDFscaledPP4(px,dummy);
636}
637
75e0cc59 638Double_t AliGenMUONlib::YJpsiBPbPb( const Double_t *px, const Double_t */*dummy*/)
17d28ba5 639{
640
641//
642// J/Psi from B->J/Psi X
643//
644//
645
646
647 Double_t c[7] = {7.37025e-02, 0., -2.94487e-03, 0., 6.07953e-06, 0., 5.39219e-07};
648
649 Double_t x = TMath::Abs(px[0]);
650 Double_t y;
651
652 if (x > 6.) {
653 y = 0.;
654 } else {
655 Int_t j;
656 y = c[j = 6];
657 while (j > 0) y = y * x + c[--j];
658 }
659
660 return y;
661}
662
663
664
fe4da5cc 665// particle composition
666//
65fb704d 667Int_t AliGenMUONlib::IpJpsi(TRandom *)
fe4da5cc 668{
d90f80fd 669// J/Psi composition
88cb7938 670 return 443;
fe4da5cc 671}
88e5db43 672Int_t AliGenMUONlib::IpPsiP(TRandom *)
673{
674// Psi prime composition
675 return 100443;
676}
0ad09590 677Int_t AliGenMUONlib::IpJpsiFamily(TRandom *)
678{
679// J/Psi composition
680 Int_t ip;
681 Float_t r = gRandom->Rndm();
682 if (r < 0.98) {
683 ip = 443;
684 } else {
685 ip = 100443;
686 }
687 return ip;
688}
689
fe4da5cc 690
88e5db43 691
fe4da5cc 692// Upsilon
693//
694//
695// pt-distribution
696//____________________________________________________________
75e0cc59 697Double_t AliGenMUONlib::PtUpsilon( const Double_t *px, const Double_t */*dummy*/ )
fe4da5cc 698{
d90f80fd 699// Upsilon pT
700 const Double_t kpt0 = 5.3;
701 const Double_t kxn = 2.5;
fe4da5cc 702 Double_t x=*px;
703 //
d90f80fd 704 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
705 return x/TMath::Power(pass1,kxn);
fe4da5cc 706}
05932df6 707
75e0cc59 708Double_t AliGenMUONlib::PtUpsilonCDFscaled( const Double_t *px, const Double_t */*dummy*/ )
0e137c25 709{
710// Upsilon pT
48416d65 711 const Double_t kpt0 = 7.753;
712 const Double_t kxn = 3.042;
0e137c25 713 Double_t x=*px;
714 //
715 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
716 return x/TMath::Power(pass1,kxn);
717}
718
75e0cc59 719Double_t AliGenMUONlib::PtUpsilonCDFscaledPP( const Double_t *px, const Double_t */*dummy*/ )
9ff13849 720{
721// Upsilon pT
722//
723// pp 14 TeV
724//
725// scaled from CDF data at 2 TeV
726
727 const Double_t kpt0 = 8.610;
728 const Double_t kxn = 3.051;
729 Double_t x=*px;
730 //
731 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
732 return x/TMath::Power(pass1,kxn);
733}
734
b33adf51 735Double_t AliGenMUONlib::PtUpsilonCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
736{
737// Upsilon pT
738//
739// pp 10 TeV
740//
741// scaled from CDF data at 2 TeV
742
743 const Double_t kpt0 = 8.235;
744 const Double_t kxn = 3.051;
745 Double_t x=*px;
746 //
747 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
748 return x/TMath::Power(pass1,kxn);
749}
750
9fd56238 751Double_t AliGenMUONlib::PtUpsilonCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
752{
753// Upsilon pT
754//
755// pp 8.8 TeV
756// scaled from CDF data at 2 TeV
757//
758 const Double_t kpt0 = 8.048;
759 const Double_t kxn = 3.051;
760 Double_t x=*px;
761 //
762 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
763 return x/TMath::Power(pass1,kxn);
764}
765
766Double_t AliGenMUONlib::PtUpsilonCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
767{
768// Upsilon pT
769//
770// pp 7 TeV
771//
772// scaled from CDF data at 2 TeV
773
774 const Double_t kpt0 = 7.817;
775 const Double_t kxn = 3.051;
776 Double_t x=*px;
777 //
778 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
779 return x/TMath::Power(pass1,kxn);
780}
781
782Double_t AliGenMUONlib::PtUpsilonCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
783{
784// Upsilon pT
785//
786// pp 3.94 TeV
787// scaled from CDF data at 2 TeV
788//
789 const Double_t kpt0 = 7.189;
790 const Double_t kxn = 3.051;
791 Double_t x=*px;
792 //
793 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
794 return x/TMath::Power(pass1,kxn);
795}
796
797Double_t AliGenMUONlib::PtUpsilonCDFscaledPPb9( const Double_t *px, const Double_t *dummy)
798{
799// Upsilon pT
800//
801// pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.90
802//
803 Double_t c[5] = {7.64952e-01, 1.12501e-04, 4.96038e-04, -3.03198e-05, 3.74035e-06};
804 Double_t x=*px;
805 Double_t y;
806 Int_t j;
807 y = c[j = 4];
808 while (j > 0) y = y * x + c[--j];
809 //
810 Double_t d = 1.+c[4]*TMath::Power(x,4);
811 return y/d * AliGenMUONlib::PtUpsilonCDFscaledPP9(px,dummy);
812}
813
814Double_t AliGenMUONlib::PtUpsilonCDFscaledPbP9( const Double_t *px, const Double_t *dummy)
815{
816// Upsilon pT
817//
818// Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.90
819//
820 Double_t c[5] = {1.09881e+00, 3.08329e-03, -2.00356e-04, 8.28991e-06, 2.52576e-06};
821 Double_t x=*px;
822 Double_t y;
823 Int_t j;
824 y = c[j = 4];
825 while (j > 0) y = y * x + c[--j];
826 //
827 Double_t d = 1.+c[4]*TMath::Power(x,4);
828 return y/d * AliGenMUONlib::PtUpsilonCDFscaledPP9(px,dummy);
829}
830
831Double_t AliGenMUONlib::PtUpsilonCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
832{
833// Upsilon pT
834//
835// PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.85
836//
837 Double_t c[5] = {8.65872e-01, 2.05465e-03, 2.56063e-04, -1.65598e-05, 2.29209e-06};
838 Double_t x=*px;
839 Double_t y;
840 Int_t j;
841 y = c[j = 4];
842 while (j > 0) y = y * x + c[--j];
843 //
844 Double_t d = 1.+c[4]*TMath::Power(x,4);
845 return y/d * AliGenMUONlib::PtUpsilonCDFscaledPP4(px,dummy);
846}
847
75e0cc59 848Double_t AliGenMUONlib::PtUpsilonFlat( const Double_t */*px*/, const Double_t */*dummy*/ )
4ca8d070 849{
850 return 1.;
851}
852
75e0cc59 853Double_t AliGenMUONlib::PtUpsilonPbPb( const Double_t *px, const Double_t */*dummy*/)
05932df6 854{
05932df6 855//
856// Upsilon pT
857//
858//
859// R. Vogt 2002
860// PbPb 5.5 TeV
861// MRST HO
862// mc = 1.4 GeV, pt-kick 1 GeV
863//
1af7144e 864 Float_t x = px[0];
865 Double_t c[8] = {
866 -1.03488e+01, 1.28065e+01, -6.60500e+00, 1.66140e+00,
867 -2.34293e-01, 1.86925e-02, -7.80708e-04, 1.30610e-05
868 };
3d905dd7 869 Double_t y;
870 if (x < 10.) {
871 Int_t j;
872 y = c[j = 7];
873 while (j > 0) y = y * x +c[--j];
874 y = x * TMath::Exp(y);
875 } else {
876 y = 0.;
877 }
1af7144e 878 return y;
05932df6 879}
880
75e0cc59 881Double_t AliGenMUONlib::PtUpsilonPP( const Double_t *px, const Double_t */*dummy*/)
bb6e81ac 882{
bb6e81ac 883//
884// Upsilon pT
885//
886//
887// R. Vogt 2002
888// pp 14 TeV
889// MRST HO
890// mc = 1.4 GeV, pt-kick 1 GeV
891//
892 Float_t x = px[0];
893 Double_t c[8] = {-7.93955e+00, 1.06306e+01, -5.21392e+00, 1.19703e+00,
894 -1.45718e-01, 8.95151e-03, -2.04806e-04, -1.13053e-06};
895
896 Double_t y;
897 if (x < 10.) {
898 Int_t j;
899 y = c[j = 7];
900 while (j > 0) y = y * x +c[--j];
901 y = x * TMath::Exp(y);
902 } else {
903 y = 0.;
904 }
905 return y;
906}
907
fe4da5cc 908//
909// y-distribution
910//
911//____________________________________________________________
75e0cc59 912Double_t AliGenMUONlib::YUpsilon(const Double_t *py, const Double_t */*dummy*/)
fe4da5cc 913{
d90f80fd 914// Upsilon y
915 const Double_t ky0 = 3.;
916 const Double_t kb=1.;
fe4da5cc 917 Double_t yu;
918 Double_t y=TMath::Abs(*py);
919 //
d90f80fd 920 if (y < ky0)
921 yu=kb;
fe4da5cc 922 else
d90f80fd 923 yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
fe4da5cc 924 return yu;
925}
05932df6 926
927
75e0cc59 928Double_t AliGenMUONlib::YUpsilonPbPb( const Double_t *px, const Double_t */*dummy*/)
05932df6 929{
930
931//
932// Upsilon y
933//
934//
935// R. Vogt 2002
936// PbPb 5.5 TeV
937// MRST HO
938// mc = 1.4 GeV, pt-kick 1 GeV
939//
940
1af7144e 941 Double_t c[7] = {3.40036e-01, -3.98882e-07, -4.48398e-03, 8.46411e-08, -6.10854e-04,
942 -2.99753e-09, 1.28895e-05};
b33adf51 943 Double_t x = TMath::Abs(px[0]);
944 if (x > 5.55) return 0.;
1af7144e 945 Int_t j;
946 Double_t y = c[j = 6];
947 while (j > 0) y = y * x +c[--j];
948 return y;
05932df6 949}
950
75e0cc59 951Double_t AliGenMUONlib::YUpsilonCDFscaled( const Double_t *px, const Double_t *dummy)
0e137c25 952{
953 // Upsilon y
954 return AliGenMUONlib::YUpsilonPbPb(px, dummy);
955
4ca8d070 956}
9ff13849 957
75e0cc59 958Double_t AliGenMUONlib::YUpsilonCDFscaledPP( const Double_t *px, const Double_t *dummy)
9ff13849 959{
960 // Upsilon y
961 return AliGenMUONlib::YUpsilonPP(px, dummy);
962
963}
964
75e0cc59 965Double_t AliGenMUONlib::YUpsilonFlat( const Double_t */*px*/, const Double_t */*dummy*/)
4ca8d070 966{
967 // Upsilon y
968 return 1.;
969
0e137c25 970}
971
b33adf51 972Double_t AliGenMUONlib::YUpsilonCDFscaledPP10( const Double_t *px, const Double_t */*dummy*/)
973{
b33adf51 974// Upsilon y
975//
976// pp 10 TeV
977// scaled from YUpsilonPP(14 TeV) using 10 TeV / 14 TeV ratio of y-spectra in LO pQCD.
978// see S.Grigoryan, PWG3 Meeting, 27th Oct 2008
979//
9fd56238 980 Double_t c[4] = {1., -2.17877e-02, -6.52830e-04, 1.40578e-05};
981 Double_t x = TMath::Abs(px[0]);
982 if (x > 6.1) return 0.;
983 Int_t j;
984 Double_t y = c[j = 3];
985 while (j > 0) y = y * x*x +c[--j];
986 return y;
987}
988
989Double_t AliGenMUONlib::YUpsilonCDFscaledPP9( const Double_t *px, const Double_t */*dummy*/)
990{
991// Upsilon y
992//
993// pp 8.8 TeV
994// rescaling of YUpsilonPP(14 TeV) using 8.8 TeV / 14 TeV ratio of y-spectra in LO QCD
995//
996 Double_t c[4] = {1., -2.37621e-02, -6.29610e-04, 1.47976e-05};
b33adf51 997 Double_t x = TMath::Abs(px[0]);
998 if (x > 6.1) return 0.;
999 Int_t j;
1000 Double_t y = c[j = 3];
1001 while (j > 0) y = y * x*x +c[--j];
1002 return y;
1003}
1004
9fd56238 1005Double_t AliGenMUONlib::YUpsilonCDFscaledPP9dummy(Double_t px)
1006{
1007 return AliGenMUONlib::YUpsilonCDFscaledPP9(&px, (Double_t*) 0);
1008}
1009
1010Double_t AliGenMUONlib::YUpsilonCDFscaledPP7( const Double_t *px, const Double_t */*dummy*/)
1011{
1012// Upsilon y
1013//
1014// pp 7 TeV
1015// scaled from YUpsilonPP(14 TeV) using 7 TeV / 14 TeV ratio of y-spectra in LO pQCD.
1016//
1017 Double_t c[4] = {1., -2.61009e-02, -6.83937e-04, 1.78451e-05};
1018 Double_t x = TMath::Abs(px[0]);
1019 if (x > 6.0) return 0.;
1020 Int_t j;
1021 Double_t y = c[j = 3];
1022 while (j > 0) y = y * x*x +c[--j];
1023 return y;
1024}
1025
1026Double_t AliGenMUONlib::YUpsilonCDFscaledPP4( const Double_t *px, const Double_t */*dummy*/)
1027{
1028// Upsilon y
1029//
1030// pp 3.94 TeV
1031// rescaling of YUpsilonPP(14 TeV) using 3.94 TeV / 14 TeV ratio of y-spectra in LO QCD
1032//
1033 Double_t c[4] = {1., -3.91924e-02, -4.26184e-04, 2.10914e-05};
1034 Double_t x = TMath::Abs(px[0]);
1035 if (x > 5.7) return 0.;
1036 Int_t j;
1037 Double_t y = c[j = 3];
1038 while (j > 0) y = y * x*x +c[--j];
1039
1040 return y;
1041}
1042
75e0cc59 1043Double_t AliGenMUONlib::YUpsilonPP( const Double_t *px, const Double_t */*dummy*/)
bb6e81ac 1044{
1045
1046//
1047// Upsilon y
1048//
1049//
1050// R. Vogt 2002
1051// p p 14. TeV
1052// MRST HO
1053// mc = 1.4 GeV, pt-kick 1 GeV
1054//
1055 Double_t c[7] = {8.91936e-01, -6.46645e-07, -1.52774e-02, 4.28677e-08, -7.01517e-04,
1056 -6.20539e-10, 1.29943e-05};
b33adf51 1057 Double_t x = TMath::Abs(px[0]);
1058 if (x > 6.2) return 0.;
bb6e81ac 1059 Int_t j;
1060 Double_t y = c[j = 6];
1061 while (j > 0) y = y * x +c[--j];
1062 return y;
1063}
1064
9fd56238 1065Double_t AliGenMUONlib::YUpsilonCDFscaledPPb9( const Double_t *px, const Double_t */*dummy*/)
1066{
1067// Upsilon y
1068//
1069// pPb 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.90
1070//
1071 Double_t c[7] = {8.71829e-01, 4.77467e-02, 8.09671e-03, 6.45294e-04, -2.15730e-04,
1072 -4.67538e-05,-2.11683e-06};
1073 Double_t y;
1074 Double_t x = px[0] + 0.47; // rapidity shift
1075 Int_t j;
1076 y = c[j = 6];
1077 while (j > 0) y = y * x + c[--j];
1078 if(y<0) y=0;
1079
1080 return y * AliGenMUONlib::YUpsilonCDFscaledPP9dummy(x);
1081}
1082
1083Double_t AliGenMUONlib::YUpsilonCDFscaledPbP9( const Double_t *px, const Double_t */*dummy*/)
1084{
1085// Upsilon y
1086//
1087// Pbp 8.8 TeV, for EKS98 with minimum bias shadowing factor 0.90
1088//
1089 Double_t c[7] = {8.71829e-01, 4.77467e-02, 8.09671e-03, 6.45294e-04, -2.15730e-04,
1090 -4.67538e-05,-2.11683e-06};
1091 Double_t y;
1092 Double_t x = -px[0] + 0.47; // rapidity shift
1093 Int_t j;
1094 y = c[j = 6];
1095 while (j > 0) y = y * x + c[--j];
1096 if(y<0) y=0;
1097
1098 return y * AliGenMUONlib::YUpsilonCDFscaledPP9dummy(x);
1099}
1100
1101Double_t AliGenMUONlib::YUpsilonCDFscaledPbPb4( const Double_t *px, const Double_t *dummy)
1102{
1103// Upsilon y
1104//
1105// PbPb 3.94 TeV, for EKS98 with minimum bias shadowing factor 0.85
1106//
1107 Double_t c[4] = {8.27837e-01, 1.70115e-02, -1.26046e-03, 1.52091e-05};
1108 Double_t x = px[0]*px[0];
1109 Double_t y;
1110 Int_t j;
1111 y = c[j = 3];
1112 while (j > 0) y = y * x + c[--j];
1113 if(y<0) y=0;
1114
1115 return y * AliGenMUONlib::YUpsilonCDFscaledPP4(px,dummy);
1116}
1117
1118
fe4da5cc 1119// particle composition
1120//
65fb704d 1121Int_t AliGenMUONlib::IpUpsilon(TRandom *)
fe4da5cc 1122{
d90f80fd 1123// y composition
88cb7938 1124 return 553;
fe4da5cc 1125}
88e5db43 1126Int_t AliGenMUONlib::IpUpsilonP(TRandom *)
1127{
1128// y composition
1129 return 100553;
1130}
1131Int_t AliGenMUONlib::IpUpsilonPP(TRandom *)
1132{
1133// y composition
1134 return 200553;
1135}
0ad09590 1136Int_t AliGenMUONlib::IpUpsilonFamily(TRandom *)
1137{
1138// y composition
1139 Int_t ip;
1140 Float_t r = gRandom->Rndm();
1141
1142 if (r < 0.712) {
1143 ip = 553;
1144 } else if (r < 0.896) {
1145 ip = 100553;
1146 } else {
1147 ip = 200553;
1148 }
1149 return ip;
1150}
fe4da5cc 1151
88e5db43 1152
fe4da5cc 1153//
1154// Phi
1155//
1156//
1157// pt-distribution (by scaling of pion distribution)
1158//____________________________________________________________
75e0cc59 1159Double_t AliGenMUONlib::PtPhi( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 1160{
d90f80fd 1161// Phi pT
0db2f441 1162 return PtScal(*px,6);
fe4da5cc 1163}
1164// y-distribution
75e0cc59 1165Double_t AliGenMUONlib::YPhi( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 1166{
d90f80fd 1167// Phi y
1168 Double_t *dum=0;
1169 return YJpsi(px,dum);
fe4da5cc 1170}
1171// particle composition
1172//
65fb704d 1173Int_t AliGenMUONlib::IpPhi(TRandom *)
fe4da5cc 1174{
d90f80fd 1175// Phi composition
89512a3b 1176 return 333;
1177}
1178
1179//
1180// omega
1181//
1182//
1183// pt-distribution (by scaling of pion distribution)
1184//____________________________________________________________
75e0cc59 1185Double_t AliGenMUONlib::PtOmega( const Double_t *px, const Double_t */*dummy*/)
89512a3b 1186{
1187// Omega pT
0db2f441 1188 return PtScal(*px,4);
89512a3b 1189}
1190// y-distribution
75e0cc59 1191Double_t AliGenMUONlib::YOmega( const Double_t *px, const Double_t */*dummy*/)
89512a3b 1192{
1193// Omega y
1194 Double_t *dum=0;
1195 return YJpsi(px,dum);
1196}
1197// particle composition
1198//
1199Int_t AliGenMUONlib::IpOmega(TRandom *)
1200{
1201// Omega composition
1202 return 223;
1203}
1204
1205
1206//
1207// Eta
1208//
1209//
1210// pt-distribution (by scaling of pion distribution)
1211//____________________________________________________________
75e0cc59 1212Double_t AliGenMUONlib::PtEta( const Double_t *px, const Double_t */*dummy*/)
89512a3b 1213{
1214// Eta pT
0db2f441 1215 return PtScal(*px,2);
89512a3b 1216}
1217// y-distribution
75e0cc59 1218Double_t AliGenMUONlib::YEta( const Double_t *px, const Double_t */*dummy*/)
89512a3b 1219{
1220// Eta y
1221 Double_t *dum=0;
1222 return YJpsi(px,dum);
1223}
1224// particle composition
1225//
1226Int_t AliGenMUONlib::IpEta(TRandom *)
1227{
1228// Eta composition
1229 return 221;
fe4da5cc 1230}
1231
1232//
1233// Charm
1234//
1235//
1236// pt-distribution
1237//____________________________________________________________
75e0cc59 1238Double_t AliGenMUONlib::PtCharm( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 1239{
d90f80fd 1240// Charm pT
bd0276a8 1241 const Double_t kpt0 = 2.25;
1242 const Double_t kxn = 3.17;
fe4da5cc 1243 Double_t x=*px;
1244 //
bd0276a8 1245 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1246 return x/TMath::Power(pass1,kxn);
1247}
1248
75e0cc59 1249Double_t AliGenMUONlib::PtCharmCentral( const Double_t *px, const Double_t */*dummy*/)
bd0276a8 1250{
1251// Charm pT
1252 const Double_t kpt0 = 2.12;
1253 const Double_t kxn = 2.78;
bd0276a8 1254 Double_t x=*px;
1255 //
1256 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
d90f80fd 1257 return x/TMath::Power(pass1,kxn);
fe4da5cc 1258}
75e0cc59 1259Double_t AliGenMUONlib::PtCharmF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1260{
1261// FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1262// PtCharmFiMjSkPP = PtCharmF0M0S0PP * (dN(i,j,k)/dpt / dN(0,0,0)/dpt)_MNR
1263// i=0,1,2; j=0,1,2; k=0,1,...,6
1264// dN(i,j,k)/dpt - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1265// http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1266// calculations for the following inputs:
1267// Peterson fragmentation function (F) with \epsilon_c = 0.02, 0.002 & 0.11
1268// for i=0,1 & 2 respectively; quark mass (M) of 1.5, 1.3 & 1.7 GeV
1269// for j=0,1 & 2 respectively;
1270// factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1271// with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2
1272// for k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1273// (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1274// June 2008, Smbat.Grigoryan@cern.ch
1275
1276// Charm pT
1277// Pythia6.214 (kCharmppMNRwmi, PDF = CTEQ5L, quark mass = 1.2 GeV, PtHard > 2.76 GeV/c)
1278// for pp collisions at 14 TeV with one c-cbar pair per event.
1279// Corresponding NLO total cross section is 5.68 mb
1280
1281
1282 const Double_t kpt0 = 2.2930;
1283 const Double_t kxn = 3.1196;
1284 Double_t c[3]={-5.2180e-01,1.8753e-01,2.8669e-02};
1285 Double_t x=*px;
1286 //
1287 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1288 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1289}
75e0cc59 1290Double_t AliGenMUONlib::PtCharmF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1291{
1292// Charm pT
1293// Corresponding NLO total cross section is 6.06 mb
1294 const Double_t kpt0 = 2.8669;
1295 const Double_t kxn = 3.1044;
1296 Double_t c[3]={-4.6714e-01,1.5005e-01,4.5003e-02};
1297 Double_t x=*px;
1298 //
1299 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1300 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1301}
75e0cc59 1302Double_t AliGenMUONlib::PtCharmF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1303{
1304// Charm pT
1305// Corresponding NLO total cross section is 6.06 mb
1306 const Double_t kpt0 = 1.8361;
1307 const Double_t kxn = 3.2966;
1308 Double_t c[3]={-6.1550e-01,2.6498e-01,1.0728e-02};
1309 Double_t x=*px;
1310 //
1311 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1312 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1313}
75e0cc59 1314Double_t AliGenMUONlib::PtCharmF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1315{
1316// Charm pT
1317// Corresponding NLO total cross section is 7.69 mb
1318 const Double_t kpt0 = 2.1280;
1319 const Double_t kxn = 3.1397;
1320 Double_t c[3]={-5.4021e-01,2.0944e-01,2.5211e-02};
1321 Double_t x=*px;
1322 //
1323 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1324 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1325}
75e0cc59 1326Double_t AliGenMUONlib::PtCharmF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1327{
1328// Charm pT
1329// Corresponding NLO total cross section is 4.81 mb
1330 const Double_t kpt0 = 2.4579;
1331 const Double_t kxn = 3.1095;
1332 Double_t c[3]={-5.1497e-01,1.7532e-01,3.2429e-02};
1333 Double_t x=*px;
1334 //
1335 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1336 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1337}
75e0cc59 1338Double_t AliGenMUONlib::PtCharmF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1339{
1340// Charm pT
1341// Corresponding NLO total cross section is 14.09 mb
1342 const Double_t kpt0 = 2.1272;
1343 const Double_t kxn = 3.1904;
1344 Double_t c[3]={-4.6088e-01,2.1918e-01,2.3055e-02};
1345 Double_t x=*px;
1346 //
1347 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1348 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1349}
75e0cc59 1350Double_t AliGenMUONlib::PtCharmF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1351{
1352// Charm pT
1353// Corresponding NLO total cross section is 1.52 mb
1354 const Double_t kpt0 = 2.8159;
1355 const Double_t kxn = 3.0857;
1356 Double_t c[3]={-6.4691e-01,2.0289e-01,2.4922e-02};
1357 Double_t x=*px;
1358 //
1359 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1360 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1361}
75e0cc59 1362Double_t AliGenMUONlib::PtCharmF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1363{
1364// Charm pT
1365// Corresponding NLO total cross section is 3.67 mb
1366 const Double_t kpt0 = 2.7297;
1367 const Double_t kxn = 3.3019;
1368 Double_t c[3]={-6.2216e-01,1.9031e-01,1.5341e-02};
1369 Double_t x=*px;
1370 //
1371 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1372 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1373}
75e0cc59 1374Double_t AliGenMUONlib::PtCharmF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1375{
1376// Charm pT
1377// Corresponding NLO total cross section is 3.38 mb
1378 const Double_t kpt0 = 2.3894;
1379 const Double_t kxn = 3.1075;
1380 Double_t c[3]={-4.9742e-01,1.7032e-01,2.5994e-02};
1381 Double_t x=*px;
1382 //
1383 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1384 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1385}
75e0cc59 1386Double_t AliGenMUONlib::PtCharmF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1387{
1388// Charm pT
1389// Corresponding NLO total cross section is 10.37 mb
1390 const Double_t kpt0 = 2.0187;
1391 const Double_t kxn = 3.3011;
1392 Double_t c[3]={-3.9869e-01,2.9248e-01,1.1763e-02};
1393 Double_t x=*px;
1394 //
1395 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1396 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1397}
75e0cc59 1398Double_t AliGenMUONlib::PtCharmF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1399{
1400// Charm pT
1401// Corresponding NLO total cross section is 7.22 mb
1402 const Double_t kpt0 = 2.1089;
1403 const Double_t kxn = 3.1848;
1404 Double_t c[3]={-4.6275e-01,1.8114e-01,2.1363e-02};
1405 Double_t x=*px;
1406 //
1407 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1408 return x/TMath::Power(pass1,kxn)*(1.+c[0]*x+c[1]*x*x)/(1.+c[2]*x*x);
1409}
1410
fe4da5cc 1411// y-distribution
75e0cc59 1412Double_t AliGenMUONlib::YCharm( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 1413{
ec772ba2 1414// Charm y :: Carrer & Dainese : ALICE-INT-2003-019 v.3 (hep-ph/0311225)
1415// Pythia tuned to reproduce the distribution given by the HVQMNR program based on NLO calculations (pQCD)
1416// shadowing + kt broadening
1417
1418 Double_t x=px[0];
1419 Double_t c[2]={-2.42985e-03,-2.31001e-04};
1420 Double_t y=1+(c[0]*TMath::Power(x,2))+(c[1]*TMath::Power(x,4));
1421 Double_t ycharm;
1422
1423 if (TMath::Abs(x)>8) {
1424 ycharm=0.;
1425 }
1426 else {
1427 ycharm=TMath::Power(y,3);
1428 }
1429
1430 return ycharm;
fe4da5cc 1431}
75e0cc59 1432Double_t AliGenMUONlib::YCharmF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1433{
1434// FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1435// YCharmFiMjSkPP = YCharmF0M0S0PP * (dN(i,j,k)/dy / dN(0,0,0)/dy)_MNR
1436// i=0,1,2; j=0,1,2; k=0,1,...,6
1437// dN(i,j,k)/dy - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1438// http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1439// calculations for the following inputs:
1440// Peterson fragmentation function (F) with \epsilon_c = 0.02, 0.002 & 0.11
1441// for i=0,1 & 2 respectively; quark mass (M) of 1.5, 1.3 & 1.7 GeV
1442// for j=0,1 & 2 respectively;
1443// factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1444// with a/b = 1/1,1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2 for
1445// k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1446// (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1447// June 2008, Smbat.Grigoryan@cern.ch
1448
1449// Charm y
1450// Pythia6.214 (kCharmppMNRwmi, PDF = CTEQ5L, quark mass = 1.2 GeV, PtHard > 2.76 GeV/c)
1451// for pp collisions at 14 TeV with one c-cbar pair per event.
1452// Corresponding NLO total cross section is 5.68 mb
1453
1454 Double_t x=px[0];
1455 Double_t c[2]={7.0909e-03,6.1967e-05};
1456 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1457 Double_t ycharm;
1458
1459 if (TMath::Abs(x)>9) {
1460 ycharm=0.;
1461 }
1462 else {
1463 ycharm=TMath::Power(y,3);
1464 }
1465
1466 return ycharm;
1467}
75e0cc59 1468Double_t AliGenMUONlib::YCharmF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1469{
1470// Charm y
1471// Corresponding NLO total cross section is 6.06 mb
1472 Double_t x=px[0];
1473 Double_t c[2]={6.9707e-03,6.0971e-05};
1474 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1475 Double_t ycharm;
1476
1477 if (TMath::Abs(x)>9) {
1478 ycharm=0.;
1479 }
1480 else {
1481 ycharm=TMath::Power(y,3);
1482 }
1483
1484 return ycharm;
1485}
75e0cc59 1486Double_t AliGenMUONlib::YCharmF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1487{
1488// Charm y
1489// Corresponding NLO total cross section is 6.06 mb
1490 Double_t x=px[0];
1491 Double_t c[2]={7.1687e-03,6.5303e-05};
1492 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1493 Double_t ycharm;
1494
1495 if (TMath::Abs(x)>9) {
1496 ycharm=0.;
1497 }
1498 else {
1499 ycharm=TMath::Power(y,3);
1500 }
1501
1502 return ycharm;
1503}
75e0cc59 1504Double_t AliGenMUONlib::YCharmF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1505{
1506// Charm y
1507// Corresponding NLO total cross section is 7.69 mb
1508 Double_t x=px[0];
1509 Double_t c[2]={5.9090e-03,7.1854e-05};
1510 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1511 Double_t ycharm;
1512
1513 if (TMath::Abs(x)>9) {
1514 ycharm=0.;
1515 }
1516 else {
1517 ycharm=TMath::Power(y,3);
1518 }
1519
1520 return ycharm;
1521}
75e0cc59 1522Double_t AliGenMUONlib::YCharmF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1523{
1524// Charm y
1525// Corresponding NLO total cross section is 4.81 mb
1526 Double_t x=px[0];
1527 Double_t c[2]={8.0882e-03,5.5872e-05};
1528 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1529 Double_t ycharm;
1530
1531 if (TMath::Abs(x)>9) {
1532 ycharm=0.;
1533 }
1534 else {
1535 ycharm=TMath::Power(y,3);
1536 }
1537
1538 return ycharm;
1539}
75e0cc59 1540Double_t AliGenMUONlib::YCharmF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1541{
1542// Charm y
1543// Corresponding NLO total cross section is 14.09 mb
1544 Double_t x=px[0];
1545 Double_t c[2]={7.2520e-03,6.2691e-05};
1546 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1547 Double_t ycharm;
1548
1549 if (TMath::Abs(x)>9) {
1550 ycharm=0.;
1551 }
1552 else {
1553 ycharm=TMath::Power(y,3);
1554 }
1555
1556 return ycharm;
1557}
75e0cc59 1558Double_t AliGenMUONlib::YCharmF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1559{
1560// Charm y
1561// Corresponding NLO total cross section is 1.52 mb
1562 Double_t x=px[0];
1563 Double_t c[2]={1.1040e-04,1.4498e-04};
1564 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1565 Double_t ycharm;
1566
1567 if (TMath::Abs(x)>9) {
1568 ycharm=0.;
1569 }
1570 else {
1571 ycharm=TMath::Power(y,3);
1572 }
1573
1574 return ycharm;
1575}
75e0cc59 1576Double_t AliGenMUONlib::YCharmF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1577{
1578// Charm y
1579// Corresponding NLO total cross section is 3.67 mb
1580 Double_t x=px[0];
1581 Double_t c[2]={-3.1328e-03,1.8270e-04};
1582 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1583 Double_t ycharm;
1584
1585 if (TMath::Abs(x)>9) {
1586 ycharm=0.;
1587 }
1588 else {
1589 ycharm=TMath::Power(y,3);
1590 }
1591
1592 return ycharm;
1593}
75e0cc59 1594Double_t AliGenMUONlib::YCharmF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1595{
1596// Charm y
1597// Corresponding NLO total cross section is 3.38 mb
1598 Double_t x=px[0];
1599 Double_t c[2]={7.0865e-03,6.2532e-05};
1600 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1601 Double_t ycharm;
1602
1603 if (TMath::Abs(x)>9) {
1604 ycharm=0.;
1605 }
1606 else {
1607 ycharm=TMath::Power(y,3);
1608 }
1609
1610 return ycharm;
1611}
75e0cc59 1612Double_t AliGenMUONlib::YCharmF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1613{
1614// Charm y
1615// Corresponding NLO total cross section is 10.37 mb
1616 Double_t x=px[0];
1617 Double_t c[2]={7.7070e-03,5.3533e-05};
1618 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1619 Double_t ycharm;
1620
1621 if (TMath::Abs(x)>9) {
1622 ycharm=0.;
1623 }
1624 else {
1625 ycharm=TMath::Power(y,3);
1626 }
1627
1628 return ycharm;
1629}
75e0cc59 1630Double_t AliGenMUONlib::YCharmF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1631{
1632// Charm y
1633// Corresponding NLO total cross section is 7.22 mb
1634 Double_t x=px[0];
1635 Double_t c[2]={7.9195e-03,5.3823e-05};
1636 Double_t y=1-(c[0]*TMath::Power(x,2))-(c[1]*TMath::Power(x,4));
1637 Double_t ycharm;
1638
1639 if (TMath::Abs(x)>9) {
1640 ycharm=0.;
1641 }
1642 else {
1643 ycharm=TMath::Power(y,3);
1644 }
1645
1646 return ycharm;
1647}
fe4da5cc 1648
ec772ba2 1649
65fb704d 1650Int_t AliGenMUONlib::IpCharm(TRandom *ran)
fe4da5cc 1651{
d90f80fd 1652// Charm composition
65fb704d 1653 Float_t random;
fe4da5cc 1654 Int_t ip;
1655// 411,421,431,4122
65fb704d 1656 random = ran->Rndm();
ec772ba2 1657// Taux de production Carrer & Dainese : ALICE-INT-2003-019 v.3
1658// >>>>> cf. tab 4 p 11
1659
1660 if (random < 0.30) {
1661 ip=421;
1662 } else if (random < 0.60) {
1663 ip=-421;
1664 } else if (random < 0.70) {
1665 ip=411;
1666 } else if (random < 0.80) {
1667 ip=-411;
1668 } else if (random < 0.86) {
1669 ip=431;
1670 } else if (random < 0.92) {
1671 ip=-431;
1672 } else if (random < 0.96) {
1673 ip=4122;
fe4da5cc 1674 } else {
ec772ba2 1675 ip=-4122;
fe4da5cc 1676 }
fe4da5cc 1677
1678 return ip;
1679}
1680
fe4da5cc 1681//
1682// Beauty
1683//
1684//
1685// pt-distribution
1686//____________________________________________________________
75e0cc59 1687Double_t AliGenMUONlib::PtBeauty( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 1688{
d90f80fd 1689// Beauty pT
bd0276a8 1690 const Double_t kpt0 = 6.53;
1691 const Double_t kxn = 3.59;
1692 Double_t x=*px;
1693 //
1694 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1695 return x/TMath::Power(pass1,kxn);
1696}
1697
75e0cc59 1698Double_t AliGenMUONlib::PtBeautyCentral( const Double_t *px, const Double_t */*dummy*/)
bd0276a8 1699{
1700// Beauty pT
1701 const Double_t kpt0 = 6.14;
1702 const Double_t kxn = 2.93;
fe4da5cc 1703 Double_t x=*px;
1704 //
d90f80fd 1705 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1706 return x/TMath::Power(pass1,kxn);
fe4da5cc 1707}
75e0cc59 1708Double_t AliGenMUONlib::PtBeautyF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1709{
1710// FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1711// PtBeautyFiMjSkPP = PtBeautyF0M0S0PP * (dN(i,j,k)/dpt / dN(0,0,0)/dpt)_MNR
1712// i=0,1,2; j=0,1,2; k=0,1,...,6
1713// dN(i,j,k)/dpt - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1714// http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1715// calculations for the following inputs:
1716// Peterson fragmentation function (F) with \epsilon_b = 0.001, 0.0002 & 0.004
1717// for i=0,1 & 2 respectively; quark mass (M) of 4.75, 4.5 & 5.0 GeV
1718// for j=0,1 & 2 respectively;
1719// factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1720// with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2 for
1721// k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1722// (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1723// June 2008, Smbat.Grigoryan@cern.ch
1724
1725// Beauty pT
1726// Pythia6.214 (kBeautyppMNRwmi, PDF = CTEQ5L, quark mass = 4.75 GeV, PtHard > 2.76 GeV/c)
1727// for pp collisions at 14 TeV with one b-bbar pair per event.
1728// Corresponding NLO total cross section is 0.494 mb
1729
1730 const Double_t kpt0 = 8.0575;
1731 const Double_t kxn = 3.1921;
1732 Double_t x=*px;
1733 //
1734 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1735 return x/TMath::Power(pass1,kxn);
1736}
75e0cc59 1737Double_t AliGenMUONlib::PtBeautyF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1738{
1739// Beauty pT
1740// Corresponding NLO total cross section is 0.445 mb
1741 const Double_t kpt0 = 8.6239;
1742 const Double_t kxn = 3.2911;
1743 Double_t x=*px;
1744 //
1745 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1746 return x/TMath::Power(pass1,kxn);
1747}
75e0cc59 1748Double_t AliGenMUONlib::PtBeautyF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1749{
1750// Beauty pT
1751// Corresponding NLO total cross section is 0.445 mb
1752 const Double_t kpt0 = 7.3367;
1753 const Double_t kxn = 3.0692;
1754 Double_t x=*px;
1755 //
1756 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1757 return x/TMath::Power(pass1,kxn);
1758}
75e0cc59 1759Double_t AliGenMUONlib::PtBeautyF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1760{
1761// Beauty pT
1762// Corresponding NLO total cross section is 0.518 mb
1763 const Double_t kpt0 = 7.6409;
1764 const Double_t kxn = 3.1364;
1765 Double_t x=*px;
1766 //
1767 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1768 return x/TMath::Power(pass1,kxn);
1769}
75e0cc59 1770Double_t AliGenMUONlib::PtBeautyF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1771{
1772// Beauty pT
1773// Corresponding NLO total cross section is 0.384 mb
1774 const Double_t kpt0 = 8.4948;
1775 const Double_t kxn = 3.2546;
1776 Double_t x=*px;
1777 //
1778 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1779 return x/TMath::Power(pass1,kxn);
1780}
75e0cc59 1781Double_t AliGenMUONlib::PtBeautyF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1782{
1783// Beauty pT
1784// Corresponding NLO total cross section is 0.648 mb
1785 const Double_t kpt0 = 7.6631;
1786 const Double_t kxn = 3.1621;
1787 Double_t x=*px;
1788 //
1789 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1790 return x/TMath::Power(pass1,kxn);
1791}
75e0cc59 1792Double_t AliGenMUONlib::PtBeautyF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1793{
1794// Beauty pT
1795// Corresponding NLO total cross section is 0.294 mb
1796 const Double_t kpt0 = 8.7245;
1797 const Double_t kxn = 3.2213;
1798 Double_t x=*px;
1799 //
1800 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1801 return x/TMath::Power(pass1,kxn);
1802}
75e0cc59 1803Double_t AliGenMUONlib::PtBeautyF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1804{
1805// Beauty pT
1806// Corresponding NLO total cross section is 0.475 mb
1807 const Double_t kpt0 = 8.5296;
1808 const Double_t kxn = 3.2187;
1809 Double_t x=*px;
1810 //
1811 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1812 return x/TMath::Power(pass1,kxn);
1813}
75e0cc59 1814Double_t AliGenMUONlib::PtBeautyF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1815{
1816// Beauty pT
1817// Corresponding NLO total cross section is 0.324 mb
1818 const Double_t kpt0 = 7.9440;
1819 const Double_t kxn = 3.1614;
1820 Double_t x=*px;
1821 //
1822 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1823 return x/TMath::Power(pass1,kxn);
1824}
75e0cc59 1825Double_t AliGenMUONlib::PtBeautyF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1826{
1827// Beauty pT
1828// Corresponding NLO total cross section is 0.536 mb
1829 const Double_t kpt0 = 8.2408;
1830 const Double_t kxn = 3.3029;
1831 Double_t x=*px;
1832 //
1833 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1834 return x/TMath::Power(pass1,kxn);
1835}
75e0cc59 1836Double_t AliGenMUONlib::PtBeautyF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1837{
1838// Beauty pT
1839// Corresponding NLO total cross section is 0.420 mb
1840 const Double_t kpt0 = 7.8041;
1841 const Double_t kxn = 3.2094;
1842 Double_t x=*px;
1843 //
1844 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
1845 return x/TMath::Power(pass1,kxn);
1846}
1847
fe4da5cc 1848// y-distribution
75e0cc59 1849Double_t AliGenMUONlib::YBeauty( const Double_t *px, const Double_t */*dummy*/)
fe4da5cc 1850{
ec772ba2 1851// Beauty y :: Carrer & Dainese : ALICE-INT-2003-019 v.3 (hep-ph/0311225)
1852// Pythia tuned to reproduce the distribution given by the HVQMNR program based on NLO calculations (pQCD)
1853// shadowing + kt broadening
1854
1855 Double_t x=px[0];
1856 Double_t c[2]={-1.27590e-02,-2.42731e-04};
1857 Double_t y=1+c[0]*TMath::Power(x,2)+c[1]*TMath::Power(x,4);
1858 Double_t ybeauty;
1859
1860 if (TMath::Abs(x)>6) {
1861 ybeauty=0.;
1862 }
1863 else {
1864 ybeauty=TMath::Power(y,3);
1865 }
1866
1867 return ybeauty;
fe4da5cc 1868}
75e0cc59 1869Double_t AliGenMUONlib::YBeautyF0M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1870{
1871// FiMjSkPP define theoretical uncertainties around F0M0S0PP as follows:
1872// YBeautyFiMjSkPP = YBeautyF0M0S0PP * (dN(i,j,k)/dy / dN(0,0,0)/dy)_MNR
1873// i=0,1,2; j=0,1,2; k=0,1,...,6
1874// dN(i,j,k)/dy - spectra obtained by A.Dainese (hep-ph/0601164, p.88;
1875// http://www-zeus.desy.de/~corradi/benchmarks) from NLO pQCD (MNR)
1876// calculations for the following inputs:
1877// Peterson fragmentation function (F) with \epsilon_b = 0.001, 0.0002 & 0.004
1878// for i=0,1 & 2 respectively; quark mass (M) of 4.75, 4.5 & 5.0 GeV
1879// for j=0,1 & 2 respectively;
1880// factorisation \mu_F = a*mt and renormalisation \mu_R = b*mt scales (S)
1881// with a/b = 1/1, 1/0.5, 0.5/1, 0.5/0.5, 1/2, 2/1 & 2/2
1882// for k = 0, 1, 2, 3, 4, 5 & 6 respectively; CTEQ6.1 PDF set
1883// (PDF uncertainty not considered since is small, see hep-ph/0601164, p.89).
1884// June 2008, Smbat.Grigoryan@cern.ch
1885
1886// Beauty y
1887// Pythia6.214 (kBeautyppMNRwmi, PDF = CTEQ5L, quark mass = 4.75 GeV, PtHard > 2.76 GeV/c)
1888// for pp collisions at 14 TeV with one b-bbar pair per event.
1889// Corresponding NLO total cross section is 0.494 mb
1890
fe4da5cc 1891
9e9ae065 1892 Double_t x=px[0];
1893 Double_t c[2]={1.2350e-02,9.2667e-05};
1894 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1895 Double_t ybeauty;
1896
1897 if (TMath::Abs(x)>7.6) {
1898 ybeauty=0.;
1899 }
1900 else {
1901 ybeauty=TMath::Power(y,3);
1902 }
1903
1904 return ybeauty;
1905}
75e0cc59 1906Double_t AliGenMUONlib::YBeautyF1M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1907{
1908// Beauty y
1909// Corresponding NLO total cross section is 0.445 mb
1910 Double_t x=px[0];
1911 Double_t c[2]={1.2292e-02,9.1847e-05};
1912 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1913 Double_t ybeauty;
1914
1915 if (TMath::Abs(x)>7.6) {
1916 ybeauty=0.;
1917 }
1918 else {
1919 ybeauty=TMath::Power(y,3);
1920 }
1921
1922 return ybeauty;
1923}
75e0cc59 1924Double_t AliGenMUONlib::YBeautyF2M0S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1925{
1926// Beauty y
1927// Corresponding NLO total cross section is 0.445 mb
1928 Double_t x=px[0];
1929 Double_t c[2]={1.2436e-02,9.3709e-05};
1930 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1931 Double_t ybeauty;
1932
1933 if (TMath::Abs(x)>7.6) {
1934 ybeauty=0.;
1935 }
1936 else {
1937 ybeauty=TMath::Power(y,3);
1938 }
1939
1940 return ybeauty;
1941}
75e0cc59 1942Double_t AliGenMUONlib::YBeautyF0M1S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1943{
1944// Beauty y
1945// Corresponding NLO total cross section is 0.518 mb
1946 Double_t x=px[0];
1947 Double_t c[2]={1.1714e-02,1.0068e-04};
1948 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1949 Double_t ybeauty;
1950
1951 if (TMath::Abs(x)>7.6) {
1952 ybeauty=0.;
1953 }
1954 else {
1955 ybeauty=TMath::Power(y,3);
1956 }
1957
1958 return ybeauty;
1959}
75e0cc59 1960Double_t AliGenMUONlib::YBeautyF0M2S0PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1961{
1962// Beauty y
1963// Corresponding NLO total cross section is 0.384 mb
1964 Double_t x=px[0];
1965 Double_t c[2]={1.2944e-02,8.5500e-05};
1966 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1967 Double_t ybeauty;
1968
1969 if (TMath::Abs(x)>7.6) {
1970 ybeauty=0.;
1971 }
1972 else {
1973 ybeauty=TMath::Power(y,3);
1974 }
1975
1976 return ybeauty;
1977}
75e0cc59 1978Double_t AliGenMUONlib::YBeautyF0M0S1PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1979{
1980// Beauty y
1981// Corresponding NLO total cross section is 0.648 mb
1982 Double_t x=px[0];
1983 Double_t c[2]={1.2455e-02,9.2713e-05};
1984 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
1985 Double_t ybeauty;
1986
1987 if (TMath::Abs(x)>7.6) {
1988 ybeauty=0.;
1989 }
1990 else {
1991 ybeauty=TMath::Power(y,3);
1992 }
1993
1994 return ybeauty;
1995}
75e0cc59 1996Double_t AliGenMUONlib::YBeautyF0M0S2PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 1997{
1998// Beauty y
1999// Corresponding NLO total cross section is 0.294 mb
2000 Double_t x=px[0];
2001 Double_t c[2]={1.0897e-02,1.1878e-04};
2002 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
2003 Double_t ybeauty;
2004
2005 if (TMath::Abs(x)>7.6) {
2006 ybeauty=0.;
2007 }
2008 else {
2009 ybeauty=TMath::Power(y,3);
2010 }
2011
2012 return ybeauty;
2013}
75e0cc59 2014Double_t AliGenMUONlib::YBeautyF0M0S3PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 2015{
2016// Beauty y
2017// Corresponding NLO total cross section is 0.475 mb
2018 Double_t x=px[0];
2019 Double_t c[2]={1.0912e-02,1.1858e-04};
2020 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
2021 Double_t ybeauty;
2022
2023 if (TMath::Abs(x)>7.6) {
2024 ybeauty=0.;
2025 }
2026 else {
2027 ybeauty=TMath::Power(y,3);
2028 }
2029
2030 return ybeauty;
2031}
75e0cc59 2032Double_t AliGenMUONlib::YBeautyF0M0S4PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 2033{
2034// Beauty y
2035// Corresponding NLO total cross section is 0.324 mb
2036 Double_t x=px[0];
2037 Double_t c[2]={1.2378e-02,9.2490e-05};
2038 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
2039 Double_t ybeauty;
2040
2041 if (TMath::Abs(x)>7.6) {
2042 ybeauty=0.;
2043 }
2044 else {
2045 ybeauty=TMath::Power(y,3);
2046 }
2047
2048 return ybeauty;
2049}
75e0cc59 2050Double_t AliGenMUONlib::YBeautyF0M0S5PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 2051{
2052// Beauty y
2053// Corresponding NLO total cross section is 0.536 mb
2054 Double_t x=px[0];
2055 Double_t c[2]={1.2886e-02,8.2912e-05};
2056 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
2057 Double_t ybeauty;
2058
2059 if (TMath::Abs(x)>7.6) {
2060 ybeauty=0.;
2061 }
2062 else {
2063 ybeauty=TMath::Power(y,3);
2064 }
2065
2066 return ybeauty;
2067}
75e0cc59 2068Double_t AliGenMUONlib::YBeautyF0M0S6PP( const Double_t *px, const Double_t */*dummy*/)
9e9ae065 2069{
2070// Beauty y
2071// Corresponding NLO total cross section is 0.420 mb
2072 Double_t x=px[0];
2073 Double_t c[2]={1.3106e-02,8.0115e-05};
2074 Double_t y=1-c[0]*TMath::Power(x,2)-c[1]*TMath::Power(x,4);
2075 Double_t ybeauty;
2076
2077 if (TMath::Abs(x)>7.6) {
2078 ybeauty=0.;
2079 }
2080 else {
2081 ybeauty=TMath::Power(y,3);
2082 }
2083
2084 return ybeauty;
2085}
ec772ba2 2086
65fb704d 2087Int_t AliGenMUONlib::IpBeauty(TRandom *ran)
fe4da5cc 2088{
d90f80fd 2089// Beauty Composition
65fb704d 2090 Float_t random;
fe4da5cc 2091 Int_t ip;
ec772ba2 2092 random = ran->Rndm();
2093
2094// Taux de production Carrer & Dainese : ALICE-INT-2003-019 v.3
2095// >>>>> cf. tab 4 p 11
2096
2097 if (random < 0.20) {
2098 ip=511;
2099 } else if (random < 0.40) {
2100 ip=-511;
2101 } else if (random < 0.605) {
2102 ip=521;
2103 } else if (random < 0.81) {
2104 ip=-521;
2105 } else if (random < 0.87) {
2106 ip=531;
2107 } else if (random < 0.93) {
2108 ip=-531;
2109 } else if (random < 0.965) {
2110 ip=5122;
fe4da5cc 2111 } else {
ec772ba2 2112 ip=-5122;
fe4da5cc 2113 }
fe4da5cc 2114
ec772ba2 2115 return ip;
fe4da5cc 2116}
2117
ec772ba2 2118
75e0cc59 2119typedef Double_t (*GenFunc) (const Double_t*, const Double_t*);
53904666 2120GenFunc AliGenMUONlib::GetPt(Int_t param, const char* tname) const
fe4da5cc 2121{
d90f80fd 2122// Return pointer to pT parameterisation
05932df6 2123 TString sname = TString(tname);
fe4da5cc 2124 GenFunc func;
753690b0 2125 switch (param)
fe4da5cc 2126 {
34f60c01 2127 case kPhi:
fe4da5cc 2128 func=PtPhi;
2129 break;
89512a3b 2130 case kOmega:
2131 func=PtOmega;
2132 break;
2133 case kEta:
2134 func=PtEta;
2135 break;
0ad09590 2136 case kJpsiFamily:
88e5db43 2137 case kPsiP:
ba8bf3f5 2138 case kChic1:
2139 case kChic2:
34f60c01 2140 case kJpsi:
bb6e81ac 2141 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 2142 func=PtJpsiPbPb;
bb6e81ac 2143 } else if (sname == "Vogt pp") {
2144 func=PtJpsiPP;
0e137c25 2145 } else if (sname == "CDF scaled") {
2146 func=PtJpsiCDFscaled;
9ff13849 2147 } else if (sname == "CDF pp") {
2148 func=PtJpsiCDFscaledPP;
9e9ae065 2149 } else if (sname == "CDF pp 10") {
2150 func=PtJpsiCDFscaledPP10;
9fd56238 2151 } else if (sname == "CDF pp 8.8") {
2152 func=PtJpsiCDFscaledPP9;
94a2994d 2153 } else if (sname == "CDF pp 7" || sname == "CDF pp 7 flat y") {
9fd56238 2154 func=PtJpsiCDFscaledPP7;
2155 } else if (sname == "CDF pp 3.94") {
2156 func=PtJpsiCDFscaledPP4;
2157 } else if (sname == "CDF pPb 8.8") {
2158 func=PtJpsiCDFscaledPPb9;
2159 } else if (sname == "CDF Pbp 8.8") {
2160 func=PtJpsiCDFscaledPbP9;
2161 } else if (sname == "CDF PbPb 3.94") {
2162 func=PtJpsiCDFscaledPbPb4;
94a2994d 2163 } else if (sname == "Flat" || sname == "CDF pp 7 flat pt") {
4ca8d070 2164 func=PtJpsiFlat;
05932df6 2165 } else {
2166 func=PtJpsi;
2167 }
fe4da5cc 2168 break;
17d28ba5 2169 case kJpsiFromB:
2170 func = PtJpsiBPbPb;
2171 break;
0ad09590 2172 case kUpsilonFamily:
88e5db43 2173 case kUpsilonP:
2174 case kUpsilonPP:
34f60c01 2175 case kUpsilon:
bb6e81ac 2176 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 2177 func=PtUpsilonPbPb;
bb6e81ac 2178 } else if (sname == "Vogt pp") {
2179 func=PtUpsilonPP;
0e137c25 2180 } else if (sname == "CDF scaled") {
2181 func=PtUpsilonCDFscaled;
9ff13849 2182 } else if (sname == "CDF pp") {
2183 func=PtUpsilonCDFscaledPP;
b33adf51 2184 } else if (sname == "CDF pp 10") {
2185 func=PtUpsilonCDFscaledPP10;
9fd56238 2186 } else if (sname == "CDF pp 8.8") {
2187 func=PtUpsilonCDFscaledPP9;
2188 } else if (sname == "CDF pp 7") {
2189 func=PtUpsilonCDFscaledPP7;
2190 } else if (sname == "CDF pp 3.94") {
2191 func=PtUpsilonCDFscaledPP4;
2192 } else if (sname == "CDF pPb 8.8") {
2193 func=PtUpsilonCDFscaledPPb9;
2194 } else if (sname == "CDF Pbp 8.8") {
2195 func=PtUpsilonCDFscaledPbP9;
2196 } else if (sname == "CDF PbPb 3.94") {
2197 func=PtUpsilonCDFscaledPbPb4;
4ca8d070 2198 } else if (sname == "Flat") {
2199 func=PtUpsilonFlat;
05932df6 2200 } else {
2201 func=PtUpsilon;
2202 }
0ad09590 2203 break;
34f60c01 2204 case kCharm:
9e9ae065 2205 if (sname == "F0M0S0 pp") {
2206 func=PtCharmF0M0S0PP;
2207 } else if (sname == "F1M0S0 pp") {
2208 func=PtCharmF1M0S0PP;
2209 } else if (sname == "F2M0S0 pp") {
2210 func=PtCharmF2M0S0PP;
2211 } else if (sname == "F0M1S0 pp") {
2212 func=PtCharmF0M1S0PP;
2213 } else if (sname == "F0M2S0 pp") {
2214 func=PtCharmF0M2S0PP;
2215 } else if (sname == "F0M0S1 pp") {
2216 func=PtCharmF0M0S1PP;
2217 } else if (sname == "F0M0S2 pp") {
2218 func=PtCharmF0M0S2PP;
2219 } else if (sname == "F0M0S3 pp") {
2220 func=PtCharmF0M0S3PP;
2221 } else if (sname == "F0M0S4 pp") {
2222 func=PtCharmF0M0S4PP;
2223 } else if (sname == "F0M0S5 pp") {
2224 func=PtCharmF0M0S5PP;
2225 } else if (sname == "F0M0S6 pp") {
2226 func=PtCharmF0M0S6PP;
2227 } else if (sname == "central") {
bd0276a8 2228 func=PtCharmCentral;
2229 } else {
2230 func=PtCharm;
2231 }
fe4da5cc 2232 break;
34f60c01 2233 case kBeauty:
9e9ae065 2234 if (sname == "F0M0S0 pp") {
2235 func=PtBeautyF0M0S0PP;
2236 } else if (sname == "F1M0S0 pp") {
2237 func=PtBeautyF1M0S0PP;
2238 } else if (sname == "F2M0S0 pp") {
2239 func=PtBeautyF2M0S0PP;
2240 } else if (sname == "F0M1S0 pp") {
2241 func=PtBeautyF0M1S0PP;
2242 } else if (sname == "F0M2S0 pp") {
2243 func=PtBeautyF0M2S0PP;
2244 } else if (sname == "F0M0S1 pp") {
2245 func=PtBeautyF0M0S1PP;
2246 } else if (sname == "F0M0S2 pp") {
2247 func=PtBeautyF0M0S2PP;
2248 } else if (sname == "F0M0S3 pp") {
2249 func=PtBeautyF0M0S3PP;
2250 } else if (sname == "F0M0S4 pp") {
2251 func=PtBeautyF0M0S4PP;
2252 } else if (sname == "F0M0S5 pp") {
2253 func=PtBeautyF0M0S5PP;
2254 } else if (sname == "F0M0S6 pp") {
2255 func=PtBeautyF0M0S6PP;
2256 } else if (sname == "central") {
bd0276a8 2257 func=PtBeautyCentral;
2258 } else {
2259 func=PtBeauty;
2260 }
fe4da5cc 2261 break;
34f60c01 2262 case kPion:
753690b0 2263 func=PtPion;
2264 break;
34f60c01 2265 case kKaon:
753690b0 2266 func=PtKaon;
2267 break;
6b450aa0 2268 case kChic0:
2269 func=PtChic0;
00935af2 2270 break;
6b450aa0 2271 case kChic:
2272 func=PtChic;
00935af2 2273 break;
119b35c7 2274 default:
2275 func=0;
2276 printf("<AliGenMUONlib::GetPt> unknown parametrisation\n");
fe4da5cc 2277 }
2278 return func;
2279}
2280
53904666 2281GenFunc AliGenMUONlib::GetY(Int_t param, const char* tname) const
fe4da5cc 2282{
ac3faee4 2283 //
2284 // Return pointer to y- parameterisation
2285 //
05932df6 2286 TString sname = TString(tname);
fe4da5cc 2287 GenFunc func;
753690b0 2288 switch (param)
fe4da5cc 2289 {
34f60c01 2290 case kPhi:
fe4da5cc 2291 func=YPhi;
2292 break;
89512a3b 2293 case kEta:
2294 func=YEta;
2295 break;
2296 case kOmega:
2297 func=YOmega;
2298 break;
0ad09590 2299 case kJpsiFamily:
88e5db43 2300 case kPsiP:
ba8bf3f5 2301 case kChic1:
2302 case kChic2:
34f60c01 2303 case kJpsi:
bb6e81ac 2304 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 2305 func=YJpsiPbPb;
bb6e81ac 2306 } else if (sname == "Vogt pp"){
2307 func=YJpsiPP;
0e137c25 2308 } else if (sname == "CDF scaled") {
2309 func=YJpsiCDFscaled;
9ff13849 2310 } else if (sname == "CDF pp") {
2311 func=YJpsiCDFscaledPP;
9e9ae065 2312 } else if (sname == "CDF pp 10") {
2313 func=YJpsiCDFscaledPP10;
9fd56238 2314 } else if (sname == "CDF pp 8.8") {
2315 func=YJpsiCDFscaledPP9;
94a2994d 2316 } else if (sname == "CDF pp 7" || sname == "CDF pp 7 flat pt") {
9fd56238 2317 func=YJpsiCDFscaledPP7;
2318 } else if (sname == "CDF pp 3.94") {
2319 func=YJpsiCDFscaledPP4;
2320 } else if (sname == "CDF pPb 8.8") {
2321 func=YJpsiCDFscaledPPb9;
2322 } else if (sname == "CDF Pbp 8.8") {
2323 func=YJpsiCDFscaledPbP9;
2324 } else if (sname == "CDF PbPb 3.94") {
2325 func=YJpsiCDFscaledPbPb4;
94a2994d 2326 } else if (sname == "Flat" || sname == "CDF pp 7 flat y") {
4ca8d070 2327 func=YJpsiFlat;
05932df6 2328 } else {
2329 func=YJpsi;
2330 }
17d28ba5 2331 break;
2332 case kJpsiFromB:
2333 func = YJpsiBPbPb;
fe4da5cc 2334 break;
0ad09590 2335 case kUpsilonFamily:
88e5db43 2336 case kUpsilonP:
2337 case kUpsilonPP:
34f60c01 2338 case kUpsilon:
bb6e81ac 2339 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 2340 func=YUpsilonPbPb;
bb6e81ac 2341 } else if (sname == "Vogt pp") {
2342 func = YUpsilonPP;
0e137c25 2343 } else if (sname == "CDF scaled") {
2344 func=YUpsilonCDFscaled;
9ff13849 2345 } else if (sname == "CDF pp") {
2346 func=YUpsilonCDFscaledPP;
b33adf51 2347 } else if (sname == "CDF pp 10") {
2348 func=YUpsilonCDFscaledPP10;
9fd56238 2349 } else if (sname == "CDF pp 8.8") {
2350 func=YUpsilonCDFscaledPP9;
2351 } else if (sname == "CDF pp 7") {
2352 func=YUpsilonCDFscaledPP7;
2353 } else if (sname == "CDF pp 3.94") {
2354 func=YUpsilonCDFscaledPP4;
2355 } else if (sname == "CDF pPb 8.8") {
2356 func=YUpsilonCDFscaledPPb9;
2357 } else if (sname == "CDF Pbp 8.8") {
2358 func=YUpsilonCDFscaledPbP9;
2359 } else if (sname == "CDF PbPb 3.94") {
2360 func=YUpsilonCDFscaledPbPb4;
4ca8d070 2361 } else if (sname == "Flat") {
2362 func=YUpsilonFlat;
05932df6 2363 } else {
2364 func=YUpsilon;
2365 }
fe4da5cc 2366 break;
34f60c01 2367 case kCharm:
9e9ae065 2368 if (sname == "F0M0S0 pp") {
2369 func=YCharmF0M0S0PP;
2370 } else if (sname == "F1M0S0 pp") {
2371 func=YCharmF1M0S0PP;
2372 } else if (sname == "F2M0S0 pp") {
2373 func=YCharmF2M0S0PP;
2374 } else if (sname == "F0M1S0 pp") {
2375 func=YCharmF0M1S0PP;
2376 } else if (sname == "F0M2S0 pp") {
2377 func=YCharmF0M2S0PP;
2378 } else if (sname == "F0M0S1 pp") {
2379 func=YCharmF0M0S1PP;
2380 } else if (sname == "F0M0S2 pp") {
2381 func=YCharmF0M0S2PP;
2382 } else if (sname == "F0M0S3 pp") {
2383 func=YCharmF0M0S3PP;
2384 } else if (sname == "F0M0S4 pp") {
2385 func=YCharmF0M0S4PP;
2386 } else if (sname == "F0M0S5 pp") {
2387 func=YCharmF0M0S5PP;
2388 } else if (sname == "F0M0S6 pp") {
2389 func=YCharmF0M0S6PP;
2390 } else {
2391 func=YCharm;
2392 }
fe4da5cc 2393 break;
34f60c01 2394 case kBeauty:
9e9ae065 2395 if (sname == "F0M0S0 pp") {
2396 func=YBeautyF0M0S0PP;
2397 } else if (sname == "F1M0S0 pp") {
2398 func=YBeautyF1M0S0PP;
2399 } else if (sname == "F2M0S0 pp") {
2400 func=YBeautyF2M0S0PP;
2401 } else if (sname == "F0M1S0 pp") {
2402 func=YBeautyF0M1S0PP;
2403 } else if (sname == "F0M2S0 pp") {
2404 func=YBeautyF0M2S0PP;
2405 } else if (sname == "F0M0S1 pp") {
2406 func=YBeautyF0M0S1PP;
2407 } else if (sname == "F0M0S2 pp") {
2408 func=YBeautyF0M0S2PP;
2409 } else if (sname == "F0M0S3 pp") {
2410 func=YBeautyF0M0S3PP;
2411 } else if (sname == "F0M0S4 pp") {
2412 func=YBeautyF0M0S4PP;
2413 } else if (sname == "F0M0S5 pp") {
2414 func=YBeautyF0M0S5PP;
2415 } else if (sname == "F0M0S6 pp") {
2416 func=YBeautyF0M0S6PP;
2417 } else {
2418 func=YBeauty;
2419 }
fe4da5cc 2420 break;
34f60c01 2421 case kPion:
753690b0 2422 func=YPion;
2423 break;
34f60c01 2424 case kKaon:
753690b0 2425 func=YKaon;
2426 break;
6b450aa0 2427 case kChic0:
2428 func=YChic0;
00935af2 2429 break;
6b450aa0 2430 case kChic:
2431 func=YChic;
00935af2 2432 break;
119b35c7 2433 default:
2434 func=0;
2435 printf("<AliGenMUONlib::GetY> unknown parametrisation\n");
fe4da5cc 2436 }
2437 return func;
2438}
00935af2 2439
2440//
2441// Chi
2442//
2443//
2444// pt-distribution
2445//____________________________________________________________
75e0cc59 2446Double_t AliGenMUONlib::PtChic0( const Double_t *px, const Double_t */*dummy*/)
00935af2 2447{
2448// Chi_c1 pT
2449 const Double_t kpt0 = 4.;
2450 const Double_t kxn = 3.6;
2451 Double_t x=*px;
2452 //
2453 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2454 return x/TMath::Power(pass1,kxn);
2455}
75e0cc59 2456Double_t AliGenMUONlib::PtChic1( const Double_t *px, const Double_t */*dummy*/)
00935af2 2457{
2458// Chi_c1 pT
2459 const Double_t kpt0 = 4.;
2460 const Double_t kxn = 3.6;
2461 Double_t x=*px;
2462 //
2463 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2464 return x/TMath::Power(pass1,kxn);
2465}
75e0cc59 2466Double_t AliGenMUONlib::PtChic2( const Double_t *px, const Double_t */*dummy*/)
00935af2 2467{
2468// Chi_c2 pT
2469 const Double_t kpt0 = 4.;
2470 const Double_t kxn = 3.6;
2471 Double_t x=*px;
2472 //
2473 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2474 return x/TMath::Power(pass1,kxn);
2475}
75e0cc59 2476Double_t AliGenMUONlib::PtChic( const Double_t *px, const Double_t */*dummy*/)
00935af2 2477{
2478// Chi_c family pT
2479 const Double_t kpt0 = 4.;
2480 const Double_t kxn = 3.6;
2481 Double_t x=*px;
2482 //
2483 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
2484 return x/TMath::Power(pass1,kxn);
2485}
2486
2487//
2488// y-distribution
2489//____________________________________________________________
75e0cc59 2490Double_t AliGenMUONlib::YChic0(const Double_t *py, const Double_t */*dummy*/)
00935af2 2491{
2492// Chi-1c y
2493 const Double_t ky0 = 4.;
2494 const Double_t kb=1.;
2495 Double_t yj;
2496 Double_t y=TMath::Abs(*py);
2497 //
2498 if (y < ky0)
2499 yj=kb;
2500 else
2501 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2502 return yj;
2503}
2504
75e0cc59 2505Double_t AliGenMUONlib::YChic1(const Double_t *py, const Double_t */*dummy*/)
00935af2 2506{
2507// Chi-1c y
2508 const Double_t ky0 = 4.;
2509 const Double_t kb=1.;
2510 Double_t yj;
2511 Double_t y=TMath::Abs(*py);
2512 //
2513 if (y < ky0)
2514 yj=kb;
2515 else
2516 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2517 return yj;
2518}
2519
75e0cc59 2520Double_t AliGenMUONlib::YChic2(const Double_t *py, const Double_t */*dummy*/)
00935af2 2521{
2522// Chi-2c y
2523 const Double_t ky0 = 4.;
2524 const Double_t kb=1.;
2525 Double_t yj;
2526 Double_t y=TMath::Abs(*py);
2527 //
2528 if (y < ky0)
2529 yj=kb;
2530 else
2531 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2532 return yj;
2533}
2534
75e0cc59 2535Double_t AliGenMUONlib::YChic(const Double_t *py, const Double_t */*dummy*/)
00935af2 2536{
2537// Chi_c family y
2538 const Double_t ky0 = 4.;
2539 const Double_t kb=1.;
2540 Double_t yj;
2541 Double_t y=TMath::Abs(*py);
2542 //
2543 if (y < ky0)
2544 yj=kb;
2545 else
2546 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
2547 return yj;
2548}
2549
2550// particle composition
2551//
6b450aa0 2552Int_t AliGenMUONlib::IpChic0(TRandom *)
00935af2 2553{
2554// Chi composition
2555 return 10441;
2556}
2557//
6b450aa0 2558Int_t AliGenMUONlib::IpChic1(TRandom *)
00935af2 2559{
2560// Chi composition
2561 return 20443;
2562}
6b450aa0 2563Int_t AliGenMUONlib::IpChic2(TRandom *)
00935af2 2564{
2565// Chi_c2 prime composition
2566 return 445;
2567}
6b450aa0 2568Int_t AliGenMUONlib::IpChic(TRandom *)
00935af2 2569{
2570// Chi composition
2571 Int_t ip;
2572 Float_t r = gRandom->Rndm();
2573 if (r < 0.001) {
2574 ip = 10441;
2575 } else if( r < 0.377 ) {
2576 ip = 20443;
2577 } else {
2578 ip = 445;
2579 }
2580 return ip;
2581}
2582
2583
2584//_____________________________________________________________
2585
65fb704d 2586typedef Int_t (*GenFuncIp) (TRandom *);
198bb1c7 2587GenFuncIp AliGenMUONlib::GetIp(Int_t param, const char* /*tname*/) const
fe4da5cc 2588{
d90f80fd 2589// Return pointer to particle type parameterisation
fe4da5cc 2590 GenFuncIp func;
753690b0 2591 switch (param)
fe4da5cc 2592 {
34f60c01 2593 case kPhi:
fe4da5cc 2594 func=IpPhi;
2595 break;
89512a3b 2596 case kEta:
2597 func=IpEta;
2598 break;
2599 case kOmega:
2600 func=IpOmega;
2601 break;
0ad09590 2602 case kJpsiFamily:
2603 func=IpJpsiFamily;
2604 break;
88e5db43 2605 case kPsiP:
2606 func=IpPsiP;
2607 break;
34f60c01 2608 case kJpsi:
17d28ba5 2609 case kJpsiFromB:
fe4da5cc 2610 func=IpJpsi;
2611 break;
34f60c01 2612 case kUpsilon:
fe4da5cc 2613 func=IpUpsilon;
2614 break;
0ad09590 2615 case kUpsilonFamily:
2616 func=IpUpsilonFamily;
2617 break;
88e5db43 2618 case kUpsilonP:
2619 func=IpUpsilonP;
2620 break;
2621 case kUpsilonPP:
2622 func=IpUpsilonPP;
2623 break;
34f60c01 2624 case kCharm:
fe4da5cc 2625 func=IpCharm;
2626 break;
34f60c01 2627 case kBeauty:
fe4da5cc 2628 func=IpBeauty;
2629 break;
34f60c01 2630 case kPion:
753690b0 2631 func=IpPion;
2632 break;
34f60c01 2633 case kKaon:
753690b0 2634 func=IpKaon;
2635 break;
6b450aa0 2636 case kChic0:
2637 func=IpChic0;
00935af2 2638 break;
6b450aa0 2639 case kChic1:
2640 func=IpChic1;
00935af2 2641 break;
6b450aa0 2642 case kChic2:
2643 func=IpChic2;
00935af2 2644 break;
6b450aa0 2645 case kChic:
2646 func=IpChic;
00935af2 2647 break;
119b35c7 2648 default:
2649 func=0;
2650 printf("<AliGenMUONlib::GetIp> unknown parametrisation\n");
fe4da5cc 2651 }
2652 return func;
2653}
2654
2655
753690b0 2656
05932df6 2657Float_t AliGenMUONlib::Interpolate(Float_t x, Float_t* y, Float_t x0,
2658 Float_t dx,
2659 Int_t n, Int_t no)
2660{
2661//
2662// Neville's alorithm for interpolation
2663//
2664// x: x-value
2665// y: Input array
2666// x0: minimum x
2667// dx: step size
2668// n: number of data points
2669// no: order of polynom
2670//
2671 Float_t* c = new Float_t[n];
2672 Float_t* d = new Float_t[n];
2673 Int_t m, i;
2674 for (i = 0; i < n; i++) {
2675 c[i] = y[i];
2676 d[i] = y[i];
2677 }
2678
2679 Int_t ns = int((x - x0)/dx);
2680
2681 Float_t y1 = y[ns];
2682 ns--;
2683 for (m = 0; m < no; m++) {
2684 for (i = 0; i < n-m; i++) {
2685 Float_t ho = x0 + Float_t(i) * dx - x;
2686 Float_t hp = x0 + Float_t(i+m+1) * dx - x;
2687 Float_t w = c[i+1] - d[i];
2688 Float_t den = ho-hp;
2689 den = w/den;
2690 d[i] = hp * den;
2691 c[i] = ho * den;
2692 }
2693 Float_t dy;
2694
2695 if (2*ns < (n-m-1)) {
2696 dy = c[ns+1];
2697 } else {
2698 dy = d[ns--];
2699 }
2700 y1 += dy;}
2701 delete[] c;
2702 delete[] d;
2703
2704 return y1;
2705}
2706
753690b0 2707