Fix bug in building local list of valid files.
[u/mrichter/AliRoot.git] / TEvtGen / AliGenEvtGen.cxx
CommitLineData
da0e9ce3 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// Class to generate decays of particles generated by a //
17// previous generator. It works as a generator, but pratically it //
18// performs only decays. It works with this scheme: first loops over //
19// particles on the stack, selects those to be decayed, decays them //
20// and then pushes the decay products on the stack. //
21// //
22// Giuseppe E. Bruno & Fiorella Fionda //
23// (Giuseppe.Bruno@ba.infn.it) (Fiorella.Fionda@ba.infn.it) //
24///////////////////////////////////////////////////////////////////////////
25
26#include "AliStack.h"
27#include "AliGenEvtGen.h"
28#include "AliRun.h"
51998ca8 29#include "AliLog.h"
da0e9ce3 30#include <TParticle.h>
31
32ClassImp(AliGenEvtGen)
33///////////////////////////////////////////////////////////////////////////
34AliGenEvtGen::AliGenEvtGen():
35 fStack(0x0),
36 fDecayer(0x0),
37 fForceDecay(kAll),
38 fSwitchOff(kBeautyPart),
39 fUserDecay(kFALSE),
40 fUserDecayTablePath(0x0)
41 {
42 //
43 // Default Construction
44 //
45 }
46///////////////////////////////////////////////////////////////////////////////////////////
47AliGenEvtGen::~AliGenEvtGen()
48 {
49 //
50 // Standard Destructor
51 //
52 if(fStack) {delete fStack;}
53 fStack = 0;
54 if(fDecayer) {delete fDecayer;}
55 fDecayer = 0;
56 if(fUserDecayTablePath) {delete fUserDecayTablePath;}
57 fUserDecayTablePath = 0;
58 }
59///////////////////////////////////////////////////////////////////////////////////////////
60
61void AliGenEvtGen::Init()
62 {
63 //
64 // Standard AliGenerator Initializer - no input
65 // 1) initialize EvtGen with default decay and particle table
66 // 2) set the decay mode to force particle
67 // 3) set a user decay table if defined
68 //
69 if(fDecayer)
70 {
71 AliWarning("AliGenEvtGen already initialized!!!");
72 return;
73 }
74 fDecayer = new AliDecayerEvtGen();
75 fDecayer->Init(); //read the default decay table DECAY.DEC and particle table
76
77 //if is set a decay mode: default decay mode is kAll
78 fDecayer->SetForceDecay(fForceDecay);
79 fDecayer->ForceDecay();
80
81 //if is defined a user decay table
82 if(fUserDecay)
83 {
84 fDecayer->SetDecayTablePath(fUserDecayTablePath);
85 fDecayer->ReadDecayTable();
86 }
87 }
88
89/////////////////////////////////////////////////////////////////////////////////////////////
90void AliGenEvtGen::Generate()
91 {
92 //
93 //Generate method - Input none - Output none
94 //For each event:
95 //1)return the stack of the previous generator and select particles to be decayed by EvtGen
96 //2)decay particles selected and put the decay products on the stack
97 //
98 //
99 Float_t polar[3]= {0,0,0}; // Polarisation of daughter particles
100 Float_t origin0[3]; // Origin of the parent particle
101 Float_t pc[3], och[3]; // Momentum and origin of the children particles from EvtGen
102 Int_t nt;
103 Float_t tof;
104 Int_t nPrimsPythia;
105 TLorentzVector *mom=new TLorentzVector();
106 static TClonesArray *particles;
107 if(!particles) particles = new TClonesArray("TParticle",1000);
108 fStack = AliRunLoader::Instance()->Stack();
109 if(!fStack) {Info("Generate","Error: No stack found!"); return;}
110 nPrimsPythia = fStack->GetNprimary();
111 AliDebug(1,Form("nPrimsPythia = %d \n",nPrimsPythia));
112 for (Int_t iTrack = 0; iTrack < nPrimsPythia; ++iTrack) {
113 TParticle *part = fStack->Particle(iTrack);
114 Int_t pdg=part->GetPdgCode();
115
116 AliDebug(1,Form("GetFlavour = %d e pdg = %d \n",GetFlavour(pdg),pdg));
117
118 switch(fSwitchOff)
119 {
120 case kAllPart:
121 break;
122 case kBeautyPart:
123 if(GetFlavour(pdg)!=5) continue;
124 break;
125 case kCharmPart:
126 if(GetFlavour(pdg)!=4) continue;
127 break;
128 }
129
130 //check if particle is already decayed by Pythia
131 if(part->GetStatusCode() != 1 || part->GetNDaughters()>0)
132 {
f2f8bbb9 133 AliDebug(1,Form("Attention: particle %d is already decayed by Pythia!",pdg));
da0e9ce3 134 continue;
135 }
136
137 part->SetStatusCode(11); //Set particle as decayed : change the status code
b05e386a 138 part->SetBit(kDoneBit);
139 part->ResetBit(kTransportBit);
140
da0e9ce3 141 mom->SetPxPyPzE(part->Px(),part->Py(),part->Pz(),part->Energy());
142 Int_t np;
143
144 do{
145 fDecayer->Decay(part->GetPdgCode(),mom);
146 np = fDecayer->ImportParticles(particles);
147 }while(np<0);
148
149 Int_t* trackIt = new Int_t[np];
150 Int_t* pParent = new Int_t[np];
151 AliDebug(1,Form("np = %d \n",np));
152
153 for (int i = 0; i < np; i++) {
154 pParent[i] = -1;
155 trackIt[i] = 0;
156 }
157 //select trackable particle
158 if (np >1) {
008d5aae 159 TParticle* iparticle = 0;//(TParticle *) particles->At(0);//parent particle
da0e9ce3 160 for (int i = 1; i<np ; i++) {
161 iparticle = (TParticle*) particles->At(i);
162 Int_t ks = iparticle->GetStatusCode();
163
164 //track last decay products
165 if(ks==1) trackIt[i]=1;
166
167 }//decay particles loop
168
169 }// if decay products
170
171 origin0[0]=part->Vx(); //[cm]
172 origin0[1]=part->Vy(); //[cm]
173 origin0[2]=part->Vz(); //[cm]
174 //
175 // Put decay products on the stack
176 //
177 for (int i = 1; i < np; i++) {
178 TParticle* iparticle = (TParticle *) particles->At(i);
179 Int_t kf = iparticle->GetPdgCode();
180 Int_t ksc = iparticle->GetStatusCode();
181 Int_t jpa = iparticle->GetFirstMother()-1; //jpa = 0 for daughters of beauty particles
182 Int_t iparent = (jpa > 0) ? pParent[jpa] : iTrack;
183
774ceaaf 184 och[0] = origin0[0]+iparticle->Vx(); //[cm]
185 och[1] = origin0[1]+iparticle->Vy(); //[cm]
186 och[2] = origin0[2]+iparticle->Vz(); //[cm]
da0e9ce3 187 pc[0] = iparticle->Px(); //[GeV/c]
188 pc[1] = iparticle->Py(); //[GeV/c]
189 pc[2] = iparticle->Pz(); //[GeV/c]
774ceaaf 190 tof = part->T()+iparticle->T();
da0e9ce3 191
192 AliDebug(1,Form("FirstMother = %d e indicePart = %d e pdg = %d \n",jpa,i,kf));
193
194 PushTrack(trackIt[i], iparent, kf, pc, och, polar,tof, kPDecay, nt, 1., ksc);
195 if(trackIt[i]==1) AliDebug(1,Form("Trackable particles: %d e pdg %d \n",i,kf));
196 pParent[i] = nt;
197 KeepTrack(nt);
198 SetHighWaterMark(nt);
199 }// Particle loop
200 particles->Clear();
201 if (trackIt) delete[] trackIt;
202 if (pParent) delete[] pParent;
203 }
204 Info("Generate","AliGenEvtGen DONE");
205 }
206
207//////////////////////////////////////////////////////////////////////////////////////////
208Int_t AliGenEvtGen::GetFlavour(Int_t pdgCode)
209 {
210 //
211 // return the flavour of a particle
212 // input: pdg code of the particle
213 // output: Int_t
214 // 3 in case of strange (open and hidden)
215 // 4 in case of charm (")
216 // 5 in case of beauty (")
217 //
218 Int_t pdg = TMath::Abs(pdgCode);
219 //Resonance
220 if (pdg > 100000) pdg %= 100000;
221 if(pdg > 10000) pdg %= 10000;
222 // meson ?
223 if(pdg > 10) pdg/=100;
224 // baryon ?
225 if(pdg > 10) pdg/=10;
226 return pdg;
227 }
228
229/////////////////////////////////////////////////////////////////////////////////////////
230Bool_t AliGenEvtGen::SetUserDecayTable(Char_t *path)
231 {
232 //
233 //Set the path of user decay table if it is defined
234 //
235 //put a comment to control if path exists
236 if(gSystem->AccessPathName(path))
237 {
238 AliWarning("Attention: This path not exist!\n");
239 return kFALSE;
240 }
241 fUserDecayTablePath = path;
242 fUserDecay = kTRUE;
243 return kTRUE;
244 }