New reader for the pedestal run and vdrift (Julian) and some bug fixing (Raphaelle)
[u/mrichter/AliRoot.git] / EVGEN / AliGenPHOSlib.cxx
CommitLineData
74c62c73 1
886b6f73 2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
88cb7938 17/* $Id$ */
886b6f73 18
19//======================================================================
20// AliGenPHOSlib class contains parameterizations of the
21// pion, kaon, eta, omega, etaprime, phi and baryon (proton,
22// antiproton, neutron and anti-neutron) particles for the
23// study of the neutral background in PHOS detector.
24// These parameterizations are used by the
25// AliGenParam class:
26// AliGenParam(npar, param, AliGenPHOSlib::GetPt(param),
27// AliGenPHOSlib::GetY(param),
28// AliGenPHOSlib::GetIp(param) )
29// param represents the particle to be simulated :
30// Pion, Kaon, Eta, Omega, Etaprime, Phi or Baryon
31// Pt distributions are calculated from the transverse mass scaling
9ae59e17 32// with Pions, using the PtScal function taken from AliGenMUONlib
886b6f73 33// version aliroot 3.01
34//
9ae59e17 35// Gines MARTINEZ. Laurent APHECETCHE and Yves SCHUTZ
36// GPS @ SUBATECH, Nantes , France (October 1999)
886b6f73 37// http://www-subatech.in2p3.fr/~photons/subatech
38// martinez@subatech.in2p3.fr
35ead6c1 39// Additional particle species simulation options has been added:
40// Charged Pion, Charged Kaons, KLong Proton, Anti-Proton, Neutron,
41// Anti-Neutron --> Changes made by Gustavo Conesa in November 2004
886b6f73 42//======================================================================
43
65fb704d 44#include "TMath.h"
45#include "TRandom.h"
46
886b6f73 47#include "AliGenPHOSlib.h"
886b6f73 48
49ClassImp(AliGenPHOSlib)
50
51//======================================================================
52// P I O N S
53// (From GetPt, GetY and GetIp as param = Pion)
54// Transverse momentum distribution" PtPion
55// Rapidity distribution YPion
56// Particle distribution IdPion 111, 211 and -211 (pi0, pi+ and pi-)
57//
75e0cc59 58 Double_t AliGenPHOSlib::PtPion(const Double_t *px, const Double_t *)
886b6f73 59{
60// Pion transverse momentum distribtuion taken
9ae59e17 61// from AliGenMUONlib class, version 3.01 of aliroot
886b6f73 62// PT-PARAMETERIZATION CDF, PRL 61(88) 1819
63// POWER LAW FOR PT > 500 MEV
64// MT SCALING BELOW (T=160 MEV)
65//
f87cfe57 66 const Double_t kp0 = 1.3;
67 const Double_t kxn = 8.28;
68 const Double_t kxlim=0.5;
69 const Double_t kt=0.160;
70 const Double_t kxmpi=0.139;
71 const Double_t kb=1.;
72 Double_t y, y1, kxmpi2, ynorm, a;
886b6f73 73 Double_t x=*px;
74 //
f87cfe57 75 y1=TMath::Power(kp0/(kp0+kxlim),kxn);
76 kxmpi2=kxmpi*kxmpi;
77 ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+kxmpi2)/kt));
886b6f73 78 a=ynorm/y1;
f87cfe57 79 if (x > kxlim)
80 y=a*TMath::Power(kp0/(kp0+x),kxn);
886b6f73 81 else
f87cfe57 82 y=kb*TMath::Exp(-sqrt(x*x+kxmpi2)/kt);
886b6f73 83 return y*x;
84}
75e0cc59 85 Double_t AliGenPHOSlib::YPion( const Double_t *py, const Double_t *)
886b6f73 86{
f87cfe57 87//
88// pion y-distribution
89//
90
91 const Double_t ka = 7000.;
92 const Double_t kdy = 4.;
886b6f73 93
94 Double_t y=TMath::Abs(*py);
95 //
f87cfe57 96 Double_t ex = y*y/(2*kdy*kdy);
97 return ka*TMath::Exp(-ex);
886b6f73 98}
f87cfe57 99
65fb704d 100 Int_t AliGenPHOSlib::IpPion(TRandom *ran)
886b6f73 101{
f87cfe57 102// particle composition pi+, pi0, pi-
103//
104
74c62c73 105 Float_t random = ran->Rndm();
106
107 if ( (3.*random) < 1. )
108 {
109 return 211 ;
110 }
111 else
112 {
113 if ( (3.*random) >= 2.)
114 {
115 return -211 ;
116 }
117 else
118 {
886b6f73 119 return 111 ;
120 }
121 }
122}
35ead6c1 123 Int_t AliGenPHOSlib::IpChargedPion(TRandom *ran)
124{
125// particle composition pi+, pi0, pi-
126//
127
128 Float_t random = ran->Rndm();
129
130 if ( (2.*random) < 1. )
131 {
132 return 211 ;
133 }
134 else
135 {
136 return -211 ;
137 }
138}
74c62c73 139
140//End Pions
141//======================================================================
142// Pi 0 Flat Distribution
143// Transverse momentum distribution PtPi0Flat
144// Rapidity distribution YPi0Flat
145// Particle distribution IdPi0Flat 111 (pi0)
146//
147
75e0cc59 148Double_t AliGenPHOSlib::PtPi0Flat(const Double_t */*px*/, const Double_t *)
74c62c73 149{
150// Pion transverse momentum flat distribution
151
152return 1;
153
154}
155
75e0cc59 156Double_t AliGenPHOSlib::YPi0Flat( const Double_t */*py*/, const Double_t *)
74c62c73 157{
158
159// pion y-distribution
160//
161 return 1.;
162}
163
164 Int_t AliGenPHOSlib::IpPi0Flat(TRandom *)
165{
166
167// particle composition pi0
168//
169 return 111 ;
170}
171// End Pi0Flat
886b6f73 172//=============================================================
173//
f87cfe57 174 Double_t AliGenPHOSlib::PtScal(Double_t pt, Int_t np)
175{
886b6f73 176// Mt-scaling
177// Fonction for the calculation of the Pt distribution for a
178// given particle np, from the pion Pt distribution using the
179// mt scaling. This function was taken from AliGenMUONlib
180// aliroot version 3.01, and was extended for baryons
181// np = 1=>Pions 2=>Kaons 3=>Etas 4=>Omegas 5=>ETA' 6=>PHI
182// 7=>BARYONS-BARYONBARS
f87cfe57 183
886b6f73 184 // SCALING EN MASSE PAR RAPPORT A PTPI
185 // MASS 1=>PI, 2=>K, 3=>ETA, 4=>OMEGA, 5=>ETA',6=>PHI
f87cfe57 186 const Double_t khm[10] = {0.1396, 0.494, 0.547, 0.782, 0.957, 1.02,
886b6f73 187 // MASS 7=>BARYON-BARYONBAR
188 0.938, 0. , 0., 0.};
189 // VALUE MESON/PI AT 5 GEV
f87cfe57 190 const Double_t kfmax[10]={1., 1., 1., 1., 1., 1., 1., 1., 1., 1.};
886b6f73 191 np--;
f87cfe57 192 Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
193 Double_t kfmax2=f5/kfmax[np];
886b6f73 194 // PIONS
195 Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
196 Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
f87cfe57 197 (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2;
886b6f73 198 return fmtscal*ptpion;
74c62c73 199
886b6f73 200}
201// End Scaling
202//============================================================================
203// K A O N S
75e0cc59 204 Double_t AliGenPHOSlib::PtKaon( const Double_t *px, const Double_t *)
886b6f73 205{
f87cfe57 206// kaon
207// pt-distribution
208//____________________________________________________________
209
886b6f73 210 return PtScal(*px,2); // 2==> Kaon in the PtScal function
211}
212
75e0cc59 213 Double_t AliGenPHOSlib::YKaon( const Double_t *py, const Double_t *)
886b6f73 214{
f87cfe57 215// y-distribution
216//____________________________________________________________
217
218 const Double_t ka = 1000.;
219 const Double_t kdy = 4.;
886b6f73 220
221
222 Double_t y=TMath::Abs(*py);
223 //
f87cfe57 224 Double_t ex = y*y/(2*kdy*kdy);
225 return ka*TMath::Exp(-ex);
886b6f73 226}
227
65fb704d 228 Int_t AliGenPHOSlib::IpKaon(TRandom *ran)
886b6f73 229{
f87cfe57 230// particle composition
231//
232
65fb704d 233 Float_t random = ran->Rndm();
234 Float_t random2 = ran->Rndm();
235 if (random2 < 0.5)
886b6f73 236 {
65fb704d 237 if (random < 0.5) {
886b6f73 238 return 321; // K+
239 } else {
240 return -321; // K-
241 }
242 }
243 else
244 {
65fb704d 245 if (random < 0.5) {
9ae59e17 246 return 130; // K^0 short
886b6f73 247 } else {
9ae59e17 248 return 310; // K^0 long
886b6f73 249 }
250 }
251}
35ead6c1 252
253 Int_t AliGenPHOSlib::IpChargedKaon(TRandom *ran)
254{
255// particle composition
256//
257
258 Float_t random = ran->Rndm();
259
260 if (random < 0.5) {
261 return 321; // K+
262 } else {
263 return -321; // K-
264 }
265
266
267}
268Int_t AliGenPHOSlib::IpKaon0L(TRandom *)
269{
270 // particle composition
271 //
272
273 return 130; // K^0 long
274}
886b6f73 275// End Kaons
276//============================================================================
277//============================================================================
278// E T A S
75e0cc59 279 Double_t AliGenPHOSlib::PtEta( const Double_t *px, const Double_t *)
886b6f73 280{
f87cfe57 281// etas
282// pt-distribution
283//____________________________________________________________
284
886b6f73 285 return PtScal(*px,3); // 3==> Eta in the PtScal function
286}
287
75e0cc59 288 Double_t AliGenPHOSlib::YEta( const Double_t *py, const Double_t *)
886b6f73 289{
f87cfe57 290// y-distribution
291//____________________________________________________________
292
293 const Double_t ka = 1000.;
294 const Double_t kdy = 4.;
886b6f73 295
296
297 Double_t y=TMath::Abs(*py);
298 //
f87cfe57 299 Double_t ex = y*y/(2*kdy*kdy);
300 return ka*TMath::Exp(-ex);
886b6f73 301}
302
65fb704d 303 Int_t AliGenPHOSlib::IpEta(TRandom *)
886b6f73 304{
f87cfe57 305// particle composition
306//
307
886b6f73 308 return 221; // eta
309}
310// End Etas
74c62c73 311
312//======================================================================
313// Eta Flat Distribution
314// Transverse momentum distribution PtEtaFlat
315// Rapidity distribution YEtaFlat
316// Particle distribution IdEtaFlat 111 (pi0)
317//
318
75e0cc59 319Double_t AliGenPHOSlib::PtEtaFlat(const Double_t */*px*/, const Double_t *)
74c62c73 320{
321// Eta transverse momentum flat distribution
322
323 return 1;
324
325}
326
75e0cc59 327Double_t AliGenPHOSlib::YEtaFlat( const Double_t */*py*/, const Double_t *)
74c62c73 328{
329//
330// pion y-distribution
331//
332 return 1.;
333}
334
335 Int_t AliGenPHOSlib::IpEtaFlat(TRandom *)
336{
337//
338// particle composition eta
339//
340 return 221 ;
341}
35ead6c1 342// End EtaFlat
886b6f73 343//============================================================================
344//============================================================================
345// O M E G A S
75e0cc59 346 Double_t AliGenPHOSlib::PtOmega( const Double_t *px, const Double_t *)
f87cfe57 347{
886b6f73 348// omegas
349// pt-distribution
350//____________________________________________________________
f87cfe57 351
886b6f73 352 return PtScal(*px,4); // 4==> Omega in the PtScal function
353}
354
75e0cc59 355 Double_t AliGenPHOSlib::YOmega( const Double_t *py, const Double_t *)
886b6f73 356{
f87cfe57 357// y-distribution
358//____________________________________________________________
359
360 const Double_t ka = 1000.;
361 const Double_t kdy = 4.;
886b6f73 362
363
364 Double_t y=TMath::Abs(*py);
365 //
f87cfe57 366 Double_t ex = y*y/(2*kdy*kdy);
367 return ka*TMath::Exp(-ex);
886b6f73 368}
369
65fb704d 370 Int_t AliGenPHOSlib::IpOmega(TRandom *)
886b6f73 371{
f87cfe57 372// particle composition
373//
374
886b6f73 375 return 223; // Omega
376}
377// End Omega
378//============================================================================
379//============================================================================
380// E T A P R I M E
75e0cc59 381 Double_t AliGenPHOSlib::PtEtaprime( const Double_t *px, const Double_t *)
f87cfe57 382{
886b6f73 383// etaprime
384// pt-distribution
385//____________________________________________________________
f87cfe57 386
886b6f73 387 return PtScal(*px,5); // 5==> Etaprime in the PtScal function
388}
389
75e0cc59 390 Double_t AliGenPHOSlib::YEtaprime( const Double_t *py, const Double_t *)
886b6f73 391{
f87cfe57 392// y-distribution
393//____________________________________________________________
394
395 const Double_t ka = 1000.;
396 const Double_t kdy = 4.;
886b6f73 397
398
399 Double_t y=TMath::Abs(*py);
400 //
f87cfe57 401 Double_t ex = y*y/(2*kdy*kdy);
402 return ka*TMath::Exp(-ex);
886b6f73 403}
404
65fb704d 405 Int_t AliGenPHOSlib::IpEtaprime(TRandom *)
886b6f73 406{
f87cfe57 407// particle composition
408//
409
886b6f73 410 return 331; // Etaprime
411}
412// End EtaPrime
413//===================================================================
414//============================================================================
415// P H I S
75e0cc59 416 Double_t AliGenPHOSlib::PtPhi( const Double_t *px, const Double_t *)
f87cfe57 417{
886b6f73 418// phi
419// pt-distribution
420//____________________________________________________________
f87cfe57 421
886b6f73 422 return PtScal(*px,6); // 6==> Phi in the PtScal function
423}
424
75e0cc59 425 Double_t AliGenPHOSlib::YPhi( const Double_t *py, const Double_t *)
886b6f73 426{
f87cfe57 427// y-distribution
428//____________________________________________________________
429
430 const Double_t ka = 1000.;
431 const Double_t kdy = 4.;
886b6f73 432
433
434 Double_t y=TMath::Abs(*py);
435 //
f87cfe57 436 Double_t ex = y*y/(2*kdy*kdy);
437 return ka*TMath::Exp(-ex);
886b6f73 438}
439
65fb704d 440 Int_t AliGenPHOSlib::IpPhi(TRandom *)
f87cfe57 441{
886b6f73 442// particle composition
443//
f87cfe57 444
886b6f73 445 return 333; // Phi
446}
447// End Phis
448//===================================================================
449//============================================================================
450// B A R Y O N S == protons, protonsbar, neutrons, and neutronsbars
75e0cc59 451 Double_t AliGenPHOSlib::PtBaryon( const Double_t *px, const Double_t *)
f87cfe57 452{
886b6f73 453// baryons
454// pt-distribution
455//____________________________________________________________
f87cfe57 456
886b6f73 457 return PtScal(*px,7); // 7==> Baryon in the PtScal function
458}
459
75e0cc59 460 Double_t AliGenPHOSlib::YBaryon( const Double_t *py, const Double_t *)
886b6f73 461{
f87cfe57 462// y-distribution
463//____________________________________________________________
464
465 const Double_t ka = 1000.;
466 const Double_t kdy = 4.;
886b6f73 467
468
469 Double_t y=TMath::Abs(*py);
470 //
f87cfe57 471 Double_t ex = y*y/(2*kdy*kdy);
472 return ka*TMath::Exp(-ex);
886b6f73 473}
474
65fb704d 475 Int_t AliGenPHOSlib::IpBaryon(TRandom *ran)
886b6f73 476{
f87cfe57 477// particle composition
478//
479
65fb704d 480 Float_t random = ran->Rndm();
481 Float_t random2 = ran->Rndm();
482 if (random2 < 0.5)
886b6f73 483 {
65fb704d 484 if (random < 0.5) {
886b6f73 485 return 2212; // p
486 } else {
487 return -2212; // pbar
488 }
489 }
490 else
491 {
65fb704d 492 if (random < 0.5) {
886b6f73 493 return 2112; // n
494 } else {
495 return -2112; // n bar
496 }
497 }
498}
35ead6c1 499
500 Int_t AliGenPHOSlib::IpProton(TRandom *)
501{
502// particle composition
503//
504 return 2212; // p
505
506}
507 Int_t AliGenPHOSlib::IpAProton(TRandom *)
508{
509// particle composition
510//
511 return -2212; // p bar
512
513}
514
515 Int_t AliGenPHOSlib::IpNeutron(TRandom *)
516{
517// particle composition
518//
519 return 2112; // n
520
521}
522 Int_t AliGenPHOSlib::IpANeutron(TRandom *)
523{
524// particle composition
525//
526 return -2112; // n
527
528}
886b6f73 529// End Baryons
530//===================================================================
531
532
35ead6c1 533
75e0cc59 534typedef Double_t (*GenFunc) (const Double_t*, const Double_t*);
198bb1c7 535GenFunc AliGenPHOSlib::GetPt(Int_t param, const char* /*tname*/) const
886b6f73 536{
f87cfe57 537// Return pinter to pT parameterisation
886b6f73 538 GenFunc func;
539
540 switch (param)
35ead6c1 541 {
542 case kPion:
886b6f73 543 func=PtPion;
544 break;
35ead6c1 545 case kPi0Flat:
74c62c73 546 func=PtPi0Flat;
547 break;
35ead6c1 548 case kKaon:
886b6f73 549 func=PtKaon;
550 break;
35ead6c1 551 case kEta:
886b6f73 552 func=PtEta;
553 break;
35ead6c1 554 case kEtaFlat:
74c62c73 555 func=PtEtaFlat;
556 break;
35ead6c1 557 case kOmega:
886b6f73 558 func=PtOmega;
559 break;
35ead6c1 560 case kEtaPrime:
886b6f73 561 func=PtEtaprime;
562 break;
35ead6c1 563 case kBaryon:
886b6f73 564 func=PtBaryon;
565 break;
35ead6c1 566 default:
886b6f73 567 func=0;
568 printf("<AliGenPHOSlib::GetPt> unknown parametrisationn");
35ead6c1 569 }
886b6f73 570 return func;
571}
572
198bb1c7 573GenFunc AliGenPHOSlib::GetY(Int_t param, const char* /*tname*/) const
886b6f73 574{
35ead6c1 575 // Return pointer to Y parameterisation
576 GenFunc func;
577 switch (param)
886b6f73 578 {
34f60c01 579 case kPion:
35ead6c1 580 func=YPion;
581 break;
74c62c73 582 case kPi0Flat:
35ead6c1 583 func=YPi0Flat;
584 break;
34f60c01 585 case kKaon:
35ead6c1 586 func=YKaon;
587 break;
34f60c01 588 case kEta:
35ead6c1 589 func=YEta;
590 break;
591 case kEtaFlat:
592 func=YEtaFlat;
593 break;
34f60c01 594 case kOmega:
35ead6c1 595 func=YOmega;
596 break;
34f60c01 597 case kEtaPrime:
35ead6c1 598 func=YEtaprime;
599 break;
34f60c01 600 case kPhi:
35ead6c1 601 func=YPhi;
602 break;
34f60c01 603 case kBaryon:
35ead6c1 604 func=YBaryon;
605 break;
886b6f73 606 default:
35ead6c1 607 func=0;
608 printf("<AliGenPHOSlib::GetY> unknown parametrisationn");
886b6f73 609 }
35ead6c1 610 return func;
886b6f73 611}
65fb704d 612typedef Int_t (*GenFuncIp) (TRandom *);
198bb1c7 613GenFuncIp AliGenPHOSlib::GetIp(Int_t param, const char* /*tname*/) const
886b6f73 614{
35ead6c1 615 // Return pointer to particle composition
616 GenFuncIp func;
617 switch (param)
886b6f73 618 {
74c62c73 619 case kPion:
35ead6c1 620 func=IpPion;
621 break;
622 case kChargedPion:
623 func=IpChargedPion;
624 break;
74c62c73 625 case kPi0Flat:
35ead6c1 626 func=IpPi0Flat;
627 break;
34f60c01 628 case kKaon:
35ead6c1 629 func=IpKaon;
630 break;
631 case kChargedKaon:
632 func=IpChargedKaon;
633 break;
634 case kKaon0L:
635 func=IpKaon0L;
636 break;
34f60c01 637 case kEta:
35ead6c1 638 func=IpEta;
639 break;
74c62c73 640 case kEtaFlat:
35ead6c1 641 func=IpEtaFlat;
642 break;
643
34f60c01 644 case kOmega:
35ead6c1 645 func=IpOmega;
646 break;
34f60c01 647 case kEtaPrime:
35ead6c1 648 func=IpEtaprime;
649 break;
34f60c01 650 case kPhi:
35ead6c1 651 func=IpPhi;
652 break;
34f60c01 653 case kBaryon:
35ead6c1 654 func=IpBaryon;
655 break;
656 case kProton:
657 func=IpProton;
658 break;
659 case kAProton:
660 func=IpAProton;
661 break;
662 case kNeutron:
663 func=IpNeutron;
664 break;
665 case kANeutron:
666 func=IpANeutron;
667 break;
668
886b6f73 669 default:
35ead6c1 670 func=0;
671 printf("<AliGenPHOSlib::GetIp> unknown parametrisationn");
886b6f73 672 }
35ead6c1 673 return func;
886b6f73 674}
675