Use AliGrayParticleModel.
[u/mrichter/AliRoot.git] / EVGEN / AliGenGrayParticles.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.3  2002/12/02 10:02:40  morsch
19 Corrections introduced by F. Silker:
20 - SetBetaSource
21 - Particle type according to charge.
22
23 Revision 1.2  2002/10/14 14:55:35  hristov
24 Merging the VirtualMC branch to the main development branch (HEAD)
25
26 Revision 1.1.2.1  2002/10/10 16:40:08  hristov
27 Updating VirtualMC to v3-09-02
28
29 Revision 1.1  2002/10/08 13:53:17  morsch
30 Gray particle generator, first commit.
31
32 */
33
34 /*
35   Generator for gray nucluons in pA interactions. 
36   Source is modelled by a relativistic Maxwell distributions.
37   Original code by  Ferenc Sikler  <sikler@rmki.kfki.hu>
38  */
39 #include "AliGenGrayParticles.h"
40 #include "AliGrayParticleModel.h"
41 #include "AliPDG.h"
42 #include <TDatabasePDG.h>
43
44  ClassImp(AliGenGrayParticles)
45     
46  AliGenGrayParticles::AliGenGrayParticles():AliGenerator(-1)
47 {
48 // Default constructor
49 }
50
51 AliGenGrayParticles::AliGenGrayParticles(Int_t npart)
52     :AliGenerator(npart)
53 {
54 // Constructor
55     fName  = "GrayParticles";
56     fTitle = "Generator for gray particles in pA collisions";
57     SetPmax();
58     SetTarget();
59     SetNominalCmsEnergy();
60     SetCharge();
61     SetTemperature();
62     SetBetaSource();
63 }
64
65 //____________________________________________________________
66 AliGenGrayParticles::~AliGenGrayParticles()
67 {
68 // Destructor
69 }
70
71
72 void AliGenGrayParticles::Init()
73 {
74   //
75   // Initialization
76   //
77     Float_t kMass  = TDatabasePDG::Instance()->GetParticle(kProton)->Mass();
78     fMomentum = fCMS/2. * fZTarget / fATarget;
79     fBeta     = fMomentum / TMath::Sqrt(kMass * kMass + fMomentum * fMomentum);
80 }
81
82
83 void AliGenGrayParticles::Generate()
84 {
85   //
86   // Generate one event
87   //
88     Float_t p[3];
89     Float_t origin[3] = {0., 0., 0.};
90     Float_t polar [3] = {0., 0., 0.};    
91     Int_t nt, i;
92     for(i = 0;i < fNpart; i++) {
93         Int_t kf;
94         if(fCharge==1) kf = kProton;
95                   else kf = kNeutron;
96         GenerateSlow(fCharge, fTemperature, fBetaSource, p);
97         
98         SetTrack(fTrackIt, -1, kf, p, origin, polar,
99                  0., kPNoProcess, nt, 1.);
100
101         KeepTrack(nt);
102     }
103 }
104
105
106
107
108 void AliGenGrayParticles::GenerateSlow(Int_t charge, Double_t T, Double_t beta, Float_t* q)
109 /* 
110    Emit a slow nucleon with "temperature" T [GeV], 
111    from a source moving with velocity beta         
112    Three-momentum [GeV/c] is given back in q[3]    
113 */
114
115 {
116  Double_t m, pmax, p, f, theta, phi;
117  
118  TDatabasePDG * pdg = TDatabasePDG::Instance();
119  const Double_t kMassProton  = pdg->GetParticle(kProton) ->Mass();
120  const Double_t kMassNeutron = pdg->GetParticle(kNeutron)->Mass();
121  
122  /* Select nucleon type */
123  if(charge == 0) m = kMassNeutron;
124  else m = kMassProton;
125
126  /* Momentum at maximum of Maxwell-distribution */
127
128  pmax = TMath::Sqrt(2*T*(T+sqrt(T*T+m*m)));
129
130  /* Try until proper momentum                                  */
131  /* for lack of primitive function of the Maxwell-distribution */
132  /* a brute force trial-accept loop, normalized at pmax        */
133
134  do
135  {
136      p = Rndm() * fPmax;
137      f = Maxwell(m, p, T) / Maxwell(m , pmax, T);
138  }
139  while(f < Rndm());
140
141  /* Spherical symmetric emission */
142  theta = TMath::ACos(2. * Rndm() - 1.);
143  phi   = 2. * TMath::Pi() * Rndm();
144
145  /* Determine momentum components in system of the moving source */
146  q[0] = p * TMath::Sin(theta) * TMath::Cos(phi);
147  q[1] = p * TMath::Sin(theta) * TMath::Sin(phi);
148  q[2] = p * TMath::Cos(theta);
149
150  /* Transform to system of the target nucleus                             */
151  /* beta is passed as negative, because the gray nucleons are slowed down */
152  Lorentz(m, -beta, q);
153
154  /* Transform to laboratory system */
155  Lorentz(m, fBeta, q);
156 }
157
158 Double_t AliGenGrayParticles::Maxwell(Double_t m, Double_t p, Double_t T)
159 {
160 /* Relativistic Maxwell-distribution */
161     Double_t ekin;
162     ekin = TMath::Sqrt(p*p+m*m)-m;
163     return (p*p * exp(-ekin/T));
164 }
165
166
167 void AliGenGrayParticles::Lorentz(Double_t m, Double_t beta, Float_t* q)
168 {
169 /* Lorentz transform in the direction of q[2] */
170  
171     Double_t gamma  = 1/sqrt(1-beta*beta); 
172     Double_t energy = sqrt(m*m + q[0]*q[0] + q[1]*q[1] + q[2]*q[2]);
173     q[2] = gamma * (q[2] + beta*energy);
174 }
175
176
177
178
179
180
181