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