Additional initialization. Do not remove fHeader, it is done in AliRunLoader
[u/mrichter/AliRoot.git] / EVGEN / AliGenCocktail.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 // Container class for AliGenerator through recursion.
19 // Container is itself an AliGenerator.
20 // What is stored are not the pointers to the generators directly but to objects of type
21 // AliGenCocktail entry.   
22 // The class provides also iterator functionality.  
23 // Author: andreas.morsch@cern.ch 
24 //
25
26 #include <TList.h>
27 #include <TObjArray.h>
28
29 #include "AliGenCocktail.h"
30 #include "AliGenCocktailEntry.h"
31 #include "AliCollisionGeometry.h"
32 #include "AliRun.h"
33 #include "AliMC.h"
34 #include "AliGenCocktailEventHeader.h"
35
36 ClassImp(AliGenCocktail)
37
38 AliGenCocktail::AliGenCocktail()
39                  :AliGenerator()
40 {
41 // Constructor
42     fName = "Cocktail";
43     fTitle= "Particle Generator using cocktail of generators";
44     flnk1 = 0;
45     flnk2 = 0;
46     fNGenerators=0;
47     fEntries = 0;
48     fRandom  = kFALSE;
49     fHeader  = 0;
50 }
51
52 AliGenCocktail::AliGenCocktail(const AliGenCocktail & cocktail):
53     AliGenerator(cocktail)
54 {
55 // Copy constructor
56     cocktail.Copy(*this);
57 }
58
59 AliGenCocktail::~AliGenCocktail()
60 {
61 // Destructor
62     delete fEntries;
63     fEntries = 0;
64     //    delete fHeader; // It is removed in AliRunLoader
65     fHeader = 0;
66 }
67
68 void AliGenCocktail::
69 AddGenerator(AliGenerator *Generator, const char* Name, Float_t RateExp)
70 {
71 //
72 // Add a generator to the list 
73 // First check that list exists
74     if (!fEntries) fEntries = new TList();
75
76 //
77 //  Forward parameters to the new generator
78     if(TestBit(kPtRange) && !(Generator->TestBit(kPtRange)) && !(Generator->TestBit(kMomentumRange))) 
79         Generator->SetPtRange(fPtMin,fPtMax);
80     if(TestBit(kMomentumRange) && !(Generator->TestBit(kPtRange)) && !(Generator->TestBit(kMomentumRange)))
81         Generator->SetMomentumRange(fPMin,fPMax);
82     
83     if (!(Generator->TestBit(kYRange)))    
84         Generator->SetYRange(fYMin,fYMax);
85     if (!(Generator->TestBit(kPhiRange)))   
86         Generator->SetPhiRange(fPhiMin*180/TMath::Pi(),fPhiMax*180/TMath::Pi());
87     if (!(Generator->TestBit(kThetaRange))) 
88         Generator->SetThetaRange(fThetaMin*180/TMath::Pi(),fThetaMax*180/TMath::Pi());
89     if (!(Generator->TestBit(kVertexRange))) 
90         Generator->SetOrigin(fOrigin[0], fOrigin[1], fOrigin[2]);
91
92     Generator->SetSigma(fOsigma[0], fOsigma[1], fOsigma[2]);
93     Generator->SetVertexSmear(fVertexSmear);
94     Generator->SetVertexSource(kContainer);
95     Generator->SetTrackingFlag(fTrackIt);
96     Generator->SetContainer(this);
97     
98         
99 //
100 //  Add generator to list   
101     char theName[256];
102     sprintf(theName, "%s_%d",Name, fNGenerators);
103     Generator->SetName(theName);
104
105     AliGenCocktailEntry *entry = 
106         new AliGenCocktailEntry(Generator, Name, RateExp);
107     
108      fEntries->Add(entry);
109      fNGenerators++;
110      flnk1 = 0;
111      flnk2 = 0;
112      fRandom  = kFALSE;
113      fHeader  = 0;
114 }
115
116   void AliGenCocktail::Init()
117 {
118 // Initialisation
119     TIter next(fEntries);
120     AliGenCocktailEntry *entry;
121     //
122     // Loop over generators and initialize
123     while((entry = (AliGenCocktailEntry*)next())) {
124         if (fStack)  entry->Generator()->SetStack(fStack);
125         entry->Generator()->Init();
126     }  
127
128     next.Reset();
129
130     if (fRandom) {
131         fProb.Set(fNGenerators);
132         next.Reset();
133         Float_t sum = 0.;
134         while((entry = (AliGenCocktailEntry*)next())) {
135             sum += entry->Rate();
136         } 
137
138         next.Reset();
139         Int_t i = 0;
140         Float_t psum = 0.;
141         while((entry = (AliGenCocktailEntry*)next())) {
142             psum +=  entry->Rate() / sum;
143             fProb[i++] = psum;
144         }
145     }
146         next.Reset();
147 }
148
149   void AliGenCocktail::FinishRun()
150 {
151 // Initialisation
152     TIter next(fEntries);
153     AliGenCocktailEntry *entry;
154     //
155     // Loop over generators and initialize
156     while((entry = (AliGenCocktailEntry*)next())) {
157         entry->Generator()->FinishRun();
158     }  
159 }
160
161  void AliGenCocktail::Generate()
162 {
163 //
164 // Generate event 
165     TIter next(fEntries);
166     AliGenCocktailEntry *entry = 0;
167     AliGenCocktailEntry *preventry = 0;
168     AliGenerator* gen = 0;
169     if (fHeader) delete fHeader;
170     fHeader = new AliGenCocktailEventHeader("Cocktail Header");
171
172     TObjArray *partArray = gAlice->GetMCApp()->Particles();
173
174 //
175 //  Generate the vertex position used by all generators
176 //    
177     if(fVertexSmear == kPerEvent) Vertex();
178
179     TArrayF eventVertex;
180     eventVertex.Set(3);
181     for (Int_t j=0; j < 3; j++) eventVertex[j] = fVertex[j];
182
183     if (!fRandom) {
184         //
185         // Loop over generators and generate events
186         Int_t igen=0;
187         
188         while((entry = (AliGenCocktailEntry*)next())) {
189             igen++;
190             if (igen ==1) {
191                 entry->SetFirst(0);
192             } else {
193                 entry->SetFirst((partArray->GetEntriesFast())+1);
194             }
195 //
196 //      Handle case in which current generator needs collision geometry from previous generator
197 //
198             gen = entry->Generator();
199             if (gen->NeedsCollisionGeometry())
200             {
201                 if (preventry && preventry->Generator()->ProvidesCollisionGeometry())
202                 {
203                     gen->SetCollisionGeometry(preventry->Generator()->CollisionGeometry());
204                 } else {
205                     Fatal("Generate()", "No Collision Geometry Provided");
206                 }
207             }
208             entry->Generator()->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2));
209             entry->Generator()->Generate();
210             entry->SetLast(partArray->GetEntriesFast());
211             preventry = entry;
212         }  
213     } else {
214         //
215         // Select a generator randomly
216         //
217         Int_t i;
218         Float_t p0 =  gRandom->Rndm();
219
220         for (i = 0; i < fNGenerators; i++) {
221             if (p0 < fProb[i]) break;
222         }
223
224         entry = (AliGenCocktailEntry*) fEntries->At(i);
225         entry->SetFirst(0);
226         gen = entry->Generator();
227         entry->Generator()->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2));
228         entry->Generator()->Generate();
229         entry->SetLast(partArray->GetEntriesFast());
230     }
231     
232     
233     next.Reset();
234
235 // Event Vertex
236     fHeader->SetPrimaryVertex(eventVertex);
237     gAlice->SetGenEventHeader(fHeader); 
238 }
239
240 void AliGenCocktail::SetVertexSmear(VertexSmear_t smear)
241 {
242 // Set vertex smearing and propagate it to the generators
243
244   AliGenerator::SetVertexSmear(smear);
245   TIter next(fEntries);
246   while (AliGenCocktailEntry* entry = (AliGenCocktailEntry*)next()) {
247     entry->Generator()->SetVertexSmear(smear);
248   }
249 }
250
251 AliGenCocktailEntry *  AliGenCocktail::FirstGenerator()
252 {
253 // Iterator over generators: Initialisation
254     flnk1 = fEntries->FirstLink();
255     if (flnk1) {
256         return (AliGenCocktailEntry*) (flnk1->GetObject());
257     } else {
258         return 0;
259     }
260 }
261
262 AliGenCocktailEntry*  AliGenCocktail::NextGenerator()
263 {
264 // Iterator over generators: Increment
265     flnk1 = flnk1->Next();
266     if (flnk1) {
267         return (AliGenCocktailEntry*) (flnk1->GetObject());
268     } else {
269         return 0;
270     }
271 }
272
273 void AliGenCocktail::
274 FirstGeneratorPair(AliGenCocktailEntry*& e1, AliGenCocktailEntry*& e2)
275 {
276 // Iterator over generator pairs: Initialisation
277     flnk2 = flnk1 = fEntries->FirstLink();
278     if (flnk1) {
279         e2 = e1 = (AliGenCocktailEntry*) (flnk1->GetObject());
280     } else {
281         e2= e1 = 0;
282     }
283 }
284
285 void AliGenCocktail::
286 NextGeneratorPair(AliGenCocktailEntry*& e1, AliGenCocktailEntry*& e2)
287 {
288 // Iterator over generators: Increment
289     flnk2 = flnk2->Next();
290     if (flnk2) {
291         e1 = (AliGenCocktailEntry*) (flnk1->GetObject());
292         e2 = (AliGenCocktailEntry*) (flnk2->GetObject());       
293     } else {
294         flnk2 = flnk1 = flnk1->Next();
295         if (flnk1) {
296             e1 = (AliGenCocktailEntry*) (flnk1->GetObject());
297             e2 = (AliGenCocktailEntry*) (flnk2->GetObject());
298         } else {
299             e1=0;
300             e2=0;
301         }
302     }
303 }
304
305 void AliGenCocktail::AddHeader(AliGenEventHeader* header)
306 {
307 // Add a header to the list 
308     if (fHeader) fHeader->AddHeader(header);
309 }
310                               
311 AliGenCocktail& AliGenCocktail::operator=(const  AliGenCocktail& rhs)
312 {
313 // Assignment operator
314     rhs.Copy(*this); 
315     return (*this);
316 }
317
318 void AliGenCocktail::Copy(TObject &) const
319 {
320     Fatal("Copy","Not implemented!\n");
321 }
322
323
324