Update master to aliroot
[u/mrichter/AliRoot.git] / EVGEN / AliGenCocktail.cxx
CommitLineData
4c039060 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
88cb7938 16/* $Id$ */
4c039060 17
675e9664 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
116cbefd 26#include <TList.h>
27#include <TObjArray.h>
7c54213e 28#include <TFormula.h>
116cbefd 29
fe4da5cc 30#include "AliGenCocktail.h"
374924b5 31#include "AliGenCocktailEntry.h"
09e2e187 32#include "AliCollisionGeometry.h"
fe4da5cc 33#include "AliRun.h"
7c54213e 34#include "AliLog.h"
5d12ce38 35#include "AliMC.h"
4f85aa78 36#include "AliGenCocktailEventHeader.h"
fe4da5cc 37
38ClassImp(AliGenCocktail)
39
40AliGenCocktail::AliGenCocktail()
1c56e311 41 :AliGenerator(),
42 fNGenerators(0),
8058ead1 43 fTotalRate(0.),
ab834258 44 fSRandom(kFALSE),
1c56e311 45 fUsePerEventRate(kFALSE),
46 fProb(0),
47 fEntries(0),
48 flnk1(0),
49 flnk2(0),
9577906a 50 fHeader(0),
51 fSeed(0)
fe4da5cc 52{
374924b5 53// Constructor
8b31bfac 54 fName = "Cocktail";
55 fTitle= "Particle Generator using cocktail of generators";
fe4da5cc 56}
57
fe4da5cc 58AliGenCocktail::~AliGenCocktail()
59{
374924b5 60// Destructor
fe4da5cc 61 delete fEntries;
373c1a55 62 fEntries = 0;
63 // delete fHeader; // It is removed in AliRunLoader
64 fHeader = 0;
fe4da5cc 65}
66
67void AliGenCocktail::
8c1efaa8 68AddGenerator(AliGenerator *Generator, const char* Name, Float_t RateExp, TFormula* formula, Int_t ntimes)
fe4da5cc 69{
70//
b1403e15 71// Add a generator to the list
72// First check that list exists
73 if (!fEntries) fEntries = new TList();
8058ead1 74 fTotalRate += RateExp;
b1403e15 75//
fe4da5cc 76// Forward parameters to the new generator
212bcf04 77 if(TestBit(kPtRange) && !(Generator->TestBit(kPtRange)) && !(Generator->TestBit(kMomentumRange)))
78 Generator->SetPtRange(fPtMin,fPtMax);
79 if(TestBit(kMomentumRange) && !(Generator->TestBit(kPtRange)) && !(Generator->TestBit(kMomentumRange)))
80 Generator->SetMomentumRange(fPMin,fPMax);
81
565cc9c0 82 if (TestBit(kYRange) && !(Generator->TestBit(kYRange)))
212bcf04 83 Generator->SetYRange(fYMin,fYMax);
565cc9c0 84 if (TestBit(kPhiRange) && !(Generator->TestBit(kPhiRange)))
212bcf04 85 Generator->SetPhiRange(fPhiMin*180/TMath::Pi(),fPhiMax*180/TMath::Pi());
565cc9c0 86 if (TestBit(kThetaRange) && !(Generator->TestBit(kThetaRange)) && !(Generator->TestBit(kEtaRange)))
212bcf04 87 Generator->SetThetaRange(fThetaMin*180/TMath::Pi(),fThetaMax*180/TMath::Pi());
2b597689 88 if (!(Generator->TestBit(kVertexRange))) {
212bcf04 89 Generator->SetOrigin(fOrigin[0], fOrigin[1], fOrigin[2]);
4a6ca8fc 90 Generator->SetSigma(fOsigma[0], fOsigma[1], fOsigma[2]);
91 Generator->SetVertexSmear(fVertexSmear);
92 Generator->SetVertexSource(kContainer);
93 }
c8f7f6f9 94 Generator->SetTrackingFlag(fTrackIt);
4f85aa78 95 Generator->SetContainer(this);
7c54213e 96
c97002cc 97
fe4da5cc 98//
99// Add generator to list
c5c3e4b0 100 char theName[256];
04a640a7 101 snprintf(theName, 256, "%s_%d",Name, fNGenerators);
c5c3e4b0 102 Generator->SetName(theName);
103
374924b5 104 AliGenCocktailEntry *entry =
fe4da5cc 105 new AliGenCocktailEntry(Generator, Name, RateExp);
8c1efaa8 106 if (formula) entry->SetFormula(formula);
107 entry->SetNTimes(ntimes);
374924b5 108 fEntries->Add(entry);
fe4da5cc 109 fNGenerators++;
373c1a55 110 flnk1 = 0;
111 flnk2 = 0;
ab834258 112 fSRandom = kFALSE;
373c1a55 113 fHeader = 0;
114}
fe4da5cc 115
116 void AliGenCocktail::Init()
374924b5 117{
118// Initialisation
119 TIter next(fEntries);
120 AliGenCocktailEntry *entry;
121 //
122 // Loop over generators and initialize
123 while((entry = (AliGenCocktailEntry*)next())) {
c97002cc 124 if (fStack) entry->Generator()->SetStack(fStack);
ab71c6e6 125 entry->Generator()->SetSeed(fSeed);
374924b5 126 entry->Generator()->Init();
127 }
1ebe6dc8 128
129 next.Reset();
130
ab834258 131 if (fSRandom) {
1ebe6dc8 132 fProb.Set(fNGenerators);
133 next.Reset();
134 Float_t sum = 0.;
135 while((entry = (AliGenCocktailEntry*)next())) {
136 sum += entry->Rate();
137 }
138
139 next.Reset();
140 Int_t i = 0;
141 Float_t psum = 0.;
142 while((entry = (AliGenCocktailEntry*)next())) {
143 psum += entry->Rate() / sum;
144 fProb[i++] = psum;
145 }
146 }
147 next.Reset();
374924b5 148}
fe4da5cc 149
09e2e187 150 void AliGenCocktail::FinishRun()
151{
152// Initialisation
153 TIter next(fEntries);
154 AliGenCocktailEntry *entry;
155 //
156 // Loop over generators and initialize
157 while((entry = (AliGenCocktailEntry*)next())) {
158 entry->Generator()->FinishRun();
159 }
160}
161
fe4da5cc 162 void AliGenCocktail::Generate()
374924b5 163{
164//
165// Generate event
166 TIter next(fEntries);
a8408367 167 AliGenCocktailEntry *entry = 0;
168 AliGenCocktailEntry *preventry = 0;
7c54213e 169 AliGenCocktailEntry *collentry = 0;
a8408367 170 AliGenerator* gen = 0;
4f85aa78 171 if (fHeader) delete fHeader;
c89fc3f1 172
173
4f85aa78 174 fHeader = new AliGenCocktailEventHeader("Cocktail Header");
09e2e187 175
d94af0c1 176 const TObjArray *partArray = gAlice->GetMCApp()->Particles();
c8f7f6f9 177
178//
179// Generate the vertex position used by all generators
180//
181 if(fVertexSmear == kPerEvent) Vertex();
182
67db3ff2 183 TArrayF eventVertex;
184 eventVertex.Set(3);
185 for (Int_t j=0; j < 3; j++) eventVertex[j] = fVertex[j];
186
ab834258 187 if (!fSRandom) {
1ebe6dc8 188 //
189 // Loop over generators and generate events
37098832 190 Int_t igen = 0;
1ebe6dc8 191 while((entry = (AliGenCocktailEntry*)next())) {
8c1efaa8 192 Int_t ntimes = entry->NTimes();
7c54213e 193 if (fUsePerEventRate && (gRandom->Rndm() > entry->Rate())) continue;
194
195 igen++;
196 if (igen ==1) {
197 entry->SetFirst(0);
198 } else {
199 entry->SetFirst((partArray->GetEntriesFast())+1);
200 }
201 gen = entry->Generator();
202 if (gen->ProvidesCollisionGeometry()) collentry = entry;
203 //
204 // Handle case in which current generator needs collision geometry from previous generator
205 //
206 if (gen->NeedsCollisionGeometry() && (entry->Formula() == 0))
a8408367 207 {
7c54213e 208 if (preventry && preventry->Generator()->ProvidesCollisionGeometry())
1ebe6dc8 209 {
7c54213e 210 gen->SetCollisionGeometry(preventry->Generator()->CollisionGeometry());
1ebe6dc8 211 } else {
7c54213e 212 Fatal("Generate()", "No Collision Geometry Provided");
213 }
214 }
215 //
216 // Number of signals is calculated from Collision Geometry
db051b0c 217 // and entry with given centrality bin is selected
7c54213e 218 //
219 if (entry->Formula() != 0)
220 {
221 if (!collentry) {
222 Fatal("Generate()", "No Collision Geometry Provided");
223 return;
224 }
225 AliCollisionGeometry* coll = (collentry->Generator())->CollisionGeometry();
226 Float_t b = coll->ImpactParameter();
227 Int_t nsig = Int_t(entry->Formula()->Eval(b));
db051b0c 228 Int_t bin = entry->Bin() - 100;
229 if (bin > 0) {
230 if (bin != nsig) continue;
231 } else {
232 if (nsig < 1) nsig = 1;
233 AliInfo(Form("Signal Events %13.3f %5d %5d\n", b, coll->HardScatters(), nsig));
234 ntimes = nsig;
235 }
a8408367 236 }
735e0707 237 gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2), fTime);
72cf8302 238
d5615b85 239 gen->GenerateN(ntimes);
7c54213e 240 entry->SetLast(partArray->GetEntriesFast());
241 preventry = entry;
242 }
ab834258 243 } else if (fSRandom) {
1ebe6dc8 244 //
245 // Select a generator randomly
246 //
247 Int_t i;
248 Float_t p0 = gRandom->Rndm();
249
250 for (i = 0; i < fNGenerators; i++) {
251 if (p0 < fProb[i]) break;
a8408367 252 }
1ebe6dc8 253
254 entry = (AliGenCocktailEntry*) fEntries->At(i);
255 entry->SetFirst(0);
256 gen = entry->Generator();
735e0707 257 gen->SetVertex(fVertex.At(0), fVertex.At(1), fVertex.At(2), fTime);
08bffa4d 258 gen->Generate();
374924b5 259 entry->SetLast(partArray->GetEntriesFast());
c89fc3f1 260 }
1ebe6dc8 261
374924b5 262 next.Reset();
4f85aa78 263
45c05df7 264 // Event Vertex
4f85aa78 265 fHeader->SetPrimaryVertex(eventVertex);
bcee5fdf 266 fHeader->CalcNProduced();
45c05df7 267 if (fContainer) {
268 fHeader->SetName(fName);
269 fContainer->AddHeader(fHeader);
270 } else {
271 gAlice->SetGenEventHeader(fHeader);
272 }
374924b5 273}
fe4da5cc 274
c44d1c08 275void AliGenCocktail::SetVertexSmear(VertexSmear_t smear)
276{
277// Set vertex smearing and propagate it to the generators
278
279 AliGenerator::SetVertexSmear(smear);
280 TIter next(fEntries);
281 while (AliGenCocktailEntry* entry = (AliGenCocktailEntry*)next()) {
282 entry->Generator()->SetVertexSmear(smear);
283 }
284}
285
fe4da5cc 286AliGenCocktailEntry * AliGenCocktail::FirstGenerator()
287{
374924b5 288// Iterator over generators: Initialisation
fe4da5cc 289 flnk1 = fEntries->FirstLink();
290 if (flnk1) {
291 return (AliGenCocktailEntry*) (flnk1->GetObject());
292 } else {
293 return 0;
294 }
295}
296
297AliGenCocktailEntry* AliGenCocktail::NextGenerator()
298{
374924b5 299// Iterator over generators: Increment
fe4da5cc 300 flnk1 = flnk1->Next();
301 if (flnk1) {
302 return (AliGenCocktailEntry*) (flnk1->GetObject());
303 } else {
304 return 0;
305 }
306}
307
308void AliGenCocktail::
309FirstGeneratorPair(AliGenCocktailEntry*& e1, AliGenCocktailEntry*& e2)
310{
374924b5 311// Iterator over generator pairs: Initialisation
fe4da5cc 312 flnk2 = flnk1 = fEntries->FirstLink();
313 if (flnk1) {
314 e2 = e1 = (AliGenCocktailEntry*) (flnk1->GetObject());
315 } else {
316 e2= e1 = 0;
317 }
318}
319
320void AliGenCocktail::
321NextGeneratorPair(AliGenCocktailEntry*& e1, AliGenCocktailEntry*& e2)
322{
374924b5 323// Iterator over generators: Increment
fe4da5cc 324 flnk2 = flnk2->Next();
325 if (flnk2) {
326 e1 = (AliGenCocktailEntry*) (flnk1->GetObject());
327 e2 = (AliGenCocktailEntry*) (flnk2->GetObject());
328 } else {
329 flnk2 = flnk1 = flnk1->Next();
330 if (flnk1) {
331 e1 = (AliGenCocktailEntry*) (flnk1->GetObject());
332 e2 = (AliGenCocktailEntry*) (flnk2->GetObject());
333 } else {
334 e1=0;
335 e2=0;
336 }
337 }
338}
339
4f85aa78 340void AliGenCocktail::AddHeader(AliGenEventHeader* header)
341{
342// Add a header to the list
343 if (fHeader) fHeader->AddHeader(header);
344}
fe4da5cc 345
198bb1c7 346
347
fe4da5cc 348