a3881924a26ed024168dca16407edb0a40854fa8
[u/mrichter/AliRoot.git] / STEER / AliSimulation.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 // class for running generation, simulation and digitization                 //
21 //                                                                           //
22 // Hits, sdigits and digits are created for all detectors by typing:         //
23 //                                                                           //
24 //   AliSimulation sim;                                                      //
25 //   sim.Run();                                                              //
26 //                                                                           //
27 // The Run method returns kTRUE in case of successful execution.             //
28 // The number of events can be given as argument to the Run method or it     //
29 // can be set by                                                             //
30 //                                                                           //
31 //   sim.SetNumberOfEvents(n);                                               //
32 //                                                                           //
33 // The name of the configuration file can be passed as argument to the       //
34 // AliSimulation constructor or can be specified by                          //
35 //                                                                           //
36 //   sim.SetConfigFile("...");                                               //
37 //                                                                           //
38 // The generation of particles and the simulation of detector hits can be    //
39 // switched on or off by                                                     //
40 //                                                                           //
41 //   sim.SetRunGeneration(kTRUE);   // generation of primary particles       //
42 //   sim.SetRunSimulation(kFALSE);  // but no tracking                       //
43 //                                                                           //
44 // For which detectors sdigits and digits will be created, can be steered    //
45 // by                                                                        //
46 //                                                                           //
47 //   sim.SetMakeSDigits("ALL");     // make sdigits for all detectors        //
48 //   sim.SetMakeDigits("ITS TPC");  // make digits only for ITS and TPC      //
49 //                                                                           //
50 // The argument is a (case sensitive) string with the names of the           //
51 // detectors separated by a space. An empty string ("") can be used to       //
52 // disable the creation of sdigits or digits. The special string "ALL"       //
53 // selects all available detectors. This is the default.                     //
54 //                                                                           //
55 // The creation of digits from hits instead of from sdigits can be selected  //
56 // by                                                                        //
57 //                                                                           //
58 //   sim.SetMakeDigitsFromHits("TRD");                                       //
59 //                                                                           //
60 // The argument is again a string with the selected detectors. Be aware that //
61 // this feature is not available for all detectors and that merging is not   //
62 // possible, when digits are created directly from hits.                     //
63 //                                                                           //
64 // Background events can be merged by calling                                //
65 //                                                                           //
66 //   sim.MergeWith("background/galice.root", 2);                             //
67 //                                                                           //
68 // The first argument is the file name of the background galice file. The    //
69 // second argument is the number of signal events per background event.      //
70 // By default this number is calculated from the number of available         //
71 // background events. MergeWith can be called several times to merge more    //
72 // than two event streams. It is assumed that the sdigits were already       //
73 // produced for the background events.                                       //
74 //                                                                           //
75 // The output of raw data can be switched on by calling                      //
76 //                                                                           //
77 //   sim.SetWriteRawData("MUON");   // write raw data for MUON               //
78 //                                                                           //
79 // The default output format of the raw data are DDL files. They are         //
80 // converted to a DATE file, if a file name is given as second argument.     //
81 // For this conversion the program "dateStream" is required. If the file     //
82 // name has the extension ".root", the DATE file is converted to a root      //
83 // file. The program "alimdc" is used for this purpose. For the conversion   //
84 // to DATE and root format the two conversion programs have to be installed. //
85 // Only the raw data in the final format is kept if the third argument is    //
86 // kTRUE.                                                                    //
87 //                                                                           //
88 // The methods RunSimulation, RunSDigitization, RunDigitization,             //
89 // RunHitsDigitization and WriteRawData can be used to run only parts of     //
90 // the full simulation chain. The creation of raw data DDL files and their   //
91 // conversion to the DATE or root format can be run directly by calling      //
92 // the methods WriteRawFiles, ConvertRawFilesToDate and ConvertDateToRoot.   //
93 //                                                                           //
94 // The default number of events per file, which is usually set in the        //
95 // config file, can be changed for individual detectors and data types       //
96 // by calling                                                                //
97 //                                                                           //
98 //   sim.SetEventsPerFile("PHOS", "Reconstructed Points", 3);                //
99 //                                                                           //
100 // The first argument is the detector, the second one the data type and the  //
101 // last one the number of events per file. Valid data types are "Hits",      //
102 // "Summable Digits", "Digits", "Reconstructed Points" and "Tracks".         //
103 // The number of events per file has to be set before the simulation of      //
104 // hits. Otherwise it has no effect.                                         //
105 //                                                                           //
106 ///////////////////////////////////////////////////////////////////////////////
107
108 #include <TVirtualMCApplication.h>
109 #include <TVirtualMC.h>
110 #include <TGeoManager.h>
111 #include <TObjString.h>
112 #include <TSystem.h>
113 #include <TFile.h>
114 #include <TROOT.h>
115
116 #include "AliCodeTimer.h"
117 #include "AliCDBStorage.h"
118 #include "AliCDBEntry.h"
119 #include "AliCDBManager.h"
120 #include "AliGeomManager.h"
121 #include "AliAlignObj.h"
122 #include "AliCentralTrigger.h"
123 #include "AliDAQ.h"
124 #include "AliDigitizer.h"
125 #include "AliGenerator.h"
126 #include "AliLog.h"
127 #include "AliModule.h"
128 #include "AliRun.h"
129 #include "AliRunDigitizer.h"
130 #include "AliRunLoader.h"
131 #include "AliSimulation.h"
132 #include "AliVertexGenFile.h"
133 #include "AliCentralTrigger.h"
134 #include "AliCTPRawData.h"
135 #include "AliRawReaderFile.h"
136 #include "AliRawReaderRoot.h"
137 #include "AliRawReaderDate.h"
138 #include "AliESD.h"
139 #include "AliHeader.h"
140 #include "AliGenEventHeader.h"
141 #include "AliMC.h"
142 #include "AliHLTSimulation.h"
143 #include "AliSysInfo.h"
144 #include "AliMagF.h"
145 #include "AliGRPObject.h"
146
147 ClassImp(AliSimulation)
148
149 AliSimulation *AliSimulation::fgInstance = 0;
150 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
151
152 //_____________________________________________________________________________
153 AliSimulation::AliSimulation(const char* configFileName,
154                              const char* name, const char* title) :
155   TNamed(name, title),
156
157   fRunGeneration(kTRUE),
158   fRunSimulation(kTRUE),
159   fLoadAlignFromCDB(kTRUE),
160   fLoadAlObjsListOfDets("ALL"),
161   fMakeSDigits("ALL"),
162   fMakeDigits("ALL"),
163   fMakeTrigger(""),
164   fMakeDigitsFromHits(""),
165   fWriteRawData(""),
166   fRawDataFileName(""),
167   fDeleteIntermediateFiles(kFALSE),
168   fWriteSelRawData(kFALSE),
169   fStopOnError(kFALSE),
170   fNEvents(1),
171   fConfigFileName(configFileName),
172   fGAliceFileName("galice.root"),
173   fEventsPerFile(),
174   fBkgrdFileNames(NULL),
175   fAlignObjArray(NULL),
176   fUseBkgrdVertex(kTRUE),
177   fRegionOfInterest(kFALSE),
178   fCDBUri(""),
179   fSpecCDBUri(),
180   fRun(-1),
181   fSeed(0),
182   fInitCDBCalled(kFALSE),
183   fInitRunNumberCalled(kFALSE),
184   fSetRunNumberFromDataCalled(kFALSE),
185   fEmbeddingFlag(kFALSE),
186   fQADetectors("ALL"),                  
187   fQATasks("ALL"),      
188   fQASteer(NULL), 
189   fRunQA(kTRUE), 
190   fEventSpecie(AliRecoParam::kDefault),
191   fRunHLT("default"),
192   fWriteGRPEntry(kTRUE)
193 {
194 // create simulation object with default parameters
195   fgInstance = this;
196   SetGAliceFile("galice.root");
197   
198 // for QA
199         fQASteer = new AliQADataMakerSteer("sim") ; 
200         fQASteer->SetActiveDetectors(fQADetectors) ; 
201         fQATasks = Form("%d %d %d", AliQA::kHITS, AliQA::kSDIGITS, AliQA::kDIGITS) ; 
202         fQASteer->SetTasks(fQATasks) ;  
203 }
204
205 //_____________________________________________________________________________
206 AliSimulation::~AliSimulation()
207 {
208 // clean up
209
210   fEventsPerFile.Delete();
211 //  if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
212 //  delete fAlignObjArray; fAlignObjArray=0;
213
214   if (fBkgrdFileNames) {
215     fBkgrdFileNames->Delete();
216     delete fBkgrdFileNames;
217   }
218
219   fSpecCDBUri.Delete();
220   if (fgInstance==this) fgInstance = 0;
221
222         delete fQASteer ; 
223         
224   AliCodeTimer::Instance()->Print();
225 }
226
227
228 //_____________________________________________________________________________
229 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
230 {
231 // set the number of events for one run
232
233   fNEvents = nEvents;
234 }
235
236 //_____________________________________________________________________________
237 void AliSimulation::InitCDB()
238 {
239 // activate a default CDB storage
240 // First check if we have any CDB storage set, because it is used 
241 // to retrieve the calibration and alignment constants
242
243   if (fInitCDBCalled) return;
244   fInitCDBCalled = kTRUE;
245
246   AliCDBManager* man = AliCDBManager::Instance();
247   if (man->IsDefaultStorageSet())
248   {
249     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
250     AliWarning("Default CDB storage has been already set !");
251     AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
252     AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
253     fCDBUri = man->GetDefaultStorage()->GetURI();
254   }
255   else {
256     if (fCDBUri.Length() > 0) 
257     {
258         AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
259         AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
260         AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
261     } else {
262         fCDBUri="local://$ALICE_ROOT";
263         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
264         AliWarning("Default CDB storage not yet set !!!!");
265         AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
266         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
267                 
268     }
269     man->SetDefaultStorage(fCDBUri);
270   }
271
272   // Now activate the detector specific CDB storage locations
273   for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
274     TObject* obj = fSpecCDBUri[i];
275     if (!obj) continue;
276     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
277     AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
278     AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
279     man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
280   }
281       
282 }
283
284 //_____________________________________________________________________________
285 void AliSimulation::InitRunNumber(){
286 // check run number. If not set, set it to 0 !!!!
287   
288   if (fInitRunNumberCalled) return;
289   fInitRunNumberCalled = kTRUE;
290   
291   AliCDBManager* man = AliCDBManager::Instance();
292   if (man->GetRun() >= 0)
293   {
294         AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
295                         "Use external variable DC_RUN or AliSimulation::SetRun()!"));
296   }
297     
298   if(fRun >= 0) {
299         AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
300         AliDebug(2, Form("Setting CDB run number to: %d",fRun));
301         AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
302   } else {
303         fRun=0;
304         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
305         AliWarning("Run number not yet set !!!!");
306         AliWarning(Form("Setting it now to: %d", fRun));
307         AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
308         
309   }
310   man->SetRun(fRun);
311
312   man->Print();
313
314 }
315
316 //_____________________________________________________________________________
317 void AliSimulation::SetCDBLock() {
318   // Set CDB lock: from now on it is forbidden to reset the run number
319   // or the default storage or to activate any further storage!
320   
321   AliCDBManager::Instance()->SetLock(1);
322 }
323
324 //_____________________________________________________________________________
325 void AliSimulation::SetDefaultStorage(const char* uri) {
326 // Store the desired default CDB storage location
327 // Activate it later within the Run() method
328
329   fCDBUri = uri;
330
331 }
332
333 //_____________________________________________________________________________
334 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
335 // Store a detector-specific CDB storage location
336 // Activate it later within the Run() method
337
338   AliCDBPath aPath(calibType);
339   if(!aPath.IsValid()){
340         AliError(Form("Not a valid path: %s", calibType));
341         return;
342   }
343
344   TObject* obj = fSpecCDBUri.FindObject(calibType);
345   if (obj) fSpecCDBUri.Remove(obj);
346   fSpecCDBUri.Add(new TNamed(calibType, uri));
347
348 }
349
350 //_____________________________________________________________________________
351 void AliSimulation::SetRunNumber(Int_t run)
352 {
353 // sets run number
354 // Activate it later within the Run() method
355
356         fRun = run;
357 }
358
359 //_____________________________________________________________________________
360 void AliSimulation::SetSeed(Int_t seed)
361 {
362 // sets seed number
363 // Activate it later within the Run() method
364
365         fSeed = seed;
366 }
367
368 //_____________________________________________________________________________
369 Bool_t AliSimulation::SetRunNumberFromData()
370 {
371   // Set the CDB manager run number
372   // The run number is retrieved from gAlice
373
374     if (fSetRunNumberFromDataCalled) return kTRUE;
375     fSetRunNumberFromDataCalled = kTRUE;    
376   
377     AliCDBManager* man = AliCDBManager::Instance();
378     Int_t runData = -1, runCDB = -1;
379   
380     AliRunLoader* runLoader = LoadRun("READ");
381     if (!runLoader) return kFALSE;
382     else {
383         runData = runLoader->GetHeader()->GetRun();
384         delete runLoader;
385     }
386   
387     runCDB = man->GetRun();
388     if(runCDB >= 0) {
389         if (runCDB != runData) {
390                 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
391                 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
392                 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
393                 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");    
394         }
395         
396     }
397       
398     man->SetRun(runData);
399     fRun = runData;
400     
401     if(man->GetRun() < 0) {
402         AliError("Run number not properly initalized!");
403         return kFALSE;
404     }
405   
406     man->Print();
407     
408     return kTRUE;
409 }
410
411 //_____________________________________________________________________________
412 void AliSimulation::SetConfigFile(const char* fileName)
413 {
414 // set the name of the config file
415
416   fConfigFileName = fileName;
417 }
418
419 //_____________________________________________________________________________
420 void AliSimulation::SetGAliceFile(const char* fileName)
421 {
422 // set the name of the galice file
423 // the path is converted to an absolute one if it is relative
424
425   fGAliceFileName = fileName;
426   if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
427     char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
428                                                 fGAliceFileName);
429     fGAliceFileName = absFileName;
430     delete[] absFileName;
431   }
432
433   AliDebug(2, Form("galice file name set to %s", fileName));
434 }
435
436 //_____________________________________________________________________________
437 void AliSimulation::SetEventsPerFile(const char* detector, const char* type, 
438                                      Int_t nEvents)
439 {
440 // set the number of events per file for the given detector and data type
441 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
442
443   TNamed* obj = new TNamed(detector, type);
444   obj->SetUniqueID(nEvents);
445   fEventsPerFile.Add(obj);
446 }
447
448 //_____________________________________________________________________________
449 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
450 {
451   // Read the alignment objects from CDB.
452   // Each detector is supposed to have the
453   // alignment objects in DET/Align/Data CDB path.
454   // All the detector objects are then collected,
455   // sorted by geometry level (starting from ALIC) and
456   // then applied to the TGeo geometry.
457   // Finally an overlaps check is performed.
458
459   if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
460     AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
461     return kFALSE;
462   }  
463   
464   // initialize CDB storage, run number, set CDB lock
465   InitCDB();
466 //  if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
467   SetCDBLock();
468     
469   Bool_t delRunLoader = kFALSE;
470   if (!runLoader) {
471     runLoader = LoadRun("READ");
472     if (!runLoader) return kFALSE;
473     delRunLoader = kTRUE;
474   }
475   
476   // Export ideal geometry 
477   if(!gAlice->IsRootGeometry()) AliGeomManager::GetGeometry()->Export("geometry.root");
478
479   // Load alignment data from CDB and apply to geometry through AliGeomManager
480   if(fLoadAlignFromCDB){
481     
482     TString detStr = fLoadAlObjsListOfDets;
483     TString loadAlObjsListOfDets = "";
484     
485     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
486     for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
487       AliModule* det = (AliModule*) detArray->At(iDet);
488       if (!det || !det->IsActive()) continue;
489       if (IsSelected(det->GetName(), detStr)) {
490         //add det to list of dets to be aligned from CDB
491         loadAlObjsListOfDets += det->GetName();
492         loadAlObjsListOfDets += " ";
493       }
494     } // end loop over detectors
495     loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
496     AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
497   }else{
498     // Check if the array with alignment objects was
499     // provided by the user. If yes, apply the objects
500     // to the present TGeo geometry
501     if (fAlignObjArray) {
502       if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
503         AliError("The misalignment of one or more volumes failed!"
504                  "Compare the list of simulated detectors and the list of detector alignment data!");
505         if (delRunLoader) delete runLoader;
506         return kFALSE;
507       }
508     }
509   }
510
511   // Update the internal geometry of modules (ITS needs it)
512   TString detStr = fLoadAlObjsListOfDets;
513   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
514   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
515
516     AliModule* det = (AliModule*) detArray->At(iDet);
517     if (!det || !det->IsActive()) continue;
518     if (IsSelected(det->GetName(), detStr)) {
519       det->UpdateInternalGeometry();
520     }
521   } // end loop over detectors
522
523
524   if (delRunLoader) delete runLoader;
525
526   return kTRUE;
527 }
528
529 //_____________________________________________________________________________
530 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
531 {
532 // add a file with background events for merging
533
534   TObjString* fileNameStr = new TObjString(fileName);
535   fileNameStr->SetUniqueID(nSignalPerBkgrd);
536   if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
537   fBkgrdFileNames->Add(fileNameStr);
538 }
539
540 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
541 {
542 // add a file with background events for embeddin
543   MergeWith(fileName, nSignalPerBkgrd);
544   fEmbeddingFlag = kTRUE;
545 }
546
547 //_____________________________________________________________________________
548 Bool_t AliSimulation::Run(Int_t nEvents)
549 {
550 // run the generation, simulation and digitization
551
552  
553   AliCodeTimerAuto("")
554   
555   // Load run number and seed from environmental vars
556   ProcessEnvironmentVars();
557
558   gRandom->SetSeed(fSeed);
559    
560   if (nEvents > 0) fNEvents = nEvents;
561
562   // generation and simulation -> hits
563   if (fRunGeneration) {
564     if (!RunSimulation()) if (fStopOnError) return kFALSE;
565   }
566            
567   // initialize CDB storage from external environment
568   // (either CDB manager or AliSimulation setters),
569   // if not already done in RunSimulation()
570   InitCDB();
571   
572   // Set run number in CDBManager from data 
573   // From this point on the run number must be always loaded from data!
574   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
575   
576   // Set CDB lock: from now on it is forbidden to reset the run number
577   // or the default storage or to activate any further storage!
578   SetCDBLock();
579
580   // If RunSimulation was not called, load the geometry and misalign it
581   if (!AliGeomManager::GetGeometry()) {
582     // Initialize the geometry manager
583     AliGeomManager::LoadGeometry("geometry.root");
584     
585 //    // Check that the consistency of symbolic names for the activated subdetectors
586 //    // in the geometry loaded by AliGeomManager
587 //    AliRunLoader* runLoader = LoadRun("READ");
588 //    if (!runLoader) return kFALSE;
589 //
590 //    TString detsToBeChecked = "";
591 //    TObjArray* detArray = runLoader->GetAliRun()->Detectors();
592 //    for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
593 //      AliModule* det = (AliModule*) detArray->At(iDet);
594 //      if (!det || !det->IsActive()) continue;
595 //      detsToBeChecked += det->GetName();
596 //      detsToBeChecked += " ";
597 //    } // end loop over detectors
598 //    if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
599     if(!AliGeomManager::CheckSymNamesLUT("ALL"))
600         AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
601         
602     if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
603     // Misalign geometry
604     if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
605   }
606
607
608   // hits -> summable digits
609   AliSysInfo::AddStamp("Start_sdigitization");
610   if (!fMakeSDigits.IsNull()) {
611     if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
612  
613   }
614   AliSysInfo::AddStamp("Stop_sdigitization");
615   
616   AliSysInfo::AddStamp("Start_digitization");  
617   // summable digits -> digits  
618   if (!fMakeDigits.IsNull()) {
619     if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
620       if (fStopOnError) return kFALSE;
621     }
622    }
623   AliSysInfo::AddStamp("Stop_digitization");
624
625   
626   
627   // hits -> digits
628   if (!fMakeDigitsFromHits.IsNull()) {
629     if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
630       AliWarning(Form("Merging and direct creation of digits from hits " 
631                  "was selected for some detectors. "
632                  "No merging will be done for the following detectors: %s",
633                  fMakeDigitsFromHits.Data()));
634     }
635     if (!RunHitsDigitization(fMakeDigitsFromHits)) {
636       if (fStopOnError) return kFALSE;
637     }
638   }
639
640   
641   
642   // digits -> trigger
643   if (!RunTrigger(fMakeTrigger,fMakeDigits)) {
644     if (fStopOnError) return kFALSE;
645   }
646
647   
648   
649   // digits -> raw data
650   if (!fWriteRawData.IsNull()) {
651     if (!WriteRawData(fWriteRawData, fRawDataFileName, 
652                       fDeleteIntermediateFiles,fWriteSelRawData)) {
653       if (fStopOnError) return kFALSE;
654     }
655   }
656
657   
658   // run HLT simulation on simulated digit data if raw data is not
659   // simulated, otherwise its called as part of WriteRawData
660   if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
661     if (!RunHLT()) {
662       if (fStopOnError) return kFALSE;
663     }
664   }
665   
666   //QA
667         if (fRunQA) {
668                 Bool_t rv = RunQA() ; 
669                 if (!rv)
670                         if (fStopOnError) 
671                                 return kFALSE ;         
672         }
673
674   // Cleanup of CDB manager: cache and active storages!
675   AliCDBManager::Instance()->ClearCache();
676
677   return kTRUE;
678 }
679
680 //_____________________________________________________________________________
681 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
682 {
683   // run the trigger
684
685   AliCodeTimerAuto("")
686
687   // initialize CDB storage from external environment
688   // (either CDB manager or AliSimulation setters),
689   // if not already done in RunSimulation()
690   InitCDB();
691   
692   // Set run number in CDBManager from data 
693   // From this point on the run number must be always loaded from data!
694   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
695   
696   // Set CDB lock: from now on it is forbidden to reset the run number
697   // or the default storage or to activate any further storage!
698   SetCDBLock();
699    
700    AliRunLoader* runLoader = LoadRun("READ");
701    if (!runLoader) return kFALSE;
702    TString trconfiguration = config;
703
704    if (trconfiguration.IsNull()) {
705      if (strcmp(gAlice->GetTriggerDescriptor(),"")) {
706        trconfiguration = gAlice->GetTriggerDescriptor();
707      }
708      else
709        AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
710    }
711
712    runLoader->MakeTree( "GG" );
713    AliCentralTrigger* aCTP = runLoader->GetTrigger();
714    // Load Configuration
715    if (!aCTP->LoadConfiguration( trconfiguration ))
716      return kFALSE;
717
718    // digits -> trigger
719    if( !aCTP->RunTrigger( runLoader , detectors ) ) {
720       if (fStopOnError) {
721         //  delete aCTP;
722         return kFALSE;
723       }
724    }
725
726    delete runLoader;
727
728    return kTRUE;
729 }
730
731 //_____________________________________________________________________________
732 Bool_t AliSimulation::WriteTriggerRawData()
733 {
734   // Writes the CTP (trigger) DDL raw data
735   // Details of the format are given in the
736   // trigger TDR - pages 134 and 135.
737   AliCTPRawData writer;
738   writer.RawData();
739
740   return kTRUE;
741 }
742
743 //_____________________________________________________________________________
744 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
745 {
746 // run the generation and simulation
747
748   AliCodeTimerAuto("")
749
750   // initialize CDB storage and run number from external environment
751   // (either CDB manager or AliSimulation setters)
752   InitCDB();
753   InitRunNumber();
754   SetCDBLock();
755   
756   if (!gAlice) {
757     AliError("no gAlice object. Restart aliroot and try again.");
758     return kFALSE;
759   }
760   if (gAlice->Modules()->GetEntries() > 0) {
761     AliError("gAlice was already run. Restart aliroot and try again.");
762     return kFALSE;
763   }
764
765   AliInfo(Form("initializing gAlice with config file %s",
766           fConfigFileName.Data()));
767   StdoutToAliInfo(StderrToAliError(
768     gAlice->Init(fConfigFileName.Data());
769   ););
770   
771   // Get the trigger descriptor string
772   // Either from AliSimulation or from
773   // gAlice
774   if (fMakeTrigger.IsNull()) {
775     if (strcmp(gAlice->GetTriggerDescriptor(),""))
776       fMakeTrigger = gAlice->GetTriggerDescriptor();
777   }
778   else
779     gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
780
781   // Set run number in CDBManager
782   AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
783
784   AliRunLoader* runLoader = AliRunLoader::GetRunLoader();
785   if (!runLoader) {
786              AliError(Form("gAlice has no run loader object. "
787                              "Check your config file: %s", fConfigFileName.Data()));
788              return kFALSE;
789   }
790   SetGAliceFile(runLoader->GetFileName());
791       
792   // Misalign geometry
793 #if ROOT_VERSION_CODE < 331527
794   AliGeomManager::SetGeometry(gGeoManager);
795   
796   // Check that the consistency of symbolic names for the activated subdetectors
797   // in the geometry loaded by AliGeomManager
798   TString detsToBeChecked = "";
799   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
800   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
801     AliModule* det = (AliModule*) detArray->At(iDet);
802     if (!det || !det->IsActive()) continue;
803     detsToBeChecked += det->GetName();
804     detsToBeChecked += " ";
805   } // end loop over detectors
806   if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
807     AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
808   MisalignGeometry(runLoader);
809 #endif
810
811 //   AliRunLoader* runLoader = AliRunLoader::GetRunLoader();
812 //   if (!runLoader) {
813 //     AliError(Form("gAlice has no run loader object. "
814 //                   "Check your config file: %s", fConfigFileName.Data()));
815 //     return kFALSE;
816 //   }
817 //   SetGAliceFile(runLoader->GetFileName());
818
819   if (!gAlice->Generator()) {
820     AliError(Form("gAlice has no generator object. "
821                   "Check your config file: %s", fConfigFileName.Data()));
822     return kFALSE;
823   }
824
825   // Write GRP entry corresponding to the setting found in Cofig.C
826   if (fWriteGRPEntry)
827     WriteGRPEntry();
828
829   if (nEvents <= 0) nEvents = fNEvents;
830
831   // get vertex from background file in case of merging
832   if (fUseBkgrdVertex &&
833       fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
834     Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
835     const char* fileName = ((TObjString*)
836                             (fBkgrdFileNames->At(0)))->GetName();
837     AliInfo(Form("The vertex will be taken from the background "
838                  "file %s with nSignalPerBackground = %d", 
839                  fileName, signalPerBkgrd));
840     AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
841     gAlice->Generator()->SetVertexGenerator(vtxGen);
842   }
843
844   if (!fRunSimulation) {
845     gAlice->Generator()->SetTrackingFlag(0);
846   }
847
848   // set the number of events per file for given detectors and data types
849   for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
850     if (!fEventsPerFile[i]) continue;
851     const char* detName = fEventsPerFile[i]->GetName();
852     const char* typeName = fEventsPerFile[i]->GetTitle();
853     TString loaderName(detName);
854     loaderName += "Loader";
855     AliLoader* loader = runLoader->GetLoader(loaderName);
856     if (!loader) {
857       AliError(Form("RunSimulation", "no loader for %s found\n"
858                     "Number of events per file not set for %s %s", 
859                     detName, typeName, detName));
860       continue;
861     }
862     AliDataLoader* dataLoader = 
863       loader->GetDataLoader(typeName);
864     if (!dataLoader) {
865       AliError(Form("no data loader for %s found\n"
866                     "Number of events per file not set for %s %s", 
867                     typeName, detName, typeName));
868       continue;
869     }
870     dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
871     AliDebug(1, Form("number of events per file set to %d for %s %s",
872                      fEventsPerFile[i]->GetUniqueID(), detName, typeName));
873   }
874
875   AliInfo("running gAlice");
876   AliSysInfo::AddStamp("Start_simulation");
877
878   gAlice->InitMC("Config.C");
879   
880   // Create the Root Tree with one branch per detector
881   //Hits moved to begin event -> now we are crating separate tree for each event
882
883   gMC->ProcessRun(nEvents);
884
885   // End of this run, close files
886   if(nEvents>0) gAlice->FinishRun();
887
888   AliSysInfo::AddStamp("Stop_simulation");
889   delete runLoader;
890
891   return kTRUE;
892 }
893
894 //_____________________________________________________________________________
895 Bool_t AliSimulation::RunSDigitization(const char* detectors)
896 {
897 // run the digitization and produce summable digits
898   static Int_t eventNr=0;
899   AliCodeTimerAuto("")
900
901   // initialize CDB storage, run number, set CDB lock
902   InitCDB();
903   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
904   SetCDBLock();
905   
906   AliRunLoader* runLoader = LoadRun();
907   if (!runLoader) return kFALSE;
908
909   TString detStr = detectors;
910   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
911   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
912     AliModule* det = (AliModule*) detArray->At(iDet);
913     if (!det || !det->IsActive()) continue;
914     if (IsSelected(det->GetName(), detStr)) {
915       AliInfo(Form("creating summable digits for %s", det->GetName()));
916       AliCodeTimerAuto(Form("creating summable digits for %s", det->GetName()));
917       det->Hits2SDigits();
918       AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
919     }
920   }
921
922   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
923     AliError(Form("the following detectors were not found: %s",
924                   detStr.Data()));
925     if (fStopOnError) return kFALSE;
926   }
927   eventNr++;
928   delete runLoader;
929
930   return kTRUE;
931 }
932
933
934 //_____________________________________________________________________________
935 Bool_t AliSimulation::RunDigitization(const char* detectors, 
936                                       const char* excludeDetectors)
937 {
938 // run the digitization and produce digits from sdigits
939
940   AliCodeTimerAuto("")
941
942   // initialize CDB storage, run number, set CDB lock
943   InitCDB();
944   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
945   SetCDBLock();
946   
947   delete AliRunLoader::GetRunLoader();
948   delete gAlice;
949   gAlice = NULL;
950
951   Int_t nStreams = 1;
952   if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
953   Int_t signalPerBkgrd = GetNSignalPerBkgrd();
954   AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
955   // manager->SetEmbeddingFlag(fEmbeddingFlag);
956   manager->SetInputStream(0, fGAliceFileName.Data());
957   for (Int_t iStream = 1; iStream < nStreams; iStream++) {
958     const char* fileName = ((TObjString*)
959                             (fBkgrdFileNames->At(iStream-1)))->GetName();
960     manager->SetInputStream(iStream, fileName);
961   }
962
963   TString detStr = detectors;
964   TString detExcl = excludeDetectors;
965   manager->GetInputStream(0)->ImportgAlice();
966   AliRunLoader* runLoader = 
967     AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
968   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
969   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
970     AliModule* det = (AliModule*) detArray->At(iDet);
971     if (!det || !det->IsActive()) continue;
972     if (IsSelected(det->GetName(), detStr) && 
973         !IsSelected(det->GetName(), detExcl)) {
974       AliDigitizer* digitizer = det->CreateDigitizer(manager);
975       
976       if (!digitizer) {
977         AliError(Form("no digitizer for %s", det->GetName()));
978         if (fStopOnError) return kFALSE;
979       } else {
980         digitizer->SetRegionOfInterest(fRegionOfInterest);
981       }
982     }
983   }
984
985   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
986     AliError(Form("the following detectors were not found: %s", 
987                   detStr.Data()));
988     if (fStopOnError) return kFALSE;
989   }
990
991   if (!manager->GetListOfTasks()->IsEmpty()) {
992     AliInfo("executing digitization");
993     manager->Exec("");
994   }
995
996   delete manager;
997
998   return kTRUE;
999 }
1000
1001 //_____________________________________________________________________________
1002 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1003 {
1004 // run the digitization and produce digits from hits
1005
1006   AliCodeTimerAuto("")
1007
1008   // initialize CDB storage, run number, set CDB lock
1009   InitCDB();
1010   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1011   SetCDBLock();
1012   
1013   AliRunLoader* runLoader = LoadRun("READ");
1014   if (!runLoader) return kFALSE;
1015
1016   TString detStr = detectors;
1017   TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1018   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1019     AliModule* det = (AliModule*) detArray->At(iDet);
1020     if (!det || !det->IsActive()) continue;
1021     if (IsSelected(det->GetName(), detStr)) {
1022       AliInfo(Form("creating digits from hits for %s", det->GetName()));
1023       det->Hits2Digits();
1024     }
1025   }
1026
1027   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1028     AliError(Form("the following detectors were not found: %s", 
1029                   detStr.Data()));
1030     if (fStopOnError) return kFALSE;
1031   }
1032
1033   return kTRUE;
1034 }
1035
1036 //_____________________________________________________________________________
1037 Bool_t AliSimulation::WriteRawData(const char* detectors, 
1038                                    const char* fileName,
1039                                    Bool_t deleteIntermediateFiles,
1040                                    Bool_t selrawdata)
1041 {
1042 // convert the digits to raw data
1043 // First DDL raw data files for the given detectors are created.
1044 // If a file name is given, the DDL files are then converted to a DATE file.
1045 // If deleteIntermediateFiles is true, the DDL raw files are deleted 
1046 // afterwards.
1047 // If the file name has the extension ".root", the DATE file is converted
1048 // to a root file.
1049 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1050 // 'selrawdata' flag can be used to enable writing of detectors raw data
1051 // accoring to the trigger cluster.
1052
1053   AliCodeTimerAuto("")
1054   
1055   TString detStr = detectors;
1056   if (!WriteRawFiles(detStr.Data())) {
1057     if (fStopOnError) return kFALSE;
1058   }
1059
1060   // run HLT simulation on simulated DDL raw files
1061   // and produce HLT ddl raw files to be included in date/root file
1062   if (IsSelected("HLT", detStr) && !fRunHLT.IsNull()) {
1063     if (!RunHLT()) {
1064       if (fStopOnError) return kFALSE;
1065     }
1066   }
1067
1068   TString dateFileName(fileName);
1069   if (!dateFileName.IsNull()) {
1070     Bool_t rootOutput = dateFileName.EndsWith(".root");
1071     if (rootOutput) dateFileName += ".date";
1072     TString selDateFileName;
1073     if (selrawdata) {
1074       selDateFileName = "selected.";
1075       selDateFileName+= dateFileName;
1076     }
1077     if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1078       if (fStopOnError) return kFALSE;
1079     }
1080     if (deleteIntermediateFiles) {
1081       AliRunLoader* runLoader = LoadRun("READ");
1082       if (runLoader) for (Int_t iEvent = 0; 
1083                           iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1084         char command[256];
1085         sprintf(command, "rm -r raw%d", iEvent);
1086         gSystem->Exec(command);
1087       }
1088     }
1089
1090     if (rootOutput) {
1091       if (!ConvertDateToRoot(dateFileName, fileName)) {
1092         if (fStopOnError) return kFALSE;
1093       }
1094       if (deleteIntermediateFiles) {
1095         gSystem->Unlink(dateFileName);
1096       }
1097       if (selrawdata) {
1098         TString selFileName = "selected.";
1099         selFileName        += fileName;
1100         if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1101           if (fStopOnError) return kFALSE;
1102         }
1103         if (deleteIntermediateFiles) {
1104           gSystem->Unlink(selDateFileName);
1105         }
1106       }
1107     }
1108   }
1109
1110   return kTRUE;
1111 }
1112
1113 //_____________________________________________________________________________
1114 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1115 {
1116 // convert the digits to raw data DDL files
1117
1118   AliCodeTimerAuto("")
1119   
1120   AliRunLoader* runLoader = LoadRun("READ");
1121   if (!runLoader) return kFALSE;
1122
1123   // write raw data to DDL files
1124   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1125     AliInfo(Form("processing event %d", iEvent));
1126     runLoader->GetEvent(iEvent);
1127     TString baseDir = gSystem->WorkingDirectory();
1128     char dirName[256];
1129     sprintf(dirName, "raw%d", iEvent);
1130     gSystem->MakeDirectory(dirName);
1131     if (!gSystem->ChangeDirectory(dirName)) {
1132       AliError(Form("couldn't change to directory %s", dirName));
1133       if (fStopOnError) return kFALSE; else continue;
1134     }
1135
1136     ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1137     runNbFile.close();
1138
1139     TString detStr = detectors;
1140     if (IsSelected("HLT", detStr)) {
1141       // Do nothing. "HLT" will be removed from detStr and HLT raw
1142       // data files are generated in RunHLT.
1143     }
1144
1145     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1146     for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1147       AliModule* det = (AliModule*) detArray->At(iDet);
1148       if (!det || !det->IsActive()) continue;
1149       if (IsSelected(det->GetName(), detStr)) {
1150         AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1151         det->Digits2Raw();
1152       }
1153     }
1154
1155     if (!WriteTriggerRawData())
1156       if (fStopOnError) return kFALSE;
1157
1158     gSystem->ChangeDirectory(baseDir);
1159     if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1160       AliError(Form("the following detectors were not found: %s", 
1161                     detStr.Data()));
1162       if (fStopOnError) return kFALSE;
1163     }
1164   }
1165
1166   delete runLoader;
1167   
1168   return kTRUE;
1169 }
1170
1171 //_____________________________________________________________________________
1172 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1173                                             const char* selDateFileName)
1174 {
1175 // convert raw data DDL files to a DATE file with the program "dateStream"
1176 // The second argument is not empty when the user decides to write
1177 // the detectors raw data according to the trigger cluster.
1178
1179   AliCodeTimerAuto("")
1180   
1181   char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1182   if (!path) {
1183     AliError("the program dateStream was not found");
1184     if (fStopOnError) return kFALSE;
1185   } else {
1186     delete[] path;
1187   }
1188
1189   AliRunLoader* runLoader = LoadRun("READ");
1190   if (!runLoader) return kFALSE;
1191
1192   AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1193   Bool_t selrawdata = kFALSE;
1194   if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1195
1196   char command[256];
1197   // Note the option -s. It is used in order to avoid
1198   // the generation of SOR/EOR events.
1199   sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d", 
1200           dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1201   FILE* pipe = gSystem->OpenPipe(command, "w");
1202
1203   if (!pipe) {
1204     AliError(Form("Cannot execute command: %s",command));
1205     return kFALSE;
1206   }
1207
1208   Int_t selEvents = 0;
1209   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1210
1211     UInt_t detectorPattern = 0;
1212     runLoader->GetEvent(iEvent);
1213     if (!runLoader->LoadTrigger()) {
1214       AliCentralTrigger *aCTP = runLoader->GetTrigger();
1215       detectorPattern = aCTP->GetClusterMask();
1216       // Check if the event was triggered by CTP
1217       if (selrawdata) {
1218         if (aCTP->GetClassMask()) selEvents++;
1219       }
1220     }
1221     else {
1222       AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1223       if (selrawdata) {
1224         AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1225         selrawdata = kFALSE;
1226       }
1227     }
1228
1229     fprintf(pipe, "GDC DetectorPattern %u\n", detectorPattern);
1230     Float_t ldc = 0;
1231     Int_t prevLDC = -1;
1232
1233     // loop over detectors and DDLs
1234     for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1235       for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1236
1237         Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1238         Int_t ldcID = Int_t(ldc + 0.0001);
1239         ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1240
1241         char rawFileName[256];
1242         sprintf(rawFileName, "raw%d/%s", 
1243                 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1244
1245         // check existence and size of raw data file
1246         FILE* file = fopen(rawFileName, "rb");
1247         if (!file) continue;
1248         fseek(file, 0, SEEK_END);
1249         unsigned long size = ftell(file);
1250         fclose(file);
1251         if (!size) continue;
1252
1253         if (ldcID != prevLDC) {
1254           fprintf(pipe, " LDC Id %d\n", ldcID);
1255           prevLDC = ldcID;
1256         }
1257         fprintf(pipe, "  Equipment Id %d Payload %s\n", ddlID, rawFileName);
1258       }
1259     }
1260   }
1261
1262   Int_t result = gSystem->ClosePipe(pipe);
1263
1264   if (!(selrawdata && selEvents > 0)) {
1265     delete runLoader;
1266     return (result == 0);
1267   }
1268
1269   AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1270   
1271   sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d", 
1272           selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1273   FILE* pipe2 = gSystem->OpenPipe(command, "w");
1274
1275   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1276
1277     // Get the trigger decision and cluster
1278     UInt_t detectorPattern = 0;
1279     TString detClust;
1280     runLoader->GetEvent(iEvent);
1281     if (!runLoader->LoadTrigger()) {
1282       AliCentralTrigger *aCTP = runLoader->GetTrigger();
1283       if (aCTP->GetClassMask() == 0) continue;
1284       detectorPattern = aCTP->GetClusterMask();
1285       detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1286       AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1287     }
1288
1289     fprintf(pipe2, "GDC DetectorPattern %u\n", detectorPattern);
1290     Float_t ldc = 0;
1291     Int_t prevLDC = -1;
1292
1293     // loop over detectors and DDLs
1294     for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1295       // Write only raw data from detectors that
1296       // are contained in the trigger cluster(s)
1297       if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1298
1299       for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1300
1301         Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1302         Int_t ldcID = Int_t(ldc + 0.0001);
1303         ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1304
1305         char rawFileName[256];
1306         sprintf(rawFileName, "raw%d/%s", 
1307                 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1308
1309         // check existence and size of raw data file
1310         FILE* file = fopen(rawFileName, "rb");
1311         if (!file) continue;
1312         fseek(file, 0, SEEK_END);
1313         unsigned long size = ftell(file);
1314         fclose(file);
1315         if (!size) continue;
1316
1317         if (ldcID != prevLDC) {
1318           fprintf(pipe2, " LDC Id %d\n", ldcID);
1319           prevLDC = ldcID;
1320         }
1321         fprintf(pipe2, "  Equipment Id %d Payload %s\n", ddlID, rawFileName);
1322       }
1323     }
1324   }
1325
1326   Int_t result2 = gSystem->ClosePipe(pipe2);
1327
1328   delete runLoader;
1329   return ((result == 0) && (result2 == 0));
1330 }
1331
1332 //_____________________________________________________________________________
1333 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1334                                         const char* rootFileName)
1335 {
1336 // convert a DATE file to a root file with the program "alimdc"
1337
1338   // ALIMDC setup
1339   const Int_t kDBSize = 2000000000;
1340   const Int_t kTagDBSize = 1000000000;
1341   const Bool_t kFilter = kFALSE;
1342   const Int_t kCompression = 1;
1343
1344   char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1345   if (!path) {
1346     AliError("the program alimdc was not found");
1347     if (fStopOnError) return kFALSE;
1348   } else {
1349     delete[] path;
1350   }
1351
1352   AliInfo(Form("converting DATE file %s to root file %s", 
1353                dateFileName, rootFileName));
1354
1355   const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1356   const char* tagDBFS    = "/tmp/mdc1/tags";
1357
1358   // User defined file system locations
1359   if (gSystem->Getenv("ALIMDC_RAWDB1")) 
1360     rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1361   if (gSystem->Getenv("ALIMDC_RAWDB2")) 
1362     rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1363   if (gSystem->Getenv("ALIMDC_TAGDB")) 
1364     tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1365
1366   gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1367   gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1368   gSystem->Exec(Form("rm -rf %s",tagDBFS));
1369
1370   gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1371   gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1372   gSystem->Exec(Form("mkdir %s",tagDBFS));
1373
1374   Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s", 
1375                                     kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1376   gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1377
1378   gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1379   gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1380   gSystem->Exec(Form("rm -rf %s",tagDBFS));
1381
1382   return (result == 0);
1383 }
1384
1385
1386 //_____________________________________________________________________________
1387 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1388 {
1389 // delete existing run loaders, open a new one and load gAlice
1390
1391   delete AliRunLoader::GetRunLoader();
1392   AliRunLoader* runLoader = 
1393     AliRunLoader::Open(fGAliceFileName.Data(), 
1394                        AliConfig::GetDefaultEventFolderName(), mode);
1395   if (!runLoader) {
1396     AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1397     return NULL;
1398   }
1399   runLoader->LoadgAlice();
1400   runLoader->LoadHeader();
1401   gAlice = runLoader->GetAliRun();
1402   if (!gAlice) {
1403     AliError(Form("no gAlice object found in file %s", 
1404                   fGAliceFileName.Data()));
1405     return NULL;
1406   }
1407   return runLoader;
1408 }
1409
1410 //_____________________________________________________________________________
1411 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1412 {
1413 // get or calculate the number of signal events per background event
1414
1415   if (!fBkgrdFileNames) return 1;
1416   Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1417   if (nBkgrdFiles == 0) return 1;
1418
1419   // get the number of signal events
1420   if (nEvents <= 0) {
1421     AliRunLoader* runLoader = 
1422         AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1423     if (!runLoader) return 1;
1424     
1425     nEvents = runLoader->GetNumberOfEvents();
1426     delete runLoader;
1427   }
1428
1429   Int_t result = 0;
1430   for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1431     // get the number of background events
1432     const char* fileName = ((TObjString*)
1433                             (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1434     AliRunLoader* runLoader =
1435       AliRunLoader::Open(fileName, "BKGRD");
1436     if (!runLoader) continue;
1437     Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1438     delete runLoader;
1439   
1440     // get or calculate the number of signal per background events
1441     Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1442     if (nSignalPerBkgrd <= 0) {
1443       nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1444     } else if (result && (result != nSignalPerBkgrd)) {
1445       AliInfo(Form("the number of signal events per background event "
1446                    "will be changed from %d to %d for stream %d", 
1447                    nSignalPerBkgrd, result, iBkgrdFile+1));
1448       nSignalPerBkgrd = result;
1449     }
1450
1451     if (!result) result = nSignalPerBkgrd;
1452     if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1453       AliWarning(Form("not enough background events (%d) for %d signal events "
1454                       "using %d signal per background events for stream %d",
1455                       nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1456     }
1457   }
1458
1459   return result;
1460 }
1461
1462 //_____________________________________________________________________________
1463 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1464 {
1465 // check whether detName is contained in detectors
1466 // if yes, it is removed from detectors
1467
1468   // check if all detectors are selected
1469   if ((detectors.CompareTo("ALL") == 0) ||
1470       detectors.BeginsWith("ALL ") ||
1471       detectors.EndsWith(" ALL") ||
1472       detectors.Contains(" ALL ")) {
1473     detectors = "ALL";
1474     return kTRUE;
1475   }
1476
1477   // search for the given detector
1478   Bool_t result = kFALSE;
1479   if ((detectors.CompareTo(detName) == 0) ||
1480       detectors.BeginsWith(detName+" ") ||
1481       detectors.EndsWith(" "+detName) ||
1482       detectors.Contains(" "+detName+" ")) {
1483     detectors.ReplaceAll(detName, "");
1484     result = kTRUE;
1485   }
1486
1487   // clean up the detectors string
1488   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
1489   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1490   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1491
1492   return result;
1493 }
1494
1495 //_____________________________________________________________________________
1496 Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName) 
1497 {
1498 //
1499 // Steering routine  to convert raw data in directory rawDirectory/ to fake SDigits. 
1500 // These can be used for embedding of MC tracks into RAW data using the standard 
1501 // merging procedure.
1502 //
1503 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1504 //
1505     if (!gAlice) {
1506         AliError("no gAlice object. Restart aliroot and try again.");
1507         return kFALSE;
1508     }
1509     if (gAlice->Modules()->GetEntries() > 0) {
1510         AliError("gAlice was already run. Restart aliroot and try again.");
1511         return kFALSE;
1512     }
1513     
1514     AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1515     StdoutToAliInfo(StderrToAliError(gAlice->Init(fConfigFileName.Data());););
1516 //
1517 //  Initialize CDB     
1518     InitCDB();
1519     //AliCDBManager* man = AliCDBManager::Instance();
1520     //man->SetRun(0); // Should this come from rawdata header ?
1521     
1522     Int_t iDet;
1523     //
1524     // Get the runloader
1525     AliRunLoader* runLoader = AliRunLoader::GetRunLoader();
1526     //
1527     // Open esd file if available
1528     TFile* esdFile = TFile::Open(esdFileName);
1529     Bool_t esdOK = (esdFile != 0);
1530     AliESD* esd = new AliESD;
1531     TTree* treeESD = 0;
1532     if (esdOK) {
1533         treeESD = (TTree*) esdFile->Get("esdTree");
1534         if (!treeESD) {
1535             AliWarning("No ESD tree found");
1536             esdOK = kFALSE;
1537         } else {
1538             treeESD->SetBranchAddress("ESD", &esd);
1539         }
1540     }
1541     //
1542     // Create the RawReader
1543     TString fileName(rawDirectory);
1544     AliRawReader* rawReader = 0x0;
1545     if (fileName.EndsWith("/")) {
1546       rawReader = new AliRawReaderFile(fileName);
1547     } else if (fileName.EndsWith(".root")) {
1548       rawReader = new AliRawReaderRoot(fileName);
1549     } else if (!fileName.IsNull()) {
1550       rawReader = new AliRawReaderDate(fileName);
1551     }
1552 //     if (!fEquipIdMap.IsNull() && fRawReader)
1553 //       fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1554     //
1555     // Get list of detectors
1556     TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1557     //
1558     // Get Header
1559     AliHeader* header = runLoader->GetHeader();
1560     //
1561     TString detStr = fMakeSDigits;
1562     // Event loop
1563     Int_t nev = 0;
1564     while(kTRUE) {
1565         if (!(rawReader->NextEvent())) break;
1566         //
1567         // Detector loop
1568         for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1569             AliModule* det = (AliModule*) detArray->At(iDet);
1570             if (!det || !det->IsActive()) continue;
1571             if (IsSelected(det->GetName(), detStr)) {
1572               AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1573               det->Raw2SDigits(rawReader);
1574               rawReader->Reset();
1575             }
1576         } // detectors
1577
1578
1579         //
1580         //  If ESD information available obtain reconstructed vertex and store in header.
1581         if (esdOK) {
1582             treeESD->GetEvent(nev);
1583             const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1584             Double_t position[3];
1585             esdVertex->GetXYZ(position);
1586             AliGenEventHeader* mcHeader = new  AliGenEventHeader("ESD");
1587             TArrayF mcV;
1588             mcV.Set(3);
1589             for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1590             mcHeader->SetPrimaryVertex(mcV);
1591             header->Reset(0,nev);
1592             header->SetGenEventHeader(mcHeader);
1593             printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1594         }
1595         nev++;
1596 //
1597 //      Finish the event
1598         runLoader->TreeE()->Fill();
1599         runLoader->SetNextEvent();
1600     } // events
1601  
1602     delete rawReader;
1603 //
1604 //  Finish the run 
1605     runLoader->CdGAFile();
1606     runLoader->WriteHeader("OVERWRITE");
1607     runLoader->WriteRunLoader();
1608
1609     return kTRUE;
1610 }
1611
1612 //_____________________________________________________________________________
1613 Int_t AliSimulation::GetDetIndex(const char* detector)
1614 {
1615   // return the detector index corresponding to detector
1616   Int_t index = -1 ; 
1617   for (index = 0; index < fgkNDetectors ; index++) {
1618     if ( strcmp(detector, fgkDetectorName[index]) == 0 )
1619           break ; 
1620   }     
1621   return index ; 
1622 }
1623
1624 //_____________________________________________________________________________
1625 Bool_t AliSimulation::RunHLT()
1626 {
1627   // Run the HLT simulation
1628   // HLT simulation is implemented in HLT/sim/AliHLTSimulation
1629   // Disabled if fRunHLT is empty, default vaule is "default".
1630   // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
1631   // The default simulation depends on the HLT component libraries and their
1632   // corresponding agents which define components and chains to run. See
1633   // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
1634   // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
1635   //
1636   // The libraries to be loaded can be specified as an option.
1637   // <pre>
1638   // AliSimulation sim;
1639   // sim.SetRunHLT("libAliHLTSample.so");
1640   // </pre>
1641   // will only load <tt>libAliHLTSample.so</tt>
1642
1643   // Other available options:
1644   // \li loglevel=<i>level</i> <br>
1645   //     logging level for this processing
1646   // \li alilog=off
1647   //     disable redirection of log messages to AliLog class
1648   // \li config=<i>macro</i>
1649   //     configuration macro
1650   // \li chains=<i>configuration</i>
1651   //     comma separated list of configurations to be run during simulation
1652   // \li rawfile=<i>file</i>
1653   //     source for the RawReader to be created, the default is <i>./</i> if
1654   //     raw data is simulated
1655
1656   int iResult=0;
1657   AliRunLoader* pRunLoader = LoadRun("READ");
1658   if (!pRunLoader) return kFALSE;
1659
1660   // initialize CDB storage, run number, set CDB lock
1661   InitCDB();
1662   if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1663   SetCDBLock();
1664   
1665   // load the library dynamically
1666   gSystem->Load(ALIHLTSIMULATION_LIBRARY);
1667
1668   // check for the library version
1669   AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
1670   if (!fctVersion) {
1671     AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
1672     return kFALSE;
1673   }
1674   if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
1675     AliError(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
1676     return kFALSE;
1677   }
1678
1679   // print compile info
1680   typedef void (*CompileInfo)( char*& date, char*& time);
1681   CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
1682   if (fctInfo) {
1683     Char_t* date=0;
1684     Char_t* time=0;
1685     (*fctInfo)(date,time);
1686     if (!date) {date=new Char_t[8]; strcpy(date,"unknown");}
1687     if (!time) {time=new Char_t[8]; strcpy(time,"unknown");}
1688     AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
1689     delete date;
1690     delete time;
1691   } else {
1692     AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
1693   }
1694
1695   // create instance of the HLT simulation
1696   AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
1697   AliHLTSimulation* pHLT=NULL;
1698   if (fctCreate==NULL || (pHLT=(fctCreate()))==NULL) {
1699     AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
1700     return kFALSE;    
1701   }
1702
1703   // init the HLT simulation
1704   TString options;
1705   if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
1706   TString detStr = fWriteRawData;
1707   if (!IsSelected("HLT", detStr)) {
1708     options+=" writerawfiles=";
1709   } else {
1710     options+=" writerawfiles=HLT";
1711   }
1712
1713   if (!detStr.IsNull() && !options.Contains("rawfile=")) {
1714     // as a matter of fact, HLT will run reconstruction and needs the RawReader
1715     // in order to get detector data. By default, RawReaderFile is used to read
1716     // the already simulated ddl files. Date and Root files from the raw data
1717     // are generated after the HLT simulation.
1718     options+=" rawfile=./";
1719   }
1720
1721   AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
1722   if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
1723     AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
1724   } else {
1725     // run the HLT simulation
1726     AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
1727     if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
1728       AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
1729     }
1730   }
1731
1732   // delete the instance
1733   AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
1734   if (fctDelete==NULL || fctDelete(pHLT)<0) {
1735     AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
1736   }
1737   pHLT=NULL;
1738
1739   return iResult>=0?kTRUE:kFALSE;
1740 }
1741
1742 //_____________________________________________________________________________
1743 Bool_t AliSimulation::RunQA()
1744 {
1745         // run the QA on summable hits, digits or digits
1746         
1747   if(!gAlice) return kFALSE;
1748         fQASteer->SetRunLoader(AliRunLoader::GetRunLoader()) ;
1749
1750         TString detectorsw("") ;  
1751         Bool_t rv = kTRUE ; 
1752   fQASteer->SetEventSpecie(fEventSpecie) ;
1753         detectorsw = fQASteer->Run(fQADetectors.Data()) ; 
1754         if ( detectorsw.IsNull() ) 
1755                 rv = kFALSE ; 
1756         return rv ; 
1757 }
1758
1759 //_____________________________________________________________________________
1760 Bool_t AliSimulation::SetRunQA(TString detAndAction) 
1761 {
1762         // Allows to run QA for a selected set of detectors
1763         // and a selected set of tasks among HITS, SDIGITS and DIGITS
1764         // all selected detectors run the same selected tasks
1765         
1766         if (!detAndAction.Contains(":")) {
1767                 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
1768                 fRunQA = kFALSE ;
1769                 return kFALSE ;                 
1770         }
1771         Int_t colon = detAndAction.Index(":") ; 
1772         fQADetectors = detAndAction(0, colon) ; 
1773         if (fQADetectors.Contains("ALL") )
1774                 fQADetectors = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ; 
1775                 fQATasks   = detAndAction(colon+1, detAndAction.Sizeof() ) ; 
1776         if (fQATasks.Contains("ALL") ) {
1777                 fQATasks = Form("%d %d %d", AliQA::kHITS, AliQA::kSDIGITS, AliQA::kDIGITS) ; 
1778         } else {
1779                 fQATasks.ToUpper() ; 
1780                 TString tempo("") ; 
1781                 if ( fQATasks.Contains("HIT") ) 
1782                         tempo = Form("%d ", AliQA::kHITS) ; 
1783                 if ( fQATasks.Contains("SDIGIT") ) 
1784                         tempo += Form("%d ", AliQA::kSDIGITS) ; 
1785                 if ( fQATasks.Contains("DIGIT") ) 
1786                         tempo += Form("%d ", AliQA::kDIGITS) ; 
1787                 fQATasks = tempo ; 
1788                 if (fQATasks.IsNull()) {
1789                         AliInfo("No QA requested\n")  ;
1790                         fRunQA = kFALSE ;
1791                         return kTRUE ; 
1792                 }
1793         }       
1794         TString tempo(fQATasks) ; 
1795     tempo.ReplaceAll(Form("%d", AliQA::kHITS), AliQA::GetTaskName(AliQA::kHITS))        ;
1796     tempo.ReplaceAll(Form("%d", AliQA::kSDIGITS), AliQA::GetTaskName(AliQA::kSDIGITS)) ;        
1797     tempo.ReplaceAll(Form("%d", AliQA::kDIGITS), AliQA::GetTaskName(AliQA::kDIGITS)) ;  
1798         AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;  
1799         fRunQA = kTRUE ;
1800         fQASteer->SetActiveDetectors(fQADetectors) ; 
1801         fQASteer->SetTasks(fQATasks) ; 
1802   for (Int_t det = 0 ; det < AliQA::kNDET ; det++) 
1803     fQASteer->SetWriteExpert(AliQA::DETECTORINDEX_t(det)) ;
1804   
1805         return kTRUE; 
1806
1807
1808 //_____________________________________________________________________________
1809 void AliSimulation::ProcessEnvironmentVars()
1810 {
1811 // Extract run number and random generator seed from env variables
1812
1813     AliInfo("Processing environment variables");
1814     
1815     // Random Number seed
1816     
1817     // first check that seed is not already set
1818     if (fSeed == 0) {
1819         if (gSystem->Getenv("CONFIG_SEED")) {
1820                 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
1821         }
1822     } else {
1823         if (gSystem->Getenv("CONFIG_SEED")) {
1824                 AliInfo(Form("Seed for random number generation already set (%d)"
1825                              ": CONFIG_SEED variable ignored!", fSeed));
1826         }
1827     }
1828    
1829     AliInfo(Form("Seed for random number generation = %d ", fSeed)); 
1830
1831     // Run Number
1832     
1833     // first check that run number is not already set
1834     if(fRun < 0) {    
1835         if (gSystem->Getenv("DC_RUN")) {
1836                 fRun = atoi(gSystem->Getenv("DC_RUN"));
1837         }
1838     } else {
1839         if (gSystem->Getenv("DC_RUN")) {
1840                 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
1841         }
1842     }
1843     
1844     AliInfo(Form("Run number = %d", fRun)); 
1845 }
1846
1847 //---------------------------------------------------------------------
1848 void AliSimulation::WriteGRPEntry()
1849 {
1850   // Get the necessary information from galice (generator, trigger etc) and
1851   // write a GRP entry corresponding to the settings in the Config.C used
1852   // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
1853
1854
1855   AliInfo("Writing global run parameters entry into the OCDB");
1856
1857   AliGRPObject* grpObj = new AliGRPObject();
1858
1859   grpObj->SetRunType("PHYSICS");
1860   grpObj->SetTimeStart(0);
1861   grpObj->SetTimeEnd(9999);
1862
1863   const AliGenerator *gen = gAlice->Generator();
1864   if (gen) {
1865     grpObj->SetBeamEnergy(gen->GetEnergyCMS());
1866     TString projectile;
1867     Int_t a,z;
1868     gen->GetProjectile(projectile,a,z);
1869     TString target;
1870     gen->GetTarget(target,a,z);
1871     TString beamType = projectile + "-" + target;
1872     beamType.ReplaceAll(" ","");
1873     if (!beamType.CompareTo("-")) {
1874       grpObj->SetBeamType("UNKNOWN");
1875     }
1876     else {
1877       grpObj->SetBeamType(beamType);
1878       // Heavy ion run, the event specie is set to kHighMult
1879       fEventSpecie = AliRecoParam::kHighMult;
1880       if ((strcmp(beamType,"p-p") == 0) ||
1881           (strcmp(beamType,"p-")  == 0) ||
1882           (strcmp(beamType,"-p")  == 0) ||
1883           (strcmp(beamType,"P-P") == 0) ||
1884           (strcmp(beamType,"P-")  == 0) ||
1885           (strcmp(beamType,"-P")  == 0)) {
1886         // Proton run, the event specie is set to kLowMult
1887         fEventSpecie = AliRecoParam::kLowMult;
1888       } 
1889     }
1890   } else {
1891     AliWarning("Unknown beam type and energy! Setting energy to 0");
1892     grpObj->SetBeamEnergy(0);
1893     grpObj->SetBeamType("UNKNOWN");
1894   }
1895
1896   UInt_t detectorPattern  = 0;
1897   Int_t nDets = 0;
1898   TObjArray *detArray = gAlice->Detectors();
1899   for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
1900     if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
1901       detectorPattern |= (1 << iDet);
1902       nDets++;
1903     }
1904   }
1905   // CTP
1906   if (!fMakeTrigger.IsNull() || strcmp(gAlice->GetTriggerDescriptor(),""))
1907     detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
1908
1909   // HLT
1910   if (!fRunHLT.IsNull())
1911     detectorPattern |= (1 << AliDAQ::kHLTId);
1912
1913   grpObj->SetNumberOfDetectors((Char_t)nDets);
1914   grpObj->SetDetectorMask((Int_t)detectorPattern);
1915   grpObj->SetLHCPeriod("LHC08c");
1916   grpObj->SetLHCState("STABLE_BEAMS");
1917   grpObj->SetLHCLuminosity(0,(AliGRPObject::Stats)0);
1918   grpObj->SetBeamIntensity(0,(AliGRPObject::Stats)0);
1919
1920   AliMagF *field = gAlice->Field();
1921   Float_t solenoidField = TMath::Abs(field->SolenoidField());
1922   Float_t factor = field->Factor();
1923   Float_t l3current = TMath::Abs(factor)*solenoidField*30000./5.;
1924   grpObj->SetL3Current(l3current,(AliGRPObject::Stats)0);
1925   
1926   if (factor > 0) {
1927     grpObj->SetL3Polarity(0);
1928     grpObj->SetDipolePolarity(0);
1929   }
1930   else {
1931     grpObj->SetL3Polarity(1);
1932     grpObj->SetDipolePolarity(1);
1933   }
1934
1935   if (TMath::Abs(factor) != 0)
1936     grpObj->SetDipoleCurrent(6000,(AliGRPObject::Stats)0);
1937   else 
1938     grpObj->SetDipoleCurrent(0,(AliGRPObject::Stats)0);
1939
1940   grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
1941   
1942   //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
1943
1944   // Now store the entry in OCDB
1945   AliCDBManager* man = AliCDBManager::Instance();
1946
1947   AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
1948   AliCDBMetaData *metadata= new AliCDBMetaData();
1949
1950   metadata->SetResponsible("alice-off@cern.ch");
1951   metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
1952  
1953   man->Put(grpObj,id,metadata);
1954 }
1955
1956