8127ab026f966d0ba307d50f7e964c18331411ab
[u/mrichter/AliRoot.git] / EVGEN / AliGenCocktailAfterBurner.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 /* $Id$ */
17
18 // 
19 // Container class for AliGenerator and AfterBurners 
20 // (which are AliGenerators as well) through recursion.
21 // The container is itself an AliGenerator a
22 // what is stored are not the pointers to the generators directly 
23 // but to objects of type
24 // AliGenCocktailAfterBurner entry.   
25 // The class provides also iterator functionality.  
26 // Author: andreas.morsch@cern.ch and piotr.skowronski@cern.ch
27 //
28 // 24.09.2001  Piotr Skowronski
29 //             debug -> gDebug,
30 //             fNEvents replaced with gAlice->GetEventsPerRun()
31 //
32
33
34 #include <Riostream.h>
35
36 #include <TList.h>
37 #include <TObjArray.h>
38 #include <TParticle.h>
39
40 #include "AliGenCocktailAfterBurner.h"
41 #include "AliGenCocktailEntry.h"
42 #include "AliStack.h"
43 #include "AliMC.h"
44
45
46 ClassImp(AliGenCocktailAfterBurner)
47 /*********************************************************************/ 
48 /*********************************************************************/ 
49
50 AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()
51 {
52 // Constructor
53     if (gDebug > 0) 
54        cout<<"AliGenCocktailAfterBurner::AliGenCocktailAfterBurner()"<<endl;
55     SetName("AliGenCocktailAfterBurner");
56     SetTitle("AliGenCocktailAfterBurner");
57     fInternalStacks =0;
58     fActiveStack = 0;
59     fCurrentGenerator = 0;
60     fCurrentEvent =0;
61     fAfterBurnerEntries = new TList();
62     fNAfterBurners = 0;
63     fGenerationDone = kFALSE;
64     
65     fActiveEvent = -1;  
66     fNBgEvents = 0;
67 }
68 /*********************************************************************/ 
69 AliGenCocktailAfterBurner::AliGenCocktailAfterBurner(const AliGenCocktailAfterBurner& cocktail):
70     AliGenCocktail(cocktail)
71
72 {
73  //Copy constructor
74     cocktail.Copy(*this);
75 }
76
77 /*********************************************************************/ 
78
79 AliGenCocktailAfterBurner::~AliGenCocktailAfterBurner()
80   {
81 //destructor
82
83     if (fInternalStacks) //delete stacks
84      { 
85        fInternalStacks->SetOwner();
86        delete fInternalStacks;
87     }
88     if (fAfterBurnerEntries) delete fAfterBurnerEntries; //delete entries
89   }
90 /*********************************************************************/ 
91 /*********************************************************************/ 
92
93 void AliGenCocktailAfterBurner::
94 AddAfterBurner(AliGenerator *AfterBurner, char* Name, Float_t RateExp)
95 {
96 //
97 //  Forward parameters to the new AfterBurner
98     
99     if (gDebug>0)cout<<"AliGenCocktailAfterBurner::AddAfterBurner  Named "<<Name<<endl;
100
101     if(TestBit(kPtRange)) 
102         AfterBurner->SetPtRange(fPtMin,fPtMax);
103     if(TestBit(kMomentumRange))
104         AfterBurner->SetMomentumRange(fPMin,fPMax);
105     
106     AfterBurner->SetYRange(fYMin,fYMax);
107     AfterBurner->SetPhiRange(fPhiMin*180/TMath::Pi(),fPhiMax*180/TMath::Pi());
108     AfterBurner->SetThetaRange(fThetaMin*180/TMath::Pi(),fThetaMax*180/TMath::Pi());
109     AfterBurner->SetOrigin(fOrigin[0], fOrigin[1], fOrigin[2]);
110     AfterBurner->SetSigma(fOsigma[0], fOsigma[1], fOsigma[2]);
111     AfterBurner->SetVertexSmear(fVertexSmear);
112     AfterBurner->SetTrackingFlag(fTrackIt);    
113 //
114 //  Add AfterBurner to list   
115     
116     AliGenCocktailEntry *entry = 
117         new AliGenCocktailEntry(AfterBurner, Name, RateExp);
118     fAfterBurnerEntries->Add(entry);
119     fNAfterBurners++;
120 //
121     
122 }
123 /*********************************************************************/ 
124 /*********************************************************************/ 
125
126 void AliGenCocktailAfterBurner::Init()
127 {
128 // Initialisation
129     fGenerationDone = kFALSE;
130     if (fInternalStacks) //delete stacks
131      { 
132        fInternalStacks->SetOwner();
133        fInternalStacks->Delete(); //clean after previous generation cycle
134      }
135
136     this->AliGenCocktail::Init(); 
137     
138     if (gDebug>0) cout<<"AliGenCocktailAfterBurner::Init"<<endl;
139     TIter next(fAfterBurnerEntries);
140     AliGenCocktailEntry *entry;
141     //
142     // Loop over generators and initialize
143     while((entry = (AliGenCocktailEntry*)next())) {
144         entry->Generator()->Init();
145     }  
146 }
147 /*********************************************************************/ 
148 /*********************************************************************/ 
149
150 void AliGenCocktailAfterBurner::Generate()
151 {
152 //
153 // Generate event
154 //  Firsts runs each generator for all events
155 //  than after burners ones for each event
156 //
157 //  It generates and processes all events during
158 //  first call only.
159 //  In next calls it just returns already generated 
160 //  and processed events to the gAlice
161
162     if (gDebug>0)
163       cout<<"#####################################"<<endl
164           <<"#AliGenCocktailAfterBurner::Generate#"<<endl
165           <<"#####################################"<<endl;
166     
167     Int_t i; //iterator
168     AliStack * stack;
169     
170     if (fGenerationDone)
171     {//if generation is done (in first call) 
172      //just copy particles from the stack to the gAlice
173       SetTracks(++fCurrentEvent);
174       cout<<"Returning event "<<fCurrentEvent<<endl;
175       return;  
176     }
177     else
178     { //Here we are in the first call of the method
179       fCurrentEvent=0;
180       Int_t numberOfEvents = gAlice->GetEventsPerRun();
181       //Create stacks
182       fInternalStacks = new TObjArray(numberOfEvents + fNBgEvents); //Create array of internal stacks
183       for(i=0;i<numberOfEvents + fNBgEvents;i++) 
184        {        
185         stack = new AliStack(10000);
186         stack->Reset();
187         fInternalStacks->Add(stack);
188        }
189 /*********************************************************************/ 
190       TIter next(fEntries);
191       AliGenCocktailEntry *entry;
192       AliGenCocktailEntry *e1;
193       AliGenCocktailEntry *e2;
194       TObjArray *partArray;
195   //
196   // Loop over generators and generate events
197       Int_t igen=0;
198       while((entry = (AliGenCocktailEntry*)next())) 
199       {
200         igen++;
201         cout<<"Generator "<<igen<<"  : "<<entry->GetName()<<endl;
202 /***********************************************/
203 //First generator for all evenets, than second for all events, etc...
204         for(i=0;i<numberOfEvents + fNBgEvents;i++) 
205           {  
206             cout<<"                  EVENT "<<i<<endl;
207             stack = GetStack(i);
208             partArray = stack->Particles();
209             fCurrentGenerator = entry->Generator();
210             fCurrentGenerator->SetStack(stack);
211             if (igen ==1) 
212               {
213                 entry->SetFirst(0);
214               } 
215             else 
216               {
217                 entry->SetFirst((partArray->GetEntriesFast())+1);
218               }
219                 fCurrentGenerator->Generate();
220                 entry->SetLast(partArray->GetEntriesFast());
221            }
222 /***********************************************/
223       }
224       next.Reset();
225       while((entry = (AliGenCocktailEntry*)next())) 
226         {
227           entry->PrintInfo();
228         }
229       for ( entry=FirstGenerator();entry;entry=NextGenerator() ) 
230         {
231           entry->PrintInfo();
232         }
233       for (FirstGeneratorPair(e1,e2); (e1&&e2); NextGeneratorPair(e1,e2) )
234         {
235           printf("\n -----------------------------");
236           e1->PrintInfo();
237           e2->PrintInfo();
238         }
239         
240         
241       /***********************************************/
242       /*******After Burners Processing****************/
243       /***********************************************/
244       TIter nextAfterBurner(fAfterBurnerEntries);
245       AliGenCocktailEntry *afterBurnerEntry;
246       Int_t iab =0; //number of current after burner / counter
247       
248       cout<<"\n\nRunning After Burners"<<endl;
249       while((afterBurnerEntry = (AliGenCocktailEntry*)nextAfterBurner()))
250         {
251           cout<<"After Burner "<<iab++<<"  :"<<afterBurnerEntry->GetName()<<endl;
252           fCurrentGenerator = afterBurnerEntry->Generator();
253           fCurrentGenerator->Generate();
254         }
255       cout<<endl<<"Finished. Processed "<<iab<<" After Burners"<<endl;
256
257       /***********************************************/
258       /***********************************************/
259       /***********************************************/       
260         
261       fGenerationDone=kTRUE; 
262       SetTracks(0); //copy event 0 to gAlice stack
263         
264 /*********************************************************************/
265         
266     }//else generated
267 }
268 /*********************************************************************/
269 /*********************************************************************/ 
270
271 AliGenCocktailAfterBurner& AliGenCocktailAfterBurner::operator=(const  AliGenCocktailAfterBurner& rhs)
272 {
273 // Assignment operator
274     rhs.Copy(*this);
275     return (*this);
276 }
277 /*********************************************************************/
278 /*********************************************************************/ 
279
280 AliStack* AliGenCocktailAfterBurner::GetStack(Int_t n) const
281 {
282 //Returns the pointer to the N'th stack (event)
283   if( ( n<0 ) || ( n>=GetNumberOfEvents() ) )
284     {
285       Fatal("AliGenCocktailAfterBurner::GetStack","Asked for non existing stack (%d)",n);
286       return 0; 
287     }
288     return ((AliStack*) fInternalStacks->At(n) );
289 }
290 /*********************************************************************/ 
291 /*********************************************************************/ 
292
293 void AliGenCocktailAfterBurner::SetActiveEventNumber(Int_t actev)
294 {
295 //Set Active Events Number and Active Stack
296 //There is only one active event number
297 //Made fo convinience of work with AfterBurners (HBT processor)
298
299     fActiveEvent = actev;
300     fActiveStack = GetStack(actev);
301 }
302 /*********************************************************************/ 
303 /*********************************************************************/ 
304
305 void AliGenCocktailAfterBurner::SetTracks(Int_t stackno)
306 {
307 //Method which copies tracks from given stack to the
308 //gAlice's stack
309     AliStack* instack = GetStack(stackno);
310     Int_t done;
311     Int_t parent; 
312     Int_t pdg;
313     Double_t px, py, pz, e, vx, vy, vz, tof, polx, poly, polz;
314     TMCProcess mech;
315     Int_t ntr;
316     Float_t weight;
317     TVector3 pol;
318     
319     TParticle * p;
320     Int_t n = instack->GetNtrack();
321     if (gDebug) 
322     {
323       cout<<"AliGenCocktailAfterBurner::SetTracks("<<stackno<<"). Number of particles is: "<<n<<"\n";
324     }
325     
326     for(Int_t i = 0; i < n; i++)
327     {
328         
329       p = instack->Particle(i);
330       done = !p->TestBit(kDoneBit);
331       parent = p->GetMother(0);
332       pdg = p->GetPdgCode();
333       px = p->Px();
334       py = p->Py();
335       pz = p->Pz();
336       e  = p->Energy();
337       vx = p->Vx();
338       vy = p->Vy();
339       vz = p->Vz();
340       tof = p->T();
341       p->GetPolarisation(pol);
342       polx = pol.X();
343       poly = pol.Y();
344       polz = pol.Z();
345       mech = AliGenCocktailAfterBurner::IntToMCProcess(p->GetUniqueID());
346       weight = p->GetWeight();
347
348       gAlice->GetMCApp()->PushTrack(done, parent, pdg, px, py, pz, e, vx, vy, vz, tof,
349                        polx, poly, polz, mech, ntr, weight);
350     }
351 }
352 /*********************************************************************/ 
353 /*********************************************************************/ 
354
355 TMCProcess AliGenCocktailAfterBurner::IntToMCProcess(Int_t no)
356 {
357  //Mothod used to convert uniqueID (integer) to TMCProcess type
358     const TMCProcess kMCprocesses[kMaxMCProcess] = 
359     {
360      kPNoProcess, kPMultipleScattering, kPEnergyLoss, kPMagneticFieldL, 
361      kPDecay, kPPair, kPCompton, kPPhotoelectric, kPBrem, kPDeltaRay,
362      kPAnnihilation, kPHadronic, kPNoProcess, kPEvaporation, kPNuclearFission,
363      kPNuclearAbsorption, kPPbarAnnihilation, kPNCapture, kPHElastic, 
364      kPHInhelastic, kPMuonNuclear, kPTOFlimit,kPPhotoFission, kPNoProcess, 
365      kPRayleigh, kPNoProcess, kPNoProcess, kPNoProcess, kPNull, kPStop
366     };
367     
368     for (Int_t i = 0;i<kMaxMCProcess;i++)
369     {
370       if (kMCprocesses[i] == no)
371         {
372           return kMCprocesses[i];
373         }
374     } 
375     return kPNoProcess;
376 }
377
378 void AliGenCocktailAfterBurner::Copy(AliGenCocktailAfterBurner &) const
379 {
380     Fatal("Copy","Not implemented!\n");
381 }