EVGEN/AliGenPMDlib.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 /* $Id$ */
  17
  18 #include <TMath.h>
  19 #include <TPDGCode.h>
  20
  21 #include "AliGenPMDlib.h"
  22
  23 ClassImp(AliGenPMDlib)
  24 //
  25 //  Neutral Pions
  26
  27 Double_t AliGenPMDlib::PtPi0(Double_t *px, Double_t */*dummy*/)
  28 {
  29 //
  30 //     PT-PARAMETERIZATION CDF, PRL 61(88) 1819
  31 //     POWER LAW FOR PT > 500 MEV
  32 //     MT SCALING BELOW (T=160 MEV)
  33 //
  34   const Double_t kp0 = 1.3;
  35   const Double_t kxn = 8.28;
  36   const Double_t kxlim=0.5;
  37   const Double_t kt=0.160;
  38   const Double_t kxmpi=0.139;
  39   const Double_t kb=1.;
  40   Double_t y, y1, xmpi2, ynorm, a;
  41   Double_t x=*px;
  42   //
  43   y1=TMath::Power(kp0/(kp0+kxlim),kxn);
  44   xmpi2=kxmpi*kxmpi;
  45   ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+xmpi2)/kt));
  46   a=ynorm/y1;
  47   if (x > kxlim)
  48     y=a*TMath::Power(kp0/(kp0+x),kxn);
  49   else
  50     y=kb*TMath::Exp(-sqrt(x*x+xmpi2)/kt);
  51   return y*x;
  52 }
  53
  54 //
  55 // y-distribution
  56 //
  57 Double_t AliGenPMDlib::YPi0( Double_t *py, Double_t */*dummy*/)
  58 {
  59   //
  60   // y parametrisation for pi0
  61   //
  62     const Double_t ka1    = 4913.;
  63     const Double_t ka2    = 1819.;
  64     const Double_t keta1  = 0.22;
  65     const Double_t keta2  = 3.66;
  66     const Double_t kdeta1 = 1.47;
  67     const Double_t kdeta2 = 1.51;
  68     Double_t y=TMath::Abs(*py);
  69     //
  70     Double_t ex1 = (y-keta1)*(y-keta1)/(2*kdeta1*kdeta1);
  71     Double_t ex2 = (y-keta2)*(y-keta2)/(2*kdeta2*kdeta2);
  72     return ka1*TMath::Exp(-ex1)+ka2*TMath::Exp(-ex2);
  73 }
  74
  75 //                 particle composition
  76 //
  77 Int_t AliGenPMDlib::IpPi0(TRandom *)
  78 {
  79 // Pi0
  80     return kPi0;
  81 }
  82
  83 //____________________________________________________________
  84 //
  85 // Mt-scaling
  86
  87 Double_t AliGenPMDlib::PtScal(Double_t pt, Int_t np)
  88 {
  89   //    SCALING EN MASSE PAR RAPPORT A PTPI
  90   //    MASS PI,K,ETA,RHO,OMEGA,ETA',PHI
  91   const Double_t khm[10] = {.13957,.493,.5488,.769,.7826,.958,1.02,0,0,0};
  92   //     VALUE MESON/PI AT 5 GEV
  93   const Double_t kfmax[10]={1.,0.3,0.55,1.0,1.0,1.0,1.0,0,0,0};
  94   np--;
  95   Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
  96   Double_t fmax2=f5/kfmax[np];
  97   // PIONS
  98   Double_t ptpion=100.*PtPi0(&pt, (Double_t*) 0);
  99   Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
 100                                  (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ fmax2;
 101   return fmtscal*ptpion;
 102 }
 103 //
 104 // kaon
 105 //
 106 //                pt-distribution
 107 //____________________________________________________________
 108
 109 Double_t AliGenPMDlib::PtEta( Double_t *px, Double_t */*dummy*/)
 110 {
 111 // Kaon pT
 112   return PtScal(*px,3);
 113 }
 114
 115 // y-distribution
 116 //____________________________________________________________
 117 Double_t AliGenPMDlib::YEta( Double_t *py, Double_t */*dummy*/)
 118 {
 119     //
 120     // y parametrisation for etas
 121     //
 122     const Double_t ka1    = 4913.;
 123     const Double_t ka2    = 1819.;
 124     const Double_t keta1  = 0.22;
 125     const Double_t keta2  = 3.66;
 126     const Double_t kdeta1 = 1.47;
 127     const Double_t kdeta2 = 1.51;
 128     Double_t y=TMath::Abs(*py);
 129     //
 130     Double_t ex1 = (y-keta1)*(y-keta1)/(2*kdeta1*kdeta1);
 131     Double_t ex2 = (y-keta2)*(y-keta2)/(2*kdeta2*kdeta2);
 132     return ka1*TMath::Exp(-ex1)+ka2*TMath::Exp(-ex2);
 133 }
 134
 135 //                 particle composition
 136 //
 137 Int_t AliGenPMDlib::IpEta(TRandom *)
 138 {
 139     return 221;
 140 }
 141
 142
 143 typedef Double_t (*GenFunc) (Double_t*,  Double_t*);
 144 GenFunc AliGenPMDlib::GetPt(Int_t param,  const char* /*tname*/) const
 145 {
 146 // Return pointer to pT parameterisation
 147     GenFunc func=NULL;
 148     switch (param)
 149     {
 150     case kPion:
 151         func=PtPi0;
 152         break;
 153     case kEta:
 154         func=PtEta;
 155         break;
 156     default:
 157         func=0;
 158         printf("<AliGenPMDlib::GetPt> unknown parametrisation\n");
 159     }
 160     return func;
 161 }
 162
 163 GenFunc AliGenPMDlib::GetY(Int_t param, const char* /*tname*/) const
 164 {
 165 // Return pointer to y- parameterisation
 166     GenFunc func=NULL;
 167     switch (param)
 168     {
 169     case kPion:
 170         func=YPi0;
 171         break;
 172     case kEta:
 173         func=YEta;
 174         break;
 175     default:
 176         func=0;
 177         printf("<AliGenPMDlib::GetY> unknown parametrisation\n");
 178     }
 179     return func;
 180
 181 }
 182 typedef Int_t (*GenFuncIp) (TRandom *);
 183 GenFuncIp AliGenPMDlib::GetIp(Int_t param,  const char* /*tname*/) const
 184 {
 185 // Return pointer to particle type parameterisation
 186     GenFuncIp func=NULL;
 187     switch (param)
 188     {
 189     case kPion:
 190         func=IpPi0;
 191         break;
 192     case kEta:
 193         func=IpEta;
 194         break;
 195     default:
 196         printf("<AliGenPMDlib::GetIp> unknown parametrisation\n");
 197     }
 198     return func;
 199 }
 200
 201
 202
 203