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