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