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