New options kPi0Flat y kEtaFlat (Gustavo Conesa)
[u/mrichter/AliRoot.git] / EVGEN / AliGenPHOSlib.cxx
1
2 /**************************************************************************
3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4  *                                                                        *
5  * Author: The ALICE Off-line Project.                                    *
6  * Contributors are mentioned in the code where appropriate.              *
7  *                                                                        *
8  * Permission to use, copy, modify and distribute this software and its   *
9  * documentation strictly for non-commercial purposes is hereby granted   *
10  * without fee, provided that the above copyright notice appears in all   *
11  * copies and that both the copyright notice and this permission notice   *
12  * appear in the supporting documentation. The authors make no claims     *
13  * about the suitability of this software for any purpose. It is          *
14  * provided "as is" without express or implied warranty.                  *
15  **************************************************************************/
16
17 /*
18 $Log$
19 Revision 1.7  2001/03/09 13:01:41  morsch
20 - enum constants for paramterisation type (particle family) moved to AliGen*lib.h
21 - use AliGenGSIlib::kUpsilon, AliGenPHOSlib::kEtaPrime to access the constants
22
23 Revision 1.6  2000/11/30 07:12:50  alibrary
24 Introducing new Rndm and QA classes
25
26 Revision 1.5  2000/06/29 21:08:27  morsch
27 All paramatrisation libraries derive from the pure virtual base class AliGenLib.
28 This allows to pass a pointer to a library directly to AliGenParam and avoids the
29 use of function pointers in Config.C.
30
31 Revision 1.4  2000/06/14 15:21:05  morsch
32 Include clean-up (IH)
33
34 Revision 1.3  2000/06/09 20:32:54  morsch
35 All coding rule violations except RS3 corrected
36
37 Revision 1.2  1999/11/04 11:30:48  fca
38 Improve comments
39
40 Revision 1.1  1999/11/03 17:43:20  fca
41 New version from G.Martinez & A.Morsch
42
43 */
44
45 //======================================================================
46 //  AliGenPHOSlib class contains parameterizations of the
47 //  pion, kaon, eta, omega, etaprime, phi and baryon (proton, 
48 //  antiproton, neutron and anti-neutron) particles for the 
49 //  study of the neutral background in PHOS detector. 
50 //  These parameterizations are used by the 
51 //  AliGenParam  class:
52 //  AliGenParam(npar, param,  AliGenPHOSlib::GetPt(param),
53 //                            AliGenPHOSlib::GetY(param),
54 //                            AliGenPHOSlib::GetIp(param) )
55 //  param represents the particle to be simulated : 
56 //  Pion, Kaon, Eta, Omega, Etaprime, Phi or Baryon    
57 //  Pt distributions are calculated from the transverse mass scaling 
58 //  with Pions, using the PtScal function taken from AliGenMUONlib 
59 //  version aliroot 3.01
60 //
61 //     Gines MARTINEZ. Laurent APHECETCHE and Yves SCHUTZ
62 //      GPS @ SUBATECH,  Nantes , France  (October 1999)
63 //     http://www-subatech.in2p3.fr/~photons/subatech
64 //     martinez@subatech.in2p3.fr
65 //======================================================================
66
67 #include "TMath.h"
68 #include "TRandom.h"
69
70 #include "AliGenPHOSlib.h"
71
72 ClassImp(AliGenPHOSlib)
73
74 //======================================================================
75 //    P  I  O  N  S
76 //    (From GetPt, GetY and GetIp as param = Pion)
77 //    Transverse momentum distribution" PtPion 
78 //    Rapidity distribution YPion
79 //    Particle distribution IdPion  111, 211 and -211 (pi0, pi+ and pi-)
80 //
81  Double_t AliGenPHOSlib::PtPion(Double_t *px, Double_t *)
82 {
83 //     Pion transverse momentum distribtuion taken 
84 //     from AliGenMUONlib class, version 3.01 of aliroot
85 //     PT-PARAMETERIZATION CDF, PRL 61(88) 1819
86 //     POWER LAW FOR PT > 500 MEV
87 //     MT SCALING BELOW (T=160 MEV)
88 //
89   const Double_t kp0 = 1.3;
90   const Double_t kxn = 8.28;
91   const Double_t kxlim=0.5;
92   const Double_t kt=0.160;
93   const Double_t kxmpi=0.139;
94   const Double_t kb=1.;
95   Double_t y, y1, kxmpi2, ynorm, a;
96   Double_t x=*px;
97   //
98   y1=TMath::Power(kp0/(kp0+kxlim),kxn);
99   kxmpi2=kxmpi*kxmpi;
100   ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+kxmpi2)/kt));
101   a=ynorm/y1;
102   if (x > kxlim)
103     y=a*TMath::Power(kp0/(kp0+x),kxn);
104   else
105     y=kb*TMath::Exp(-sqrt(x*x+kxmpi2)/kt);
106   return y*x;
107 }
108  Double_t AliGenPHOSlib::YPion( Double_t *py, Double_t *)
109 {
110 //
111 // pion y-distribution
112 //
113
114   const Double_t ka    = 7000.;   
115   const Double_t kdy   = 4.;
116
117   Double_t y=TMath::Abs(*py);
118   //
119   Double_t ex = y*y/(2*kdy*kdy);
120   return ka*TMath::Exp(-ex);
121 }
122
123  Int_t AliGenPHOSlib::IpPion(TRandom *ran)
124 {
125 //                 particle composition pi+, pi0, pi-
126 //
127
128      Float_t random = ran->Rndm();
129
130      if ( (3.*random)  < 1. ) 
131      {
132            return 211 ;
133      } 
134      else
135      {  
136        if ( (3.*random) >= 2.)
137        {
138           return -211 ;
139        }
140        else 
141        {
142         return 111  ;
143       }
144     }
145 }
146
147 //End Pions
148 //======================================================================
149 //    Pi 0 Flat Distribution
150 //    Transverse momentum distribution PtPi0Flat
151 //    Rapidity distribution YPi0Flat
152 //    Particle distribution IdPi0Flat  111 (pi0)
153 //
154
155 Double_t AliGenPHOSlib::PtPi0Flat(Double_t *px, Double_t *)
156 {
157 //     Pion transverse momentum flat distribution 
158
159 return 1;
160
161 }
162
163 Double_t AliGenPHOSlib::YPi0Flat( Double_t *py, Double_t *)
164 {
165
166 // pion y-distribution
167 //
168   return 1.;
169 }
170
171  Int_t AliGenPHOSlib::IpPi0Flat(TRandom *)
172 {
173
174 //                 particle composition pi0
175 //
176         return 111 ;
177 }
178 // End Pi0Flat
179 //============================================================= 
180 //
181  Double_t AliGenPHOSlib::PtScal(Double_t pt, Int_t np)
182 {
183 // Mt-scaling
184 // Fonction for the calculation of the Pt distribution for a 
185 // given particle np, from the pion Pt distribution using the 
186 // mt scaling. This function was taken from AliGenMUONlib 
187 // aliroot version 3.01, and was extended for baryons
188 // np = 1=>Pions 2=>Kaons 3=>Etas 4=>Omegas 5=>ETA' 6=>PHI
189 //      7=>BARYONS-BARYONBARS
190
191   //    SCALING EN MASSE PAR RAPPORT A PTPI
192   //    MASS                1=>PI,  2=>K, 3=>ETA, 4=>OMEGA,  5=>ETA',6=>PHI
193   const Double_t khm[10] = {0.1396, 0.494,  0.547,    0.782,   0.957,   1.02, 
194   //    MASS               7=>BARYON-BARYONBAR  
195                                          0.938, 0. , 0., 0.};
196   //     VALUE MESON/PI AT 5 GEV
197   const Double_t kfmax[10]={1., 1., 1., 1., 1., 1., 1., 1., 1., 1.};
198   np--;
199   Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
200   Double_t kfmax2=f5/kfmax[np];
201   // PIONS
202   Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
203   Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
204                                  (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2;
205   return fmtscal*ptpion;
206
207 }
208 // End Scaling
209 //============================================================================
210 //    K  A  O  N  S
211  Double_t AliGenPHOSlib::PtKaon( Double_t *px, Double_t *)
212 {
213 //                kaon
214 //                pt-distribution
215 //____________________________________________________________
216
217   return PtScal(*px,2);  //  2==> Kaon in the PtScal function
218 }
219
220  Double_t AliGenPHOSlib::YKaon( Double_t *py, Double_t *)
221 {
222 // y-distribution
223 //____________________________________________________________
224
225   const Double_t ka    = 1000.;
226   const Double_t kdy   = 4.;
227
228
229   Double_t y=TMath::Abs(*py);
230   //
231   Double_t ex = y*y/(2*kdy*kdy);
232   return ka*TMath::Exp(-ex);
233 }
234
235  Int_t AliGenPHOSlib::IpKaon(TRandom *ran)
236 {
237 //                 particle composition
238 //
239
240     Float_t random = ran->Rndm();
241     Float_t random2 = ran->Rndm();
242     if (random2 < 0.5) 
243     {
244       if (random < 0.5) {       
245         return  321;   //   K+
246       } else {
247         return -321;   // K-
248       }
249     }
250     else
251     {  
252       if (random < 0.5) {       
253         return  130;   // K^0 short
254       } else {  
255         return  310;   // K^0 long
256       }
257     }
258 }
259 // End Kaons
260 //============================================================================
261 //============================================================================
262 //   E  T  A  S
263  Double_t AliGenPHOSlib::PtEta( Double_t *px, Double_t *)
264 {
265 //                etas
266 //                pt-distribution
267 //____________________________________________________________
268
269   return PtScal(*px,3);  //  3==> Eta in the PtScal function
270 }
271
272  Double_t AliGenPHOSlib::YEta( Double_t *py, Double_t *)
273 {
274 // y-distribution
275 //____________________________________________________________
276
277   const Double_t ka    = 1000.;
278   const Double_t kdy   = 4.;
279
280
281   Double_t y=TMath::Abs(*py);
282   //
283   Double_t ex = y*y/(2*kdy*kdy);
284   return ka*TMath::Exp(-ex);
285 }
286
287  Int_t AliGenPHOSlib::IpEta(TRandom *)
288 {
289 //                 particle composition
290 //
291
292         return  221;   //   eta
293 }
294 // End Etas
295
296 //======================================================================
297 //    Eta Flat Distribution
298 //    Transverse momentum distribution PtEtaFlat
299 //    Rapidity distribution YEtaFlat
300 //    Particle distribution IdEtaFlat  111 (pi0)
301 //
302
303 Double_t AliGenPHOSlib::PtEtaFlat(Double_t *px, Double_t *)
304 {
305 //     Eta transverse momentum flat distribution 
306
307   return 1;
308
309 }
310
311 Double_t AliGenPHOSlib::YEtaFlat( Double_t *py, Double_t *)
312 {
313 //
314 // pion y-distribution
315 //
316   return 1.;
317 }
318
319  Int_t AliGenPHOSlib::IpEtaFlat(TRandom *)
320 {
321 //
322 //                 particle composition eta
323 //
324         return 221 ;
325 }
326 // End Pi0Flat
327 //============================================================================
328 //============================================================================
329 //    O  M  E  G  A  S
330  Double_t AliGenPHOSlib::PtOmega( Double_t *px, Double_t *)
331 {
332 // omegas
333 //                pt-distribution
334 //____________________________________________________________
335
336   return PtScal(*px,4);  //  4==> Omega in the PtScal function
337 }
338
339  Double_t AliGenPHOSlib::YOmega( Double_t *py, Double_t *)
340 {
341 // y-distribution
342 //____________________________________________________________
343
344   const Double_t ka    = 1000.;
345   const Double_t kdy   = 4.;
346
347
348   Double_t y=TMath::Abs(*py);
349   //
350   Double_t ex = y*y/(2*kdy*kdy);
351   return ka*TMath::Exp(-ex);
352 }
353
354  Int_t AliGenPHOSlib::IpOmega(TRandom *)
355 {
356 //                 particle composition
357 //
358
359         return  223;   // Omega
360 }
361 // End Omega
362 //============================================================================
363 //============================================================================
364 //    E  T  A  P  R  I  M  E
365  Double_t AliGenPHOSlib::PtEtaprime( Double_t *px, Double_t *)
366 {
367 // etaprime
368 //                pt-distribution
369 //____________________________________________________________
370
371   return PtScal(*px,5);  //  5==> Etaprime in the PtScal function
372 }
373
374  Double_t AliGenPHOSlib::YEtaprime( Double_t *py, Double_t *)
375 {
376 // y-distribution
377 //____________________________________________________________
378
379   const Double_t ka    = 1000.;
380   const Double_t kdy   = 4.;
381
382
383   Double_t y=TMath::Abs(*py);
384   //
385   Double_t ex = y*y/(2*kdy*kdy);
386   return ka*TMath::Exp(-ex);
387 }
388
389  Int_t AliGenPHOSlib::IpEtaprime(TRandom *)
390 {
391 //                 particle composition
392 //
393
394         return  331;   //   Etaprime
395 }
396 // End EtaPrime
397 //===================================================================
398 //============================================================================
399 //    P  H  I   S
400  Double_t AliGenPHOSlib::PtPhi( Double_t *px, Double_t *)
401 {
402 // phi
403 //                pt-distribution
404 //____________________________________________________________
405
406   return PtScal(*px,6);  //  6==> Phi in the PtScal function
407 }
408
409  Double_t AliGenPHOSlib::YPhi( Double_t *py, Double_t *)
410 {
411 // y-distribution
412 //____________________________________________________________
413
414   const Double_t ka    = 1000.;
415   const Double_t kdy   = 4.;
416
417
418   Double_t y=TMath::Abs(*py);
419   //
420   Double_t ex = y*y/(2*kdy*kdy);
421   return ka*TMath::Exp(-ex);
422 }
423
424  Int_t AliGenPHOSlib::IpPhi(TRandom *)
425 {
426 //                 particle composition
427 //
428     
429         return  333;   //   Phi      
430 }
431 // End Phis
432 //===================================================================
433 //============================================================================
434 //    B  A  R  Y  O  N  S  == protons, protonsbar, neutrons, and neutronsbars
435  Double_t AliGenPHOSlib::PtBaryon( Double_t *px, Double_t *)
436 {
437 // baryons
438 //                pt-distribution
439 //____________________________________________________________
440
441   return PtScal(*px,7);  //  7==> Baryon in the PtScal function
442 }
443
444  Double_t AliGenPHOSlib::YBaryon( Double_t *py, Double_t *)
445 {
446 // y-distribution
447 //____________________________________________________________
448
449   const Double_t ka    = 1000.;
450   const Double_t kdy   = 4.;
451
452
453   Double_t y=TMath::Abs(*py);
454   //
455   Double_t ex = y*y/(2*kdy*kdy);
456   return ka*TMath::Exp(-ex);
457 }
458
459  Int_t AliGenPHOSlib::IpBaryon(TRandom *ran)
460 {
461 //                 particle composition
462 //
463
464     Float_t random = ran->Rndm();
465     Float_t random2 = ran->Rndm();
466     if (random2 < 0.5) 
467     {
468       if (random < 0.5) {       
469         return  2212;   //   p
470       } else {
471         return -2212;   // pbar
472       }
473     }
474     else
475     {  
476       if (random < 0.5) {       
477         return  2112;   // n
478       } else {  
479         return -2112;   // n bar
480       }
481     }
482 }
483 // End Baryons
484 //===================================================================
485
486
487 typedef Double_t (*GenFunc) (Double_t*,  Double_t*);
488  GenFunc AliGenPHOSlib::GetPt(Int_t param, const char* tname)
489 {
490 // Return pinter to pT parameterisation
491     GenFunc func;
492     
493     switch (param)
494     {
495     case kPion:     
496         func=PtPion;
497         break;
498     case kPi0Flat:     
499         func=PtPi0Flat;
500         break;
501     case kKaon:
502         func=PtKaon;
503         break;
504     case kEta:
505         func=PtEta;
506         break;
507     case kEtaFlat:
508         func=PtEtaFlat;
509         break;
510     case kOmega:
511         func=PtOmega;
512         break;
513     case kEtaPrime:
514         func=PtEtaprime;
515         break;
516     case kBaryon:
517         func=PtBaryon;
518         break;
519     default:
520         func=0;
521         printf("<AliGenPHOSlib::GetPt> unknown parametrisationn");
522     }
523     return func;
524 }
525
526  GenFunc AliGenPHOSlib::GetY(Int_t param, const char* tname)
527 {
528 // Return pointer to Y parameterisation
529     GenFunc func;
530     switch (param)
531     {
532     case kPion:
533         func=YPion;
534         break;
535     case kPi0Flat:
536         func=YPi0Flat;
537         break;
538     case kKaon:
539         func=YKaon;
540         break;
541     case kEta:
542         func=YEta;
543         break;
544    case kEtaFlat:
545         func=YEtaFlat;
546         break;
547     case kOmega:
548         func=YOmega;
549         break;
550     case kEtaPrime:
551         func=YEtaprime;
552         break;
553     case kPhi:
554         func=YPhi;
555         break;
556     case kBaryon:
557         func=YBaryon;
558         break;
559     default:
560         func=0;
561         printf("<AliGenPHOSlib::GetY> unknown parametrisationn");
562     }
563     return func;
564 }
565 typedef Int_t (*GenFuncIp) (TRandom *);
566  GenFuncIp AliGenPHOSlib::GetIp(Int_t param,  const char* tname)
567 {
568 // Return pointer to particle composition
569     GenFuncIp func;
570     switch (param)
571     {
572     case kPion:       
573         func=IpPion;
574         break;
575     case kPi0Flat:       
576         func=IpPi0Flat;
577         break;
578     case kKaon:
579         func=IpKaon;
580         break;
581     case kEta:
582         func=IpEta;
583         break;
584     case kEtaFlat:
585         func=IpEtaFlat;
586         break;
587         
588     case kOmega:
589         func=IpOmega;
590         break;
591     case kEtaPrime:
592         func=IpEtaprime;
593         break;
594     case kPhi:
595         func=IpPhi;
596         break;
597     case kBaryon:
598         func=IpBaryon;
599         break;
600     default:
601         func=0;
602         printf("<AliGenPHOSlib::GetIp> unknown parametrisationn");
603     }
604     return func;
605 }
606