1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 Revision 1.2 1999/11/04 11:30:48 fca
21 Revision 1.1 1999/11/03 17:43:20 fca
22 New version from G.Martinez & A.Morsch
26 //======================================================================
27 // AliGenPHOSlib class contains parameterizations of the
28 // pion, kaon, eta, omega, etaprime, phi and baryon (proton,
29 // antiproton, neutron and anti-neutron) particles for the
30 // study of the neutral background in PHOS detector.
31 // These parameterizations are used by the
33 // AliGenParam(npar, param, AliGenPHOSlib::GetPt(param),
34 // AliGenPHOSlib::GetY(param),
35 // AliGenPHOSlib::GetIp(param) )
36 // param represents the particle to be simulated :
37 // Pion, Kaon, Eta, Omega, Etaprime, Phi or Baryon
38 // Pt distributions are calculated from the transverse mass scaling
39 // with Pions, using the PtScal function taken from AliGenMUONlib
40 // version aliroot 3.01
42 // Gines MARTINEZ. Laurent APHECETCHE and Yves SCHUTZ
43 // GPS @ SUBATECH, Nantes , France (October 1999)
44 // http://www-subatech.in2p3.fr/~photons/subatech
45 // martinez@subatech.in2p3.fr
46 //======================================================================
48 #include "AliGenPHOSlib.h"
52 ClassImp(AliGenPHOSlib)
54 //======================================================================
56 // (From GetPt, GetY and GetIp as param = Pion)
57 // Transverse momentum distribution" PtPion
58 // Rapidity distribution YPion
59 // Particle distribution IdPion 111, 211 and -211 (pi0, pi+ and pi-)
61 Double_t AliGenPHOSlib::PtPion(Double_t *px, Double_t *)
63 // Pion transverse momentum distribtuion taken
64 // from AliGenMUONlib class, version 3.01 of aliroot
65 // PT-PARAMETERIZATION CDF, PRL 61(88) 1819
66 // POWER LAW FOR PT > 500 MEV
67 // MT SCALING BELOW (T=160 MEV)
69 const Double_t kp0 = 1.3;
70 const Double_t kxn = 8.28;
71 const Double_t kxlim=0.5;
72 const Double_t kt=0.160;
73 const Double_t kxmpi=0.139;
75 Double_t y, y1, kxmpi2, ynorm, a;
78 y1=TMath::Power(kp0/(kp0+kxlim),kxn);
80 ynorm=kb*(TMath::Exp(-sqrt(kxlim*kxlim+kxmpi2)/kt));
83 y=a*TMath::Power(kp0/(kp0+x),kxn);
85 y=kb*TMath::Exp(-sqrt(x*x+kxmpi2)/kt);
88 Double_t AliGenPHOSlib::YPion( Double_t *py, Double_t *)
91 // pion y-distribution
94 const Double_t ka = 7000.;
95 const Double_t kdy = 4.;
97 Double_t y=TMath::Abs(*py);
99 Double_t ex = y*y/(2*kdy*kdy);
100 return ka*TMath::Exp(-ex);
103 Int_t AliGenPHOSlib::IpPion()
105 // particle composition pi+, pi0, pi-
111 if ( (3.*random[0]) < 1. )
117 if ( (3.*random[0]) >= 2.)
128 //=============================================================
130 Double_t AliGenPHOSlib::PtScal(Double_t pt, Int_t np)
133 // Fonction for the calculation of the Pt distribution for a
134 // given particle np, from the pion Pt distribution using the
135 // mt scaling. This function was taken from AliGenMUONlib
136 // aliroot version 3.01, and was extended for baryons
137 // np = 1=>Pions 2=>Kaons 3=>Etas 4=>Omegas 5=>ETA' 6=>PHI
138 // 7=>BARYONS-BARYONBARS
140 // SCALING EN MASSE PAR RAPPORT A PTPI
141 // MASS 1=>PI, 2=>K, 3=>ETA, 4=>OMEGA, 5=>ETA',6=>PHI
142 const Double_t khm[10] = {0.1396, 0.494, 0.547, 0.782, 0.957, 1.02,
143 // MASS 7=>BARYON-BARYONBAR
145 // VALUE MESON/PI AT 5 GEV
146 const Double_t kfmax[10]={1., 1., 1., 1., 1., 1., 1., 1., 1., 1.};
148 Double_t f5=TMath::Power(((sqrt(100.018215)+2.)/(sqrt(100.+khm[np]*khm[np])+2.0)),12.3);
149 Double_t kfmax2=f5/kfmax[np];
151 Double_t ptpion=100.*PtPion(&pt, (Double_t*) 0);
152 Double_t fmtscal=TMath::Power(((sqrt(pt*pt+0.018215)+2.)/
153 (sqrt(pt*pt+khm[np]*khm[np])+2.0)),12.3)/ kfmax2;
154 return fmtscal*ptpion;
157 //============================================================================
159 Double_t AliGenPHOSlib::PtKaon( Double_t *px, Double_t *)
163 //____________________________________________________________
165 return PtScal(*px,2); // 2==> Kaon in the PtScal function
168 Double_t AliGenPHOSlib::YKaon( Double_t *py, Double_t *)
171 //____________________________________________________________
173 const Double_t ka = 1000.;
174 const Double_t kdy = 4.;
177 Double_t y=TMath::Abs(*py);
179 Double_t ex = y*y/(2*kdy*kdy);
180 return ka*TMath::Exp(-ex);
183 Int_t AliGenPHOSlib::IpKaon()
185 // particle composition
188 Float_t random[1],random2[1];
190 gMC->Rndm(random2,1);
191 if (random2[0] < 0.5)
193 if (random[0] < 0.5) {
201 if (random[0] < 0.5) {
202 return 130; // K^0 short
204 return 310; // K^0 long
209 //============================================================================
210 //============================================================================
212 Double_t AliGenPHOSlib::PtEta( Double_t *px, Double_t *)
216 //____________________________________________________________
218 return PtScal(*px,3); // 3==> Eta in the PtScal function
221 Double_t AliGenPHOSlib::YEta( Double_t *py, Double_t *)
224 //____________________________________________________________
226 const Double_t ka = 1000.;
227 const Double_t kdy = 4.;
230 Double_t y=TMath::Abs(*py);
232 Double_t ex = y*y/(2*kdy*kdy);
233 return ka*TMath::Exp(-ex);
236 Int_t AliGenPHOSlib::IpEta()
238 // particle composition
244 //============================================================================
245 //============================================================================
247 Double_t AliGenPHOSlib::PtOmega( Double_t *px, Double_t *)
251 //____________________________________________________________
253 return PtScal(*px,4); // 4==> Omega in the PtScal function
256 Double_t AliGenPHOSlib::YOmega( Double_t *py, Double_t *)
259 //____________________________________________________________
261 const Double_t ka = 1000.;
262 const Double_t kdy = 4.;
265 Double_t y=TMath::Abs(*py);
267 Double_t ex = y*y/(2*kdy*kdy);
268 return ka*TMath::Exp(-ex);
271 Int_t AliGenPHOSlib::IpOmega()
273 // particle composition
279 //============================================================================
280 //============================================================================
282 Double_t AliGenPHOSlib::PtEtaprime( Double_t *px, Double_t *)
286 //____________________________________________________________
288 return PtScal(*px,5); // 5==> Etaprime in the PtScal function
291 Double_t AliGenPHOSlib::YEtaprime( Double_t *py, Double_t *)
294 //____________________________________________________________
296 const Double_t ka = 1000.;
297 const Double_t kdy = 4.;
300 Double_t y=TMath::Abs(*py);
302 Double_t ex = y*y/(2*kdy*kdy);
303 return ka*TMath::Exp(-ex);
306 Int_t AliGenPHOSlib::IpEtaprime()
308 // particle composition
311 return 331; // Etaprime
314 //===================================================================
315 //============================================================================
317 Double_t AliGenPHOSlib::PtPhi( Double_t *px, Double_t *)
321 //____________________________________________________________
323 return PtScal(*px,6); // 6==> Phi in the PtScal function
326 Double_t AliGenPHOSlib::YPhi( Double_t *py, Double_t *)
329 //____________________________________________________________
331 const Double_t ka = 1000.;
332 const Double_t kdy = 4.;
335 Double_t y=TMath::Abs(*py);
337 Double_t ex = y*y/(2*kdy*kdy);
338 return ka*TMath::Exp(-ex);
341 Int_t AliGenPHOSlib::IpPhi()
343 // particle composition
349 //===================================================================
350 //============================================================================
351 // B A R Y O N S == protons, protonsbar, neutrons, and neutronsbars
352 Double_t AliGenPHOSlib::PtBaryon( Double_t *px, Double_t *)
356 //____________________________________________________________
358 return PtScal(*px,7); // 7==> Baryon in the PtScal function
361 Double_t AliGenPHOSlib::YBaryon( Double_t *py, Double_t *)
364 //____________________________________________________________
366 const Double_t ka = 1000.;
367 const Double_t kdy = 4.;
370 Double_t y=TMath::Abs(*py);
372 Double_t ex = y*y/(2*kdy*kdy);
373 return ka*TMath::Exp(-ex);
376 Int_t AliGenPHOSlib::IpBaryon()
378 // particle composition
381 Float_t random[1],random2[1];
383 gMC->Rndm(random2,1);
384 if (random2[0] < 0.5)
386 if (random[0] < 0.5) {
389 return -2212; // pbar
394 if (random[0] < 0.5) {
397 return -2112; // n bar
402 //===================================================================
405 typedef Double_t (*GenFunc) (Double_t*, Double_t*);
406 GenFunc AliGenPHOSlib::GetPt(Param_t param)
408 // Return pinter to pT parameterisation
433 printf("<AliGenPHOSlib::GetPt> unknown parametrisationn");
438 GenFunc AliGenPHOSlib::GetY(Param_t param)
440 // Return pointer to Y parameterisation
467 printf("<AliGenPHOSlib::GetY> unknown parametrisationn");
471 typedef Int_t (*GenFuncIp) ();
472 GenFuncIp AliGenPHOSlib::GetIp(Param_t param)
474 // Return pointer to particle composition
502 printf("<AliGenPHOSlib::GetIp> unknown parametrisationn");