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.8  2000/06/09 20:32:11  morsch
  19 All coding rule violations except RS3 corrected
  20
  21 Revision 1.7  2000/05/02 08:12:13  morsch
  22 Coding rule violations corrected.
  23
  24 Revision 1.6  1999/09/29 09:24:14  fca
  25 Introduction of the Copyright and cvs Log
  26
  27 */
  28
  29 #include "AliGenMUONlib.h"
  30 #include "AliMC.h"
  31
  32 ClassImp(AliGenMUONlib)
  33 //
  34 //  Pions
  35 Double_t AliGenMUONlib::PtPion(Double_t *px, Double_t *dummy)
  36 {
  37 //
  38 //     PT-PARAMETERIZATION CDF, PRL 61(88) 1819
  39 //     POWER LAW FOR PT > 500 MEV
  40 //     MT SCALING BELOW (T=160 MEV)
  41 //
  42   const Double_t kp0 = 1.3;
  43   const Double_t kxn = 8.28;
  44   const Double_t kxlim=0.5;
  45   const Double_t kt=0.160;
  46   const Double_t kxmpi=0.139;
  47   const Double_t kb=1.;
  48   Double_t y, y1, xmpi2, ynorm, a;
  49   Double_t x=*px;
  50   //
  51   y1=TMath::Power(kp0/(kp0+kxlim),kxn);
  52   xmpi2=kxmpi*kxmpi;
  53   ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
  54   a=ynorm/y1;
  55   if (x > kxlim)
  56     y=a*TMath::Power(kp0/(kp0+x),kxn);
  57   else
  58     y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
  59   return y*x;
  60 }
  61 //
  62 // y-distribution
  63 //
  64 Double_t AliGenMUONlib::YPion( Double_t *py, Double_t *dummy)
  65 {
  66 // Pion y
  67   const Double_t ka    = 7000.;
  68   const Double_t kdy   = 4.;
  69
  70   Double_t y=TMath::Abs(*py);
  71   //
  72   Double_t ex = y*y/(2*kdy*kdy);
  73   return ka*TMath::Exp(-ex);
  74 }
  75 //                 particle composition
  76 //
  77 Int_t AliGenMUONlib::IpPion()
  78 {
  79 // Pion composition
  80     Float_t random[1];
  81     gMC->Rndm(random,1);
  82     if (random[0] < 0.5) {
  83         return  211;
  84     } else {
  85         return -211;
  86     }
  87 }
  88
  89 //____________________________________________________________
  90 //
  91 // Mt-scaling
  92
  93 Double_t AliGenMUONlib::PtScal(Double_t pt, Int_t np)
  94 {
  95   //    SCALING EN MASSE PAR RAPPORT A PTPI
  96   //    MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
  97   const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
  98   //     VALUE MESON/PI AT 5 GEV
  99   const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
 100   np--;
 101   Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
 102   Double_t fmax2=f5/kfmax[np];
 103   // PIONS
 104   Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
 105   Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
 106                                  (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
 107   return fmtscal*ptpion;
 108 }
 109 //
 110 // kaon
 111 //
 112 //                pt-distribution
 113 //____________________________________________________________
 114 Double_t AliGenMUONlib::PtKaon( Double_t *px, Double_t *dummy)
 115 {
 116 // Kaon pT
 117   return PtScal(*px,2);
 118 }
 119
 120 // y-distribution
 121 //____________________________________________________________
 122 Double_t AliGenMUONlib::YKaon( Double_t *py, Double_t *dummy)
 123 {
 124 // Kaon y
 125   const Double_t ka    = 1000.;
 126   const Double_t kdy   = 4.;
 127
 128
 129   Double_t y=TMath::Abs(*py);
 130   //
 131   Double_t ex = y*y/(2*kdy*kdy);
 132   return ka*TMath::Exp(-ex);
 133 }
 134
 135 //                 particle composition
 136 //
 137 Int_t AliGenMUONlib::IpKaon()
 138 {
 139 // Kaon composition
 140     Float_t random[1];
 141     gMC->Rndm(random,1);
 142     if (random[0] < 0.5) {
 143         return  321;
 144     } else {
 145         return -321;
 146     }
 147 }
 148
 149 //                    J/Psi
 150 //
 151 //
 152 //                pt-distribution
 153 //____________________________________________________________
 154 Double_t AliGenMUONlib::PtJpsi( Double_t *px, Double_t *dummy)
 155 {
 156 // J/Psi pT
 157   const Double_t kpt0 = 4.;
 158   const Double_t kxn  = 3.6;
 159   Double_t x=*px;
 160   //
 161   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 162   return x/TMath::Power(pass1,kxn);
 163 }
 164 //
 165 //               y-distribution
 166 //____________________________________________________________
 167 Double_t AliGenMUONlib::YJpsi(Double_t *py, Double_t *dummy)
 168 {
 169 // J/psi y
 170   const Double_t ky0 = 4.;
 171   const Double_t kb=1.;
 172   Double_t yj;
 173   Double_t y=TMath::Abs(*py);
 174   //
 175   if (y < ky0)
 176     yj=kb;
 177   else
 178     yj=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
 179   return yj;
 180 }
 181 //                 particle composition
 182 //
 183 Int_t AliGenMUONlib::IpJpsi()
 184 {
 185 // J/Psi composition
 186     return 443;
 187 }
 188
 189 //                      Upsilon
 190 //
 191 //
 192 //                  pt-distribution
 193 //____________________________________________________________
 194 Double_t AliGenMUONlib::PtUpsilon( Double_t *px, Double_t *dummy )
 195 {
 196 // Upsilon pT
 197   const Double_t kpt0 = 5.3;
 198   const Double_t kxn  = 2.5;
 199   Double_t x=*px;
 200   //
 201   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 202   return x/TMath::Power(pass1,kxn);
 203 }
 204 //
 205 //                    y-distribution
 206 //
 207 //____________________________________________________________
 208 Double_t AliGenMUONlib::YUpsilon(Double_t *py, Double_t *dummy)
 209 {
 210 // Upsilon y
 211   const Double_t ky0 = 3.;
 212   const Double_t kb=1.;
 213   Double_t yu;
 214   Double_t y=TMath::Abs(*py);
 215   //
 216   if (y < ky0)
 217     yu=kb;
 218   else
 219     yu=kb*TMath::Exp(-(y-ky0)*(y-ky0)/2);
 220   return yu;
 221 }
 222 //                 particle composition
 223 //
 224 Int_t AliGenMUONlib::IpUpsilon()
 225 {
 226 // y composition
 227     return 553;
 228 }
 229
 230 //
 231 //                        Phi
 232 //
 233 //
 234 //    pt-distribution (by scaling of pion distribution)
 235 //____________________________________________________________
 236 Double_t AliGenMUONlib::PtPhi( Double_t *px, Double_t *dummy)
 237 {
 238 // Phi pT
 239   return PtScal(*px,7);
 240 }
 241 //    y-distribution
 242 Double_t AliGenMUONlib::YPhi( Double_t *px, Double_t *dummy)
 243 {
 244 // Phi y
 245     Double_t *dum=0;
 246     return YJpsi(px,dum);
 247 }
 248 //                 particle composition
 249 //
 250 Int_t AliGenMUONlib::IpPhi()
 251 {
 252 // Phi composition
 253     return 41;
 254 }
 255
 256 //
 257 //                        Charm
 258 //
 259 //
 260 //                    pt-distribution
 261 //____________________________________________________________
 262 Double_t AliGenMUONlib::PtCharm( Double_t *px, Double_t *dummy)
 263 {
 264 // Charm pT
 265   const Double_t kpt0 = 4.08;
 266   const Double_t kxn  = 9.40;
 267   Double_t x=*px;
 268   //
 269   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 270   return x/TMath::Power(pass1,kxn);
 271 }
 272 //                  y-distribution
 273 Double_t AliGenMUONlib::YCharm( Double_t *px, Double_t *dummy)
 274 {
 275 // Charm y
 276     Double_t *dum=0;
 277     return YJpsi(px,dum);
 278 }
 279
 280 Int_t AliGenMUONlib::IpCharm()
 281 {
 282 // Charm composition
 283     Float_t random[2];
 284     Int_t ip;
 285 //    411,421,431,4122
 286     gMC->Rndm(random,2);
 287     if (random[0] < 0.5) {
 288         ip=411;
 289     } else if (random[0] < 0.75) {
 290         ip=421;
 291     } else if (random[0] < 0.90) {
 292         ip=431;
 293     } else {
 294         ip=4122;
 295     }
 296     if (random[1] < 0.5) {ip=-ip;}
 297
 298     return ip;
 299 }
 300
 301
 302 //
 303 //                        Beauty
 304 //
 305 //
 306 //                    pt-distribution
 307 //____________________________________________________________
 308 Double_t AliGenMUONlib::PtBeauty( Double_t *px, Double_t *dummy)
 309 {
 310 // Beauty pT
 311   const Double_t kpt0 = 4.;
 312   const Double_t kxn  = 3.6;
 313   Double_t x=*px;
 314   //
 315   Double_t pass1 = 1.+(x/kpt0)*(x/kpt0);
 316   return x/TMath::Power(pass1,kxn);
 317 }
 318 //                     y-distribution
 319 Double_t AliGenMUONlib::YBeauty( Double_t *px, Double_t *dummy)
 320 {
 321 // Beauty y
 322     Double_t *dum=0;
 323     return YJpsi(px,dum);
 324 }
 325
 326 Int_t AliGenMUONlib::IpBeauty()
 327 {
 328 // Beauty Composition
 329     Float_t random[2];
 330     Int_t ip;
 331     gMC->Rndm(random,2);
 332     if (random[0] < 0.5) {
 333         ip=511;
 334     } else if (random[0] < 0.75) {
 335         ip=521;
 336     } else if (random[0] < 0.90) {
 337         ip=531;
 338     } else {
 339         ip=5122;
 340     }
 341     if (random[1] < 0.5) {ip=-ip;}
 342
 343     return ip;
 344 }
 345
 346 typedef Double_t (*GenFunc) (Double_t*,  Double_t*);
 347 GenFunc AliGenMUONlib::GetPt(Param_t param)
 348 {
 349 // Return pointer to pT parameterisation
 350     GenFunc func;
 351     switch (param)
 352     {
 353     case phi_p:
 354         func=PtPhi;
 355         break;
 356     case jpsi_p:
 357         func=PtJpsi;
 358         break;
 359     case upsilon_p:
 360         func=PtUpsilon;
 361         break;
 362     case charm_p:
 363         func=PtCharm;
 364         break;
 365     case beauty_p:
 366         func=PtBeauty;
 367         break;
 368     case pion_p:
 369         func=PtPion;
 370         break;
 371     case kaon_p:
 372         func=PtKaon;
 373         break;
 374     default:
 375         func=0;
 376         printf("<AliGenMUONlib::GetPt> unknown parametrisation\n");
 377     }
 378     return func;
 379 }
 380
 381 GenFunc AliGenMUONlib::GetY(Param_t param)
 382 {
 383 // Return pointer to y- parameterisation
 384     GenFunc func;
 385     switch (param)
 386     {
 387     case phi_p:
 388         func=YPhi;
 389         break;
 390     case jpsi_p:
 391         func=YJpsi;
 392         break;
 393     case upsilon_p:
 394         func=YUpsilon;
 395         break;
 396     case charm_p:
 397         func=YCharm;
 398         break;
 399     case beauty_p:
 400         func=YBeauty;
 401         break;
 402     case pion_p:
 403         func=YPion;
 404         break;
 405     case kaon_p:
 406         func=YKaon;
 407         break;
 408     default:
 409         func=0;
 410         printf("<AliGenMUONlib::GetY> unknown parametrisation\n");
 411     }
 412     return func;
 413 }
 414 typedef Int_t (*GenFuncIp) ();
 415 GenFuncIp AliGenMUONlib::GetIp(Param_t param)
 416 {
 417 // Return pointer to particle type parameterisation
 418     GenFuncIp func;
 419     switch (param)
 420     {
 421     case phi_p:
 422         func=IpPhi;
 423         break;
 424     case jpsi_p:
 425         func=IpJpsi;
 426         break;
 427     case upsilon_p:
 428         func=IpUpsilon;
 429         break;
 430     case charm_p:
 431         func=IpCharm;
 432         break;
 433     case beauty_p:
 434         func=IpBeauty;
 435         break;
 436     case pion_p:
 437         func=IpPion;
 438         break;
 439     case kaon_p:
 440         func=IpKaon;
 441         break;
 442     default:
 443         func=0;
 444         printf("<AliGenMUONlib::GetIp> unknown parametrisation\n");
 445     }
 446     return func;
 447 }
 448
 449
 450
 451