EVGEN/AliGenPHOSlib.cxx

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