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