]> git.uio.no Git - u/mrichter/AliRoot.git/blob - STEER/AliSimulation.cxx
Added the writing of the trigger DDL raw data. New raw data header class for simulati...
[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 // Background 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 // By default this number is calculated from the number of available         //
71 // background events. MergeWith can be called several times to merge more    //
72 // than two event streams. It is assumed that the sdigits were already       //
73 // produced for the background events.                                       //
74 //                                                                           //
75 // The output of raw data can be switched on by calling                      //
76 //                                                                           //
77 //   sim.SetWriteRawData("MUON");   // write raw data for MUON               //
78 //                                                                           //
79 // The default output format of the raw data are DDL files. They are         //
80 // converted to a DATE file, if a file name is given as second argument.     //
81 // For this conversion the program "dateStream" is required. If the file     //
82 // name has the extension ".root", the DATE file is converted to a root      //
83 // file. The program "alimdc" is used for this purpose. For the conversion   //
84 // to DATE and root format the two conversion programs have to be installed. //
85 // Only the raw data in the final format is kept if the third argument is    //
86 // kTRUE.                                                                    //
87 //                                                                           //
88 // The methods RunSimulation, RunSDigitization, RunDigitization,             //
89 // RunHitsDigitization and WriteRawData can be used to run only parts of     //
90 // the full simulation chain. The creation of raw data DDL files and their   //
91 // conversion to the DATE or root format can be run directly by calling      //
92 // the methods WriteRawFiles, ConvertRawFilesToDate and ConvertDateToRoot.   //
93 //                                                                           //
94 // The default number of events per file, which is usually set in the        //
95 // config file, can be changed for individual detectors and data types       //
96 // by calling                                                                //
97 //                                                                           //
98 //   sim.SetEventsPerFile("PHOS", "Reconstructed Points", 3);                //
99 //                                                                           //
100 // The first argument is the detector, the second one the data type and the  //
101 // last one the number of events per file. Valid data types are "Hits",      //
102 // "Summable Digits", "Digits", "Reconstructed Points" and "Tracks".         //
103 // The number of events per file has to be set before the simulation of      //
104 // hits. Otherwise it has no effect.                                         //
105 //                                                                           //
106 ///////////////////////////////////////////////////////////////////////////////
107
108 #include <TGeoManager.h>
109 #include <TObjString.h>
110 #include <TStopwatch.h>
111 #include <TSystem.h>
112 #include <TFile.h>
113
114 #include "AliCDBStorage.h"
115 #include "AliCDBEntry.h"
116 #include "AliCDBManager.h"
117 #include "AliAlignObj.h"
118 #include "AliCentralTrigger.h"
119 #include "AliDAQConfig.h"
120 #include "AliDigitizer.h"
121 #include "AliGenerator.h"
122 #include "AliLog.h"
123 #include "AliModule.h"
124 #include "AliRun.h"
125 #include "AliRunDigitizer.h"
126 #include "AliRunLoader.h"
127 #include "AliSimulation.h"
128 #include "AliVertexGenFile.h"
129 #include "AliCentralTrigger.h"
130 #include "AliCTPRawData.h"
131
132 ClassImp(AliSimulation)
133
134
135 //_____________________________________________________________________________
136 AliSimulation::AliSimulation(const char* configFileName, const char* cdbUri,
137                              const char* name, const char* title) :
138   TNamed(name, title),
139
140   fRunGeneration(kTRUE),
141   fRunSimulation(kTRUE),
142   fLoadAlignFromCDB(kTRUE),
143   fLoadAlignData("ALL"),
144   fMakeSDigits("ALL"),
145   fMakeDigits("ALL"),
146   fMakeTrigger(""),
147   fMakeDigitsFromHits(""),
148   fWriteRawData(""),
149   fRawDataFileName(""),
150   fDeleteIntermediateFiles(kFALSE),
151   fStopOnError(kFALSE),
152
153   fNEvents(1),
154   fConfigFileName(configFileName),
155   fGAliceFileName("galice.root"),
156   fEventsPerFile(),
157   fBkgrdFileNames(NULL),
158   fAlignObjArray(NULL),
159   fUseBkgrdVertex(kTRUE),
160   fRegionOfInterest(kFALSE),
161   fCDBUri(cdbUri)
162 {
163 // create simulation object with default parameters
164
165   SetGAliceFile("galice.root");
166 }
167
168 //_____________________________________________________________________________
169 AliSimulation::AliSimulation(const AliSimulation& sim) :
170   TNamed(sim),
171
172   fRunGeneration(sim.fRunGeneration),
173   fRunSimulation(sim.fRunSimulation),
174   fLoadAlignFromCDB(sim.fLoadAlignFromCDB),
175   fLoadAlignData(sim.fLoadAlignData),
176   fMakeSDigits(sim.fMakeSDigits),
177   fMakeDigits(sim.fMakeDigits),
178   fMakeTrigger(sim.fMakeTrigger),
179   fMakeDigitsFromHits(sim.fMakeDigitsFromHits),
180   fWriteRawData(sim.fWriteRawData),
181   fRawDataFileName(""),
182   fDeleteIntermediateFiles(kFALSE),
183   fStopOnError(sim.fStopOnError),
184
185   fNEvents(sim.fNEvents),
186   fConfigFileName(sim.fConfigFileName),
187   fGAliceFileName(sim.fGAliceFileName),
188   fEventsPerFile(),
189   fBkgrdFileNames(NULL),
190   fAlignObjArray(NULL),
191   fUseBkgrdVertex(sim.fUseBkgrdVertex),
192   fRegionOfInterest(sim.fRegionOfInterest),
193   fCDBUri(sim.fCDBUri)
194 {
195 // copy constructor
196
197   for (Int_t i = 0; i < sim.fEventsPerFile.GetEntriesFast(); i++) {
198     if (!sim.fEventsPerFile[i]) continue;
199     fEventsPerFile.Add(sim.fEventsPerFile[i]->Clone());
200   }
201
202   fBkgrdFileNames = new TObjArray;
203   for (Int_t i = 0; i < sim.fBkgrdFileNames->GetEntriesFast(); i++) {
204     if (!sim.fBkgrdFileNames->At(i)) continue;
205     fBkgrdFileNames->Add(sim.fBkgrdFileNames->At(i)->Clone());
206   }
207 }
208
209 //_____________________________________________________________________________
210 AliSimulation& AliSimulation::operator = (const AliSimulation& sim)
211 {
212 // assignment operator
213
214   this->~AliSimulation();
215   new(this) AliSimulation(sim);
216   return *this;
217 }
218
219 //_____________________________________________________________________________
220 AliSimulation::~AliSimulation()
221 {
222 // clean up
223
224   fEventsPerFile.Delete();
225 //  if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
226 //  delete fAlignObjArray; fAlignObjArray=0;
227
228   if (fBkgrdFileNames) {
229     fBkgrdFileNames->Delete();
230     delete fBkgrdFileNames;
231   }
232 }
233
234
235 //_____________________________________________________________________________
236 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
237 {
238 // set the number of events for one run
239
240   fNEvents = nEvents;
241 }
242
243 //_____________________________________________________________________________
244 void AliSimulation::InitCDBStorage(const char* uri)
245 {
246 // activate a default CDB storage
247 // First check if we have any CDB storage set, because it is used 
248 // to retrieve the calibration and alignment constants
249
250   AliCDBManager* man = AliCDBManager::Instance();
251   if (!man->IsDefaultStorageSet())
252   {
253     AliWarningClass("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
254     AliWarningClass("Default CDB storage not yet set");
255     AliWarningClass(Form("Using default storage declared in AliSimulation: %s",uri));
256     AliWarningClass("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
257     SetDefaultStorage(uri);
258   }  
259   
260 }
261
262 //_____________________________________________________________________________
263 void AliSimulation::SetDefaultStorage(const char* uri) {
264 // activate a default CDB storage 
265
266    AliCDBManager::Instance()->SetDefaultStorage(uri);
267
268 }
269
270 //_____________________________________________________________________________
271 void AliSimulation::SetSpecificStorage(const char* detName, const char* uri) {
272 // activate a detector-specific CDB storage 
273
274    AliCDBManager::Instance()->SetSpecificStorage(detName, uri);
275
276 }
277
278 //_____________________________________________________________________________
279 void AliSimulation::SetConfigFile(const char* fileName)
280 {
281 // set the name of the config file
282
283   fConfigFileName = fileName;
284 }
285
286 //_____________________________________________________________________________
287 void AliSimulation::SetGAliceFile(const char* fileName)
288 {
289 // set the name of the galice file
290 // the path is converted to an absolute one if it is relative
291
292   fGAliceFileName = fileName;
293   if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
294     char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
295                                                 fGAliceFileName);
296     fGAliceFileName = absFileName;
297     delete[] absFileName;
298   }
299
300   AliDebug(2, Form("galice file name set to %s", fileName));
301 }
302
303 //_____________________________________________________________________________
304 void AliSimulation::SetEventsPerFile(const char* detector, const char* type, 
305                                      Int_t nEvents)
306 {
307 // set the number of events per file for the given detector and data type
308 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
309
310   TNamed* obj = new TNamed(detector, type);
311   obj->SetUniqueID(nEvents);
312   fEventsPerFile.Add(obj);
313 }
314
315 //_____________________________________________________________________________
316 Bool_t AliSimulation::ApplyAlignObjsToGeom(TObjArray* alObjArray)
317 {
318   // Read collection of alignment objects (AliAlignObj derived) saved
319   // in the TClonesArray ClArrayName and apply them to the geometry
320   // manager singleton.
321   //
322   alObjArray->Sort();
323   Int_t nvols = alObjArray->GetEntriesFast();
324
325   Bool_t flag = kTRUE;
326
327   for(Int_t j=0; j<nvols; j++)
328     {
329       AliAlignObj* alobj = (AliAlignObj*) alObjArray->UncheckedAt(j);
330       if (alobj->ApplyToGeometry() == kFALSE) flag = kFALSE;
331     }
332
333   if (AliDebugLevelClass() >= 1) {
334     gGeoManager->CheckOverlaps(20);
335     TObjArray* ovexlist = gGeoManager->GetListOfOverlaps();
336     if(ovexlist->GetEntriesFast()){  
337       AliError("The application of alignment objects to the geometry caused huge overlaps/extrusions!");
338    }
339   }
340
341   return flag;
342
343 }
344
345 //_____________________________________________________________________________
346 Bool_t AliSimulation::ApplyAlignObjsToGeom(const char* fileName, const char* clArrayName)
347 {
348   // read collection of alignment objects (AliAlignObj derived) saved
349   // in the TClonesArray ClArrayName in the file fileName and apply
350   // them to the TGeo geometry passed as argument
351   //
352
353   TFile* inFile = TFile::Open(fileName,"READ");
354   if (!inFile || !inFile->IsOpen()) {
355     AliErrorClass(Form("Could not open file %s !",fileName));
356     return kFALSE;
357   }
358
359   TClonesArray* alObjArray = ((TClonesArray*) inFile->Get(clArrayName));
360   inFile->Close();
361   if (!alObjArray) {
362     AliErrorClass(Form("Could not get array (%s) from file (%s) !",clArrayName,fileName));
363     return kFALSE;
364   }
365
366   return ApplyAlignObjsToGeom(alObjArray);
367
368 }
369
370 //_____________________________________________________________________________
371 Bool_t AliSimulation::ApplyAlignObjsToGeom(AliCDBParam* param, AliCDBId& Id)
372 {
373   // read collection of alignment objects (AliAlignObj derived) saved
374   // in the TClonesArray ClArrayName in the AliCDBEntry identified by
375   // param (to get the AliCDBStorage) and Id; apply the alignment objects
376   // to the TGeo geometry passed as argument
377   //
378
379   AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage(param);
380   AliCDBEntry* entry = storage->Get(Id);
381   TClonesArray* AlObjArray = ((TClonesArray*) entry->GetObject());
382
383   return ApplyAlignObjsToGeom(AlObjArray);
384
385 }
386
387 //_____________________________________________________________________________
388 Bool_t AliSimulation::ApplyAlignObjsToGeom(const char* uri, const char* path, Int_t runnum, Int_t version, Int_t sversion)
389 {
390   // read collection of alignment objects (AliAlignObj derived) saved
391   // in the TClonesArray ClArrayName in the AliCDBEntry identified by
392   // param (to get the AliCDBStorage) and Id; apply the alignment objects
393   // to the TGeo geometry passed as argument
394   //
395
396   AliCDBParam* param = AliCDBManager::Instance()->CreateParameter(uri);
397   AliCDBId id(path, runnum, runnum, version, sversion);
398
399   return ApplyAlignObjsToGeom(param, id);
400
401 }
402
403 //_____________________________________________________________________________
404 Bool_t AliSimulation::ApplyAlignObjsToGeom(const char* detName, Int_t runnum, Int_t version, Int_t sversion)
405 {
406   // read collection of alignment objects (AliAlignObj derived) saved
407   // in the TClonesArray ClArrayName in the AliCDBEntry identified by
408   // param (to get the AliCDBStorage) and Id; apply the alignment objects
409   // to the TGeo geometry passed as argument
410   //
411
412   InitCDBStorage("local://$ALICE_ROOT");
413   AliCDBPath path(detName,"Align","Data");
414   AliCDBEntry* entry = AliCDBManager::Instance()->Get(path.GetPath(),runnum,version,sversion);
415
416   if(!entry) return kFALSE;
417   TClonesArray* AlObjArray = ((TClonesArray*) entry->GetObject());
418
419   return ApplyAlignObjsToGeom(AlObjArray);
420 }
421
422 //_____________________________________________________________________________
423 Bool_t AliSimulation::SetAlignObjArraySingleDet(const char* detName)
424 {
425   // Fills array of single detector's alignable objects from CDB
426   
427   AliDebug(2, Form("Loading alignment data for detector: %s",detName));
428   
429   AliCDBEntry *entry;
430         
431   AliCDBPath path(detName,"Align","Data");
432         
433   entry=AliCDBManager::Instance()->Get(path.GetPath());
434   if(!entry){ 
435         AliDebug(2,Form("Couldn't load alignment data for detector %s",detName));
436         return kFALSE;
437   }
438   entry->SetOwner(1);
439   TClonesArray *alignArray = (TClonesArray*) entry->GetObject();        
440   alignArray->SetOwner(0);
441   AliDebug(2,Form("Found %d alignment objects for %s",
442                         alignArray->GetEntries(),detName));
443
444   AliAlignObj *alignObj=0;
445   TIter iter(alignArray);
446         
447   // loop over align objects in detector
448   while( ( alignObj=(AliAlignObj *) iter.Next() ) ){
449         fAlignObjArray->Add(alignObj);
450   }
451   // delete entry --- Don't delete, it is cached!
452         
453   AliDebug(2, Form("fAlignObjArray entries: %d",fAlignObjArray->GetEntries() ));
454   return kTRUE;
455
456 }
457
458 //_____________________________________________________________________________
459 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
460 {
461   // Read the alignment objects from CDB.
462   // Each detector is supposed to have the
463   // alignment objects in DET/Align/Data CDB path.
464   // All the detector objects are then collected,
465   // sorted by geometry level (starting from ALIC) and
466   // then applied to the TGeo geometry.
467   // Finally an overlaps check is performed.
468
469   Bool_t delRunLoader = kFALSE;
470   if (!runLoader) {
471     runLoader = LoadRun("READ");
472     if (!runLoader) return kFALSE;
473     delRunLoader = kTRUE;
474   }
475
476   // Load alignment data from CDB and fill fAlignObjArray 
477   if(fLoadAlignFromCDB){
478         if(!fAlignObjArray) fAlignObjArray = new TObjArray();
479         
480         //fAlignObjArray->RemoveAll(); 
481         fAlignObjArray->Clear();        
482         fAlignObjArray->SetOwner(0);
483  
484         TString detStr = fLoadAlignData;
485         TString dataNotLoaded="";
486         TString dataLoaded="";
487   
488         TObjArray* detArray = runLoader->GetAliRun()->Detectors();
489         for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
490                 AliModule* det = (AliModule*) detArray->At(iDet);
491                 if (!det || !det->IsActive()) continue;
492                 if (IsSelected(det->GetName(), detStr)) {
493                         if(!SetAlignObjArraySingleDet(det->GetName())){
494                                 dataNotLoaded += det->GetName();
495                                 dataNotLoaded += " ";
496                         } else {
497                                 dataLoaded += det->GetName();
498                                 dataLoaded += " ";
499                         }
500                 }
501         } // end loop over detectors
502   
503         if ((detStr.CompareTo("ALL") == 0)) detStr = "";
504         dataNotLoaded += detStr;
505         AliInfo(Form("Alignment data loaded for: %s",
506                           dataLoaded.Data()));
507         AliInfo(Form("Didn't/couldn't load alignment data for: %s",
508                           dataNotLoaded.Data()));
509   } // fLoadAlignFromCDB flag
510  
511   // Check if the array with alignment objects was
512   // provided by the user. If yes, apply the objects
513   // to the present TGeo geometry
514   if (fAlignObjArray) {
515     if (gGeoManager && gGeoManager->IsClosed()) {
516       if (ApplyAlignObjsToGeom(fAlignObjArray) == kFALSE) {
517         AliError("The misalignment of one or more volumes failed!"
518                  "Compare the list of simulated detectors and the list of detector alignment data!");
519         if (delRunLoader) delete runLoader;
520         return kFALSE;
521       }
522     }
523     else {
524       AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
525       if (delRunLoader) delete runLoader;
526       return kFALSE;
527     }
528   }
529
530   if (delRunLoader) delete runLoader;
531
532   return kTRUE;
533 }
534
535
536 //_____________________________________________________________________________
537 Bool_t AliSimulation::SetRunNumber()
538 {
539   // Set the CDB manager run number
540   // The run number is retrieved from gAlice
541
542   if(AliCDBManager::Instance()->GetRun() < 0) {
543     AliRunLoader* runLoader = LoadRun("READ");
544     if (!runLoader) return kFALSE;
545     else {
546       AliCDBManager::Instance()->SetRun(runLoader->GetAliRun()->GetRunNumber());
547       AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
548       delete runLoader;
549     }
550   }
551   return kTRUE;
552 }
553
554 //_____________________________________________________________________________
555 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
556 {
557 // add a file with background events for merging
558
559   TObjString* fileNameStr = new TObjString(fileName);
560   fileNameStr->SetUniqueID(nSignalPerBkgrd);
561   if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
562   fBkgrdFileNames->Add(fileNameStr);
563 }
564
565
566 //_____________________________________________________________________________
567 Bool_t AliSimulation::Run(Int_t nEvents)
568 {
569 // run the generation, simulation and digitization
570
571   InitCDBStorage(fCDBUri);
572
573   if (nEvents > 0) fNEvents = nEvents;
574
575   // generation and simulation -> hits
576   if (fRunGeneration) {
577     if (!RunSimulation()) if (fStopOnError) return kFALSE;
578   }
579
580   // Set run number in CDBManager (if it is not already set in RunSimulation)
581   if (!SetRunNumber()) if (fStopOnError) return kFALSE;
582
583   // Load and misalign the geometry
584   if (!gGeoManager) {
585     TGeoManager::Import("geometry.root");
586     if (!gGeoManager) if (fStopOnError) return kFALSE;
587     if (!MisalignGeometry()) if (fStopOnError) return kFALSE;
588   }
589
590   // hits -> summable digits
591   if (!fMakeSDigits.IsNull()) {
592     if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
593   }
594
595   // summable digits -> digits
596   if (!fMakeDigits.IsNull()) {
597     if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
598       if (fStopOnError) return kFALSE;
599     }
600   }
601
602   // hits -> digits
603   if (!fMakeDigitsFromHits.IsNull()) {
604     if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
605       AliWarning(Form("Merging and direct creation of digits from hits " 
606                  "was selected for some detectors. "
607                  "No merging will be done for the following detectors: %s",
608                  fMakeDigitsFromHits.Data()));
609     }
610     if (!RunHitsDigitization(fMakeDigitsFromHits)) {
611       if (fStopOnError) return kFALSE;
612     }
613   }
614
615   // digits -> trigger
616   if (!fMakeTrigger.IsNull()) {
617     if (!RunTrigger(fMakeTrigger)) {
618       if (fStopOnError) return kFALSE;
619     }
620   }
621
622   // digits -> raw data
623   if (!fWriteRawData.IsNull()) {
624     if (!WriteRawData(fWriteRawData, fRawDataFileName, 
625                       fDeleteIntermediateFiles)) {
626       if (fStopOnError) return kFALSE;
627     }
628   }
629
630   return kTRUE;
631 }
632
633 //_____________________________________________________________________________
634 Bool_t AliSimulation::RunTrigger(const char* descriptors)
635 {
636   // run the trigger
637
638    TStopwatch stopwatch;
639    stopwatch.Start();
640
641    AliRunLoader* runLoader = LoadRun("READ");
642    if (!runLoader) return kFALSE;
643    TString des = descriptors;
644
645    runLoader->MakeTree( "CT" );
646    AliCentralTrigger* aCTP = runLoader->GetTrigger();
647   // Load Descriptors
648    aCTP->LoadDescriptor( des );
649
650   // digits -> trigger
651    if( !aCTP->RunTrigger( runLoader ) ) {
652       if (fStopOnError) {
653     //  delete aCTP;
654          return kFALSE;
655       }
656    }
657
658    AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
659            stopwatch.RealTime(),stopwatch.CpuTime()));
660
661    delete runLoader;
662
663    return kTRUE;
664 }
665
666 //_____________________________________________________________________________
667 Bool_t AliSimulation::WriteTriggerRawData()
668 {
669   // Writes the CTP (trigger) DDL raw data
670   // Details of the format are given in the
671   // trigger TDR - pages 134 and 135.
672   AliCTPRawData writer;
673   writer.RawData();
674
675   return kTRUE;
676 }
677
678 //_____________________________________________________________________________
679 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
680 {
681 // run the generation and simulation
682
683   TStopwatch stopwatch;
684   stopwatch.Start();
685
686   if (!gAlice) {
687     AliError("no gAlice object. Restart aliroot and try again.");
688     return kFALSE;
689   }
690   if (gAlice->Modules()->GetEntries() > 0) {
691     AliError("gAlice was already run. Restart aliroot and try again.");
692     return kFALSE;
693   }
694
695   AliInfo(Form("initializing gAlice with config file %s",
696           fConfigFileName.Data()));
697   StdoutToAliInfo(StderrToAliError(
698     gAlice->Init(fConfigFileName.Data());
699   ););
700
701   // Set run number in CDBManager
702   AliCDBManager::Instance()->SetRun(gAlice->GetRunNumber());
703   AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
704
705   AliRunLoader* runLoader = gAlice->GetRunLoader();
706   if (!runLoader) {
707              AliError(Form("gAlice has no run loader object. "
708                              "Check your config file: %s", fConfigFileName.Data()));
709              return kFALSE;
710   }
711   SetGAliceFile(runLoader->GetFileName());
712  
713   // Export ideal geometry 
714   if (gGeoManager) gGeoManager->Export("geometry.root");
715
716   // Misalign geometry
717   if (!MisalignGeometry(runLoader)) return kFALSE;
718
719   // Temporary fix by A.Gheata
720   // Could be removed with the next Root version (>5.11)
721   if (gGeoManager) {
722     TIter next(gGeoManager->GetListOfVolumes());
723     TGeoVolume *vol;
724     while ((vol = (TGeoVolume *)next())) {
725       if (vol->GetVoxels()) {
726         if (vol->GetVoxels()->NeedRebuild()) {
727           vol->GetVoxels()->Voxelize();
728           vol->FindOverlaps();
729         }
730       }
731     }
732   }
733
734   // Export (mis)aligned geometry 
735   if (gGeoManager) gGeoManager->Export("misaligned_geometry.root");
736
737 //   AliRunLoader* runLoader = gAlice->GetRunLoader();
738 //   if (!runLoader) {
739 //     AliError(Form("gAlice has no run loader object. "
740 //                   "Check your config file: %s", fConfigFileName.Data()));
741 //     return kFALSE;
742 //   }
743 //   SetGAliceFile(runLoader->GetFileName());
744
745   if (!gAlice->Generator()) {
746     AliError(Form("gAlice has no generator object. "
747                   "Check your config file: %s", fConfigFileName.Data()));
748     return kFALSE;
749   }
750   if (nEvents <= 0) nEvents = fNEvents;
751
752   // get vertex from background file in case of merging
753   if (fUseBkgrdVertex &&
754       fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
755     Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
756     const char* fileName = ((TObjString*)
757                             (fBkgrdFileNames->At(0)))->GetName();
758     AliInfo(Form("The vertex will be taken from the background "
759                  "file %s with nSignalPerBackground = %d", 
760                  fileName, signalPerBkgrd));
761     AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
762     gAlice->Generator()->SetVertexGenerator(vtxGen);
763   }
764
765   if (!fRunSimulation) {
766     gAlice->Generator()->SetTrackingFlag(0);
767   }
768
769   // set the number of events per file for given detectors and data types
770   for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
771     if (!fEventsPerFile[i]) continue;
772     const char* detName = fEventsPerFile[i]->GetName();
773     const char* typeName = fEventsPerFile[i]->GetTitle();
774     TString loaderName(detName);
775     loaderName += "Loader";
776     AliLoader* loader = runLoader->GetLoader(loaderName);
777     if (!loader) {
778       AliError(Form("RunSimulation", "no loader for %s found\n"
779                     "Number of events per file not set for %s %s", 
780                     detName, typeName, detName));
781       continue;
782     }
783     AliDataLoader* dataLoader = 
784       loader->GetDataLoader(typeName);
785     if (!dataLoader) {
786       AliError(Form("no data loader for %s found\n"
787                     "Number of events per file not set for %s %s", 
788                     typeName, detName, typeName));
789       continue;
790     }
791     dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
792     AliDebug(1, Form("number of events per file set to %d for %s %s",
793                      fEventsPerFile[i]->GetUniqueID(), detName, typeName));
794   }
795
796   AliInfo("running gAlice");
797   StdoutToAliInfo(StderrToAliError(
798     gAlice->Run(nEvents);
799   ););
800
801   delete runLoader;
802
803   AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
804                stopwatch.RealTime(),stopwatch.CpuTime()));
805
806   return kTRUE;
807 }
808
809 //_____________________________________________________________________________
810 Bool_t AliSimulation::RunSDigitization(const char* detectors)
811 {
812 // run the digitization and produce summable digits
813
814   TStopwatch stopwatch;
815   stopwatch.Start();
816
817   AliRunLoader* runLoader = LoadRun();
818   if (!runLoader) return kFALSE;
819
820   TString detStr = detectors;
821   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
822   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
823     AliModule* det = (AliModule*) detArray->At(iDet);
824     if (!det || !det->IsActive()) continue;
825     if (IsSelected(det->GetName(), detStr)) {
826       AliInfo(Form("creating summable digits for %s", det->GetName()));
827       TStopwatch stopwatchDet;
828       stopwatchDet.Start();
829       det->Hits2SDigits();
830       AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
831            det->GetName(),stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
832     }
833   }
834
835   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
836     AliError(Form("the following detectors were not found: %s",
837                   detStr.Data()));
838     if (fStopOnError) return kFALSE;
839   }
840
841   delete runLoader;
842
843   AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
844            stopwatch.RealTime(),stopwatch.CpuTime()));
845
846   return kTRUE;
847 }
848
849
850 //_____________________________________________________________________________
851 Bool_t AliSimulation::RunDigitization(const char* detectors, 
852                                       const char* excludeDetectors)
853 {
854 // run the digitization and produce digits from sdigits
855
856   TStopwatch stopwatch;
857   stopwatch.Start();
858
859   while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
860   if (gAlice) delete gAlice;
861   gAlice = NULL;
862
863   Int_t nStreams = 1;
864   if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
865   Int_t signalPerBkgrd = GetNSignalPerBkgrd();
866   AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
867   manager->SetInputStream(0, fGAliceFileName.Data());
868   for (Int_t iStream = 1; iStream < nStreams; iStream++) {
869     const char* fileName = ((TObjString*)
870                             (fBkgrdFileNames->At(iStream-1)))->GetName();
871     manager->SetInputStream(iStream, fileName);
872   }
873
874   TString detStr = detectors;
875   TString detExcl = excludeDetectors;
876   manager->GetInputStream(0)->ImportgAlice();
877   AliRunLoader* runLoader = 
878     AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
879   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
880   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
881     AliModule* det = (AliModule*) detArray->At(iDet);
882     if (!det || !det->IsActive()) continue;
883     if (IsSelected(det->GetName(), detStr) && 
884         !IsSelected(det->GetName(), detExcl)) {
885       AliDigitizer* digitizer = det->CreateDigitizer(manager);
886       if (!digitizer) {
887         AliError(Form("no digitizer for %s", det->GetName()));
888         if (fStopOnError) return kFALSE;
889       } else {
890         digitizer->SetRegionOfInterest(fRegionOfInterest);
891       }
892     }
893   }
894
895   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
896     AliError(Form("the following detectors were not found: %s", 
897                   detStr.Data()));
898     if (fStopOnError) return kFALSE;
899   }
900
901   if (!manager->GetListOfTasks()->IsEmpty()) {
902     AliInfo("executing digitization");
903     manager->Exec("");
904   }
905
906   delete manager;
907
908   AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
909                stopwatch.RealTime(),stopwatch.CpuTime()));
910   
911   return kTRUE;
912 }
913
914 //_____________________________________________________________________________
915 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
916 {
917 // run the digitization and produce digits from hits
918
919   TStopwatch stopwatch;
920   stopwatch.Start();
921
922   AliRunLoader* runLoader = LoadRun("READ");
923   if (!runLoader) return kFALSE;
924
925   TString detStr = detectors;
926   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
927   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
928     AliModule* det = (AliModule*) detArray->At(iDet);
929     if (!det || !det->IsActive()) continue;
930     if (IsSelected(det->GetName(), detStr)) {
931       AliInfo(Form("creating digits from hits for %s", det->GetName()));
932       det->Hits2Digits();
933     }
934   }
935
936   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
937     AliError(Form("the following detectors were not found: %s", 
938                   detStr.Data()));
939     if (fStopOnError) return kFALSE;
940   }
941
942   delete runLoader;
943   //PH Temporary fix to avoid interference with the PHOS loder/getter
944   //PH The problem has to be solved in more general way 09/06/05
945
946   AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
947                stopwatch.RealTime(),stopwatch.CpuTime()));
948
949   return kTRUE;
950 }
951
952 //_____________________________________________________________________________
953 Bool_t AliSimulation::WriteRawData(const char* detectors, 
954                                    const char* fileName,
955                                    Bool_t deleteIntermediateFiles)
956 {
957 // convert the digits to raw data
958 // First DDL raw data files for the given detectors are created.
959 // If a file name is given, the DDL files are then converted to a DATE file.
960 // If deleteIntermediateFiles is true, the DDL raw files are deleted 
961 // afterwards.
962 // If the file name has the extension ".root", the DATE file is converted
963 // to a root file.
964 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
965
966   TStopwatch stopwatch;
967   stopwatch.Start();
968
969   if (!WriteRawFiles(detectors)) {
970     if (fStopOnError) return kFALSE;
971   }
972
973   TString dateFileName(fileName);
974   if (!dateFileName.IsNull()) {
975     Bool_t rootOutput = dateFileName.EndsWith(".root");
976     if (rootOutput) dateFileName += ".date";
977     if (!ConvertRawFilesToDate(dateFileName)) {
978       if (fStopOnError) return kFALSE;
979     }
980     if (deleteIntermediateFiles) {
981       AliRunLoader* runLoader = LoadRun("READ");
982       if (runLoader) for (Int_t iEvent = 0; 
983                           iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
984         char command[256];
985         sprintf(command, "rm -r raw%d", iEvent);
986         gSystem->Exec(command);
987       }
988     }
989
990     if (rootOutput) {
991       if (!ConvertDateToRoot(dateFileName, fileName)) {
992         if (fStopOnError) return kFALSE;
993       }
994       if (deleteIntermediateFiles) {
995         gSystem->Unlink(dateFileName);
996       }
997     }
998   }
999
1000   AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
1001                stopwatch.RealTime(),stopwatch.CpuTime()));
1002
1003   return kTRUE;
1004 }
1005
1006 //_____________________________________________________________________________
1007 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1008 {
1009 // convert the digits to raw data DDL files
1010
1011   AliRunLoader* runLoader = LoadRun("READ");
1012   if (!runLoader) return kFALSE;
1013
1014   // write raw data to DDL files
1015   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1016     AliInfo(Form("processing event %d", iEvent));
1017     runLoader->GetEvent(iEvent);
1018     TString baseDir = gSystem->WorkingDirectory();
1019     char dirName[256];
1020     sprintf(dirName, "raw%d", iEvent);
1021     gSystem->MakeDirectory(dirName);
1022     if (!gSystem->ChangeDirectory(dirName)) {
1023       AliError(Form("couldn't change to directory %s", dirName));
1024       if (fStopOnError) return kFALSE; else continue;
1025     }
1026
1027     TString detStr = detectors;
1028     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1029     for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1030       AliModule* det = (AliModule*) detArray->At(iDet);
1031       if (!det || !det->IsActive()) continue;
1032       if (IsSelected(det->GetName(), detStr)) {
1033         AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1034         det->Digits2Raw();
1035       }
1036     }
1037
1038     if (!WriteTriggerRawData())
1039       if (fStopOnError) return kFALSE;
1040
1041     gSystem->ChangeDirectory(baseDir);
1042     if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1043       AliError(Form("the following detectors were not found: %s", 
1044                     detStr.Data()));
1045       if (fStopOnError) return kFALSE;
1046     }
1047   }
1048
1049   delete runLoader;
1050   return kTRUE;
1051 }
1052
1053 //_____________________________________________________________________________
1054 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName)
1055 {
1056 // convert raw data DDL files to a DATE file with the program "dateStream"
1057
1058   char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1059   if (!path) {
1060     AliError("the program dateStream was not found");
1061     if (fStopOnError) return kFALSE;
1062   } else {
1063     delete[] path;
1064   }
1065
1066   AliRunLoader* runLoader = LoadRun("READ");
1067   if (!runLoader) return kFALSE;
1068
1069   AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1070   char command[256];
1071   sprintf(command, "dateStream -o %s -# %d -C", 
1072           dateFileName, runLoader->GetNumberOfEvents());
1073   FILE* pipe = gSystem->OpenPipe(command, "w");
1074
1075   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1076     fprintf(pipe, "GDC\n");
1077     Float_t ldc = 0;
1078     Int_t prevLDC = -1;
1079
1080     // loop over detectors and DDLs
1081     for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1082       for (Int_t iDDL = 0; iDDL < kDetectorDDLs[iDet]; iDDL++) {
1083
1084         Int_t ddlID = 0x100*iDet + iDDL;
1085         Int_t ldcID = Int_t(ldc + 0.0001);
1086         ldc += kDetectorLDCs[iDet] / kDetectorDDLs[iDet];
1087
1088         char rawFileName[256];
1089         sprintf(rawFileName, "raw%d/%s_%d.ddl", 
1090                 iEvent, kDetectors[iDet], ddlID);
1091
1092         // check existence and size of raw data file
1093         FILE* file = fopen(rawFileName, "rb");
1094         if (!file) continue;
1095         fseek(file, 0, SEEK_END);
1096         unsigned long size = ftell(file);
1097         fclose(file);
1098         if (!size) continue;
1099
1100         if (ldcID != prevLDC) {
1101           fprintf(pipe, " LDC Id %d\n", ldcID);
1102           prevLDC = ldcID;
1103         }
1104         fprintf(pipe, "  Equipment Id %d Payload %s\n", ddlID, rawFileName);
1105       }
1106     }
1107   }
1108
1109   Int_t result = gSystem->ClosePipe(pipe);
1110
1111   delete runLoader;
1112   return (result == 0);
1113 }
1114
1115 //_____________________________________________________________________________
1116 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1117                                         const char* rootFileName)
1118 {
1119 // convert a DATE file to a root file with the program "alimdc"
1120
1121   // ALIMDC setup
1122   const Int_t kDBSize = 1000000000;
1123   const Int_t kTagDBSize = 1000000000;
1124   const Bool_t kFilter = kFALSE;
1125   const Int_t kCompression = 1;
1126
1127   char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1128   if (!path) {
1129     AliError("the program alimdc was not found");
1130     if (fStopOnError) return kFALSE;
1131   } else {
1132     delete[] path;
1133   }
1134
1135   AliInfo(Form("converting DATE file %s to root file %s", 
1136                dateFileName, rootFileName));
1137
1138   const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1139   const char* tagDBFS    = "/tmp/mdc1/tags";
1140   const char* runDBFS    = "/tmp/mdc1/meta";
1141
1142   // User defined file system locations
1143   if (gSystem->Getenv("ALIMDC_RAWDB1")) 
1144     rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1145   if (gSystem->Getenv("ALIMDC_RAWDB2")) 
1146     rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1147   if (gSystem->Getenv("ALIMDC_TAGDB")) 
1148     tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1149   if (gSystem->Getenv("ALIMDC_RUNDB")) 
1150     runDBFS = gSystem->Getenv("ALIMDC_RUNDB");
1151
1152   gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1153   gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1154   gSystem->Exec(Form("rm -rf %s",tagDBFS));
1155   gSystem->Exec(Form("rm -rf %s",runDBFS));
1156
1157   gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1158   gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1159   gSystem->Exec(Form("mkdir %s",tagDBFS));
1160   gSystem->Exec(Form("mkdir %s",runDBFS));
1161
1162   Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s", 
1163                                     kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1164   gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1165
1166   gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1167   gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1168   gSystem->Exec(Form("rm -rf %s",tagDBFS));
1169   gSystem->Exec(Form("rm -rf %s",runDBFS));
1170
1171   return (result == 0);
1172 }
1173
1174
1175 //_____________________________________________________________________________
1176 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1177 {
1178 // delete existing run loaders, open a new one and load gAlice
1179
1180   while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
1181   AliRunLoader* runLoader = 
1182     AliRunLoader::Open(fGAliceFileName.Data(), 
1183                        AliConfig::GetDefaultEventFolderName(), mode);
1184   if (!runLoader) {
1185     AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1186     return NULL;
1187   }
1188   runLoader->LoadgAlice();
1189   gAlice = runLoader->GetAliRun();
1190   if (!gAlice) {
1191     AliError(Form("no gAlice object found in file %s", 
1192                   fGAliceFileName.Data()));
1193     return NULL;
1194   }
1195   return runLoader;
1196 }
1197
1198 //_____________________________________________________________________________
1199 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1200 {
1201 // get or calculate the number of signal events per background event
1202
1203   if (!fBkgrdFileNames) return 1;
1204   Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1205   if (nBkgrdFiles == 0) return 1;
1206
1207   // get the number of signal events
1208   if (nEvents <= 0) {
1209     AliRunLoader* runLoader = 
1210       AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1211     if (!runLoader) return 1;
1212     nEvents = runLoader->GetNumberOfEvents();
1213     delete runLoader;
1214   }
1215
1216   Int_t result = 0;
1217   for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1218     // get the number of background events
1219     const char* fileName = ((TObjString*)
1220                             (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1221     AliRunLoader* runLoader = 
1222       AliRunLoader::Open(fileName, "BKGRD");
1223     if (!runLoader) continue;
1224     Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1225     delete runLoader;
1226
1227     // get or calculate the number of signal per background events
1228     Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1229     if (nSignalPerBkgrd <= 0) {
1230       nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1231     } else if (result && (result != nSignalPerBkgrd)) {
1232       AliInfo(Form("the number of signal events per background event "
1233                    "will be changed from %d to %d for stream %d", 
1234                    nSignalPerBkgrd, result, iBkgrdFile+1));
1235       nSignalPerBkgrd = result;
1236     }
1237
1238     if (!result) result = nSignalPerBkgrd;
1239     if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1240       AliWarning(Form("not enough background events (%d) for %d signal events "
1241                       "using %d signal per background events for stream %d",
1242                       nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1243     }
1244   }
1245
1246   return result;
1247 }
1248
1249 //_____________________________________________________________________________
1250 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1251 {
1252 // check whether detName is contained in detectors
1253 // if yes, it is removed from detectors
1254
1255   // check if all detectors are selected
1256   if ((detectors.CompareTo("ALL") == 0) ||
1257       detectors.BeginsWith("ALL ") ||
1258       detectors.EndsWith(" ALL") ||
1259       detectors.Contains(" ALL ")) {
1260     detectors = "ALL";
1261     return kTRUE;
1262   }
1263
1264   // search for the given detector
1265   Bool_t result = kFALSE;
1266   if ((detectors.CompareTo(detName) == 0) ||
1267       detectors.BeginsWith(detName+" ") ||
1268       detectors.EndsWith(" "+detName) ||
1269       detectors.Contains(" "+detName+" ")) {
1270     detectors.ReplaceAll(detName, "");
1271     result = kTRUE;
1272   }
1273
1274   // clean up the detectors string
1275   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
1276   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1277   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1278
1279   return result;
1280 }