EVGEN/AliGenPHOSlib.cxx

   1
   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
  17 /*
  18 $Log$
  19 Revision 1.9  2002/04/23 12:54:29  morsch
  20 New options kPi0Flat y kEtaFlat (Gustavo Conesa)
  21
  22 Revision 1.7  2001/03/09 13:01:41  morsch
  23 - enum constants for paramterisation type (particle family) moved to AliGen*lib.h
  24 - use AliGenGSIlib::kUpsilon, AliGenPHOSlib::kEtaPrime to access the constants
  25
  26 Revision 1.6  2000/11/30 07:12:50  alibrary
  27 Introducing new Rndm and QA classes
  28
  29 Revision 1.5  2000/06/29 21:08:27  morsch
  30 All paramatrisation libraries derive from the pure virtual base class AliGenLib.
  31 This allows to pass a pointer to a library directly to AliGenParam and avoids the
  32 use of function pointers in Config.C.
  33
  34 Revision 1.4  2000/06/14 15:21:05  morsch
  35 Include clean-up (IH)
  36
  37 Revision 1.3  2000/06/09 20:32:54  morsch
  38 All coding rule violations except RS3 corrected
  39
  40 Revision 1.2  1999/11/04 11:30:48  fca
  41 Improve comments
  42
  43 Revision 1.1  1999/11/03 17:43:20  fca
  44 New version from G.Martinez & A.Morsch
  45
  46 */
  47
  48 //======================================================================
  49 //  AliGenPHOSlib class contains parameterizations of the
  50 //  pion, kaon, eta, omega, etaprime, phi and baryon (proton,
  51 //  antiproton, neutron and anti-neutron) particles for the
  52 //  study of the neutral background in PHOS detector.
  53 //  These parameterizations are used by the
  54 //  AliGenParam  class:
  55 //  AliGenParam(npar, param,  AliGenPHOSlib::GetPt(param),
  56 //                            AliGenPHOSlib::GetY(param),
  57 //                            AliGenPHOSlib::GetIp(param) )
  58 //  param represents the particle to be simulated :
  59 //  Pion, Kaon, Eta, Omega, Etaprime, Phi or Baryon
  60 //  Pt distributions are calculated from the transverse mass scaling
  61 //  with Pions, using the PtScal function taken from AliGenMUONlib
  62 //  version aliroot 3.01
  63 //
  64 //     Gines MARTINEZ. Laurent APHECETCHE and Yves SCHUTZ
  65 //      GPS @ SUBATECH,  Nantes , France  (October 1999)
  66 //     http://www-subatech.in2p3.fr/~photons/subatech
  67 //     martinez@subatech.in2p3.fr
  68 //======================================================================
  69
  70 #include "TMath.h"
  71 #include "TRandom.h"
  72
  73 #include "AliGenPHOSlib.h"
  74
  75 ClassImp(AliGenPHOSlib)
  76
  77 //======================================================================
  78 //    P  I  O  N  S
  79 //    (From GetPt, GetY and GetIp as param = Pion)
  80 //    Transverse momentum distribution" PtPion
  81 //    Rapidity distribution YPion
  82 //    Particle distribution IdPion  111, 211 and -211 (pi0, pi+ and pi-)
  83 //
  84  Double_t AliGenPHOSlib::PtPion(Double_t *px, Double_t *)
  85 {
  86 //     Pion transverse momentum distribtuion taken
  87 //     from AliGenMUONlib class, version 3.01 of aliroot
  88 //     PT-PARAMETERIZATION CDF, PRL 61(88) 1819
  89 //     POWER LAW FOR PT > 500 MEV
  90 //     MT SCALING BELOW (T=160 MEV)
  91 //
  92   const Double_t kp0 = 1.3;
  93   const Double_t kxn = 8.28;
  94   const Double_t kxlim=0.5;
  95   const Double_t kt=0.160;
  96   const Double_t kxmpi=0.139;
  97   const Double_t kb=1.;
  98   Double_t y, y1, kxmpi2, ynorm, a;
  99   Double_t x=*px;
 100   //
 101   y1=TMath::Power(kp0/(kp0+kxlim),kxn);
 102   kxmpi2=kxmpi*kxmpi;
 103   ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+kxmpi2)/kt));
 104   a=ynorm/y1;
 105   if (x > kxlim)
 106     y=a*TMath::Power(kp0/(kp0+x),kxn);
 107   else
 108     y=kb*TMath::Exp(-sqrt(x*x+kxmpi2)/kt);
 109   return y*x;
 110 }
 111  Double_t AliGenPHOSlib::YPion( Double_t *py, Double_t *)
 112 {
 113 //
 114 // pion y-distribution
 115 //
 116
 117   const Double_t ka    = 7000.;
 118   const Double_t kdy   = 4.;
 119
 120   Double_t y=TMath::Abs(*py);
 121   //
 122   Double_t ex = y*y/(2*kdy*kdy);
 123   return ka*TMath::Exp(-ex);
 124 }
 125
 126  Int_t AliGenPHOSlib::IpPion(TRandom *ran)
 127 {
 128 //                 particle composition pi+, pi0, pi-
 129 //
 130
 131      Float_t random = ran->Rndm();
 132
 133      if ( (3.*random)  < 1. )
 134      {
 135            return 211 ;
 136      }
 137      else
 138      {
 139        if ( (3.*random) >= 2.)
 140        {
 141           return -211 ;
 142        }
 143        else
 144        {
 145         return 111  ;
 146       }
 147     }
 148 }
 149
 150 //End Pions
 151 //======================================================================
 152 //    Pi 0 Flat Distribution
 153 //    Transverse momentum distribution PtPi0Flat
 154 //    Rapidity distribution YPi0Flat
 155 //    Particle distribution IdPi0Flat  111 (pi0)
 156 //
 157
 158 Double_t AliGenPHOSlib::PtPi0Flat(Double_t *px, Double_t *)
 159 {
 160 //     Pion transverse momentum flat distribution
 161
 162 return 1;
 163
 164 }
 165
 166 Double_t AliGenPHOSlib::YPi0Flat( Double_t *py, Double_t *)
 167 {
 168
 169 // pion y-distribution
 170 //
 171   return 1.;
 172 }
 173
 174  Int_t AliGenPHOSlib::IpPi0Flat(TRandom *)
 175 {
 176
 177 //                 particle composition pi0
 178 //
 179         return 111 ;
 180 }
 181 // End Pi0Flat
 182 //=============================================================
 183 //
 184  Double_t AliGenPHOSlib::PtScal(Double_t pt, Int_t np)
 185 {
 186 // Mt-scaling
 187 // Fonction for the calculation of the Pt distribution for a
 188 // given particle np, from the pion Pt distribution using the
 189 // mt scaling. This function was taken from AliGenMUONlib
 190 // aliroot version 3.01, and was extended for baryons
 191 // np = 1=>Pions 2=>Kaons 3=>Etas 4=>Omegas 5=>ETA' 6=>PHI
 192 //      7=>BARYONS-BARYONBARS
 193
 194   //    SCALING EN MASSE PAR RAPPORT A PTPI
 195   //    MASS                1=>PI,  2=>K, 3=>ETA, 4=>OMEGA,  5=>ETA',6=>PHI
 196   const Double_t khm[10] = {0.1396, 0.494,  0.547,    0.782,   0.957,   1.02,
 197   //    MASS               7=>BARYON-BARYONBAR
 198                                          0.938, 0. , 0., 0.};
 199   //     VALUE MESON/PI AT 5 GEV
 200   const Double_t kfmax[10]={1., 1., 1., 1., 1., 1., 1., 1., 1., 1.};
 201   np--;
 202   Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
 203   Double_t kfmax2=f5/kfmax[np];
 204   // PIONS
 205   Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
 206   Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
 207                                  (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2;
 208   return fmtscal*ptpion;
 209
 210 }
 211 // End Scaling
 212 //============================================================================
 213 //    K  A  O  N  S
 214  Double_t AliGenPHOSlib::PtKaon( Double_t *px, Double_t *)
 215 {
 216 //                kaon
 217 //                pt-distribution
 218 //____________________________________________________________
 219
 220   return PtScal(*px,2);  //  2==> Kaon in the PtScal function
 221 }
 222
 223  Double_t AliGenPHOSlib::YKaon( Double_t *py, Double_t *)
 224 {
 225 // y-distribution
 226 //____________________________________________________________
 227
 228   const Double_t ka    = 1000.;
 229   const Double_t kdy   = 4.;
 230
 231
 232   Double_t y=TMath::Abs(*py);
 233   //
 234   Double_t ex = y*y/(2*kdy*kdy);
 235   return ka*TMath::Exp(-ex);
 236 }
 237
 238  Int_t AliGenPHOSlib::IpKaon(TRandom *ran)
 239 {
 240 //                 particle composition
 241 //
 242
 243     Float_t random = ran->Rndm();
 244     Float_t random2 = ran->Rndm();
 245     if (random2 < 0.5)
 246     {
 247       if (random < 0.5) {
 248         return  321;   //   K+
 249       } else {
 250         return -321;   // K-
 251       }
 252     }
 253     else
 254     {
 255       if (random < 0.5) {
 256         return  130;   // K^0 short
 257       } else {
 258         return  310;   // K^0 long
 259       }
 260     }
 261 }
 262 // End Kaons
 263 //============================================================================
 264 //============================================================================
 265 //   E  T  A  S
 266  Double_t AliGenPHOSlib::PtEta( Double_t *px, Double_t *)
 267 {
 268 //                etas
 269 //                pt-distribution
 270 //____________________________________________________________
 271
 272   return PtScal(*px,3);  //  3==> Eta in the PtScal function
 273 }
 274
 275  Double_t AliGenPHOSlib::YEta( Double_t *py, Double_t *)
 276 {
 277 // y-distribution
 278 //____________________________________________________________
 279
 280   const Double_t ka    = 1000.;
 281   const Double_t kdy   = 4.;
 282
 283
 284   Double_t y=TMath::Abs(*py);
 285   //
 286   Double_t ex = y*y/(2*kdy*kdy);
 287   return ka*TMath::Exp(-ex);
 288 }
 289
 290  Int_t AliGenPHOSlib::IpEta(TRandom *)
 291 {
 292 //                 particle composition
 293 //
 294
 295         return  221;   //   eta
 296 }
 297 // End Etas
 298
 299 //======================================================================
 300 //    Eta Flat Distribution
 301 //    Transverse momentum distribution PtEtaFlat
 302 //    Rapidity distribution YEtaFlat
 303 //    Particle distribution IdEtaFlat  111 (pi0)
 304 //
 305
 306 Double_t AliGenPHOSlib::PtEtaFlat(Double_t *px, Double_t *)
 307 {
 308 //     Eta transverse momentum flat distribution
 309
 310   return 1;
 311
 312 }
 313
 314 Double_t AliGenPHOSlib::YEtaFlat( Double_t *py, Double_t *)
 315 {
 316 //
 317 // pion y-distribution
 318 //
 319   return 1.;
 320 }
 321
 322  Int_t AliGenPHOSlib::IpEtaFlat(TRandom *)
 323 {
 324 //
 325 //                 particle composition eta
 326 //
 327         return 221 ;
 328 }
 329 // End Pi0Flat
 330 //============================================================================
 331 //============================================================================
 332 //    O  M  E  G  A  S
 333  Double_t AliGenPHOSlib::PtOmega( Double_t *px, Double_t *)
 334 {
 335 // omegas
 336 //                pt-distribution
 337 //____________________________________________________________
 338
 339   return PtScal(*px,4);  //  4==> Omega in the PtScal function
 340 }
 341
 342  Double_t AliGenPHOSlib::YOmega( Double_t *py, Double_t *)
 343 {
 344 // y-distribution
 345 //____________________________________________________________
 346
 347   const Double_t ka    = 1000.;
 348   const Double_t kdy   = 4.;
 349
 350
 351   Double_t y=TMath::Abs(*py);
 352   //
 353   Double_t ex = y*y/(2*kdy*kdy);
 354   return ka*TMath::Exp(-ex);
 355 }
 356
 357  Int_t AliGenPHOSlib::IpOmega(TRandom *)
 358 {
 359 //                 particle composition
 360 //
 361
 362         return  223;   // Omega
 363 }
 364 // End Omega
 365 //============================================================================
 366 //============================================================================
 367 //    E  T  A  P  R  I  M  E
 368  Double_t AliGenPHOSlib::PtEtaprime( Double_t *px, Double_t *)
 369 {
 370 // etaprime
 371 //                pt-distribution
 372 //____________________________________________________________
 373
 374   return PtScal(*px,5);  //  5==> Etaprime in the PtScal function
 375 }
 376
 377  Double_t AliGenPHOSlib::YEtaprime( Double_t *py, Double_t *)
 378 {
 379 // y-distribution
 380 //____________________________________________________________
 381
 382   const Double_t ka    = 1000.;
 383   const Double_t kdy   = 4.;
 384
 385
 386   Double_t y=TMath::Abs(*py);
 387   //
 388   Double_t ex = y*y/(2*kdy*kdy);
 389   return ka*TMath::Exp(-ex);
 390 }
 391
 392  Int_t AliGenPHOSlib::IpEtaprime(TRandom *)
 393 {
 394 //                 particle composition
 395 //
 396
 397         return  331;   //   Etaprime
 398 }
 399 // End EtaPrime
 400 //===================================================================
 401 //============================================================================
 402 //    P  H  I   S
 403  Double_t AliGenPHOSlib::PtPhi( Double_t *px, Double_t *)
 404 {
 405 // phi
 406 //                pt-distribution
 407 //____________________________________________________________
 408
 409   return PtScal(*px,6);  //  6==> Phi in the PtScal function
 410 }
 411
 412  Double_t AliGenPHOSlib::YPhi( Double_t *py, Double_t *)
 413 {
 414 // y-distribution
 415 //____________________________________________________________
 416
 417   const Double_t ka    = 1000.;
 418   const Double_t kdy   = 4.;
 419
 420
 421   Double_t y=TMath::Abs(*py);
 422   //
 423   Double_t ex = y*y/(2*kdy*kdy);
 424   return ka*TMath::Exp(-ex);
 425 }
 426
 427  Int_t AliGenPHOSlib::IpPhi(TRandom *)
 428 {
 429 //                 particle composition
 430 //
 431
 432         return  333;   //   Phi
 433 }
 434 // End Phis
 435 //===================================================================
 436 //============================================================================
 437 //    B  A  R  Y  O  N  S  == protons, protonsbar, neutrons, and neutronsbars
 438  Double_t AliGenPHOSlib::PtBaryon( Double_t *px, Double_t *)
 439 {
 440 // baryons
 441 //                pt-distribution
 442 //____________________________________________________________
 443
 444   return PtScal(*px,7);  //  7==> Baryon in the PtScal function
 445 }
 446
 447  Double_t AliGenPHOSlib::YBaryon( Double_t *py, Double_t *)
 448 {
 449 // y-distribution
 450 //____________________________________________________________
 451
 452   const Double_t ka    = 1000.;
 453   const Double_t kdy   = 4.;
 454
 455
 456   Double_t y=TMath::Abs(*py);
 457   //
 458   Double_t ex = y*y/(2*kdy*kdy);
 459   return ka*TMath::Exp(-ex);
 460 }
 461
 462  Int_t AliGenPHOSlib::IpBaryon(TRandom *ran)
 463 {
 464 //                 particle composition
 465 //
 466
 467     Float_t random = ran->Rndm();
 468     Float_t random2 = ran->Rndm();
 469     if (random2 < 0.5)
 470     {
 471       if (random < 0.5) {
 472         return  2212;   //   p
 473       } else {
 474         return -2212;   // pbar
 475       }
 476     }
 477     else
 478     {
 479       if (random < 0.5) {
 480         return  2112;   // n
 481       } else {
 482         return -2112;   // n bar
 483       }
 484     }
 485 }
 486 // End Baryons
 487 //===================================================================
 488
 489
 490 typedef Double_t (*GenFunc) (Double_t*,  Double_t*);
 491  GenFunc AliGenPHOSlib::GetPt(Int_t param, const char* tname) const
 492 {
 493 // Return pinter to pT parameterisation
 494     GenFunc func;
 495
 496     switch (param)
 497     {
 498     case kPion:
 499         func=PtPion;
 500         break;
 501     case kPi0Flat:
 502         func=PtPi0Flat;
 503         break;
 504     case kKaon:
 505         func=PtKaon;
 506         break;
 507     case kEta:
 508         func=PtEta;
 509         break;
 510     case kEtaFlat:
 511         func=PtEtaFlat;
 512         break;
 513     case kOmega:
 514         func=PtOmega;
 515         break;
 516     case kEtaPrime:
 517         func=PtEtaprime;
 518         break;
 519     case kBaryon:
 520         func=PtBaryon;
 521         break;
 522     default:
 523         func=0;
 524         printf("<AliGenPHOSlib::GetPt> unknown parametrisationn");
 525     }
 526     return func;
 527 }
 528
 529  GenFunc AliGenPHOSlib::GetY(Int_t param, const char* tname) const
 530 {
 531 // Return pointer to Y parameterisation
 532     GenFunc func;
 533     switch (param)
 534     {
 535     case kPion:
 536         func=YPion;
 537         break;
 538     case kPi0Flat:
 539         func=YPi0Flat;
 540         break;
 541     case kKaon:
 542         func=YKaon;
 543         break;
 544     case kEta:
 545         func=YEta;
 546         break;
 547    case kEtaFlat:
 548         func=YEtaFlat;
 549         break;
 550     case kOmega:
 551         func=YOmega;
 552         break;
 553     case kEtaPrime:
 554         func=YEtaprime;
 555         break;
 556     case kPhi:
 557         func=YPhi;
 558         break;
 559     case kBaryon:
 560         func=YBaryon;
 561         break;
 562     default:
 563         func=0;
 564         printf("<AliGenPHOSlib::GetY> unknown parametrisationn");
 565     }
 566     return func;
 567 }
 568 typedef Int_t (*GenFuncIp) (TRandom *);
 569  GenFuncIp AliGenPHOSlib::GetIp(Int_t param,  const char* tname) const
 570 {
 571 // Return pointer to particle composition
 572     GenFuncIp func;
 573     switch (param)
 574     {
 575     case kPion:
 576         func=IpPion;
 577         break;
 578     case kPi0Flat:
 579         func=IpPi0Flat;
 580         break;
 581     case kKaon:
 582         func=IpKaon;
 583         break;
 584     case kEta:
 585         func=IpEta;
 586         break;
 587     case kEtaFlat:
 588         func=IpEtaFlat;
 589         break;
 590
 591     case kOmega:
 592         func=IpOmega;
 593         break;
 594     case kEtaPrime:
 595         func=IpEtaprime;
 596         break;
 597     case kPhi:
 598         func=IpPhi;
 599         break;
 600     case kBaryon:
 601         func=IpBaryon;
 602         break;
 603     default:
 604         func=0;
 605         printf("<AliGenPHOSlib::GetIp> unknown parametrisationn");
 606     }
 607     return func;
 608 }
 609