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