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