Quarkonia for pp @ 14 TeV added.
[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
16/*
17$Log$
bb6e81ac 18Revision 1.20 2003/03/13 11:54:39 morsch
19Limited pT range for parameterized Upsilon and J/Psi pT distributions.
20
3d905dd7 21Revision 1.19 2003/02/24 16:46:11 morsch
22New parameterisation for Psi and Upsilon (PbPb)
23
1af7144e 24Revision 1.18 2002/11/07 09:13:42 morsch
25Use "Vogt" to label new distributions.
26
113049ac 27Revision 1.17 2002/11/07 09:06:10 morsch
28J/Psi and Upsilon pt and y distributions from R. Vogt 2002 added.
29
05932df6 30Revision 1.16 2002/10/14 14:55:35 hristov
31Merging the VirtualMC branch to the main development branch (HEAD)
32
b9d0a01d 33Revision 1.14.6.1 2002/06/10 14:57:41 hristov
34Merged with v3-08-02
35
36Revision 1.15 2002/04/17 10:11:51 morsch
37Coding Rule violations corrected.
38
53904666 39Revision 1.14 2002/02/22 17:26:43 morsch
40Eta and omega added.
41
89512a3b 42Revision 1.13 2001/03/27 11:01:04 morsch
43Charm pt-distribution corrected. More realistic y-distribution for pi and K.
44
2280e6af 45Revision 1.12 2001/03/09 13:01:41 morsch
46- enum constants for paramterisation type (particle family) moved to AliGen*lib.h
47- use AliGenGSIlib::kUpsilon, AliGenPHOSlib::kEtaPrime to access the constants
48
34f60c01 49Revision 1.11 2000/11/30 07:12:50 alibrary
50Introducing new Rndm and QA classes
51
65fb704d 52Revision 1.10 2000/06/29 21:08:27 morsch
53All paramatrisation libraries derive from the pure virtual base class AliGenLib.
54This allows to pass a pointer to a library directly to AliGenParam and avoids the
55use of function pointers in Config.C.
56
b22ee262 57Revision 1.9 2000/06/14 15:20:56 morsch
58Include clean-up (IH)
59
5c3fd7ea 60Revision 1.8 2000/06/09 20:32:11 morsch
61All coding rule violations except RS3 corrected
62
f87cfe57 63Revision 1.7 2000/05/02 08:12:13 morsch
64Coding rule violations corrected.
65
d90f80fd 66Revision 1.6 1999/09/29 09:24:14 fca
67Introduction of the Copyright and cvs Log
68
4c039060 69*/
70
53904666 71// Library class for particle pt and y distributions used for
72// muon spectrometer simulations.
73// To be used with AliGenParam.
74// The following particle typed can be simulated:
75// pi, K, phi, omega, eta, J/Psi, Upsilon, charm and beauty mesons.
76//
77// andreas.morsch@cern.ch
78//
79
65fb704d 80#include "TMath.h"
81#include "TRandom.h"
82
fe4da5cc 83#include "AliGenMUONlib.h"
5c3fd7ea 84
fe4da5cc 85ClassImp(AliGenMUONlib)
86//
87// Pions
d90f80fd 88Double_t AliGenMUONlib::PtPion(Double_t *px, Double_t *dummy)
fe4da5cc 89{
90//
91// PT-PARAMETERIZATION CDF, PRL 61(88) 1819
92// POWER LAW FOR PT > 500 MEV
93// MT SCALING BELOW (T=160 MEV)
94//
d90f80fd 95 const Double_t kp0 = 1.3;
96 const Double_t kxn = 8.28;
97 const Double_t kxlim=0.5;
98 const Double_t kt=0.160;
99 const Double_t kxmpi=0.139;
100 const Double_t kb=1.;
fe4da5cc 101 Double_t y, y1, xmpi2, ynorm, a;
102 Double_t x=*px;
103 //
d90f80fd 104 y1=TMath::Power(kp0/(kp0+kxlim),kxn);
105 xmpi2=kxmpi*kxmpi;
106 ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
fe4da5cc 107 a=ynorm/y1;
d90f80fd 108 if (x > kxlim)
109 y=a*TMath::Power(kp0/(kp0+x),kxn);
fe4da5cc 110 else
d90f80fd 111 y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
fe4da5cc 112 return y*x;
113}
753690b0 114//
115// y-distribution
116//
d90f80fd 117Double_t AliGenMUONlib::YPion( Double_t *py, Double_t *dummy)
753690b0 118{
d90f80fd 119// Pion y
2280e6af 120 Double_t y=TMath::Abs(*py);
121/*
d90f80fd 122 const Double_t ka = 7000.;
123 const Double_t kdy = 4.;
d90f80fd 124 Double_t ex = y*y/(2*kdy*kdy);
125 return ka*TMath::Exp(-ex);
2280e6af 126*/
127 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
128
753690b0 129}
130// particle composition
131//
65fb704d 132Int_t AliGenMUONlib::IpPion(TRandom *ran)
753690b0 133{
d90f80fd 134// Pion composition
65fb704d 135 if (ran->Rndm() < 0.5) {
753690b0 136 return 211;
137 } else {
138 return -211;
139 }
140}
fe4da5cc 141
142//____________________________________________________________
143//
144// Mt-scaling
145
146Double_t AliGenMUONlib::PtScal(Double_t pt, Int_t np)
147{
148 // SCALING EN MASSE PAR RAPPORT A PTPI
149 // MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
d90f80fd 150 const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
fe4da5cc 151 // VALUE MESON/PI AT 5 GEV
d90f80fd 152 const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
fe4da5cc 153 np--;
d90f80fd 154 Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
155 Double_t fmax2=f5/kfmax[np];
fe4da5cc 156 // PIONS
157 Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
158 Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
d90f80fd 159 (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
fe4da5cc 160 return fmtscal*ptpion;
161}
162//
753690b0 163// kaon
164//
165// pt-distribution
166//____________________________________________________________
d90f80fd 167Double_t AliGenMUONlib::PtKaon( Double_t *px, Double_t *dummy)
753690b0 168{
d90f80fd 169// Kaon pT
753690b0 170 return PtScal(*px,2);
171}
172
173// y-distribution
fe4da5cc 174//____________________________________________________________
d90f80fd 175Double_t AliGenMUONlib::YKaon( Double_t *py, Double_t *dummy)
fe4da5cc 176{
d90f80fd 177// Kaon y
2280e6af 178 Double_t y=TMath::Abs(*py);
179/*
d90f80fd 180 const Double_t ka = 1000.;
181 const Double_t kdy = 4.;
fe4da5cc 182 //
d90f80fd 183 Double_t ex = y*y/(2*kdy*kdy);
184 return ka*TMath::Exp(-ex);
2280e6af 185*/
186
187 return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
753690b0 188}
189
190// particle composition
191//
65fb704d 192Int_t AliGenMUONlib::IpKaon(TRandom *ran)
753690b0 193{
d90f80fd 194// Kaon composition
65fb704d 195 if (ran->Rndm() < 0.5) {
753690b0 196 return 321;
197 } else {
198 return -321;
199 }
fe4da5cc 200}
753690b0 201
fe4da5cc 202// J/Psi
203//
204//
205// pt-distribution
206//____________________________________________________________
d90f80fd 207Double_t AliGenMUONlib::PtJpsi( Double_t *px, Double_t *dummy)
fe4da5cc 208{
d90f80fd 209// J/Psi pT
210 const Double_t kpt0 = 4.;
211 const Double_t kxn = 3.6;
fe4da5cc 212 Double_t x=*px;
213 //
d90f80fd 214 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
215 return x/TMath::Power(pass1,kxn);
fe4da5cc 216}
05932df6 217
218Double_t AliGenMUONlib::PtJpsiPbPb( Double_t *px, Double_t *dummy)
219{
1af7144e 220// J/Psi pT spectrum
05932df6 221//
222// R. Vogt 2002
223// PbPb 5.5 TeV
224// MRST HO
225// mc = 1.4 GeV, pt-kick 1 GeV
226//
1af7144e 227 Float_t x = px[0];
228 Float_t c[8] = {
229 -2.13098e+00, 9.46552e+00, -5.06799e+00, 1.27260e+00,
230 -1.83806e-01, 1.55853e-02, -7.23241e-04, 1.42105e-05
05932df6 231 };
1af7144e 232
3d905dd7 233 Double_t y;
234 if (x < 10.) {
235 Int_t j;
236 y = c[j = 7];
237 while (j > 0) y = y * x +c[--j];
238 y = x * TMath::Exp(y);
239 } else {
240 y = 0.;
241 }
1af7144e 242 return y;
05932df6 243}
bb6e81ac 244Double_t AliGenMUONlib::PtJpsiPP( Double_t *px, Double_t *dummy)
245{
246// J/Psi pT spectrum
247//
248// R. Vogt 2002
249// pp 14 TeV
250// MRST HO
251// mc = 1.4 GeV, pt-kick 1 GeV
252//
253 Float_t x = px[0];
254 Float_t c[4] = {8.47471e+00, -1.93567e+00, 1.50271e-01, -5.51212e-03};
255
256 Double_t y;
257 if (x < 10.) {
258 Int_t j;
259 y = c[j = 3];
260 while (j > 0) y = y * x +c[--j];
261 y = x * TMath::Exp(y);
262 } else {
263 y = 0.;
264 }
265 return y;
266}
267
fe4da5cc 268//
269// y-distribution
270//____________________________________________________________
d90f80fd 271Double_t AliGenMUONlib::YJpsi(Double_t *py, Double_t *dummy)
fe4da5cc 272{
d90f80fd 273// J/psi y
274 const Double_t ky0 = 4.;
275 const Double_t kb=1.;
fe4da5cc 276 Double_t yj;
277 Double_t y=TMath::Abs(*py);
278 //
d90f80fd 279 if (y < ky0)
280 yj=kb;
fe4da5cc 281 else
d90f80fd 282 yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
fe4da5cc 283 return yj;
284}
05932df6 285
286
287Double_t AliGenMUONlib::YJpsiPbPb( Double_t *px, Double_t *dummy)
288{
289
290//
291// J/Psi y
292//
293//
294// R. Vogt 2002
295// PbPb 5.5 TeV
296// MRST HO
297// mc = 1.4 GeV, pt-kick 1 GeV
298//
1af7144e 299 Double_t c[5] = {-6.03425e+02, 4.98257e+02, -1.38794e+02, 1.62209e+01, -6.85955e-01};
300 Double_t x = TMath::Abs(px[0]);
301 Double_t y;
302
303 if (x < 4.) {
304 y = 31.754;
305 } else if (x < 6) {
306 Int_t j;
307 y = c[j = 4];
308 while (j > 0) y = y * x + c[--j];
309 } else {
310 y =0.;
311 }
312
313 return y;
05932df6 314}
315
bb6e81ac 316Double_t AliGenMUONlib::YJpsiPP( Double_t *px, Double_t *dummy)
317{
318
319//
320// J/Psi y
321//
322//
323// R. Vogt 2002
324// pp 14 TeV
325// MRST HO
326// mc = 1.4 GeV, pt-kick 1 GeV
327//
328
329 Double_t c[5] = {1.38532e+00, 1.00596e+02, -3.46378e+01, 3.94172e+00, -1.48319e-01};
330 Double_t x = TMath::Abs(px[0]);
331 Double_t y;
332
333 if (x < 2.5) {
334 y = 96.455 - 0.8483 * x * x;
335 } else if (x < 7.9) {
336 Int_t j;
337 y = c[j = 4];
338 while (j > 0) y = y * x + c[--j];
339 } else {
340 y =0.;
341 }
342
343 return y;
344}
345
fe4da5cc 346// particle composition
347//
65fb704d 348Int_t AliGenMUONlib::IpJpsi(TRandom *)
fe4da5cc 349{
d90f80fd 350// J/Psi composition
fe4da5cc 351 return 443;
352}
353
354// Upsilon
355//
356//
357// pt-distribution
358//____________________________________________________________
d90f80fd 359Double_t AliGenMUONlib::PtUpsilon( Double_t *px, Double_t *dummy )
fe4da5cc 360{
d90f80fd 361// Upsilon pT
362 const Double_t kpt0 = 5.3;
363 const Double_t kxn = 2.5;
fe4da5cc 364 Double_t x=*px;
365 //
d90f80fd 366 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
367 return x/TMath::Power(pass1,kxn);
fe4da5cc 368}
05932df6 369
370Double_t AliGenMUONlib::PtUpsilonPbPb( Double_t *px, Double_t *dummy)
371{
372
373//
374// Upsilon pT
375//
376//
377// R. Vogt 2002
378// PbPb 5.5 TeV
379// MRST HO
380// mc = 1.4 GeV, pt-kick 1 GeV
381//
1af7144e 382 Float_t x = px[0];
383 Double_t c[8] = {
384 -1.03488e+01, 1.28065e+01, -6.60500e+00, 1.66140e+00,
385 -2.34293e-01, 1.86925e-02, -7.80708e-04, 1.30610e-05
386 };
3d905dd7 387 Double_t y;
388 if (x < 10.) {
389 Int_t j;
390 y = c[j = 7];
391 while (j > 0) y = y * x +c[--j];
392 y = x * TMath::Exp(y);
393 } else {
394 y = 0.;
395 }
1af7144e 396 return y;
05932df6 397}
398
bb6e81ac 399Double_t AliGenMUONlib::PtUpsilonPP( Double_t *px, Double_t *dummy)
400{
401
402//
403// Upsilon pT
404//
405//
406// R. Vogt 2002
407// pp 14 TeV
408// MRST HO
409// mc = 1.4 GeV, pt-kick 1 GeV
410//
411 Float_t x = px[0];
412 Double_t c[8] = {-7.93955e+00, 1.06306e+01, -5.21392e+00, 1.19703e+00,
413 -1.45718e-01, 8.95151e-03, -2.04806e-04, -1.13053e-06};
414
415 Double_t y;
416 if (x < 10.) {
417 Int_t j;
418 y = c[j = 7];
419 while (j > 0) y = y * x +c[--j];
420 y = x * TMath::Exp(y);
421 } else {
422 y = 0.;
423 }
424 return y;
425}
426
fe4da5cc 427//
428// y-distribution
429//
430//____________________________________________________________
d90f80fd 431Double_t AliGenMUONlib::YUpsilon(Double_t *py, Double_t *dummy)
fe4da5cc 432{
d90f80fd 433// Upsilon y
434 const Double_t ky0 = 3.;
435 const Double_t kb=1.;
fe4da5cc 436 Double_t yu;
437 Double_t y=TMath::Abs(*py);
438 //
d90f80fd 439 if (y < ky0)
440 yu=kb;
fe4da5cc 441 else
d90f80fd 442 yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
fe4da5cc 443 return yu;
444}
05932df6 445
446
447Double_t AliGenMUONlib::YUpsilonPbPb( Double_t *px, Double_t *dummy)
448{
449
450//
451// Upsilon y
452//
453//
454// R. Vogt 2002
455// PbPb 5.5 TeV
456// MRST HO
457// mc = 1.4 GeV, pt-kick 1 GeV
458//
459
1af7144e 460 Double_t c[7] = {3.40036e-01, -3.98882e-07, -4.48398e-03, 8.46411e-08, -6.10854e-04,
461 -2.99753e-09, 1.28895e-05};
462
463 Double_t x = px[0];
464 if (TMath::Abs(x) > 5.55) return 0.;
465 Int_t j;
466 Double_t y = c[j = 6];
467 while (j > 0) y = y * x +c[--j];
468 return y;
05932df6 469}
470
bb6e81ac 471Double_t AliGenMUONlib::YUpsilonPP( Double_t *px, Double_t *dummy)
472{
473
474//
475// Upsilon y
476//
477//
478// R. Vogt 2002
479// p p 14. TeV
480// MRST HO
481// mc = 1.4 GeV, pt-kick 1 GeV
482//
483 Double_t c[7] = {8.91936e-01, -6.46645e-07, -1.52774e-02, 4.28677e-08, -7.01517e-04,
484 -6.20539e-10, 1.29943e-05};
485
486 Double_t x = px[0];
487 if (TMath::Abs(x) > 6.2) return 0.;
488 Int_t j;
489 Double_t y = c[j = 6];
490 while (j > 0) y = y * x +c[--j];
491 return y;
492}
493
fe4da5cc 494// particle composition
495//
65fb704d 496Int_t AliGenMUONlib::IpUpsilon(TRandom *)
fe4da5cc 497{
d90f80fd 498// y composition
fe4da5cc 499 return 553;
500}
501
502//
503// Phi
504//
505//
506// pt-distribution (by scaling of pion distribution)
507//____________________________________________________________
d90f80fd 508Double_t AliGenMUONlib::PtPhi( Double_t *px, Double_t *dummy)
fe4da5cc 509{
d90f80fd 510// Phi pT
fe4da5cc 511 return PtScal(*px,7);
512}
513// y-distribution
d90f80fd 514Double_t AliGenMUONlib::YPhi( Double_t *px, Double_t *dummy)
fe4da5cc 515{
d90f80fd 516// Phi y
517 Double_t *dum=0;
518 return YJpsi(px,dum);
fe4da5cc 519}
520// particle composition
521//
65fb704d 522Int_t AliGenMUONlib::IpPhi(TRandom *)
fe4da5cc 523{
d90f80fd 524// Phi composition
89512a3b 525 return 333;
526}
527
528//
529// omega
530//
531//
532// pt-distribution (by scaling of pion distribution)
533//____________________________________________________________
534Double_t AliGenMUONlib::PtOmega( Double_t *px, Double_t *dummy)
535{
536// Omega pT
537 return PtScal(*px,5);
538}
539// y-distribution
540Double_t AliGenMUONlib::YOmega( Double_t *px, Double_t *dummy)
541{
542// Omega y
543 Double_t *dum=0;
544 return YJpsi(px,dum);
545}
546// particle composition
547//
548Int_t AliGenMUONlib::IpOmega(TRandom *)
549{
550// Omega composition
551 return 223;
552}
553
554
555//
556// Eta
557//
558//
559// pt-distribution (by scaling of pion distribution)
560//____________________________________________________________
561Double_t AliGenMUONlib::PtEta( Double_t *px, Double_t *dummy)
562{
563// Eta pT
564 return PtScal(*px,3);
565}
566// y-distribution
567Double_t AliGenMUONlib::YEta( Double_t *px, Double_t *dummy)
568{
569// Eta y
570 Double_t *dum=0;
571 return YJpsi(px,dum);
572}
573// particle composition
574//
575Int_t AliGenMUONlib::IpEta(TRandom *)
576{
577// Eta composition
578 return 221;
fe4da5cc 579}
580
581//
582// Charm
583//
584//
585// pt-distribution
586//____________________________________________________________
d90f80fd 587Double_t AliGenMUONlib::PtCharm( Double_t *px, Double_t *dummy)
fe4da5cc 588{
d90f80fd 589// Charm pT
590 const Double_t kpt0 = 4.08;
591 const Double_t kxn = 9.40;
2280e6af 592
fe4da5cc 593 Double_t x=*px;
594 //
2280e6af 595 Double_t pass1 = 1.+(x/kpt0);
d90f80fd 596 return x/TMath::Power(pass1,kxn);
fe4da5cc 597}
598// y-distribution
d90f80fd 599Double_t AliGenMUONlib::YCharm( Double_t *px, Double_t *dummy)
fe4da5cc 600{
d90f80fd 601// Charm y
602 Double_t *dum=0;
603 return YJpsi(px,dum);
fe4da5cc 604}
605
65fb704d 606Int_t AliGenMUONlib::IpCharm(TRandom *ran)
fe4da5cc 607{
d90f80fd 608// Charm composition
65fb704d 609 Float_t random;
fe4da5cc 610 Int_t ip;
611// 411,421,431,4122
65fb704d 612 random = ran->Rndm();
613 if (random < 0.5) {
fe4da5cc 614 ip=411;
65fb704d 615 } else if (random < 0.75) {
fe4da5cc 616 ip=421;
65fb704d 617 } else if (random < 0.90) {
fe4da5cc 618 ip=431;
619 } else {
620 ip=4122;
621 }
65fb704d 622 if (ran->Rndm() < 0.5) {ip=-ip;}
fe4da5cc 623
624 return ip;
625}
626
627
628//
629// Beauty
630//
631//
632// pt-distribution
633//____________________________________________________________
d90f80fd 634Double_t AliGenMUONlib::PtBeauty( Double_t *px, Double_t *dummy)
fe4da5cc 635{
d90f80fd 636// Beauty pT
637 const Double_t kpt0 = 4.;
638 const Double_t kxn = 3.6;
fe4da5cc 639 Double_t x=*px;
640 //
d90f80fd 641 Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
642 return x/TMath::Power(pass1,kxn);
fe4da5cc 643}
644// y-distribution
d90f80fd 645Double_t AliGenMUONlib::YBeauty( Double_t *px, Double_t *dummy)
fe4da5cc 646{
d90f80fd 647// Beauty y
648 Double_t *dum=0;
649 return YJpsi(px,dum);
fe4da5cc 650}
651
65fb704d 652Int_t AliGenMUONlib::IpBeauty(TRandom *ran)
fe4da5cc 653{
d90f80fd 654// Beauty Composition
65fb704d 655 Float_t random;
fe4da5cc 656 Int_t ip;
65fb704d 657 random = ran->Rndm();
658 if (random < 0.5) {
fe4da5cc 659 ip=511;
65fb704d 660 } else if (random < 0.75) {
fe4da5cc 661 ip=521;
65fb704d 662 } else if (random < 0.90) {
fe4da5cc 663 ip=531;
664 } else {
665 ip=5122;
666 }
65fb704d 667 if (ran->Rndm() < 0.5) {ip=-ip;}
fe4da5cc 668
669 return ip;
670}
671
672typedef Double_t (*GenFunc) (Double_t*, Double_t*);
53904666 673GenFunc AliGenMUONlib::GetPt(Int_t param, const char* tname) const
fe4da5cc 674{
d90f80fd 675// Return pointer to pT parameterisation
05932df6 676 TString sname = TString(tname);
fe4da5cc 677 GenFunc func;
753690b0 678 switch (param)
fe4da5cc 679 {
34f60c01 680 case kPhi:
fe4da5cc 681 func=PtPhi;
682 break;
89512a3b 683 case kOmega:
684 func=PtOmega;
685 break;
686 case kEta:
687 func=PtEta;
688 break;
34f60c01 689 case kJpsi:
bb6e81ac 690 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 691 func=PtJpsiPbPb;
bb6e81ac 692 } else if (sname == "Vogt pp") {
693 func=PtJpsiPP;
05932df6 694 } else {
695 func=PtJpsi;
696 }
fe4da5cc 697 break;
34f60c01 698 case kUpsilon:
bb6e81ac 699 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 700 func=PtUpsilonPbPb;
bb6e81ac 701 } else if (sname == "Vogt pp") {
702 func=PtUpsilonPP;
05932df6 703 } else {
704 func=PtUpsilon;
705 }
fe4da5cc 706 break;
34f60c01 707 case kCharm:
fe4da5cc 708 func=PtCharm;
709 break;
34f60c01 710 case kBeauty:
fe4da5cc 711 func=PtBeauty;
712 break;
34f60c01 713 case kPion:
753690b0 714 func=PtPion;
715 break;
34f60c01 716 case kKaon:
753690b0 717 func=PtKaon;
718 break;
119b35c7 719 default:
720 func=0;
721 printf("<AliGenMUONlib::GetPt> unknown parametrisation\n");
fe4da5cc 722 }
723 return func;
724}
725
53904666 726GenFunc AliGenMUONlib::GetY(Int_t param, const char* tname) const
fe4da5cc 727{
05932df6 728 TString sname = TString(tname);
729
d90f80fd 730// Return pointer to y- parameterisation
fe4da5cc 731 GenFunc func;
753690b0 732 switch (param)
fe4da5cc 733 {
34f60c01 734 case kPhi:
fe4da5cc 735 func=YPhi;
736 break;
89512a3b 737 case kEta:
738 func=YEta;
739 break;
740 case kOmega:
741 func=YOmega;
742 break;
34f60c01 743 case kJpsi:
bb6e81ac 744 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 745 func=YJpsiPbPb;
bb6e81ac 746 } else if (sname == "Vogt pp"){
747 func=YJpsiPP;
05932df6 748 } else {
749 func=YJpsi;
750 }
bb6e81ac 751
fe4da5cc 752 break;
34f60c01 753 case kUpsilon:
bb6e81ac 754 if (sname == "Vogt" || sname == "Vogt PbPb") {
05932df6 755 func=YUpsilonPbPb;
bb6e81ac 756 } else if (sname == "Vogt pp") {
757 func = YUpsilonPP;
05932df6 758 } else {
759 func=YUpsilon;
760 }
fe4da5cc 761 break;
34f60c01 762 case kCharm:
fe4da5cc 763 func=YCharm;
764 break;
34f60c01 765 case kBeauty:
fe4da5cc 766 func=YBeauty;
767 break;
34f60c01 768 case kPion:
753690b0 769 func=YPion;
770 break;
34f60c01 771 case kKaon:
753690b0 772 func=YKaon;
773 break;
119b35c7 774 default:
775 func=0;
776 printf("<AliGenMUONlib::GetY> unknown parametrisation\n");
fe4da5cc 777 }
778 return func;
779}
65fb704d 780typedef Int_t (*GenFuncIp) (TRandom *);
53904666 781GenFuncIp AliGenMUONlib::GetIp(Int_t param, const char* tname) const
fe4da5cc 782{
d90f80fd 783// Return pointer to particle type parameterisation
fe4da5cc 784 GenFuncIp func;
753690b0 785 switch (param)
fe4da5cc 786 {
34f60c01 787 case kPhi:
fe4da5cc 788 func=IpPhi;
789 break;
89512a3b 790 case kEta:
791 func=IpEta;
792 break;
793 case kOmega:
794 func=IpOmega;
795 break;
34f60c01 796 case kJpsi:
fe4da5cc 797 func=IpJpsi;
798 break;
34f60c01 799 case kUpsilon:
fe4da5cc 800 func=IpUpsilon;
801 break;
34f60c01 802 case kCharm:
fe4da5cc 803 func=IpCharm;
804 break;
34f60c01 805 case kBeauty:
fe4da5cc 806 func=IpBeauty;
807 break;
34f60c01 808 case kPion:
753690b0 809 func=IpPion;
810 break;
34f60c01 811 case kKaon:
753690b0 812 func=IpKaon;
813 break;
119b35c7 814 default:
815 func=0;
816 printf("<AliGenMUONlib::GetIp> unknown parametrisation\n");
fe4da5cc 817 }
818 return func;
819}
820
821
753690b0 822
05932df6 823Float_t AliGenMUONlib::Interpolate(Float_t x, Float_t* y, Float_t x0,
824 Float_t dx,
825 Int_t n, Int_t no)
826{
827//
828// Neville's alorithm for interpolation
829//
830// x: x-value
831// y: Input array
832// x0: minimum x
833// dx: step size
834// n: number of data points
835// no: order of polynom
836//
837 Float_t* c = new Float_t[n];
838 Float_t* d = new Float_t[n];
839 Int_t m, i;
840 for (i = 0; i < n; i++) {
841 c[i] = y[i];
842 d[i] = y[i];
843 }
844
845 Int_t ns = int((x - x0)/dx);
846
847 Float_t y1 = y[ns];
848 ns--;
849 for (m = 0; m < no; m++) {
850 for (i = 0; i < n-m; i++) {
851 Float_t ho = x0 + Float_t(i) * dx - x;
852 Float_t hp = x0 + Float_t(i+m+1) * dx - x;
853 Float_t w = c[i+1] - d[i];
854 Float_t den = ho-hp;
855 den = w/den;
856 d[i] = hp * den;
857 c[i] = ho * den;
858 }
859 Float_t dy;
860
861 if (2*ns < (n-m-1)) {
862 dy = c[ns+1];
863 } else {
864 dy = d[ns--];
865 }
866 y1 += dy;}
867 delete[] c;
868 delete[] d;
869
870 return y1;
871}
872
753690b0 873