Using vertex from underlying event in AliSimulation (T.Kuhr)
[u/mrichter/AliRoot.git] / STEER / AliSimulation.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 //                                                                           //
20 // class for running generation, simulation and digitization                 //
21 //                                                                           //
22 // Hits, sdigits and digits are created for all detectors by typing:         //
23 //                                                                           //
24 //   AliSimulation sim;                                                      //
25 //   sim.Run();                                                              //
26 //                                                                           //
27 // The Run method returns kTRUE in case of successful execution.             //
28 // The number of events can be given as argument to the Run method or it     //
29 // can be set by                                                             //
30 //                                                                           //
31 //   sim.SetNumberOfEvents(n);                                               //
32 //                                                                           //
33 // The name of the configuration file can be passed as argument to the       //
34 // AliSimulation constructor or can be specified by                          //
35 //                                                                           //
36 //   sim.SetConfigFile("...");                                               //
37 //                                                                           //
38 // The generation of particles and the simulation of detector hits can be    //
39 // switched on or off by                                                     //
40 //                                                                           //
41 //   sim.SetRunGeneration(kTRUE);   // generation of primary particles       //
42 //   sim.SetRunSimulation(kFALSE);  // but no tracking                       //
43 //                                                                           //
44 // For which detectors sdigits and digits will be created, can be steered    //
45 // by                                                                        //
46 //                                                                           //
47 //   sim.SetMakeSDigits("ALL");     // make sdigits for all detectors        //
48 //   sim.SetMakeDigits("ITS TPC");  // make digits only for ITS and TPC      //
49 //                                                                           //
50 // The argument is a (case sensitive) string with the names of the           //
51 // detectors separated by a space. An empty string ("") can be used to       //
52 // disable the creation of sdigits or digits. The special string "ALL"       //
53 // selects all available detectors. This is the default.                     //
54 //                                                                           //
55 // The creation of digits from hits instead of from sdigits can be selected  //
56 // by                                                                        //
57 //                                                                           //
58 //   sim.SetMakeDigitsFromHits("TRD");                                       //
59 //                                                                           //
60 // The argument is again a string with the selected detectors. Be aware that //
61 // this feature is not available for all detectors and that merging is not   //
62 // possible, when digits are created directly from hits.                     //
63 //                                                                           //
64 // Backgound events can be merged by calling                                 //
65 //                                                                           //
66 //   sim.MergeWith("background/galice.root", 2);                             //
67 //                                                                           //
68 // The first argument is the file name of the background galice file. The    //
69 // second argument is the number of signal events per background event.      //
70 // The default value for this is 1. MergeWith can be called several times    //
71 // to merge more than two event streams. It is assumed that the sdigits      //
72 // were already produced for the background events.                          //
73 //                                                                           //
74 // The methods RunSimulation, RunSDigitization, RunDigitization and          //
75 // RunHitsDigitization can be used to run only parts of the full simulation  //
76 // chain.                                                                    //
77 //                                                                           //
78 ///////////////////////////////////////////////////////////////////////////////
79
80
81 #include "AliSimulation.h"
82 #include "AliRunLoader.h"
83 #include "AliRun.h"
84 #include "AliModule.h"
85 #include "AliGenerator.h"
86 #include "AliVertexGenFile.h"
87 #include "AliRunDigitizer.h"
88 #include "AliDigitizer.h"
89 #include <TObjString.h>
90
91
92 ClassImp(AliSimulation)
93
94
95 //_____________________________________________________________________________
96 AliSimulation::AliSimulation(const char* configFileName,
97                              const char* name, const char* title) :
98   TNamed(name, title),
99
100   fRunGeneration(kTRUE),
101   fRunSimulation(kTRUE),
102   fMakeSDigits("ALL"),
103   fMakeDigits("ALL"),
104   fMakeDigitsFromHits(""),
105   fStopOnError(kFALSE),
106
107   fNEvents(1),
108   fConfigFileName(configFileName),
109   fGAliceFileName("galice.root"),
110   fBkgrdFileNames(NULL),
111   fUseBkgrdVertex(kTRUE),
112   fRegionOfInterest(kTRUE)
113 {
114 // create simulation object with default parameters
115
116 }
117
118 //_____________________________________________________________________________
119 AliSimulation::AliSimulation(const AliSimulation& sim) :
120   TNamed(sim),
121
122   fRunGeneration(sim.fRunGeneration),
123   fRunSimulation(sim.fRunSimulation),
124   fMakeSDigits(sim.fMakeSDigits),
125   fMakeDigits(sim.fMakeDigits),
126   fMakeDigitsFromHits(sim.fMakeDigitsFromHits),
127   fStopOnError(sim.fStopOnError),
128
129   fNEvents(sim.fNEvents),
130   fConfigFileName(sim.fConfigFileName),
131   fGAliceFileName(sim.fGAliceFileName),
132   fBkgrdFileNames(NULL),
133   fUseBkgrdVertex(sim.fUseBkgrdVertex),
134   fRegionOfInterest(sim.fRegionOfInterest)
135 {
136 // copy constructor
137
138   fBkgrdFileNames = new TObjArray;
139   for (Int_t i = 0; i < sim.fBkgrdFileNames->GetEntriesFast(); i++) {
140     if (!sim.fBkgrdFileNames->At(i)) continue;
141     fBkgrdFileNames->Add(sim.fBkgrdFileNames->At(i)->Clone());
142   }
143 }
144
145 //_____________________________________________________________________________
146 AliSimulation& AliSimulation::operator = (const AliSimulation& sim)
147 {
148 // assignment operator
149
150   this->~AliSimulation();
151   new(this) AliSimulation(sim);
152   return *this;
153 }
154
155 //_____________________________________________________________________________
156 AliSimulation::~AliSimulation()
157 {
158 // clean up
159
160   if (fBkgrdFileNames) {
161     fBkgrdFileNames->Delete();
162     delete fBkgrdFileNames;
163   }
164 }
165
166
167 //_____________________________________________________________________________
168 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
169 {
170 // set the number of events for one run
171
172   fNEvents = nEvents;
173 }
174
175 //_____________________________________________________________________________
176 void AliSimulation::SetConfigFile(const char* fileName)
177 {
178 // set the name of the config file
179
180   fConfigFileName = fileName;
181 }
182
183 //_____________________________________________________________________________
184 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
185 {
186 // add a file with background events for merging
187
188   TObjString* fileNameStr = new TObjString(fileName);
189   fileNameStr->SetUniqueID(nSignalPerBkgrd);
190   if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
191   fBkgrdFileNames->Add(fileNameStr);
192 }
193
194
195 //_____________________________________________________________________________
196 Bool_t AliSimulation::Run(Int_t nEvents)
197 {
198 // run the generation, simulation and digitization
199
200   if (nEvents > 0) fNEvents = nEvents;
201
202   // generation and simulation -> hits
203   if (fRunGeneration) {
204     if (!RunSimulation()) if (fStopOnError) return kFALSE;
205   }
206
207   // hits -> summable digits
208   if (!fMakeSDigits.IsNull()) {
209     if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
210   }
211
212   // summable digits -> digits
213   if (!fMakeDigits.IsNull()) {
214     if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
215       if (fStopOnError) return kFALSE;
216     }
217   }
218
219   // hits -> digits
220   if (!fMakeDigitsFromHits.IsNull()) {
221     if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
222       Warning("Run", "Merging and direct creation of digits from hits " 
223               "was selected for some detectors. "
224               "No merging will be done for the following detectors: %s",
225               fMakeDigitsFromHits.Data());
226     }
227     if (!RunHitsDigitization(fMakeDigitsFromHits)) {
228       if (fStopOnError) return kFALSE;
229     }
230   }
231
232   return kTRUE;
233 }
234
235 //_____________________________________________________________________________
236 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
237 {
238 // run the generation and simulation
239
240   TStopwatch stopwatch;
241   stopwatch.Start();
242
243   if (!gAlice) {
244     Error("RunSimulation", "no gAlice object. Restart aliroot and try again.");
245     return kFALSE;
246   }
247   if (gAlice->Modules()->GetEntries() > 0) {
248     Error("RunSimulation", 
249           "gAlice was already run. Restart aliroot and try again.");
250     return kFALSE;
251   }
252
253   Info("RunSimulation", "initializing gAlice with config file %s",
254        fConfigFileName.Data());
255   gAlice->Init(fConfigFileName.Data());
256   AliRunLoader* runLoader = gAlice->GetRunLoader();
257   if (!runLoader) {
258     Error("RunSimulation", "gAlice has no run loader object. "
259           "Check your config file: %s", fConfigFileName.Data());
260     return kFALSE;
261   }
262   fGAliceFileName = runLoader->GetFileName();
263
264   if (!gAlice->Generator()) {
265     Error("RunSimulation", "gAlice has no generator object. "
266           "Check your config file: %s", fConfigFileName.Data());
267     return kFALSE;
268   }
269
270   // get vertex from background file in case of merging
271   if (fUseBkgrdVertex &&
272       fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
273     Int_t signalPerBkgrd = fBkgrdFileNames->At(0)->GetUniqueID();
274     const char* fileName = ((TObjString*)
275                             (fBkgrdFileNames->At(0)))->GetName();
276     Info("RunSimulation", "The vertex will be taken from the background "
277          "file %s with nSignalPerBackground = %d", 
278          fileName, signalPerBkgrd);
279     AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
280     gAlice->Generator()->SetVertexGenerator(vtxGen);
281   }
282
283   if (!fRunSimulation) {
284     gAlice->Generator()->SetTrackingFlag(0);
285   }
286
287   Info("RunSimulation", "running gAlice");
288   if (nEvents <= 0) nEvents = fNEvents;
289   gAlice->Run(nEvents);
290
291   delete runLoader;
292
293   Info("RunSimulation", "execution time:");
294   stopwatch.Print();
295
296   return kTRUE;
297 }
298
299 //_____________________________________________________________________________
300 Bool_t AliSimulation::RunSDigitization(const char* detectors)
301 {
302 // run the digitization and produce summable digits
303
304   TStopwatch stopwatch;
305   stopwatch.Start();
306
307   AliRunLoader* runLoader = LoadRun();
308   if (!runLoader) return kFALSE;
309
310   TString detStr = detectors;
311   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
312   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
313     AliModule* det = (AliModule*) detArray->At(iDet);
314     if (!det || !det->IsActive()) continue;
315     if (IsSelected(det->GetName(), detStr)) {
316       Info("RunSDigitization", "creating summable digits for %s", 
317            det->GetName());
318       det->Hits2SDigits();
319     }
320   }
321
322   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
323     Error("RunSDigitization", "the following detectors were not found: %s", 
324           detStr.Data());
325     if (fStopOnError) return kFALSE;
326   }
327
328   delete runLoader;
329
330   Info("RunSDigitization", "execution time:");
331   stopwatch.Print();
332
333   return kTRUE;
334 }
335
336
337 //_____________________________________________________________________________
338 Bool_t AliSimulation::RunDigitization(const char* detectors, 
339                                       const char* excludeDetectors)
340 {
341 // run the digitization and produce digits from sdigits
342
343   TStopwatch stopwatch;
344   stopwatch.Start();
345
346   while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
347   if (gAlice) delete gAlice;
348   gAlice = NULL;
349
350   Int_t nStreams = 1;
351   Int_t signalPerBkgrd = 1;
352   if (fBkgrdFileNames) {
353     nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
354     if (nStreams > 1) signalPerBkgrd = fBkgrdFileNames->At(0)->GetUniqueID();
355   }
356   AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
357   manager->SetInputStream(0, fGAliceFileName.Data());
358   for (Int_t iStream = 1; iStream < nStreams; iStream++) {
359     const char* fileName = ((TObjString*)
360                             (fBkgrdFileNames->At(iStream-1)))->GetName();
361     manager->SetInputStream(iStream, fileName);
362   }
363
364   TString detStr = detectors;
365   TString detExcl = excludeDetectors;
366   manager->GetInputStream(0)->ImportgAlice();
367   AliRunLoader* runLoader = 
368     AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
369   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
370   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
371     AliModule* det = (AliModule*) detArray->At(iDet);
372     if (!det || !det->IsActive()) continue;
373     if (IsSelected(det->GetName(), detStr) && 
374         !IsSelected(det->GetName(), detExcl)) {
375       AliDigitizer* digitizer = det->CreateDigitizer(manager);
376       if (!digitizer) {
377         Error("RunDigitization", "no digitizer for %s", det->GetName());
378         if (fStopOnError) return kFALSE;
379       } else {
380         digitizer->SetRegionOfInterest(fRegionOfInterest);
381       }
382     }
383   }
384
385   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
386     Error("RunDigitization", "the following detectors were not found: %s", 
387           detStr.Data());
388     if (fStopOnError) return kFALSE;
389   }
390
391   if (!manager->GetListOfTasks()->IsEmpty()) {
392     Info("RunDigitization", "executing digitization");
393     manager->Exec("");
394   }
395
396   delete manager;
397
398   Info("RunDigitization", "execution time:");
399   stopwatch.Print();
400
401   return kTRUE;
402 }
403
404 //_____________________________________________________________________________
405 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
406 {
407 // run the digitization and produce digits from hits
408
409   TStopwatch stopwatch;
410   stopwatch.Start();
411
412   AliRunLoader* runLoader = LoadRun();
413   if (!runLoader) return kFALSE;
414
415   TString detStr = detectors;
416   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
417   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
418     AliModule* det = (AliModule*) detArray->At(iDet);
419     if (!det || !det->IsActive()) continue;
420     if (IsSelected(det->GetName(), detStr)) {
421       Info("RunHitsDigitization", "creating digits from hits for %s", 
422            det->GetName());
423       det->Hits2Digits();
424     }
425   }
426
427   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
428     Error("RunHitsDigitization", "the following detectors were not found: %s", 
429           detStr.Data());
430     if (fStopOnError) return kFALSE;
431   }
432
433   delete runLoader;
434
435   Info("RunHitsDigitization", "execution time:");
436   stopwatch.Print();
437
438   return kTRUE;
439 }
440
441
442 //_____________________________________________________________________________
443 AliRunLoader* AliSimulation::LoadRun() const
444 {
445 // delete existing run loaders, open a new one and load gAlice
446
447   while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
448   AliRunLoader* runLoader = 
449     AliRunLoader::Open(fGAliceFileName.Data(), 
450                        AliConfig::fgkDefaultEventFolderName, "UPDATE");
451   if (!runLoader) {
452     Error("LoadRun", "no run loader found in file %s", 
453           fGAliceFileName.Data());
454     return NULL;
455   }
456   runLoader->LoadgAlice();
457   gAlice = runLoader->GetAliRun();
458   if (!gAlice) {
459     Error("LoadRun", "no gAlice object found in file %s", 
460           fGAliceFileName.Data());
461     return NULL;
462   }
463   return runLoader;
464 }
465
466 //_____________________________________________________________________________
467 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
468 {
469 // check whether detName is contained in detectors
470 // if yes, it is removed from detectors
471
472   // check if all detectors are selected
473   if ((detectors.CompareTo("ALL") == 0) ||
474       detectors.BeginsWith("ALL ") ||
475       detectors.EndsWith(" ALL") ||
476       detectors.Contains(" ALL ")) {
477     detectors = "ALL";
478     return kTRUE;
479   }
480
481   // search for the given detector
482   Bool_t result = kFALSE;
483   if ((detectors.CompareTo(detName) == 0) ||
484       detectors.BeginsWith(detName+" ") ||
485       detectors.EndsWith(" "+detName) ||
486       detectors.Contains(" "+detName+" ")) {
487     detectors.ReplaceAll(detName, "");
488     result = kTRUE;
489   }
490
491   // clean up the detectors string
492   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
493   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
494   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
495
496   return result;
497 }