Some warning going to error and viceversa:
[u/mrichter/AliRoot.git] / STEER / 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 <TCint.h>
109 #include <TFile.h>
110 #include <TGeoGlobalMagField.h>
111 #include <TGeoManager.h>
112 #include <TObjString.h>
113 #include <TROOT.h>
114 #include <TSystem.h>
115 #include <TVirtualMC.h>
116 #include <TVirtualMCApplication.h>
117 #include <TDatime.h>
118 #include <TInterpreter.h>
119
120 #include "AliAlignObj.h"
121 #include "AliCDBEntry.h"
122 #include "AliCDBManager.h"
123 #include "AliGRPManager.h"
124 #include "AliCDBStorage.h"
125 #include "AliCTPRawData.h"
126 #include "AliCentralTrigger.h"
127 #include "AliCentralTrigger.h"
128 #include "AliCodeTimer.h"
129 #include "AliDAQ.h"
130 #include "AliDigitizer.h"
131 #include "AliESDEvent.h"
132 #include "AliGRPObject.h"
133 #include "AliGenEventHeader.h"
134 #include "AliGenerator.h"
135 #include "AliGeomManager.h"
136 #include "AliHLTSimulation.h"
137 #include "AliHeader.h"
138 #include "AliLego.h"
139 #include "AliLegoGenerator.h"
140 #include "AliLog.h"
141 #include "AliMC.h"
142 #include "AliMagF.h"
143 #include "AliModule.h"
144 #include "AliPDG.h"
145 #include "AliRawReaderDate.h"
146 #include "AliRawReaderFile.h"
147 #include "AliRawReaderRoot.h"
148 #include "AliRun.h"
149 #include "AliDigitizationInput.h"
150 #include "AliRunLoader.h"
151 #include "AliSimulation.h"
152 #include "AliSysInfo.h"
153 #include "AliVertexGenFile.h"
154
155 using std::ofstream;
156 ClassImp(AliSimulation)
157
158 AliSimulation *AliSimulation::fgInstance = 0;
159 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE"
160 // #ifdef MFT_UPGRADE
161 //                                                                             ,"MFT"
162 // #endif 
163                                                                             ,"MFT"    // AU
164                                                                             ,"HLT"
165 };
166
167 //_____________________________________________________________________________
168 AliSimulation::AliSimulation(const char* configFileName,
169                              const char* name, const char* title) :
170   TNamed(name, title),
171
172   fRunGeneration(kTRUE),
173   fRunSimulation(kTRUE),
174   fLoadAlignFromCDB(kTRUE),
175   fLoadAlObjsListOfDets("ALL"),
176   fMakeSDigits("ALL"),
177   fMakeDigits("ALL"),
178   fTriggerConfig(""),
179   fMakeDigitsFromHits(""),
180   fWriteRawData(""),
181   fRawDataFileName(""),
182   fDeleteIntermediateFiles(kFALSE),
183   fWriteSelRawData(kFALSE),
184   fStopOnError(kFALSE),
185   fUseMonitoring(kFALSE),
186   fNEvents(1),
187   fConfigFileName(configFileName),
188   fGAliceFileName("galice.root"),
189   fEventsPerFile(),
190   fBkgrdFileNames(NULL),
191   fAlignObjArray(NULL),
192   fUseBkgrdVertex(kTRUE),
193   fRegionOfInterest(kFALSE),
194   fCDBUri(""),
195   fQARefUri(""), 
196   fSpecCDBUri(),
197   fRun(-1),
198   fSeed(0),
199   fInitCDBCalled(kFALSE),
200   fInitRunNumberCalled(kFALSE),
201   fSetRunNumberFromDataCalled(kFALSE),
202   fEmbeddingFlag(kFALSE),
203   fLego(NULL),
204   fKey(0),
205   fUseVertexFromCDB(0),
206   fUseMagFieldFromGRP(0),
207   fGRPWriteLocation(Form("local://%s", gSystem->pwd())),
208   fUseTimeStampFromCDB(0),
209   fTimeStart(0),
210   fTimeEnd(0),
211   fQADetectors("ALL"),                  
212   fQATasks("ALL"),      
213   fRunQA(kTRUE), 
214   fEventSpecie(AliRecoParam::kDefault),
215   fWriteQAExpertData(kTRUE), 
216   fGeometryFile(),
217   fRunHLT("default"),
218   fpHLT(NULL),
219   fWriteGRPEntry(kTRUE)
220 {
221 // create simulation object with default parameters
222   fgInstance = this;
223   SetGAliceFile("galice.root");
224   
225 // for QA
226         AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ; 
227         qam->SetActiveDetectors(fQADetectors) ; 
228         fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ; 
229         qam->SetTasks(fQATasks) ;       
230 }
231
232 //_____________________________________________________________________________
233 AliSimulation::~AliSimulation()
234 {
235 // clean up
236
237   fEventsPerFile.Delete();
238 //  if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
239 //  delete fAlignObjArray; fAlignObjArray=0;
240
241   if (fBkgrdFileNames) {
242     fBkgrdFileNames->Delete();
243     delete fBkgrdFileNames;
244   }
245
246   fSpecCDBUri.Delete();
247   if (fgInstance==this) fgInstance = 0;
248
249   AliQAManager::QAManager()->ShowQA() ; 
250   AliQAManager::Destroy() ;     
251   AliCodeTimer::Instance()->Print();
252 }
253
254
255 //_____________________________________________________________________________
256 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
257 {
258 // set the number of events for one run
259
260   fNEvents = nEvents;
261 }
262
263 //_____________________________________________________________________________
264 void AliSimulation::InitQA()
265 {
266   // activate a default CDB storage
267   // First check if we have any CDB storage set, because it is used 
268   // to retrieve the calibration and alignment constants
269   
270   if (fInitCDBCalled) return;
271   fInitCDBCalled = kTRUE;
272
273   AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ; 
274   qam->SetActiveDetectors(fQADetectors) ; 
275   fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ; 
276   qam->SetTasks(fQATasks) ;
277         if (fWriteQAExpertData)
278     qam->SetWriteExpert() ; 
279   
280   if (qam->IsDefaultStorageSet()) {
281     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
282     AliWarning("Default QA reference storage has been already set !");
283     AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
284     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
285     fQARefUri = qam->GetDefaultStorage()->GetURI();
286   } else {
287       if (fQARefUri.Length() > 0) {
288         AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
289         AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
290         AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
291       } else {
292         fQARefUri="local://$ALICE_ROOT/QARef";
293         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
294         AliWarning("Default QA reference storage not yet set !!!!");
295         AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
296         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
297       }
298     qam->SetDefaultStorage(fQARefUri);
299   }
300 }
301
302 //_____________________________________________________________________________
303 void AliSimulation::InitCDB()
304 {
305 // activate a default CDB storage
306 // First check if we have any CDB storage set, because it is used 
307 // to retrieve the calibration and alignment constants
308
309   if (fInitCDBCalled) return;
310   fInitCDBCalled = kTRUE;
311
312   AliCDBManager* man = AliCDBManager::Instance();
313   if (man->IsDefaultStorageSet())
314   {
315     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
316     AliWarning("Default CDB storage has been already set !");
317     AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
318     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
319     fCDBUri = man->GetDefaultStorage()->GetURI();
320   }
321   else {
322     if (fCDBUri.Length() > 0) 
323     {
324         AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
325         AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
326         AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
327     } else {
328         fCDBUri="local://$ALICE_ROOT/OCDB";
329         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
330         AliWarning("Default CDB storage not yet set !!!!");
331         AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
332         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
333                 
334     }
335     man->SetDefaultStorage(fCDBUri);
336   }
337
338   // Now activate the detector specific CDB storage locations
339   for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
340     TObject* obj = fSpecCDBUri[i];
341     if (!obj) continue;
342     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
343     AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
344     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
345     man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
346   }
347       
348 }
349
350 //_____________________________________________________________________________
351 void AliSimulation::InitRunNumber(){
352 // check run number. If not set, set it to 0 !!!!
353   
354   if (fInitRunNumberCalled) return;
355   fInitRunNumberCalled = kTRUE;
356   
357   AliCDBManager* man = AliCDBManager::Instance();
358   if (man->GetRun() >= 0)
359   {
360         AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
361                         "Use external variable DC_RUN or AliSimulation::SetRun()!"));
362   }
363     
364   if(fRun >= 0) {
365         AliDebug(2, Form("Setting CDB run number to: %d",fRun));
366   } else {
367         fRun=0;
368         AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
369                         fRun));
370   }
371   man->SetRun(fRun);
372
373   man->Print();
374
375 }
376
377 //_____________________________________________________________________________
378 void AliSimulation::SetCDBLock() {
379   // Set CDB lock: from now on it is forbidden to reset the run number
380   // or the default storage or to activate any further storage!
381   
382   ULong64_t key = AliCDBManager::Instance()->SetLock(1);
383   if (key) fKey = key;
384 }
385
386 //_____________________________________________________________________________
387 void AliSimulation::SetDefaultStorage(const char* uri) {
388   // Store the desired default CDB storage location
389   // Activate it later within the Run() method
390   
391   fCDBUri = uri;
392   
393 }
394
395 //_____________________________________________________________________________
396 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
397   // Store the desired default CDB storage location
398   // Activate it later within the Run() method
399   
400   fQARefUri = uri;
401   AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
402 }
403
404 //_____________________________________________________________________________
405 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
406 // Store a detector-specific CDB storage location
407 // Activate it later within the Run() method
408
409   AliCDBPath aPath(calibType);
410   if(!aPath.IsValid()){
411         AliError(Form("Not a valid path: %s", calibType));
412         return;
413   }
414
415   TObject* obj = fSpecCDBUri.FindObject(calibType);
416   if (obj) fSpecCDBUri.Remove(obj);
417   fSpecCDBUri.Add(new TNamed(calibType, uri));
418
419 }
420
421 //_____________________________________________________________________________
422 void AliSimulation::SetRunNumber(Int_t run)
423 {
424 // sets run number
425 // Activate it later within the Run() method
426
427         fRun = run;
428 }
429
430 //_____________________________________________________________________________
431 void AliSimulation::SetSeed(Int_t seed)
432 {
433 // sets seed number
434 // Activate it later within the Run() method
435
436         fSeed = seed;
437 }
438
439 //_____________________________________________________________________________
440 Bool_t AliSimulation::SetRunNumberFromData()
441 {
442   // Set the CDB manager run number
443   // The run number is retrieved from gAlice
444
445     if (fSetRunNumberFromDataCalled) return kTRUE;
446     fSetRunNumberFromDataCalled = kTRUE;    
447   
448     AliCDBManager* man = AliCDBManager::Instance();
449     Int_t runData = -1, runCDB = -1;
450   
451     AliRunLoader* runLoader = LoadRun("READ");
452     if (!runLoader) return kFALSE;
453     else {
454         runData = runLoader->GetHeader()->GetRun();
455         delete runLoader;
456     }
457   
458     runCDB = man->GetRun();
459     if(runCDB >= 0) {
460         if (runCDB != runData) {
461                 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
462                 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
463                 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
464                 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");    
465         }
466         
467     }
468       
469     man->SetRun(runData);
470     fRun = runData;
471     
472     if(man->GetRun() < 0) {
473         AliError("Run number not properly initalized!");
474         return kFALSE;
475     }
476   
477     man->Print();
478     
479     return kTRUE;
480 }
481
482 //_____________________________________________________________________________
483 void AliSimulation::SetConfigFile(const char* fileName)
484 {
485 // set the name of the config file
486
487   fConfigFileName = fileName;
488 }
489
490 //_____________________________________________________________________________
491 void AliSimulation::SetGAliceFile(const char* fileName)
492 {
493 // set the name of the galice file
494 // the path is converted to an absolute one if it is relative
495
496   fGAliceFileName = fileName;
497   if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
498     char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
499                                                 fGAliceFileName);
500     fGAliceFileName = absFileName;
501     delete[] absFileName;
502   }
503
504   AliDebug(2, Form("galice file name set to %s", fileName));
505 }
506
507 //_____________________________________________________________________________
508 void AliSimulation::SetEventsPerFile(const char* detector, const char* type, 
509                                      Int_t nEvents)
510 {
511 // set the number of events per file for the given detector and data type
512 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
513
514   TNamed* obj = new TNamed(detector, type);
515   obj->SetUniqueID(nEvents);
516   fEventsPerFile.Add(obj);
517 }
518
519 //_____________________________________________________________________________
520 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
521 {
522   // Read the alignment objects from CDB.
523   // Each detector is supposed to have the
524   // alignment objects in DET/Align/Data CDB path.
525   // All the detector objects are then collected,
526   // sorted by geometry level (starting from ALIC) and
527   // then applied to the TGeo geometry.
528   // Finally an overlaps check is performed.
529
530   if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
531     AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
532     return kFALSE;
533   }  
534   
535   // initialize CDB storage, run number, set CDB lock
536   InitCDB();
537 //  if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
538   SetCDBLock();
539     
540   Bool_t delRunLoader = kFALSE;
541   if (!runLoader) {
542     runLoader = LoadRun("READ");
543     if (!runLoader) return kFALSE;
544     delRunLoader = kTRUE;
545   }
546   
547   // Export ideal geometry 
548   if(!IsGeometryFromFile()) AliGeomManager::GetGeometry()->Export("geometry.root");
549
550   // Load alignment data from CDB and apply to geometry through AliGeomManager
551   if(fLoadAlignFromCDB){
552     
553     TString detStr = fLoadAlObjsListOfDets;
554     TString loadAlObjsListOfDets = "";
555     
556     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
557     for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
558       AliModule* det = (AliModule*) detArray->At(iDet);
559       if (!det || !det->IsActive()) continue;
560       if (IsSelected(det->GetName(), detStr)) {
561         //add det to list of dets to be aligned from CDB
562         loadAlObjsListOfDets += det->GetName();
563         loadAlObjsListOfDets += " ";
564       }
565     } // end loop over detectors
566     loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
567     AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
568   }else{
569     // Check if the array with alignment objects was
570     // provided by the user. If yes, apply the objects
571     // to the present TGeo geometry
572     if (fAlignObjArray) {
573       if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
574         AliError("The misalignment of one or more volumes failed!"
575                  "Compare the list of simulated detectors and the list of detector alignment data!");
576         if (delRunLoader) delete runLoader;
577         return kFALSE;
578       }
579     }
580   }
581
582   // Update the internal geometry of modules (ITS needs it)
583   TString detStr = fLoadAlObjsListOfDets;
584   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
585   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
586
587     AliModule* det = (AliModule*) detArray->At(iDet);
588     if (!det || !det->IsActive()) continue;
589     if (IsSelected(det->GetName(), detStr)) {
590       det->UpdateInternalGeometry();
591     }
592   } // end loop over detectors
593
594
595   if (delRunLoader) delete runLoader;
596
597   return kTRUE;
598 }
599
600 //_____________________________________________________________________________
601 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
602 {
603 // add a file with background events for merging
604
605   TObjString* fileNameStr = new TObjString(fileName);
606   fileNameStr->SetUniqueID(nSignalPerBkgrd);
607   if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
608   fBkgrdFileNames->Add(fileNameStr);
609 }
610
611 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
612 {
613 // add a file with background events for embeddin
614   MergeWith(fileName, nSignalPerBkgrd);
615   fEmbeddingFlag = kTRUE;
616 }
617
618 //_____________________________________________________________________________
619 Bool_t AliSimulation::Run(Int_t nEvents)
620 {
621 // run the generation, simulation and digitization
622
623  
624   AliCodeTimerAuto("",0)
625   AliSysInfo::AddStamp("Start_Run");
626   
627   // Load run number and seed from environmental vars
628   ProcessEnvironmentVars();
629   AliSysInfo::AddStamp("ProcessEnvironmentVars");
630
631   gRandom->SetSeed(fSeed);
632    
633   if (nEvents > 0) fNEvents = nEvents;
634
635   // create and setup the HLT instance
636   if (!fRunHLT.IsNull() && !CreateHLT()) {
637     if (fStopOnError) return kFALSE;
638     // disable HLT
639     fRunHLT="";
640   }
641   
642   // generation and simulation -> hits
643   if (fRunGeneration) {
644     if (!RunSimulation()) if (fStopOnError) return kFALSE;
645   }
646   AliSysInfo::AddStamp("RunSimulation");
647            
648   // initialize CDB storage from external environment
649   // (either CDB manager or AliSimulation setters),
650   // if not already done in RunSimulation()
651   InitCDB();
652   AliSysInfo::AddStamp("InitCDB");
653   
654   // Set run number in CDBManager from data 
655   // From this point on the run number must be always loaded from data!
656   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
657   
658   // Set CDB lock: from now on it is forbidden to reset the run number
659   // or the default storage or to activate any further storage!
660   SetCDBLock();
661
662   // If RunSimulation was not called, load the geometry and misalign it
663   if (!AliGeomManager::GetGeometry()) {
664     // Initialize the geometry manager
665     AliGeomManager::LoadGeometry("geometry.root");
666     AliSysInfo::AddStamp("GetGeometry");
667 //    // Check that the consistency of symbolic names for the activated subdetectors
668 //    // in the geometry loaded by AliGeomManager
669 //    AliRunLoader* runLoader = LoadRun("READ");
670 //    if (!runLoader) return kFALSE;
671 //
672 //    TString detsToBeChecked = "";
673 //    TObjArray* detArray = runLoader->GetAliRun()->Detectors();
674 //    for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
675 //      AliModule* det = (AliModule*) detArray->At(iDet);
676 //      if (!det || !det->IsActive()) continue;
677 //      detsToBeChecked += det->GetName();
678 //      detsToBeChecked += " ";
679 //    } // end loop over detectors
680 //    if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
681     if(!AliGeomManager::CheckSymNamesLUT("ALL"))
682         AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
683         
684     if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
685     // Misalign geometry
686     if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
687   }
688   AliSysInfo::AddStamp("MissalignGeometry");
689
690
691   // hits -> summable digits
692   AliSysInfo::AddStamp("Start_sdigitization");
693   if (!fMakeSDigits.IsNull()) {
694     if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
695  
696   }
697   AliSysInfo::AddStamp("Stop_sdigitization");
698   
699   AliSysInfo::AddStamp("Start_digitization");  
700   // summable digits -> digits  
701   if (!fMakeDigits.IsNull()) {
702     if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
703       if (fStopOnError) return kFALSE;
704     }
705    }
706   AliSysInfo::AddStamp("Stop_digitization");
707
708   
709   
710   // hits -> digits
711   if (!fMakeDigitsFromHits.IsNull()) {
712     if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
713       AliWarning(Form("Merging and direct creation of digits from hits " 
714                  "was selected for some detectors. "
715                  "No merging will be done for the following detectors: %s",
716                  fMakeDigitsFromHits.Data()));
717     }
718     if (!RunHitsDigitization(fMakeDigitsFromHits)) {
719       if (fStopOnError) return kFALSE;
720     }
721   }
722
723   AliSysInfo::AddStamp("Hits2Digits");
724   
725   
726   // digits -> trigger
727   if (!fTriggerConfig.IsNull() && !RunTrigger(fTriggerConfig,fMakeDigits)) {
728     if (fStopOnError) return kFALSE;
729   }
730
731   AliSysInfo::AddStamp("RunTrigger");
732   
733   
734   // digits -> raw data
735   if (!fWriteRawData.IsNull()) {
736     if (!WriteRawData(fWriteRawData, fRawDataFileName, 
737                       fDeleteIntermediateFiles,fWriteSelRawData)) {
738       if (fStopOnError) return kFALSE;
739     }
740   }
741
742   AliSysInfo::AddStamp("WriteRaw");
743   
744   // run HLT simulation on simulated digit data if raw data is not
745   // simulated, otherwise its called as part of WriteRawData
746   if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
747     if (!RunHLT()) {
748       if (fStopOnError) return kFALSE;
749     }
750   }
751
752   AliSysInfo::AddStamp("RunHLT");
753   
754   //QA
755   if (fRunQA) {
756       Bool_t rv = RunQA() ; 
757       if (!rv)
758           if (fStopOnError) 
759               return kFALSE ;
760   }
761
762   AliSysInfo::AddStamp("RunQA");
763
764   TString snapshotFileOut("");
765   if(TString(gSystem->Getenv("OCDB_SNAPSHOT_CREATE")) == TString("kTRUE")){ 
766       AliInfo(" ******** Creating the snapshot! *********");
767       TString snapshotFile(gSystem->Getenv("OCDB_SNAPSHOT_FILENAME")); 
768       if(!(snapshotFile.IsNull() || snapshotFile.IsWhitespace())) 
769           snapshotFileOut = snapshotFile;
770       else 
771           snapshotFileOut="OCDB.root"; 
772       AliCDBManager::Instance()->DumpToSnapshotFile(snapshotFileOut.Data(),kFALSE); 
773   }
774
775   // Cleanup of CDB manager: cache and active storages!
776   AliCDBManager::Instance()->ClearCache();
777
778   return kTRUE;
779 }
780
781 //_______________________________________________________________________
782 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
783                      Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
784                      Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
785 {
786   //
787   // Generates lego plots of:
788   //    - radiation length map phi vs theta
789   //    - radiation length map phi vs eta
790   //    - interaction length map
791   //    - g/cm2 length map
792   //
793   //  ntheta    bins in theta, eta
794   //  themin    minimum angle in theta (degrees)
795   //  themax    maximum angle in theta (degrees)
796   //  nphi      bins in phi
797   //  phimin    minimum angle in phi (degrees)
798   //  phimax    maximum angle in phi (degrees)
799   //  rmin      minimum radius
800   //  rmax      maximum radius
801   //  
802   //
803   //  The number of events generated = ntheta*nphi
804   //  run input parameters in macro setup (default="Config.C")
805   //
806   //  Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
807   //Begin_Html
808   /*
809     <img src="picts/AliRunLego1.gif">
810   */
811   //End_Html
812   //Begin_Html
813   /*
814     <img src="picts/AliRunLego2.gif">
815   */
816   //End_Html
817   //Begin_Html
818   /*
819     <img src="picts/AliRunLego3.gif">
820   */
821   //End_Html
822   //
823
824 // run the generation and simulation
825
826   AliCodeTimerAuto("",0)
827
828   // initialize CDB storage and run number from external environment
829   // (either CDB manager or AliSimulation setters)
830   InitCDB();
831   InitRunNumber();
832   SetCDBLock();
833   
834   if (!gAlice) {
835     AliError("no gAlice object. Restart aliroot and try again.");
836     return kFALSE;
837   }
838   if (gAlice->Modules()->GetEntries() > 0) {
839     AliError("gAlice was already run. Restart aliroot and try again.");
840     return kFALSE;
841   }
842
843   AliInfo(Form("initializing gAlice with config file %s",
844           fConfigFileName.Data()));
845
846   // Number of events 
847     if (nev == -1) nev  = nc1 * nc2;
848     
849   // check if initialisation has been done
850   // If runloader has been initialized, set the number of events per file to nc1 * nc2
851     
852   // Set new generator
853   if (!gener) gener  = new AliLegoGenerator();
854   //
855   // Configure Generator
856
857   gener->SetRadiusRange(rmin, rmax);
858   gener->SetZMax(zmax);
859   gener->SetCoor1Range(nc1, c1min, c1max);
860   gener->SetCoor2Range(nc2, c2min, c2max);
861   
862   
863   //Create Lego object  
864   fLego = new AliLego("lego",gener);
865
866   //__________________________________________________________________________
867
868   gAlice->Announce();
869
870   gROOT->LoadMacro(setup);
871   gInterpreter->ProcessLine(gAlice->GetConfigFunction());
872
873   if(AliCDBManager::Instance()->GetRun() >= 0) { 
874     SetRunNumber(AliCDBManager::Instance()->GetRun());
875   } else {
876     AliWarning("Run number not initialized!!");
877   }
878
879   AliRunLoader::Instance()->CdGAFile();
880   
881   AliPDG::AddParticlesToPdgDataBase();  
882   
883   gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
884
885   gAlice->GetMCApp()->Init();
886   
887   
888   //Must be here because some MCs (G4) adds detectors here and not in Config.C
889   gAlice->InitLoaders();
890   AliRunLoader::Instance()->MakeTree("E");
891   
892   //
893   // Save stuff at the beginning of the file to avoid file corruption
894   AliRunLoader::Instance()->CdGAFile();
895   gAlice->Write();
896
897   //Save current generator
898   AliGenerator *gen=gAlice->GetMCApp()->Generator();
899   gAlice->GetMCApp()->ResetGenerator(gener);
900   //Prepare MC for Lego Run
901   gMC->InitLego();
902   
903   //Run Lego Object
904   
905   
906   AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
907   gMC->ProcessRun(nev);
908   
909   // End of this run, close files
910   FinishRun();
911   // Restore current generator
912   gAlice->GetMCApp()->ResetGenerator(gen);
913   // Delete Lego Object
914   delete fLego;
915
916   return kTRUE;
917 }
918
919 //_____________________________________________________________________________
920 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
921 {
922   // run the trigger
923
924   AliCodeTimerAuto("",0)
925
926   // initialize CDB storage from external environment
927   // (either CDB manager or AliSimulation setters),
928   // if not already done in RunSimulation()
929   InitCDB();
930   
931   // Set run number in CDBManager from data 
932   // From this point on the run number must be always loaded from data!
933   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
934   
935   // Set CDB lock: from now on it is forbidden to reset the run number
936   // or the default storage or to activate any further storage!
937   SetCDBLock();
938    
939    AliRunLoader* runLoader = LoadRun("READ");
940    if (!runLoader) return kFALSE;
941    TString trconfiguration = config;
942
943    if (trconfiguration.IsNull()) {
944      if(!fTriggerConfig.IsNull()) {
945        trconfiguration = fTriggerConfig;
946      }
947      else
948        AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
949    }
950
951    runLoader->MakeTree( "GG" );
952    AliCentralTrigger* aCTP = runLoader->GetTrigger();
953    // Load Configuration
954    if (!aCTP->LoadConfiguration( trconfiguration ))
955      return kFALSE;
956
957    // digits -> trigger
958    if( !aCTP->RunTrigger( runLoader , detectors ) ) {
959       if (fStopOnError) {
960         //  delete aCTP;
961         return kFALSE;
962       }
963    }
964
965    delete runLoader;
966
967    return kTRUE;
968 }
969
970 //_____________________________________________________________________________
971 Bool_t AliSimulation::WriteTriggerRawData()
972 {
973   // Writes the CTP (trigger) DDL raw data
974   // Details of the format are given in the
975   // trigger TDR - pages 134 and 135.
976   AliCTPRawData writer;
977   writer.RawData();
978
979   return kTRUE;
980 }
981
982 //_____________________________________________________________________________
983 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
984 {
985 // run the generation and simulation
986
987   AliCodeTimerAuto("",0)
988
989   // initialize CDB storage and run number from external environment
990   // (either CDB manager or AliSimulation setters)
991   AliSysInfo::AddStamp("RunSimulation_Begin");
992   InitCDB();
993   AliSysInfo::AddStamp("RunSimulation_InitCDB");
994   InitRunNumber();
995
996   SetCDBLock();
997   AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
998   
999   if (!gAlice) {
1000     AliError("no gAlice object. Restart aliroot and try again.");
1001     return kFALSE;
1002   }
1003   if (gAlice->Modules()->GetEntries() > 0) {
1004     AliError("gAlice was already run. Restart aliroot and try again.");
1005     return kFALSE;
1006   }
1007   
1008   // Setup monitoring if requested
1009   gAlice->GetMCApp()->SetUseMonitoring(fUseMonitoring);
1010
1011   AliInfo(Form("initializing gAlice with config file %s",
1012           fConfigFileName.Data()));
1013
1014   //
1015   // Initialize ALICE Simulation run
1016   //
1017   gAlice->Announce();
1018
1019   //
1020   // If requested set the mag. field from the GRP entry.
1021   // After this the field is loccked and cannot be changed by Config.C
1022   if (fUseMagFieldFromGRP) {
1023       AliGRPManager grpM;
1024       grpM.ReadGRPEntry();
1025       grpM.SetMagField();
1026       AliInfo("Field is locked now. It cannot be changed in Config.C");
1027   }
1028 //
1029 // Execute Config.C
1030   TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1031   gROOT->LoadMacro(fConfigFileName.Data());
1032   Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1033   if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1034     AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1035   }
1036   AliSysInfo::AddStamp("RunSimulation_Config");
1037
1038 //
1039 // If requested obtain the vertex position and vertex sigma_z from the CDB
1040 // This overwrites the settings from the Config.C  
1041   if (fUseVertexFromCDB) {
1042       Double_t vtxPos[3] = {0., 0., 0.}; 
1043       Double_t vtxSig[3] = {0., 0., 0.};
1044       AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1045       if (entry) {
1046           AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1047           if (vertex) {
1048               if(vertex->GetXRes()>2.8) { // > pipe radius --> it's a dummy object, don't use it 
1049                   entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1050                   if (entry) vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1051               }
1052           }
1053           if (vertex) {
1054               vertex->GetXYZ(vtxPos);
1055               vertex->GetSigmaXYZ(vtxSig);
1056               AliInfo("Overwriting Config.C vertex settings !");
1057               AliInfo(Form("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
1058                            vtxPos[0], vtxPos[1], vtxPos[2], vtxSig[2]));
1059               
1060               AliGenerator *gen = gAlice->GetMCApp()->Generator();
1061               gen->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]);   // vertex position
1062               gen->SetSigmaZ(vtxSig[2]);
1063           }
1064       }
1065   }
1066
1067   // If requested we take the SOR and EOR time-stamps from the GRP and use them
1068   // in order to generate the event time-stamps
1069   if (fUseTimeStampFromCDB) {
1070     AliGRPManager grpM;
1071     grpM.ReadGRPEntry();
1072     const AliGRPObject *grpObj = grpM.GetGRPData();
1073     if (!grpObj || (grpObj->GetTimeEnd() <= grpObj->GetTimeStart())) {
1074       AliError("Missing GRP or bad SOR/EOR time-stamps! Switching off the time-stamp generation from GRP!");
1075       fTimeStart = fTimeEnd = 0;
1076       fUseTimeStampFromCDB = kFALSE;
1077     }
1078     else {
1079       fTimeStart = grpObj->GetTimeStart();
1080       fTimeEnd = grpObj->GetTimeEnd();
1081     }
1082   }
1083   
1084   if(AliCDBManager::Instance()->GetRun() >= 0) { 
1085     AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
1086     AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
1087   } else {
1088         AliWarning("Run number not initialized!!");
1089   }
1090   
1091    AliRunLoader::Instance()->CdGAFile();
1092    
1093
1094    AliPDG::AddParticlesToPdgDataBase();  
1095
1096    gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1097    AliSysInfo::AddStamp("RunSimulation_GetField");
1098    gAlice->GetMCApp()->Init();
1099    AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1100
1101    //Must be here because some MCs (G4) adds detectors here and not in Config.C
1102    gAlice->InitLoaders();
1103    AliRunLoader::Instance()->MakeTree("E");
1104    AliRunLoader::Instance()->LoadKinematics("RECREATE");
1105    AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1106    AliRunLoader::Instance()->LoadHits("all","RECREATE");
1107    //
1108    // Save stuff at the beginning of the file to avoid file corruption
1109    AliRunLoader::Instance()->CdGAFile();
1110    gAlice->Write();
1111    gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1112    AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1113   //___________________________________________________________________________________________
1114   
1115   AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1116
1117   // Set run number in CDBManager
1118   AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1119
1120   AliRunLoader* runLoader = AliRunLoader::Instance();
1121   if (!runLoader) {
1122              AliError(Form("gAlice has no run loader object. "
1123                              "Check your config file: %s", fConfigFileName.Data()));
1124              return kFALSE;
1125   }
1126   SetGAliceFile(runLoader->GetFileName());
1127       
1128   // Misalign geometry
1129 #if ROOT_VERSION_CODE < 331527
1130   AliGeomManager::SetGeometry(gGeoManager);
1131   
1132   // Check that the consistency of symbolic names for the activated subdetectors
1133   // in the geometry loaded by AliGeomManager
1134   TString detsToBeChecked = "";
1135   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1136   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1137     AliModule* det = (AliModule*) detArray->At(iDet);
1138     if (!det || !det->IsActive()) continue;
1139     detsToBeChecked += det->GetName();
1140     detsToBeChecked += " ";
1141   } // end loop over detectors
1142   if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1143     AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1144   MisalignGeometry(runLoader);
1145   AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1146 #endif
1147
1148 //   AliRunLoader* runLoader = AliRunLoader::Instance();
1149 //   if (!runLoader) {
1150 //     AliError(Form("gAlice has no run loader object. "
1151 //                   "Check your config file: %s", fConfigFileName.Data()));
1152 //     return kFALSE;
1153 //   }
1154 //   SetGAliceFile(runLoader->GetFileName());
1155
1156   if (!gAlice->GetMCApp()->Generator()) {
1157     AliError(Form("gAlice has no generator object. "
1158                   "Check your config file: %s", fConfigFileName.Data()));
1159     return kFALSE;
1160   }
1161
1162   // Write GRP entry corresponding to the setting found in Cofig.C
1163   if (fWriteGRPEntry)
1164     WriteGRPEntry();
1165   AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1166
1167   if (nEvents <= 0) nEvents = fNEvents;
1168
1169   // get vertex from background file in case of merging
1170   if (fUseBkgrdVertex &&
1171       fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1172     Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1173     const char* fileName = ((TObjString*)
1174                             (fBkgrdFileNames->At(0)))->GetName();
1175     AliInfo(Form("The vertex will be taken from the background "
1176                  "file %s with nSignalPerBackground = %d", 
1177                  fileName, signalPerBkgrd));
1178     AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1179     gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1180   }
1181
1182   if (!fRunSimulation) {
1183     gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1184   }
1185
1186   // set the number of events per file for given detectors and data types
1187   for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1188     if (!fEventsPerFile[i]) continue;
1189     const char* detName = fEventsPerFile[i]->GetName();
1190     const char* typeName = fEventsPerFile[i]->GetTitle();
1191     TString loaderName(detName);
1192     loaderName += "Loader";
1193     AliLoader* loader = runLoader->GetLoader(loaderName);
1194     if (!loader) {
1195       AliError(Form("RunSimulation no loader for %s found\n Number of events per file not set for %s %s", 
1196                     detName, typeName, detName));
1197       continue;
1198     }
1199     AliDataLoader* dataLoader = 
1200       loader->GetDataLoader(typeName);
1201     if (!dataLoader) {
1202       AliError(Form("no data loader for %s found\n"
1203                     "Number of events per file not set for %s %s", 
1204                     typeName, detName, typeName));
1205       continue;
1206     }
1207     dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1208     AliDebug(1, Form("number of events per file set to %d for %s %s",
1209                      fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1210   }
1211
1212   AliInfo("running gAlice");
1213   AliSysInfo::AddStamp("Start_ProcessRun");
1214
1215   // Create the Root Tree with one branch per detector
1216   //Hits moved to begin event -> now we are crating separate tree for each event
1217   gMC->ProcessRun(nEvents);
1218
1219   // End of this run, close files
1220   if(nEvents>0) FinishRun();
1221
1222   AliSysInfo::AddStamp("Stop_ProcessRun");
1223   delete runLoader;
1224
1225   return kTRUE;
1226 }
1227
1228 //_____________________________________________________________________________
1229 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1230 {
1231 // run the digitization and produce summable digits
1232   static Int_t eventNr=0;
1233   AliCodeTimerAuto("",0) ;
1234
1235   // initialize CDB storage, run number, set CDB lock
1236   InitCDB();
1237   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1238   SetCDBLock();
1239   
1240   AliRunLoader* runLoader = LoadRun();
1241   if (!runLoader) return kFALSE;
1242
1243   TString detStr = detectors;
1244   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1245   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1246     AliModule* det = (AliModule*) detArray->At(iDet);
1247     if (!det || !det->IsActive()) continue;
1248     if (IsSelected(det->GetName(), detStr)) {
1249       AliInfo(Form("creating summable digits for %s", det->GetName()));
1250       AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1251       det->Hits2SDigits();
1252       AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1253       AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1254     }
1255   }
1256
1257   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1258     AliError(Form("the following detectors were not found: %s",
1259                   detStr.Data()));
1260     if (fStopOnError) return kFALSE;
1261   }
1262   eventNr++;
1263   delete runLoader;
1264
1265   return kTRUE;
1266 }
1267
1268
1269 //_____________________________________________________________________________
1270 Bool_t AliSimulation::RunDigitization(const char* detectors, 
1271                                       const char* excludeDetectors)
1272 {
1273 // run the digitization and produce digits from sdigits
1274
1275   AliCodeTimerAuto("",0)
1276
1277   // initialize CDB storage, run number, set CDB lock
1278   InitCDB();
1279   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1280   SetCDBLock();
1281   
1282   delete AliRunLoader::Instance();
1283   delete gAlice;
1284   gAlice = NULL;
1285
1286   Int_t nStreams = 1;
1287   if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1288   Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1289   AliDigitizationInput digInp(nStreams, signalPerBkgrd);
1290   // digInp.SetEmbeddingFlag(fEmbeddingFlag);
1291   digInp.SetRegionOfInterest(fRegionOfInterest);
1292   digInp.SetInputStream(0, fGAliceFileName.Data());
1293   for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1294     const char* fileName = ((TObjString*)(fBkgrdFileNames->At(iStream-1)))->GetName();
1295     digInp.SetInputStream(iStream, fileName);
1296   }
1297   TObjArray detArr;
1298   detArr.SetOwner(kTRUE);
1299   TString detStr = detectors;
1300   TString detExcl = excludeDetectors;
1301   if (!static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice()) {
1302     AliError("Error occured while getting gAlice from Input 0");
1303     return kFALSE;
1304   }
1305   AliRunLoader* runLoader = AliRunLoader::GetRunLoader(digInp.GetInputStream(0)->GetFolderName());
1306   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1307   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1308     AliModule* det = (AliModule*) detArray->At(iDet);
1309     if (!det || !det->IsActive()) continue;
1310     if (!IsSelected(det->GetName(), detStr) || IsSelected(det->GetName(), detExcl)) continue;
1311     AliDigitizer* digitizer = det->CreateDigitizer(&digInp);
1312     if (!digitizer || !digitizer->Init()) {
1313       AliError(Form("no digitizer for %s", det->GetName()));
1314       if (fStopOnError) return kFALSE;
1315       else continue;
1316     }
1317     detArr.AddLast(digitizer);    
1318     AliInfo(Form("Created digitizer from SDigits -> Digits for %s", det->GetName()));    
1319   }
1320   //
1321   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1322     AliError(Form("the following detectors were not found: %s", detStr.Data()));
1323     if (fStopOnError) return kFALSE;
1324   }
1325   //
1326   Int_t ndigs = detArr.GetEntriesFast();
1327   Int_t eventsCreated = 0;
1328   AliRunLoader* outRl =  digInp.GetOutRunLoader();
1329   while ((eventsCreated++ < fNEvents) || (fNEvents < 0)) {
1330     if (!digInp.ConnectInputTrees()) break;
1331     digInp.InitEvent(); //this must be after call of Connect Input tress.
1332     if (outRl) outRl->SetEventNumber(eventsCreated-1);
1333     static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice(); // use gAlice of the first input stream
1334     for (int id=0;id<ndigs;id++) ((AliDigitizer*)detArr[id])->Digitize("");
1335     digInp.FinishEvent();
1336   };
1337   digInp.FinishGlobal();
1338   //
1339   return kTRUE;
1340 }
1341
1342 //_____________________________________________________________________________
1343 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1344 {
1345 // run the digitization and produce digits from hits
1346
1347   AliCodeTimerAuto("",0)
1348
1349   // initialize CDB storage, run number, set CDB lock
1350   InitCDB();
1351   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1352   SetCDBLock();
1353   
1354   AliRunLoader* runLoader = LoadRun("READ");
1355   if (!runLoader) return kFALSE;
1356
1357   TString detStr = detectors;
1358   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1359   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1360     AliModule* det = (AliModule*) detArray->At(iDet);
1361     if (!det || !det->IsActive()) continue;
1362     if (IsSelected(det->GetName(), detStr)) {
1363       AliInfo(Form("creating digits from hits for %s", det->GetName()));
1364       det->Hits2Digits();
1365     }
1366   }
1367
1368   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1369     AliError(Form("the following detectors were not found: %s", 
1370                   detStr.Data()));
1371     if (fStopOnError) return kFALSE;
1372   }
1373
1374   return kTRUE;
1375 }
1376
1377 //_____________________________________________________________________________
1378 Bool_t AliSimulation::WriteRawData(const char* detectors, 
1379                                    const char* fileName,
1380                                    Bool_t deleteIntermediateFiles,
1381                                    Bool_t selrawdata)
1382 {
1383 // convert the digits to raw data
1384 // First DDL raw data files for the given detectors are created.
1385 // If a file name is given, the DDL files are then converted to a DATE file.
1386 // If deleteIntermediateFiles is true, the DDL raw files are deleted 
1387 // afterwards.
1388 // If the file name has the extension ".root", the DATE file is converted
1389 // to a root file.
1390 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1391 // 'selrawdata' flag can be used to enable writing of detectors raw data
1392 // accoring to the trigger cluster.
1393
1394   AliCodeTimerAuto("",0)
1395   AliSysInfo::AddStamp("WriteRawData_Start");
1396   
1397   TString detStr = detectors;
1398   if (!WriteRawFiles(detStr.Data())) {
1399     if (fStopOnError) return kFALSE;
1400   }
1401   AliSysInfo::AddStamp("WriteRawFiles");
1402
1403   // run HLT simulation on simulated DDL raw files
1404   // and produce HLT ddl raw files to be included in date/root file
1405   // bugfix 2009-06-26: the decision whether to write HLT raw data
1406   // is taken in RunHLT. Here HLT always needs to be run in order to
1407   // create HLT digits, unless its switched off. This is due to the
1408   // special placement of the HLT between the generation of DDL files
1409   // and conversion to DATE/Root file.
1410   detStr.ReplaceAll("HLT", "");
1411   if (!fRunHLT.IsNull()) {
1412     if (!RunHLT()) {
1413       if (fStopOnError) return kFALSE;
1414     }
1415   }
1416   AliSysInfo::AddStamp("WriteRawData_RunHLT");
1417
1418   TString dateFileName(fileName);
1419   if (!dateFileName.IsNull()) {
1420     Bool_t rootOutput = dateFileName.EndsWith(".root");
1421     if (rootOutput) dateFileName += ".date";
1422     TString selDateFileName;
1423     if (selrawdata) {
1424       selDateFileName = "selected.";
1425       selDateFileName+= dateFileName;
1426     }
1427     if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1428       if (fStopOnError) return kFALSE;
1429     }
1430     AliSysInfo::AddStamp("ConvertRawFilesToDate");
1431     if (deleteIntermediateFiles) {
1432       AliRunLoader* runLoader = LoadRun("READ");
1433       if (runLoader) for (Int_t iEvent = 0; 
1434                           iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1435         char command[256];
1436         snprintf(command, 256, "rm -r raw%d", iEvent);
1437         gSystem->Exec(command);
1438       }
1439       delete runLoader;
1440     }
1441
1442     if (rootOutput) {
1443       if (!ConvertDateToRoot(dateFileName, fileName)) {
1444         if (fStopOnError) return kFALSE;
1445       }
1446       AliSysInfo::AddStamp("ConvertDateToRoot");
1447       if (deleteIntermediateFiles) {
1448         gSystem->Unlink(dateFileName);
1449       }
1450       if (selrawdata) {
1451         TString selFileName = "selected.";
1452         selFileName        += fileName;
1453         if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1454           if (fStopOnError) return kFALSE;
1455         }
1456         if (deleteIntermediateFiles) {
1457           gSystem->Unlink(selDateFileName);
1458         }
1459       }
1460     }
1461   }
1462
1463   return kTRUE;
1464 }
1465
1466 //_____________________________________________________________________________
1467 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1468 {
1469 // convert the digits to raw data DDL files
1470
1471   AliCodeTimerAuto("",0)
1472   
1473   AliRunLoader* runLoader = LoadRun("READ");
1474   if (!runLoader) return kFALSE;
1475
1476   // write raw data to DDL files
1477   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1478     AliInfo(Form("processing event %d", iEvent));
1479     runLoader->GetEvent(iEvent);
1480     TString baseDir = gSystem->WorkingDirectory();
1481     char dirName[256];
1482     snprintf(dirName, 256, "raw%d", iEvent);
1483     gSystem->MakeDirectory(dirName);
1484     if (!gSystem->ChangeDirectory(dirName)) {
1485       AliError(Form("couldn't change to directory %s", dirName));
1486       if (fStopOnError) return kFALSE; else continue;
1487     }
1488
1489     ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1490     runNbFile.close();
1491
1492     TString detStr = detectors;
1493     if (IsSelected("HLT", detStr)) {
1494       // Do nothing. "HLT" will be removed from detStr and HLT raw
1495       // data files are generated in RunHLT.
1496     }
1497
1498     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1499     for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1500       AliModule* det = (AliModule*) detArray->At(iDet);
1501       if (!det || !det->IsActive()) continue;
1502       if (IsSelected(det->GetName(), detStr)) {
1503         AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1504         det->Digits2Raw();
1505       }
1506     }
1507
1508     if (!WriteTriggerRawData())
1509       if (fStopOnError) return kFALSE;
1510
1511     gSystem->ChangeDirectory(baseDir);
1512     if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1513       AliError(Form("the following detectors were not found: %s", 
1514                     detStr.Data()));
1515       if (fStopOnError) return kFALSE;
1516     }
1517   }
1518
1519   delete runLoader;
1520   
1521   return kTRUE;
1522 }
1523
1524 //_____________________________________________________________________________
1525 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1526                                             const char* selDateFileName)
1527 {
1528 // convert raw data DDL files to a DATE file with the program "dateStream"
1529 // The second argument is not empty when the user decides to write
1530 // the detectors raw data according to the trigger cluster.
1531
1532   AliCodeTimerAuto("",0)
1533   
1534   char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1535   if (!path) {
1536     AliError("the program dateStream was not found");
1537     if (fStopOnError) return kFALSE;
1538   } else {
1539     delete[] path;
1540   }
1541
1542   AliRunLoader* runLoader = LoadRun("READ");
1543   if (!runLoader) return kFALSE;
1544
1545   AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1546   Bool_t selrawdata = kFALSE;
1547   if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1548
1549   char command[256];
1550   // Note the option -s. It is used in order to avoid
1551   // the generation of SOR/EOR events.
1552   snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d", 
1553           dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1554   FILE* pipe = gSystem->OpenPipe(command, "w");
1555
1556   if (!pipe) {
1557     AliError(Form("Cannot execute command: %s",command));
1558     return kFALSE;
1559   }
1560
1561   Int_t selEvents = 0;
1562   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1563
1564     UInt_t detectorPattern = 0;
1565     runLoader->GetEvent(iEvent);
1566     if (!runLoader->LoadTrigger()) {
1567       AliCentralTrigger *aCTP = runLoader->GetTrigger();
1568       detectorPattern = aCTP->GetClusterMask();
1569       // Check if the event was triggered by CTP
1570       if (selrawdata) {
1571         if (aCTP->GetClassMask()) selEvents++;
1572       }
1573     }
1574     else {
1575       AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1576       if (selrawdata) {
1577         AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1578         selrawdata = kFALSE;
1579       }
1580     }
1581
1582     fprintf(pipe, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1583     Float_t ldc = 0;
1584     Int_t prevLDC = -1;
1585
1586     // loop over detectors and DDLs
1587     for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1588       for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1589
1590         Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1591         Int_t ldcID = Int_t(ldc + 0.0001);
1592         ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1593
1594         char rawFileName[256];
1595         snprintf(rawFileName, 256, "raw%d/%s", 
1596                 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1597
1598         // check existence and size of raw data file
1599         FILE* file = fopen(rawFileName, "rb");
1600         if (!file) continue;
1601         fseek(file, 0, SEEK_END);
1602         unsigned long size = ftell(file);
1603         fclose(file);
1604         if (!size) continue;
1605
1606         if (ldcID != prevLDC) {
1607           fprintf(pipe, " LDC Id %d\n", ldcID);
1608           prevLDC = ldcID;
1609         }
1610         fprintf(pipe, "  Equipment Id %d Payload %s\n", ddlID, rawFileName);
1611       }
1612     }
1613   }
1614
1615   Int_t result = gSystem->ClosePipe(pipe);
1616
1617   if (!(selrawdata && selEvents > 0)) {
1618     delete runLoader;
1619     return (result == 0);
1620   }
1621
1622   AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1623   
1624   snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d", 
1625           selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1626   FILE* pipe2 = gSystem->OpenPipe(command, "w");
1627
1628   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1629
1630     // Get the trigger decision and cluster
1631     UInt_t detectorPattern = 0;
1632     TString detClust;
1633     runLoader->GetEvent(iEvent);
1634     if (!runLoader->LoadTrigger()) {
1635       AliCentralTrigger *aCTP = runLoader->GetTrigger();
1636       if (aCTP->GetClassMask() == 0) continue;
1637       detectorPattern = aCTP->GetClusterMask();
1638       detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1639       AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1640     }
1641
1642     fprintf(pipe2, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1643     Float_t ldc = 0;
1644     Int_t prevLDC = -1;
1645
1646     // loop over detectors and DDLs
1647     for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1648       // Write only raw data from detectors that
1649       // are contained in the trigger cluster(s)
1650       if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1651
1652       for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1653
1654         Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1655         Int_t ldcID = Int_t(ldc + 0.0001);
1656         ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1657
1658         char rawFileName[256];
1659         snprintf(rawFileName, 256, "raw%d/%s", 
1660                 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1661
1662         // check existence and size of raw data file
1663         FILE* file = fopen(rawFileName, "rb");
1664         if (!file) continue;
1665         fseek(file, 0, SEEK_END);
1666         unsigned long size = ftell(file);
1667         fclose(file);
1668         if (!size) continue;
1669
1670         if (ldcID != prevLDC) {
1671           fprintf(pipe2, " LDC Id %d\n", ldcID);
1672           prevLDC = ldcID;
1673         }
1674         fprintf(pipe2, "  Equipment Id %d Payload %s\n", ddlID, rawFileName);
1675       }
1676     }
1677   }
1678
1679   Int_t result2 = gSystem->ClosePipe(pipe2);
1680
1681   delete runLoader;
1682   return ((result == 0) && (result2 == 0));
1683 }
1684
1685 //_____________________________________________________________________________
1686 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1687                                         const char* rootFileName)
1688 {
1689 // convert a DATE file to a root file with the program "alimdc"
1690
1691   // ALIMDC setup
1692   const Int_t kDBSize = 2000000000;
1693   const Int_t kTagDBSize = 1000000000;
1694   const Bool_t kFilter = kFALSE;
1695   const Int_t kCompression = 1;
1696
1697   char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1698   if (!path) {
1699     AliError("the program alimdc was not found");
1700     if (fStopOnError) return kFALSE;
1701   } else {
1702     delete[] path;
1703   }
1704
1705   AliInfo(Form("converting DATE file %s to root file %s", 
1706                dateFileName, rootFileName));
1707
1708   const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1709   const char* tagDBFS    = "/tmp/mdc1/tags";
1710
1711   // User defined file system locations
1712   if (gSystem->Getenv("ALIMDC_RAWDB1")) 
1713     rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1714   if (gSystem->Getenv("ALIMDC_RAWDB2")) 
1715     rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1716   if (gSystem->Getenv("ALIMDC_TAGDB")) 
1717     tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1718
1719   gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1720   gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1721   gSystem->Exec(Form("rm -rf %s",tagDBFS));
1722
1723   gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1724   gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1725   gSystem->Exec(Form("mkdir %s",tagDBFS));
1726
1727   Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s", 
1728                                     kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1729   gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1730
1731   gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1732   gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1733   gSystem->Exec(Form("rm -rf %s",tagDBFS));
1734
1735   return (result == 0);
1736 }
1737
1738
1739 //_____________________________________________________________________________
1740 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1741 {
1742 // delete existing run loaders, open a new one and load gAlice
1743
1744   delete AliRunLoader::Instance();
1745   AliRunLoader* runLoader = 
1746     AliRunLoader::Open(fGAliceFileName.Data(), 
1747                        AliConfig::GetDefaultEventFolderName(), mode);
1748   if (!runLoader) {
1749     AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1750     return NULL;
1751   }
1752   runLoader->LoadgAlice();
1753   runLoader->LoadHeader();
1754   gAlice = runLoader->GetAliRun();
1755   if (!gAlice) {
1756     AliError(Form("no gAlice object found in file %s", 
1757                   fGAliceFileName.Data()));
1758     return NULL;
1759   }
1760   return runLoader;
1761 }
1762
1763 //_____________________________________________________________________________
1764 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1765 {
1766 // get or calculate the number of signal events per background event
1767
1768   if (!fBkgrdFileNames) return 1;
1769   Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1770   if (nBkgrdFiles == 0) return 1;
1771
1772   // get the number of signal events
1773   if (nEvents <= 0) {
1774     AliRunLoader* runLoader = 
1775         AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1776     if (!runLoader) return 1;
1777     
1778     nEvents = runLoader->GetNumberOfEvents();
1779     delete runLoader;
1780   }
1781
1782   Int_t result = 0;
1783   for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1784     // get the number of background events
1785     const char* fileName = ((TObjString*)
1786                             (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1787     AliRunLoader* runLoader =
1788       AliRunLoader::Open(fileName, "BKGRD");
1789     if (!runLoader) continue;
1790     Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1791     delete runLoader;
1792   
1793     // get or calculate the number of signal per background events
1794     Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1795     if (nSignalPerBkgrd <= 0) {
1796       nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1797     } else if (result && (result != nSignalPerBkgrd)) {
1798       AliInfo(Form("the number of signal events per background event "
1799                    "will be changed from %d to %d for stream %d", 
1800                    nSignalPerBkgrd, result, iBkgrdFile+1));
1801       nSignalPerBkgrd = result;
1802     }
1803
1804     if (!result) result = nSignalPerBkgrd;
1805     if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1806       AliWarning(Form("not enough background events (%d) for %d signal events "
1807                       "using %d signal per background events for stream %d",
1808                       nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1809     }
1810   }
1811
1812   return result;
1813 }
1814
1815 //_____________________________________________________________________________
1816 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1817 {
1818 // check whether detName is contained in detectors
1819 // if yes, it is removed from detectors
1820
1821   // check if all detectors are selected
1822   if ((detectors.CompareTo("ALL") == 0) ||
1823       detectors.BeginsWith("ALL ") ||
1824       detectors.EndsWith(" ALL") ||
1825       detectors.Contains(" ALL ")) {
1826     detectors = "ALL";
1827     return kTRUE;
1828   }
1829
1830   // search for the given detector
1831   Bool_t result = kFALSE;
1832   if ((detectors.CompareTo(detName) == 0) ||
1833       detectors.BeginsWith(detName+" ") ||
1834       detectors.EndsWith(" "+detName) ||
1835       detectors.Contains(" "+detName+" ")) {
1836     detectors.ReplaceAll(detName, "");
1837     result = kTRUE;
1838   }
1839
1840   // clean up the detectors string
1841   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
1842   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1843   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1844
1845   return result;
1846 }
1847
1848 //_____________________________________________________________________________
1849 Int_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName, Int_t N, Int_t nSkip)
1850 {
1851 //
1852 // Steering routine  to convert raw data in directory rawDirectory/ to fake SDigits. 
1853 // These can be used for embedding of MC tracks into RAW data using the standard 
1854 // merging procedure.
1855 //
1856 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1857 //
1858   if (!gAlice) {
1859     AliError("no gAlice object. Restart aliroot and try again.");
1860     return kFALSE;
1861   }
1862   if (gAlice->Modules()->GetEntries() > 0) {
1863     AliError("gAlice was already run. Restart aliroot and try again.");
1864     return kFALSE;
1865   }
1866   
1867   AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1868   
1869   gAlice->Announce();
1870   
1871   gROOT->LoadMacro(fConfigFileName.Data());
1872   gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1873   
1874   if(AliCDBManager::Instance()->GetRun() >= 0) { 
1875         SetRunNumber(AliCDBManager::Instance()->GetRun());
1876   } else {
1877         AliWarning("Run number not initialized!!");
1878   }
1879   
1880    AliRunLoader::Instance()->CdGAFile();
1881     
1882    AliPDG::AddParticlesToPdgDataBase();  
1883
1884    gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1885    
1886    gAlice->GetMCApp()->Init();
1887
1888    //Must be here because some MCs (G4) adds detectors here and not in Config.C
1889    gAlice->InitLoaders();
1890    AliRunLoader::Instance()->MakeTree("E");
1891    AliRunLoader::Instance()->LoadKinematics("RECREATE");
1892    AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1893    AliRunLoader::Instance()->LoadHits("all","RECREATE");
1894
1895    //
1896    // Save stuff at the beginning of the file to avoid file corruption
1897    AliRunLoader::Instance()->CdGAFile();
1898    gAlice->Write();
1899 //
1900 //  Initialize CDB     
1901     InitCDB();
1902     //AliCDBManager* man = AliCDBManager::Instance();
1903     //man->SetRun(0); // Should this come from rawdata header ?
1904     
1905     Int_t iDet;
1906     //
1907     // Get the runloader
1908     AliRunLoader* runLoader = AliRunLoader::Instance();
1909     //
1910     // Open esd file if available
1911     TFile* esdFile = 0;
1912     TTree* treeESD = 0;
1913     AliESDEvent* esd = 0;
1914     if (esdFileName && (strlen(esdFileName)>0)) {
1915       esdFile = TFile::Open(esdFileName);
1916       if (esdFile) {
1917         esd = new AliESDEvent();
1918         esdFile->GetObject("esdTree", treeESD);
1919                   if (treeESD) {
1920                           esd->ReadFromTree(treeESD);
1921                           if (nSkip>0) {
1922                                   AliInfo(Form("Asking to skip first %d ESDs events",nSkip));
1923                           } else {
1924                                   nSkip=0;
1925                           }
1926                   }
1927       }
1928     }
1929
1930     //
1931     // Create the RawReader
1932     TString fileName(rawDirectory);
1933     AliRawReader* rawReader = AliRawReader::Create(fileName.Data());
1934     if (!rawReader) return (kFALSE);
1935     
1936 //     if (!fEquipIdMap.IsNull() && fRawReader)
1937 //       fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1938     //
1939     // Get list of detectors
1940     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1941     //
1942     // Get Header
1943     AliHeader* header = runLoader->GetHeader();
1944     // Event loop
1945     Int_t nev = 0;
1946     while(kTRUE) {
1947         if (!(rawReader->NextEvent())) break;
1948         runLoader->SetEventNumber(nev);
1949         runLoader->GetHeader()->Reset(rawReader->GetRunNumber(), 
1950                                       nev, nev);
1951         runLoader->GetEvent(nev);
1952         AliInfo(Form("We are at event %d",nev));
1953         //
1954         // Detector loop
1955         TString detStr = fMakeSDigits;
1956         for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1957             AliModule* det = (AliModule*) detArray->At(iDet);
1958             if (!det || !det->IsActive()) continue;
1959             if (IsSelected(det->GetName(), detStr)) {
1960               AliInfo(Form("Calling Raw2SDigits for %s", det->GetName()));
1961               det->Raw2SDigits(rawReader);
1962               rawReader->Reset();
1963             }
1964         } // detectors
1965         
1966
1967         //
1968         //  If ESD information available obtain reconstructed vertex and store in header.
1969         if (treeESD) {
1970                 AliInfo(Form("Selected event %d correspond to event %d in raw and to %d in esd",nev,rawReader->GetEventIndex(),nSkip+rawReader->GetEventIndex()));
1971             treeESD->GetEvent(nSkip+rawReader->GetEventIndex());
1972             const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1973             Double_t position[3];
1974             esdVertex->GetXYZ(position);
1975             AliGenEventHeader* mcHeader = new  AliGenEventHeader("ESD");
1976             TArrayF mcV;
1977             mcV.Set(3);
1978             for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1979             mcHeader->SetPrimaryVertex(mcV);
1980             header->Reset(0,nev);
1981             header->SetGenEventHeader(mcHeader);
1982             AliInfo(Form("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]));
1983         }
1984 //
1985 //      Finish the event
1986         runLoader->TreeE()->Fill();
1987         AliInfo(Form("Finished event %d",nev));
1988         nev++;
1989         if (N>0&&nev>=N)
1990           break;
1991     } // events
1992  
1993     delete rawReader;
1994 //
1995 //  Finish the run 
1996     runLoader->CdGAFile();
1997     runLoader->WriteHeader("OVERWRITE");
1998     runLoader->WriteRunLoader();
1999
2000     return nev;
2001 }
2002
2003 //_____________________________________________________________________________
2004 void AliSimulation::FinishRun()
2005 {
2006   //
2007   // Called at the end of the run.
2008   //
2009
2010   if(IsLegoRun()) 
2011    {
2012     AliDebug(1, "Finish Lego");
2013     AliRunLoader::Instance()->CdGAFile();
2014     fLego->FinishRun();
2015    }
2016   
2017   // Clean detector information
2018   TIter next(gAlice->Modules());
2019   AliModule *detector;
2020   while((detector = dynamic_cast<AliModule*>(next()))) {
2021     AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
2022     detector->FinishRun();
2023   }
2024   
2025   AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
2026   AliRunLoader::Instance()->WriteHeader("OVERWRITE");
2027
2028   // Write AliRun info and all detectors parameters
2029   AliRunLoader::Instance()->CdGAFile();
2030   gAlice->Write(0,TObject::kOverwrite);//write AliRun
2031   AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
2032   
2033   if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();  
2034   AliRunLoader::Instance()->Synchronize();
2035 }
2036
2037 //_____________________________________________________________________________
2038 Int_t AliSimulation::GetDetIndex(const char* detector)
2039 {
2040   // return the detector index corresponding to detector
2041   Int_t index = -1 ; 
2042   for (index = 0; index < fgkNDetectors ; index++) {
2043     if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2044           break ; 
2045   }     
2046   return index ; 
2047 }
2048
2049 //_____________________________________________________________________________
2050 Bool_t AliSimulation::CreateHLT()
2051 {
2052   // Init the HLT simulation.
2053   // The function  loads the library and creates the instance of AliHLTSimulation.
2054   // the main reason for the decoupled creation is to set the transient OCDB
2055   // objects before the OCDB is locked
2056
2057   // load the library dynamically
2058   gSystem->Load(ALIHLTSIMULATION_LIBRARY);
2059
2060   // check for the library version
2061   AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
2062   if (!fctVersion) {
2063     AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
2064     return kFALSE;
2065   }
2066   if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
2067     AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
2068   }
2069
2070   // print compile info
2071   typedef void (*CompileInfo)( const char*& date, const char*& time);
2072   CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
2073   if (fctInfo) {
2074     const char* date="";
2075     const char* time="";
2076     (*fctInfo)(date, time);
2077     if (!date) date="unknown";
2078     if (!time) time="unknown";
2079     AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
2080   } else {
2081     AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
2082   }
2083
2084   // create instance of the HLT simulation
2085   AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
2086   if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
2087     AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
2088     return kFALSE;    
2089   }
2090
2091   TString specObjects;
2092   for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
2093     if (specObjects.Length()>0) specObjects+=" ";
2094     specObjects+=fSpecCDBUri[i]->GetName();
2095   }
2096
2097   AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
2098   if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
2099     AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2100   }
2101
2102   return kTRUE;
2103 }
2104
2105 //_____________________________________________________________________________
2106 Bool_t AliSimulation::RunHLT()
2107 {
2108   // Run the HLT simulation
2109   // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2110   // Disabled if fRunHLT is empty, default vaule is "default".
2111   // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2112   // The default simulation depends on the HLT component libraries and their
2113   // corresponding agents which define components and chains to run. See
2114   // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2115   // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2116   //
2117   // The libraries to be loaded can be specified as an option.
2118   // <pre>
2119   // AliSimulation sim;
2120   // sim.SetRunHLT("libAliHLTSample.so");
2121   // </pre>
2122   // will only load <tt>libAliHLTSample.so</tt>
2123
2124   // Other available options:
2125   // \li loglevel=<i>level</i> <br>
2126   //     logging level for this processing
2127   // \li alilog=off
2128   //     disable redirection of log messages to AliLog class
2129   // \li config=<i>macro</i>
2130   //     configuration macro
2131   // \li chains=<i>configuration</i>
2132   //     comma separated list of configurations to be run during simulation
2133   // \li rawfile=<i>file</i>
2134   //     source for the RawReader to be created, the default is <i>./</i> if
2135   //     raw data is simulated
2136
2137   int iResult=0;
2138
2139   if (!fpHLT && !CreateHLT()) {
2140     return kFALSE;
2141   }
2142   AliHLTSimulation* pHLT=fpHLT;
2143
2144   AliRunLoader* pRunLoader = LoadRun("READ");
2145   if (!pRunLoader) return kFALSE;
2146
2147   // initialize CDB storage, run number, set CDB lock
2148   // thats for the case of running HLT simulation without all the other steps
2149   // multiple calls are handled by the function, so we can just call
2150   InitCDB();
2151   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2152   SetCDBLock();
2153   
2154   // init the HLT simulation
2155   TString options;
2156   if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2157   TString detStr = fWriteRawData;
2158   if (!IsSelected("HLT", detStr)) {
2159     options+=" writerawfiles=";
2160   } else {
2161     options+=" writerawfiles=HLT";
2162   }
2163
2164   if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2165     // as a matter of fact, HLT will run reconstruction and needs the RawReader
2166     // in order to get detector data. By default, RawReaderFile is used to read
2167     // the already simulated ddl files. Date and Root files from the raw data
2168     // are generated after the HLT simulation.
2169     options+=" rawfile=./";
2170   }
2171
2172   AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2173   if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2174     AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2175   } else {
2176     // run the HLT simulation
2177     AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2178     if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2179       AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2180     }
2181   }
2182
2183   // delete the instance
2184   AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2185   if (fctDelete==NULL || fctDelete(pHLT)<0) {
2186     AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2187   }
2188   pHLT=NULL;
2189
2190   return iResult>=0?kTRUE:kFALSE;
2191 }
2192
2193 //_____________________________________________________________________________
2194 Bool_t AliSimulation::RunQA()
2195 {
2196         // run the QA on summable hits, digits or digits
2197         
2198     //if(!gAlice) return kFALSE;
2199         AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2200
2201         TString detectorsw("") ;  
2202         Bool_t rv = kTRUE ; 
2203   AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2204         detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ; 
2205         if ( detectorsw.IsNull() ) 
2206                 rv = kFALSE ; 
2207         return rv ; 
2208 }
2209
2210 //_____________________________________________________________________________
2211 Bool_t AliSimulation::SetRunQA(TString detAndAction) 
2212 {
2213         // Allows to run QA for a selected set of detectors
2214         // and a selected set of tasks among HITS, SDIGITS and DIGITS
2215         // all selected detectors run the same selected tasks
2216         
2217         if (!detAndAction.Contains(":")) {
2218                 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2219                 fRunQA = kFALSE ;
2220                 return kFALSE ;                 
2221         }
2222         Int_t colon = detAndAction.Index(":") ; 
2223         fQADetectors = detAndAction(0, colon) ; 
2224         if (fQADetectors.Contains("ALL") ){
2225     TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ; 
2226     Int_t minus = fQADetectors.Last('-') ; 
2227     TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ; 
2228     TString toRemove("") ;
2229     while (minus >= 0) {
2230       toRemove = fQADetectors(minus+1, fQADetectors.Length()) ; 
2231       toRemove = toRemove.Strip() ; 
2232       toKeep.ReplaceAll(toRemove, "") ; 
2233       fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ; 
2234       minus = fQADetectors.Last('-') ; 
2235     }
2236     fQADetectors = toKeep ; 
2237   }
2238   fQATasks   = detAndAction(colon+1, detAndAction.Sizeof() ) ; 
2239         if (fQATasks.Contains("ALL") ) {
2240                 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ; 
2241         } else {
2242                 fQATasks.ToUpper() ; 
2243                 TString tempo("") ; 
2244                 if ( fQATasks.Contains("HIT") ) 
2245                         tempo = Form("%d ", AliQAv1::kHITS) ; 
2246                 if ( fQATasks.Contains("SDIGIT") ) 
2247                         tempo += Form("%d ", AliQAv1::kSDIGITS) ; 
2248                 if ( fQATasks.Contains("DIGIT") ) 
2249                         tempo += Form("%d ", AliQAv1::kDIGITS) ; 
2250                 fQATasks = tempo ; 
2251                 if (fQATasks.IsNull()) {
2252                         AliInfo("No QA requested\n")  ;
2253                         fRunQA = kFALSE ;
2254                         return kTRUE ; 
2255                 }
2256         }       
2257         TString tempo(fQATasks) ; 
2258     tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS))  ;
2259     tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;  
2260     tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;    
2261         AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;  
2262         fRunQA = kTRUE ;
2263         AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ; 
2264         AliQAManager::QAManager()->SetTasks(fQATasks) ; 
2265   for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++) 
2266     AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2267   
2268         return kTRUE; 
2269
2270
2271 //_____________________________________________________________________________
2272 void AliSimulation::ProcessEnvironmentVars()
2273 {
2274 // Extract run number and random generator seed from env variables
2275
2276     AliInfo("Processing environment variables");
2277     
2278     // Random Number seed
2279     
2280     // first check that seed is not already set
2281     if (fSeed == 0) {
2282         if (gSystem->Getenv("CONFIG_SEED")) {
2283                 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2284         }
2285     } else {
2286         if (gSystem->Getenv("CONFIG_SEED")) {
2287                 AliInfo(Form("Seed for random number generation already set (%d)"
2288                              ": CONFIG_SEED variable ignored!", fSeed));
2289         }
2290     }
2291    
2292     AliInfo(Form("Seed for random number generation = %d ", fSeed)); 
2293
2294     // Run Number
2295     
2296     // first check that run number is not already set
2297     if(fRun < 0) {    
2298         if (gSystem->Getenv("DC_RUN")) {
2299                 fRun = atoi(gSystem->Getenv("DC_RUN"));
2300         }
2301     } else {
2302         if (gSystem->Getenv("DC_RUN")) {
2303                 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2304         }
2305     }
2306     
2307     AliInfo(Form("Run number = %d", fRun)); 
2308 }
2309
2310 //---------------------------------------------------------------------
2311 void AliSimulation::WriteGRPEntry()
2312 {
2313   // Get the necessary information from galice (generator, trigger etc) and
2314   // write a GRP entry corresponding to the settings in the Config.C used
2315   // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2316
2317
2318   AliInfo("Writing global run parameters entry into the OCDB");
2319
2320   AliGRPObject* grpObj = new AliGRPObject();
2321
2322   grpObj->SetRunType("PHYSICS");
2323   grpObj->SetTimeStart(fTimeStart);
2324   grpObj->SetTimeEnd(fTimeEnd); 
2325   grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2326
2327   const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2328   Int_t a = 0;
2329   Int_t z = 0;
2330
2331   if (gen) {
2332     TString projectile;
2333     gen->GetProjectile(projectile,a,z);
2334     TString target;
2335     gen->GetTarget(target,a,z);
2336     TString beamType = projectile + "-" + target;
2337     beamType.ReplaceAll(" ","");
2338     if (!beamType.CompareTo("-")) {
2339       grpObj->SetBeamType("UNKNOWN");
2340       grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2341     }
2342     else {
2343       grpObj->SetBeamType(beamType);
2344       if (z != 0) {
2345           grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 * a / z);
2346       } else {
2347           grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 );
2348       }
2349       // Heavy ion run, the event specie is set to kHighMult
2350       fEventSpecie = AliRecoParam::kHighMult;
2351       if ((strcmp(beamType,"p-p") == 0) ||
2352           (strcmp(beamType,"p-")  == 0) ||
2353           (strcmp(beamType,"-p")  == 0) ||
2354           (strcmp(beamType,"P-P") == 0) ||
2355           (strcmp(beamType,"P-")  == 0) ||
2356           (strcmp(beamType,"-P")  == 0)) {
2357         // Proton run, the event specie is set to kLowMult
2358         fEventSpecie = AliRecoParam::kLowMult;
2359       } 
2360     }
2361   } else {
2362     AliWarning("Unknown beam type and energy! Setting energy to 0");
2363     grpObj->SetBeamEnergy(0);
2364     grpObj->SetBeamType("UNKNOWN");
2365   }
2366
2367   UInt_t detectorPattern  = 0;
2368   Int_t nDets = 0;
2369   TObjArray *detArray = gAlice->Detectors();
2370   for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2371     if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2372       AliDebug(1, Form("Detector #%d found: %s", iDet, AliDAQ::OfflineModuleName(iDet)));
2373       detectorPattern |= (1 << iDet);
2374       nDets++;
2375     }
2376   }
2377   // CTP
2378   if (!fTriggerConfig.IsNull())
2379     detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2380
2381   // HLT
2382   if (!fRunHLT.IsNull())
2383     detectorPattern |= (1 << AliDAQ::kHLTId);
2384
2385   grpObj->SetNumberOfDetectors((Char_t)nDets);
2386   grpObj->SetDetectorMask((Int_t)detectorPattern);
2387   grpObj->SetLHCPeriod("LHC08c");
2388   grpObj->SetLHCState("STABLE_BEAMS");
2389   //
2390   AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2391   Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2392
2393   Float_t factorSol     = field ? field->GetFactorSol() : 0;
2394   Float_t currentSol    = TMath::Abs(factorSol)>1E-6 ? 
2395     TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2396   //
2397   Float_t factorDip     = field ? field->GetFactorDip() : 0;
2398   Float_t currentDip    = 6000.*TMath::Abs(factorDip);
2399   //
2400   grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2401   grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);  
2402   grpObj->SetL3Polarity(factorSol>0 ? 0:1);  
2403   grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2404   if (field) grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2405   grpObj->SetPolarityConventionLHC();                    // LHC convention +/+ current -> -/- field main components
2406   //
2407   grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2408   
2409   //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2410
2411   // Now store the entry in OCDB
2412   AliCDBManager* man = AliCDBManager::Instance();
2413   
2414   man->SetLock(0, fKey);
2415   
2416   AliCDBStorage* sto = man->GetStorage(fGRPWriteLocation.Data());
2417   
2418
2419   AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2420   AliCDBMetaData *metadata= new AliCDBMetaData();
2421
2422   metadata->SetResponsible("alice-off@cern.ch");
2423   metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2424  
2425   sto->Put(grpObj,id,metadata);
2426   man->SetLock(1, fKey);
2427 }
2428
2429 //_____________________________________________________________________________
2430 time_t AliSimulation::GenerateTimeStamp() const
2431 {
2432   // Generate event time-stamp according to
2433   // SOR/EOR time from GRP
2434   if (fUseTimeStampFromCDB)
2435     return fTimeStart + gRandom->Integer(fTimeEnd-fTimeStart);
2436   else
2437     return 0;
2438 }