EVGEN/AliGenMUONlib.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.14  2002/02/22 17:26:43  morsch
  19 Eta and omega added.
  20
  21 Revision 1.13  2001/03/27 11:01:04  morsch
  22 Charm pt-distribution corrected. More realistic y-distribution for pi and K.
  23
  24 Revision 1.12  2001/03/09 13:01:41  morsch
  25 - enum constants for paramterisation type (particle family) moved to AliGen*lib.h
  26 - use AliGenGSIlib::kUpsilon, AliGenPHOSlib::kEtaPrime to access the constants
  27
  28 Revision 1.11  2000/11/30 07:12:50  alibrary
  29 Introducing new Rndm and QA classes
  30
  31 Revision 1.10  2000/06/29 21:08:27  morsch
  32 All paramatrisation libraries derive from the pure virtual base class AliGenLib.
  33 This allows to pass a pointer to a library directly to AliGenParam and avoids the
  34 use of function pointers in Config.C.
  35
  36 Revision 1.9  2000/06/14 15:20:56  morsch
  37 Include clean-up (IH)
  38
  39 Revision 1.8  2000/06/09 20:32:11  morsch
  40 All coding rule violations except RS3 corrected
  41
  42 Revision 1.7  2000/05/02 08:12:13  morsch
  43 Coding rule violations corrected.
  44
  45 Revision 1.6  1999/09/29 09:24:14  fca
  46 Introduction of the Copyright and cvs Log
  47
  48 */
  49
  50 // Library class for particle pt and y distributions used for
  51 // muon spectrometer simulations.
  52 // To be used with AliGenParam.
  53 // The following particle typed can be simulated:
  54 // pi, K, phi, omega, eta, J/Psi, Upsilon, charm and beauty mesons.
  55 //
  56 // andreas.morsch@cern.ch
  57 //
  58
  59 #include "TMath.h"
  60 #include "TRandom.h"
  61
  62 #include "AliGenMUONlib.h"
  63
  64 ClassImp(AliGenMUONlib)
  65 //
  66 //  Pions
  67 Double_t AliGenMUONlib::PtPion(Double_t *px, Double_t *dummy)
  68 {
  69 //
  70 //     PT-PARAMETERIZATION CDF, PRL 61(88) 1819
  71 //     POWER LAW FOR PT > 500 MEV
  72 //     MT SCALING BELOW (T=160 MEV)
  73 //
  74   const Double_t kp0 = 1.3;
  75   const Double_t kxn = 8.28;
  76   const Double_t kxlim=0.5;
  77   const Double_t kt=0.160;
  78   const Double_t kxmpi=0.139;
  79   const Double_t kb=1.;
  80   Double_t y, y1, xmpi2, ynorm, a;
  81   Double_t x=*px;
  82   //
  83   y1=TMath::Power(kp0/(kp0+kxlim),kxn);
  84   xmpi2=kxmpi*kxmpi;
  85   ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
  86   a=ynorm/y1;
  87   if (x > kxlim)
  88     y=a*TMath::Power(kp0/(kp0+x),kxn);
  89   else
  90     y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
  91   return y*x;
  92 }
  93 //
  94 // y-distribution
  95 //
  96 Double_t AliGenMUONlib::YPion( Double_t *py, Double_t *dummy)
  97 {
  98 // Pion y
  99   Double_t y=TMath::Abs(*py);
 100 /*
 101   const Double_t ka    = 7000.;
 102   const Double_t kdy   = 4.;
 103   Double_t ex = y*y/(2*kdy*kdy);
 104   return ka*TMath::Exp(-ex);
 105 */
 106   return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
 107
 108 }
 109 //                 particle composition
 110 //
 111 Int_t AliGenMUONlib::IpPion(TRandom *ran)
 112 {
 113 // Pion composition
 114     if (ran->Rndm() < 0.5) {
 115         return  211;
 116     } else {
 117         return -211;
 118     }
 119 }
 120
 121 //____________________________________________________________
 122 //
 123 // Mt-scaling
 124
 125 Double_t AliGenMUONlib::PtScal(Double_t pt, Int_t np)
 126 {
 127   //    SCALING EN MASSE PAR RAPPORT A PTPI
 128   //    MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
 129   const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
 130   //     VALUE MESON/PI AT 5 GEV
 131   const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
 132   np--;
 133   Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
 134   Double_t fmax2=f5/kfmax[np];
 135   // PIONS
 136   Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
 137   Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
 138                                  (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
 139   return fmtscal*ptpion;
 140 }
 141 //
 142 // kaon
 143 //
 144 //                pt-distribution
 145 //____________________________________________________________
 146 Double_t AliGenMUONlib::PtKaon( Double_t *px, Double_t *dummy)
 147 {
 148 // Kaon pT
 149   return PtScal(*px,2);
 150 }
 151
 152 // y-distribution
 153 //____________________________________________________________
 154 Double_t AliGenMUONlib::YKaon( Double_t *py, Double_t *dummy)
 155 {
 156 // Kaon y
 157   Double_t y=TMath::Abs(*py);
 158 /*
 159   const Double_t ka    = 1000.;
 160   const Double_t kdy   = 4.;
 161   //
 162   Double_t ex = y*y/(2*kdy*kdy);
 163   return ka*TMath::Exp(-ex);
 164 */
 165
 166   return 1.16526e+04+y*-3.79886e+03+y*y*4.31130e+02;
 167 }
 168
 169 //                 particle composition
 170 //
 171 Int_t AliGenMUONlib::IpKaon(TRandom *ran)
 172 {
 173 // Kaon composition
 174     if (ran->Rndm() < 0.5) {
 175         return  321;
 176     } else {
 177         return -321;
 178     }
 179 }
 180
 181 //                    J/Psi
 182 //
 183 //
 184 //                pt-distribution
 185 //____________________________________________________________
 186 Double_t AliGenMUONlib::PtJpsi( Double_t *px, Double_t *dummy)
 187 {
 188 // J/Psi pT
 189   const Double_t kpt0 = 4.;
 190   const Double_t kxn  = 3.6;
 191   Double_t x=*px;
 192   //
 193   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 194   return x/TMath::Power(pass1,kxn);
 195 }
 196 //
 197 //               y-distribution
 198 //____________________________________________________________
 199 Double_t AliGenMUONlib::YJpsi(Double_t *py, Double_t *dummy)
 200 {
 201 // J/psi y
 202   const Double_t ky0 = 4.;
 203   const Double_t kb=1.;
 204   Double_t yj;
 205   Double_t y=TMath::Abs(*py);
 206   //
 207   if (y < ky0)
 208     yj=kb;
 209   else
 210     yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
 211   return yj;
 212 }
 213 //                 particle composition
 214 //
 215 Int_t AliGenMUONlib::IpJpsi(TRandom *)
 216 {
 217 // J/Psi composition
 218     return 443;
 219 }
 220
 221 //                      Upsilon
 222 //
 223 //
 224 //                  pt-distribution
 225 //____________________________________________________________
 226 Double_t AliGenMUONlib::PtUpsilon( Double_t *px, Double_t *dummy )
 227 {
 228 // Upsilon pT
 229   const Double_t kpt0 = 5.3;
 230   const Double_t kxn  = 2.5;
 231   Double_t x=*px;
 232   //
 233   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 234   return x/TMath::Power(pass1,kxn);
 235 }
 236 //
 237 //                    y-distribution
 238 //
 239 //____________________________________________________________
 240 Double_t AliGenMUONlib::YUpsilon(Double_t *py, Double_t *dummy)
 241 {
 242 // Upsilon y
 243   const Double_t ky0 = 3.;
 244   const Double_t kb=1.;
 245   Double_t yu;
 246   Double_t y=TMath::Abs(*py);
 247   //
 248   if (y < ky0)
 249     yu=kb;
 250   else
 251     yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
 252   return yu;
 253 }
 254 //                 particle composition
 255 //
 256 Int_t AliGenMUONlib::IpUpsilon(TRandom *)
 257 {
 258 // y composition
 259     return 553;
 260 }
 261
 262 //
 263 //                        Phi
 264 //
 265 //
 266 //    pt-distribution (by scaling of pion distribution)
 267 //____________________________________________________________
 268 Double_t AliGenMUONlib::PtPhi( Double_t *px, Double_t *dummy)
 269 {
 270 // Phi pT
 271   return PtScal(*px,7);
 272 }
 273 //    y-distribution
 274 Double_t AliGenMUONlib::YPhi( Double_t *px, Double_t *dummy)
 275 {
 276 // Phi y
 277     Double_t *dum=0;
 278     return YJpsi(px,dum);
 279 }
 280 //                 particle composition
 281 //
 282 Int_t AliGenMUONlib::IpPhi(TRandom *)
 283 {
 284 // Phi composition
 285     return 333;
 286 }
 287
 288 //
 289 //                        omega
 290 //
 291 //
 292 //    pt-distribution (by scaling of pion distribution)
 293 //____________________________________________________________
 294 Double_t AliGenMUONlib::PtOmega( Double_t *px, Double_t *dummy)
 295 {
 296 // Omega pT
 297   return PtScal(*px,5);
 298 }
 299 //    y-distribution
 300 Double_t AliGenMUONlib::YOmega( Double_t *px, Double_t *dummy)
 301 {
 302 // Omega y
 303     Double_t *dum=0;
 304     return YJpsi(px,dum);
 305 }
 306 //                 particle composition
 307 //
 308 Int_t AliGenMUONlib::IpOmega(TRandom *)
 309 {
 310 // Omega composition
 311     return 223;
 312 }
 313
 314
 315 //
 316 //                        Eta
 317 //
 318 //
 319 //    pt-distribution (by scaling of pion distribution)
 320 //____________________________________________________________
 321 Double_t AliGenMUONlib::PtEta( Double_t *px, Double_t *dummy)
 322 {
 323 // Eta pT
 324   return PtScal(*px,3);
 325 }
 326 //    y-distribution
 327 Double_t AliGenMUONlib::YEta( Double_t *px, Double_t *dummy)
 328 {
 329 // Eta y
 330     Double_t *dum=0;
 331     return YJpsi(px,dum);
 332 }
 333 //                 particle composition
 334 //
 335 Int_t AliGenMUONlib::IpEta(TRandom *)
 336 {
 337 // Eta composition
 338     return 221;
 339 }
 340
 341 //
 342 //                        Charm
 343 //
 344 //
 345 //                    pt-distribution
 346 //____________________________________________________________
 347 Double_t AliGenMUONlib::PtCharm( Double_t *px, Double_t *dummy)
 348 {
 349 // Charm pT
 350   const Double_t kpt0 = 4.08;
 351   const Double_t kxn  = 9.40;
 352
 353   Double_t x=*px;
 354   //
 355   Double_t pass1 = 1.+(x/kpt0);
 356   return x/TMath::Power(pass1,kxn);
 357 }
 358 //                  y-distribution
 359 Double_t AliGenMUONlib::YCharm( Double_t *px, Double_t *dummy)
 360 {
 361 // Charm y
 362     Double_t *dum=0;
 363     return YJpsi(px,dum);
 364 }
 365
 366 Int_t AliGenMUONlib::IpCharm(TRandom *ran)
 367 {
 368 // Charm composition
 369     Float_t random;
 370     Int_t ip;
 371 //    411,421,431,4122
 372     random = ran->Rndm();
 373     if (random < 0.5) {
 374         ip=411;
 375     } else if (random < 0.75) {
 376         ip=421;
 377     } else if (random < 0.90) {
 378         ip=431;
 379     } else {
 380         ip=4122;
 381     }
 382     if (ran->Rndm() < 0.5) {ip=-ip;}
 383
 384     return ip;
 385 }
 386
 387
 388 //
 389 //                        Beauty
 390 //
 391 //
 392 //                    pt-distribution
 393 //____________________________________________________________
 394 Double_t AliGenMUONlib::PtBeauty( Double_t *px, Double_t *dummy)
 395 {
 396 // Beauty pT
 397   const Double_t kpt0 = 4.;
 398   const Double_t kxn  = 3.6;
 399   Double_t x=*px;
 400   //
 401   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 402   return x/TMath::Power(pass1,kxn);
 403 }
 404 //                     y-distribution
 405 Double_t AliGenMUONlib::YBeauty( Double_t *px, Double_t *dummy)
 406 {
 407 // Beauty y
 408     Double_t *dum=0;
 409     return YJpsi(px,dum);
 410 }
 411
 412 Int_t AliGenMUONlib::IpBeauty(TRandom *ran)
 413 {
 414 // Beauty Composition
 415     Float_t random;
 416     Int_t ip;
 417     random = ran->Rndm();
 418     if (random < 0.5) {
 419         ip=511;
 420     } else if (random < 0.75) {
 421         ip=521;
 422     } else if (random < 0.90) {
 423         ip=531;
 424     } else {
 425         ip=5122;
 426     }
 427     if (ran->Rndm() < 0.5) {ip=-ip;}
 428
 429     return ip;
 430 }
 431
 432 typedef Double_t (*GenFunc) (Double_t*,  Double_t*);
 433 GenFunc AliGenMUONlib::GetPt(Int_t param,  const char* tname) const
 434 {
 435 // Return pointer to pT parameterisation
 436     GenFunc func;
 437     switch (param)
 438     {
 439     case kPhi:
 440         func=PtPhi;
 441         break;
 442     case kOmega:
 443         func=PtOmega;
 444         break;
 445     case kEta:
 446         func=PtEta;
 447         break;
 448     case kJpsi:
 449         func=PtJpsi;
 450         break;
 451     case kUpsilon:
 452         func=PtUpsilon;
 453         break;
 454     case kCharm:
 455         func=PtCharm;
 456         break;
 457     case kBeauty:
 458         func=PtBeauty;
 459         break;
 460     case kPion:
 461         func=PtPion;
 462         break;
 463     case kKaon:
 464         func=PtKaon;
 465         break;
 466     default:
 467         func=0;
 468         printf("<AliGenMUONlib::GetPt> unknown parametrisation\n");
 469     }
 470     return func;
 471 }
 472
 473 GenFunc AliGenMUONlib::GetY(Int_t param, const char* tname) const
 474 {
 475 // Return pointer to y- parameterisation
 476     GenFunc func;
 477     switch (param)
 478     {
 479     case kPhi:
 480         func=YPhi;
 481         break;
 482     case kEta:
 483         func=YEta;
 484         break;
 485     case kOmega:
 486         func=YOmega;
 487         break;
 488     case kJpsi:
 489         func=YJpsi;
 490         break;
 491     case kUpsilon:
 492         func=YUpsilon;
 493         break;
 494     case kCharm:
 495         func=YCharm;
 496         break;
 497     case kBeauty:
 498         func=YBeauty;
 499         break;
 500     case kPion:
 501         func=YPion;
 502         break;
 503     case kKaon:
 504         func=YKaon;
 505         break;
 506     default:
 507         func=0;
 508         printf("<AliGenMUONlib::GetY> unknown parametrisation\n");
 509     }
 510     return func;
 511 }
 512 typedef Int_t (*GenFuncIp) (TRandom *);
 513 GenFuncIp AliGenMUONlib::GetIp(Int_t param,  const char* tname) const
 514 {
 515 // Return pointer to particle type parameterisation
 516     GenFuncIp func;
 517     switch (param)
 518     {
 519     case kPhi:
 520         func=IpPhi;
 521         break;
 522     case kEta:
 523         func=IpEta;
 524         break;
 525     case kOmega:
 526         func=IpOmega;
 527         break;
 528     case kJpsi:
 529         func=IpJpsi;
 530         break;
 531     case kUpsilon:
 532         func=IpUpsilon;
 533         break;
 534     case kCharm:
 535         func=IpCharm;
 536         break;
 537     case kBeauty:
 538         func=IpBeauty;
 539         break;
 540     case kPion:
 541         func=IpPion;
 542         break;
 543     case kKaon:
 544         func=IpKaon;
 545         break;
 546     default:
 547         func=0;
 548         printf("<AliGenMUONlib::GetIp> unknown parametrisation\n");
 549     }
 550     return func;
 551 }
 552
 553
 554
 555