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 ***************************************************************************/
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. //
22 // Giuseppe E. Bruno & Fiorella Fionda //
23 // (Giuseppe.Bruno@ba.infn.it) (Fiorella.Fionda@ba.infn.it) //
24 ///////////////////////////////////////////////////////////////////////////
27 #include "AliGenEvtGen.h"
30 #include <TParticle.h>
32 ClassImp(AliGenEvtGen)
33 ///////////////////////////////////////////////////////////////////////////
34 AliGenEvtGen::AliGenEvtGen():
38 fSwitchOff(kBeautyPart),
40 fUserDecayTablePath(0x0)
43 // Default Construction
46 ///////////////////////////////////////////////////////////////////////////////////////////
47 AliGenEvtGen::~AliGenEvtGen()
50 // Standard Destructor
52 if(fStack) {delete fStack;}
54 if(fDecayer) {delete fDecayer;}
56 if(fUserDecayTablePath) {delete fUserDecayTablePath;}
57 fUserDecayTablePath = 0;
59 ///////////////////////////////////////////////////////////////////////////////////////////
61 void AliGenEvtGen::Init()
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
71 AliWarning("AliGenEvtGen already initialized!!!");
74 fDecayer = new AliDecayerEvtGen();
75 fDecayer->Init(); //read the default decay table DECAY.DEC and particle table
77 //if is set a decay mode: default decay mode is kAll
78 fDecayer->SetForceDecay(fForceDecay);
79 fDecayer->ForceDecay();
81 //if is defined a user decay table
84 fDecayer->SetDecayTablePath(fUserDecayTablePath);
85 fDecayer->ReadDecayTable();
89 /////////////////////////////////////////////////////////////////////////////////////////////
90 void AliGenEvtGen::Generate()
93 //Generate method - Input none - Output none
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
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
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();
116 AliDebug(1,Form("GetFlavour = %d e pdg = %d \n",GetFlavour(pdg),pdg));
123 if(GetFlavour(pdg)!=5) continue;
126 if(GetFlavour(pdg)!=4) continue;
130 //check if particle is already decayed by Pythia
131 if(part->GetStatusCode() != 1 || part->GetNDaughters()>0)
133 Info("AliGenEvtGen","Attention: particle %d is already decayed by Pythia!",pdg);
137 part->SetStatusCode(11); //Set particle as decayed : change the status code
139 mom->SetPxPyPzE(part->Px(),part->Py(),part->Pz(),part->Energy());
143 fDecayer->Decay(part->GetPdgCode(),mom);
144 np = fDecayer->ImportParticles(particles);
147 Int_t* trackIt = new Int_t[np];
148 Int_t* pParent = new Int_t[np];
149 AliDebug(1,Form("np = %d \n",np));
151 for (int i = 0; i < np; i++) {
155 //select trackable particle
157 TParticle* iparticle = (TParticle *) particles->At(0);//parent particle
158 for (int i = 1; i<np ; i++) {
159 iparticle = (TParticle*) particles->At(i);
160 Int_t ks = iparticle->GetStatusCode();
162 //track last decay products
163 if(ks==1) trackIt[i]=1;
165 }//decay particles loop
167 }// if decay products
169 origin0[0]=part->Vx(); //[cm]
170 origin0[1]=part->Vy(); //[cm]
171 origin0[2]=part->Vz(); //[cm]
173 // Put decay products on the stack
175 for (int i = 1; i < np; i++) {
176 TParticle* iparticle = (TParticle *) particles->At(i);
177 Int_t kf = iparticle->GetPdgCode();
178 Int_t ksc = iparticle->GetStatusCode();
179 Int_t jpa = iparticle->GetFirstMother()-1; //jpa = 0 for daughters of beauty particles
180 Int_t iparent = (jpa > 0) ? pParent[jpa] : iTrack;
182 och[0] = origin0[0]+iparticle->Vx()/10; //[cm]
183 och[1] = origin0[1]+iparticle->Vy()/10; //[cm]
184 och[2] = origin0[2]+iparticle->Vz()/10; //[cm]
185 pc[0] = iparticle->Px(); //[GeV/c]
186 pc[1] = iparticle->Py(); //[GeV/c]
187 pc[2] = iparticle->Pz(); //[GeV/c]
188 tof = part->T()+kconv*iparticle->T();
190 AliDebug(1,Form("FirstMother = %d e indicePart = %d e pdg = %d \n",jpa,i,kf));
192 PushTrack(trackIt[i], iparent, kf, pc, och, polar,tof, kPDecay, nt, 1., ksc);
193 if(trackIt[i]==1) AliDebug(1,Form("Trackable particles: %d e pdg %d \n",i,kf));
196 SetHighWaterMark(nt);
199 if (trackIt) delete[] trackIt;
200 if (pParent) delete[] pParent;
202 Info("Generate","AliGenEvtGen DONE");
205 //////////////////////////////////////////////////////////////////////////////////////////
206 Int_t AliGenEvtGen::GetFlavour(Int_t pdgCode)
209 // return the flavour of a particle
210 // input: pdg code of the particle
212 // 3 in case of strange (open and hidden)
213 // 4 in case of charm (")
214 // 5 in case of beauty (")
216 Int_t pdg = TMath::Abs(pdgCode);
218 if (pdg > 100000) pdg %= 100000;
219 if(pdg > 10000) pdg %= 10000;
221 if(pdg > 10) pdg/=100;
223 if(pdg > 10) pdg/=10;
227 /////////////////////////////////////////////////////////////////////////////////////////
228 Bool_t AliGenEvtGen::SetUserDecayTable(Char_t *path)
231 //Set the path of user decay table if it is defined
233 //put a comment to control if path exists
234 if(gSystem->AccessPathName(path))
236 AliWarning("Attention: This path not exist!\n");
239 fUserDecayTablePath = path;