pt and y-parameterisations for PMD physics simulation.
[u/mrichter/AliRoot.git] / 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