1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for running generation, simulation and digitization //
22 // Hits, sdigits and digits are created for all detectors by typing: //
24 // AliSimulation sim; //
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 //
31 // sim.SetNumberOfEvents(n); //
33 // The name of the configuration file can be passed as argument to the //
34 // AliSimulation constructor or can be specified by //
36 // sim.SetConfigFile("..."); //
38 // The generation of particles and the simulation of detector hits can be //
39 // switched on or off by //
41 // sim.SetRunGeneration(kTRUE); // generation of primary particles //
42 // sim.SetRunSimulation(kFALSE); // but no tracking //
44 // For which detectors sdigits and digits will be created, can be steered //
47 // sim.SetMakeSDigits("ALL"); // make sdigits for all detectors //
48 // sim.SetMakeDigits("ITS TPC"); // make digits only for ITS and TPC //
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. //
55 // The creation of digits from hits instead of from sdigits can be selected //
58 // sim.SetMakeDigitsFromHits("TRD"); //
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. //
64 // Background events can be merged by calling //
66 // sim.MergeWith("background/galice.root", 2); //
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. //
75 // The output of raw data can be switched on by calling //
77 // sim.SetWriteRawData("MUON"); // write raw data for MUON //
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 //
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. //
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 //
98 // sim.SetEventsPerFile("PHOS", "Reconstructed Points", 3); //
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. //
106 ///////////////////////////////////////////////////////////////////////////////
110 #include <TGeoGlobalMagField.h>
111 #include <TGeoManager.h>
112 #include <TObjString.h>
115 #include <TVirtualMC.h>
116 #include <TVirtualMCApplication.h>
119 #include "AliAlignObj.h"
120 #include "AliCDBEntry.h"
121 #include "AliCDBManager.h"
122 #include "AliCDBStorage.h"
123 #include "AliCTPRawData.h"
124 #include "AliCentralTrigger.h"
125 #include "AliCentralTrigger.h"
126 #include "AliCodeTimer.h"
128 #include "AliDigitizer.h"
129 #include "AliESDEvent.h"
130 #include "AliGRPObject.h"
131 #include "AliGenEventHeader.h"
132 #include "AliGenerator.h"
133 #include "AliGeomManager.h"
134 #include "AliHLTSimulation.h"
135 #include "AliHeader.h"
137 #include "AliLegoGenerator.h"
141 #include "AliModule.h"
143 #include "AliRawReaderDate.h"
144 #include "AliRawReaderFile.h"
145 #include "AliRawReaderRoot.h"
147 #include "AliRunDigitizer.h"
148 #include "AliRunLoader.h"
149 #include "AliSimulation.h"
150 #include "AliSysInfo.h"
151 #include "AliVertexGenFile.h"
153 ClassImp(AliSimulation)
155 AliSimulation *AliSimulation::fgInstance = 0;
156 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
158 //_____________________________________________________________________________
159 AliSimulation::AliSimulation(const char* configFileName,
160 const char* name, const char* title) :
163 fRunGeneration(kTRUE),
164 fRunSimulation(kTRUE),
165 fLoadAlignFromCDB(kTRUE),
166 fLoadAlObjsListOfDets("ALL"),
170 fMakeDigitsFromHits(""),
172 fRawDataFileName(""),
173 fDeleteIntermediateFiles(kFALSE),
174 fWriteSelRawData(kFALSE),
175 fStopOnError(kFALSE),
177 fConfigFileName(configFileName),
178 fGAliceFileName("galice.root"),
180 fBkgrdFileNames(NULL),
181 fAlignObjArray(NULL),
182 fUseBkgrdVertex(kTRUE),
183 fRegionOfInterest(kFALSE),
189 fInitCDBCalled(kFALSE),
190 fInitRunNumberCalled(kFALSE),
191 fSetRunNumberFromDataCalled(kFALSE),
192 fEmbeddingFlag(kFALSE),
197 fEventSpecie(AliRecoParam::kDefault),
198 fWriteQAExpertData(kTRUE),
201 fWriteGRPEntry(kTRUE)
203 // create simulation object with default parameters
205 SetGAliceFile("galice.root");
208 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
209 qam->SetActiveDetectors(fQADetectors) ;
210 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
211 qam->SetTasks(fQATasks) ;
214 //_____________________________________________________________________________
215 AliSimulation::~AliSimulation()
219 fEventsPerFile.Delete();
220 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
221 // delete fAlignObjArray; fAlignObjArray=0;
223 if (fBkgrdFileNames) {
224 fBkgrdFileNames->Delete();
225 delete fBkgrdFileNames;
228 fSpecCDBUri.Delete();
229 if (fgInstance==this) fgInstance = 0;
231 AliQAManager::QAManager()->ShowQA() ;
232 AliQAManager::Destroy() ;
233 AliCodeTimer::Instance()->Print();
237 //_____________________________________________________________________________
238 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
240 // set the number of events for one run
245 //_____________________________________________________________________________
246 void AliSimulation::InitQA()
248 // activate a default CDB storage
249 // First check if we have any CDB storage set, because it is used
250 // to retrieve the calibration and alignment constants
252 if (fInitCDBCalled) return;
253 fInitCDBCalled = kTRUE;
255 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
256 qam->SetActiveDetectors(fQADetectors) ;
257 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
258 qam->SetTasks(fQATasks) ;
259 if (fWriteQAExpertData)
260 qam->SetWriteExpert() ;
262 if (qam->IsDefaultStorageSet()) {
263 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
264 AliWarning("Default QA reference storage has been already set !");
265 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
266 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
267 fQARefUri = qam->GetDefaultStorage()->GetURI();
269 if (fQARefUri.Length() > 0) {
270 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
271 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
272 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
274 fQARefUri="local://$ALICE_ROOT/QARef";
275 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
276 AliWarning("Default QA reference storage not yet set !!!!");
277 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
278 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
280 qam->SetDefaultStorage(fQARefUri);
284 //_____________________________________________________________________________
285 void AliSimulation::InitCDB()
287 // activate a default CDB storage
288 // First check if we have any CDB storage set, because it is used
289 // to retrieve the calibration and alignment constants
291 if (fInitCDBCalled) return;
292 fInitCDBCalled = kTRUE;
294 AliCDBManager* man = AliCDBManager::Instance();
295 if (man->IsDefaultStorageSet())
297 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
298 AliWarning("Default CDB storage has been already set !");
299 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
300 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
301 fCDBUri = man->GetDefaultStorage()->GetURI();
304 if (fCDBUri.Length() > 0)
306 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
307 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
308 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
310 fCDBUri="local://$ALICE_ROOT/OCDB";
311 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
312 AliWarning("Default CDB storage not yet set !!!!");
313 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
314 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
317 man->SetDefaultStorage(fCDBUri);
320 // Now activate the detector specific CDB storage locations
321 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
322 TObject* obj = fSpecCDBUri[i];
324 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
325 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
326 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
327 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
332 //_____________________________________________________________________________
333 void AliSimulation::InitRunNumber(){
334 // check run number. If not set, set it to 0 !!!!
336 if (fInitRunNumberCalled) return;
337 fInitRunNumberCalled = kTRUE;
339 AliCDBManager* man = AliCDBManager::Instance();
340 if (man->GetRun() >= 0)
342 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
343 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
347 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
350 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
359 //_____________________________________________________________________________
360 void AliSimulation::SetCDBLock() {
361 // Set CDB lock: from now on it is forbidden to reset the run number
362 // or the default storage or to activate any further storage!
364 AliCDBManager::Instance()->SetLock(1);
367 //_____________________________________________________________________________
368 void AliSimulation::SetDefaultStorage(const char* uri) {
369 // Store the desired default CDB storage location
370 // Activate it later within the Run() method
376 //_____________________________________________________________________________
377 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
378 // Store the desired default CDB storage location
379 // Activate it later within the Run() method
382 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
385 //_____________________________________________________________________________
386 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
387 // Store a detector-specific CDB storage location
388 // Activate it later within the Run() method
390 AliCDBPath aPath(calibType);
391 if(!aPath.IsValid()){
392 AliError(Form("Not a valid path: %s", calibType));
396 TObject* obj = fSpecCDBUri.FindObject(calibType);
397 if (obj) fSpecCDBUri.Remove(obj);
398 fSpecCDBUri.Add(new TNamed(calibType, uri));
402 //_____________________________________________________________________________
403 void AliSimulation::SetRunNumber(Int_t run)
406 // Activate it later within the Run() method
411 //_____________________________________________________________________________
412 void AliSimulation::SetSeed(Int_t seed)
415 // Activate it later within the Run() method
420 //_____________________________________________________________________________
421 Bool_t AliSimulation::SetRunNumberFromData()
423 // Set the CDB manager run number
424 // The run number is retrieved from gAlice
426 if (fSetRunNumberFromDataCalled) return kTRUE;
427 fSetRunNumberFromDataCalled = kTRUE;
429 AliCDBManager* man = AliCDBManager::Instance();
430 Int_t runData = -1, runCDB = -1;
432 AliRunLoader* runLoader = LoadRun("READ");
433 if (!runLoader) return kFALSE;
435 runData = runLoader->GetHeader()->GetRun();
439 runCDB = man->GetRun();
441 if (runCDB != runData) {
442 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
443 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
444 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
445 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
450 man->SetRun(runData);
453 if(man->GetRun() < 0) {
454 AliError("Run number not properly initalized!");
463 //_____________________________________________________________________________
464 void AliSimulation::SetConfigFile(const char* fileName)
466 // set the name of the config file
468 fConfigFileName = fileName;
471 //_____________________________________________________________________________
472 void AliSimulation::SetGAliceFile(const char* fileName)
474 // set the name of the galice file
475 // the path is converted to an absolute one if it is relative
477 fGAliceFileName = fileName;
478 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
479 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
481 fGAliceFileName = absFileName;
482 delete[] absFileName;
485 AliDebug(2, Form("galice file name set to %s", fileName));
488 //_____________________________________________________________________________
489 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
492 // set the number of events per file for the given detector and data type
493 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
495 TNamed* obj = new TNamed(detector, type);
496 obj->SetUniqueID(nEvents);
497 fEventsPerFile.Add(obj);
500 //_____________________________________________________________________________
501 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
503 // Read the alignment objects from CDB.
504 // Each detector is supposed to have the
505 // alignment objects in DET/Align/Data CDB path.
506 // All the detector objects are then collected,
507 // sorted by geometry level (starting from ALIC) and
508 // then applied to the TGeo geometry.
509 // Finally an overlaps check is performed.
511 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
512 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
516 // initialize CDB storage, run number, set CDB lock
518 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
521 Bool_t delRunLoader = kFALSE;
523 runLoader = LoadRun("READ");
524 if (!runLoader) return kFALSE;
525 delRunLoader = kTRUE;
528 // Export ideal geometry
529 if(!gAlice->IsRootGeometry()) AliGeomManager::GetGeometry()->Export("geometry.root");
531 // Load alignment data from CDB and apply to geometry through AliGeomManager
532 if(fLoadAlignFromCDB){
534 TString detStr = fLoadAlObjsListOfDets;
535 TString loadAlObjsListOfDets = "";
537 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
538 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
539 AliModule* det = (AliModule*) detArray->At(iDet);
540 if (!det || !det->IsActive()) continue;
541 if (IsSelected(det->GetName(), detStr)) {
542 //add det to list of dets to be aligned from CDB
543 loadAlObjsListOfDets += det->GetName();
544 loadAlObjsListOfDets += " ";
546 } // end loop over detectors
547 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
548 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
550 // Check if the array with alignment objects was
551 // provided by the user. If yes, apply the objects
552 // to the present TGeo geometry
553 if (fAlignObjArray) {
554 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
555 AliError("The misalignment of one or more volumes failed!"
556 "Compare the list of simulated detectors and the list of detector alignment data!");
557 if (delRunLoader) delete runLoader;
563 // Update the internal geometry of modules (ITS needs it)
564 TString detStr = fLoadAlObjsListOfDets;
565 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
566 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
568 AliModule* det = (AliModule*) detArray->At(iDet);
569 if (!det || !det->IsActive()) continue;
570 if (IsSelected(det->GetName(), detStr)) {
571 det->UpdateInternalGeometry();
573 } // end loop over detectors
576 if (delRunLoader) delete runLoader;
581 //_____________________________________________________________________________
582 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
584 // add a file with background events for merging
586 TObjString* fileNameStr = new TObjString(fileName);
587 fileNameStr->SetUniqueID(nSignalPerBkgrd);
588 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
589 fBkgrdFileNames->Add(fileNameStr);
592 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
594 // add a file with background events for embeddin
595 MergeWith(fileName, nSignalPerBkgrd);
596 fEmbeddingFlag = kTRUE;
599 //_____________________________________________________________________________
600 Bool_t AliSimulation::Run(Int_t nEvents)
602 // run the generation, simulation and digitization
605 AliCodeTimerAuto("",0)
606 AliSysInfo::AddStamp("Start_Run");
608 // Load run number and seed from environmental vars
609 ProcessEnvironmentVars();
610 AliSysInfo::AddStamp("ProcessEnvironmentVars");
612 gRandom->SetSeed(fSeed);
614 if (nEvents > 0) fNEvents = nEvents;
616 // create and setup the HLT instance
617 if (!fRunHLT.IsNull() && !CreateHLT()) {
618 if (fStopOnError) return kFALSE;
623 // generation and simulation -> hits
624 if (fRunGeneration) {
625 if (!RunSimulation()) if (fStopOnError) return kFALSE;
627 AliSysInfo::AddStamp("RunSimulation");
629 // initialize CDB storage from external environment
630 // (either CDB manager or AliSimulation setters),
631 // if not already done in RunSimulation()
633 AliSysInfo::AddStamp("InitCDB");
635 // Set run number in CDBManager from data
636 // From this point on the run number must be always loaded from data!
637 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
639 // Set CDB lock: from now on it is forbidden to reset the run number
640 // or the default storage or to activate any further storage!
643 // If RunSimulation was not called, load the geometry and misalign it
644 if (!AliGeomManager::GetGeometry()) {
645 // Initialize the geometry manager
646 AliGeomManager::LoadGeometry("geometry.root");
647 AliSysInfo::AddStamp("GetGeometry");
650 // // Check that the consistency of symbolic names for the activated subdetectors
651 // // in the geometry loaded by AliGeomManager
652 // AliRunLoader* runLoader = LoadRun("READ");
653 // if (!runLoader) return kFALSE;
655 // TString detsToBeChecked = "";
656 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
657 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
658 // AliModule* det = (AliModule*) detArray->At(iDet);
659 // if (!det || !det->IsActive()) continue;
660 // detsToBeChecked += det->GetName();
661 // detsToBeChecked += " ";
662 // } // end loop over detectors
663 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
664 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
665 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
667 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
669 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
671 AliSysInfo::AddStamp("MissalignGeometry");
674 // hits -> summable digits
675 AliSysInfo::AddStamp("Start_sdigitization");
676 if (!fMakeSDigits.IsNull()) {
677 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
680 AliSysInfo::AddStamp("Stop_sdigitization");
682 AliSysInfo::AddStamp("Start_digitization");
683 // summable digits -> digits
684 if (!fMakeDigits.IsNull()) {
685 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
686 if (fStopOnError) return kFALSE;
689 AliSysInfo::AddStamp("Stop_digitization");
694 if (!fMakeDigitsFromHits.IsNull()) {
695 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
696 AliWarning(Form("Merging and direct creation of digits from hits "
697 "was selected for some detectors. "
698 "No merging will be done for the following detectors: %s",
699 fMakeDigitsFromHits.Data()));
701 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
702 if (fStopOnError) return kFALSE;
706 AliSysInfo::AddStamp("Hits2Digits");
710 if (!RunTrigger(fMakeTrigger,fMakeDigits)) {
711 if (fStopOnError) return kFALSE;
714 AliSysInfo::AddStamp("RunTrigger");
717 // digits -> raw data
718 if (!fWriteRawData.IsNull()) {
719 if (!WriteRawData(fWriteRawData, fRawDataFileName,
720 fDeleteIntermediateFiles,fWriteSelRawData)) {
721 if (fStopOnError) return kFALSE;
725 AliSysInfo::AddStamp("WriteRaw");
727 // run HLT simulation on simulated digit data if raw data is not
728 // simulated, otherwise its called as part of WriteRawData
729 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
731 if (fStopOnError) return kFALSE;
735 AliSysInfo::AddStamp("RunHLT");
739 Bool_t rv = RunQA() ;
745 AliSysInfo::AddStamp("RunQA");
747 // Cleanup of CDB manager: cache and active storages!
748 AliCDBManager::Instance()->ClearCache();
753 //_______________________________________________________________________
754 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
755 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
756 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
759 // Generates lego plots of:
760 // - radiation length map phi vs theta
761 // - radiation length map phi vs eta
762 // - interaction length map
763 // - g/cm2 length map
765 // ntheta bins in theta, eta
766 // themin minimum angle in theta (degrees)
767 // themax maximum angle in theta (degrees)
769 // phimin minimum angle in phi (degrees)
770 // phimax maximum angle in phi (degrees)
771 // rmin minimum radius
772 // rmax maximum radius
775 // The number of events generated = ntheta*nphi
776 // run input parameters in macro setup (default="Config.C")
778 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
781 <img src="picts/AliRunLego1.gif">
786 <img src="picts/AliRunLego2.gif">
791 <img src="picts/AliRunLego3.gif">
796 // run the generation and simulation
798 AliCodeTimerAuto("",0)
800 // initialize CDB storage and run number from external environment
801 // (either CDB manager or AliSimulation setters)
807 AliError("no gAlice object. Restart aliroot and try again.");
810 if (gAlice->Modules()->GetEntries() > 0) {
811 AliError("gAlice was already run. Restart aliroot and try again.");
815 AliInfo(Form("initializing gAlice with config file %s",
816 fConfigFileName.Data()));
819 if (nev == -1) nev = nc1 * nc2;
821 // check if initialisation has been done
822 // If runloader has been initialized, set the number of events per file to nc1 * nc2
825 if (!gener) gener = new AliLegoGenerator();
827 // Configure Generator
829 gener->SetRadiusRange(rmin, rmax);
830 gener->SetZMax(zmax);
831 gener->SetCoor1Range(nc1, c1min, c1max);
832 gener->SetCoor2Range(nc2, c2min, c2max);
836 fLego = new AliLego("lego",gener);
838 //__________________________________________________________________________
842 gROOT->LoadMacro(setup);
843 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
845 if(AliCDBManager::Instance()->GetRun() >= 0) {
846 SetRunNumber(AliCDBManager::Instance()->GetRun());
848 AliWarning("Run number not initialized!!");
851 AliRunLoader::Instance()->CdGAFile();
853 AliPDG::AddParticlesToPdgDataBase();
855 gAlice->GetMCApp()->Init();
857 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
860 //Must be here because some MCs (G4) adds detectors here and not in Config.C
861 gAlice->InitLoaders();
862 AliRunLoader::Instance()->MakeTree("E");
865 // Save stuff at the beginning of the file to avoid file corruption
866 AliRunLoader::Instance()->CdGAFile();
869 //Save current generator
870 AliGenerator *gen=gAlice->GetMCApp()->Generator();
871 gAlice->GetMCApp()->ResetGenerator(gener);
872 //Prepare MC for Lego Run
878 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
879 gMC->ProcessRun(nev);
881 // End of this run, close files
883 // Restore current generator
884 gAlice->GetMCApp()->ResetGenerator(gen);
885 // Delete Lego Object
891 //_____________________________________________________________________________
892 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
896 AliCodeTimerAuto("",0)
898 // initialize CDB storage from external environment
899 // (either CDB manager or AliSimulation setters),
900 // if not already done in RunSimulation()
903 // Set run number in CDBManager from data
904 // From this point on the run number must be always loaded from data!
905 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
907 // Set CDB lock: from now on it is forbidden to reset the run number
908 // or the default storage or to activate any further storage!
911 AliRunLoader* runLoader = LoadRun("READ");
912 if (!runLoader) return kFALSE;
913 TString trconfiguration = config;
915 if (trconfiguration.IsNull()) {
916 if (strcmp(gAlice->GetTriggerDescriptor(),"")) {
917 trconfiguration = gAlice->GetTriggerDescriptor();
920 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
923 runLoader->MakeTree( "GG" );
924 AliCentralTrigger* aCTP = runLoader->GetTrigger();
925 // Load Configuration
926 if (!aCTP->LoadConfiguration( trconfiguration ))
930 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
942 //_____________________________________________________________________________
943 Bool_t AliSimulation::WriteTriggerRawData()
945 // Writes the CTP (trigger) DDL raw data
946 // Details of the format are given in the
947 // trigger TDR - pages 134 and 135.
948 AliCTPRawData writer;
954 //_____________________________________________________________________________
955 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
957 // run the generation and simulation
959 AliCodeTimerAuto("",0)
961 // initialize CDB storage and run number from external environment
962 // (either CDB manager or AliSimulation setters)
963 AliSysInfo::AddStamp("RunSimulation_Begin");
965 AliSysInfo::AddStamp("RunSimulation_InitCDB");
968 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
971 AliError("no gAlice object. Restart aliroot and try again.");
974 if (gAlice->Modules()->GetEntries() > 0) {
975 AliError("gAlice was already run. Restart aliroot and try again.");
979 AliInfo(Form("initializing gAlice with config file %s",
980 fConfigFileName.Data()));
983 // Initialize ALICE Simulation run
988 gROOT->LoadMacro(fConfigFileName.Data());
989 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
990 AliSysInfo::AddStamp("RunSimulation_Config");
992 if(AliCDBManager::Instance()->GetRun() >= 0) {
993 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
994 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
996 AliWarning("Run number not initialized!!");
999 AliRunLoader::Instance()->CdGAFile();
1001 AliPDG::AddParticlesToPdgDataBase();
1003 gAlice->GetMCApp()->Init();
1004 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1006 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1007 AliSysInfo::AddStamp("RunSimulation_GetField");
1009 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1010 gAlice->InitLoaders();
1011 AliRunLoader::Instance()->MakeTree("E");
1012 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1013 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1014 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1016 // Save stuff at the beginning of the file to avoid file corruption
1017 AliRunLoader::Instance()->CdGAFile();
1019 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1020 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1021 //___________________________________________________________________________________________
1023 // Get the trigger descriptor string
1024 // Either from AliSimulation or from
1026 if (fMakeTrigger.IsNull()) {
1027 if (strcmp(gAlice->GetTriggerDescriptor(),""))
1028 fMakeTrigger = gAlice->GetTriggerDescriptor();
1031 gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
1032 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1034 // Set run number in CDBManager
1035 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1037 AliRunLoader* runLoader = AliRunLoader::Instance();
1039 AliError(Form("gAlice has no run loader object. "
1040 "Check your config file: %s", fConfigFileName.Data()));
1043 SetGAliceFile(runLoader->GetFileName());
1045 // Misalign geometry
1046 #if ROOT_VERSION_CODE < 331527
1047 AliGeomManager::SetGeometry(gGeoManager);
1049 // Check that the consistency of symbolic names for the activated subdetectors
1050 // in the geometry loaded by AliGeomManager
1051 TString detsToBeChecked = "";
1052 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1053 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1054 AliModule* det = (AliModule*) detArray->At(iDet);
1055 if (!det || !det->IsActive()) continue;
1056 detsToBeChecked += det->GetName();
1057 detsToBeChecked += " ";
1058 } // end loop over detectors
1059 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1060 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1061 MisalignGeometry(runLoader);
1062 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1065 // AliRunLoader* runLoader = AliRunLoader::Instance();
1066 // if (!runLoader) {
1067 // AliError(Form("gAlice has no run loader object. "
1068 // "Check your config file: %s", fConfigFileName.Data()));
1071 // SetGAliceFile(runLoader->GetFileName());
1073 if (!gAlice->GetMCApp()->Generator()) {
1074 AliError(Form("gAlice has no generator object. "
1075 "Check your config file: %s", fConfigFileName.Data()));
1079 // Write GRP entry corresponding to the setting found in Cofig.C
1082 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1084 if (nEvents <= 0) nEvents = fNEvents;
1086 // get vertex from background file in case of merging
1087 if (fUseBkgrdVertex &&
1088 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1089 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1090 const char* fileName = ((TObjString*)
1091 (fBkgrdFileNames->At(0)))->GetName();
1092 AliInfo(Form("The vertex will be taken from the background "
1093 "file %s with nSignalPerBackground = %d",
1094 fileName, signalPerBkgrd));
1095 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1096 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1099 if (!fRunSimulation) {
1100 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1103 // set the number of events per file for given detectors and data types
1104 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1105 if (!fEventsPerFile[i]) continue;
1106 const char* detName = fEventsPerFile[i]->GetName();
1107 const char* typeName = fEventsPerFile[i]->GetTitle();
1108 TString loaderName(detName);
1109 loaderName += "Loader";
1110 AliLoader* loader = runLoader->GetLoader(loaderName);
1112 AliError(Form("RunSimulation", "no loader for %s found\n"
1113 "Number of events per file not set for %s %s",
1114 detName, typeName, detName));
1117 AliDataLoader* dataLoader =
1118 loader->GetDataLoader(typeName);
1120 AliError(Form("no data loader for %s found\n"
1121 "Number of events per file not set for %s %s",
1122 typeName, detName, typeName));
1125 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1126 AliDebug(1, Form("number of events per file set to %d for %s %s",
1127 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1130 AliInfo("running gAlice");
1131 AliSysInfo::AddStamp("Start_ProcessRun");
1133 // Create the Root Tree with one branch per detector
1134 //Hits moved to begin event -> now we are crating separate tree for each event
1136 gMC->ProcessRun(nEvents);
1138 // End of this run, close files
1139 if(nEvents>0) FinishRun();
1141 AliSysInfo::AddStamp("Stop_ProcessRun");
1147 //_____________________________________________________________________________
1148 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1150 // run the digitization and produce summable digits
1151 static Int_t eventNr=0;
1152 AliCodeTimerAuto("",0) ;
1154 // initialize CDB storage, run number, set CDB lock
1156 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1159 AliRunLoader* runLoader = LoadRun();
1160 if (!runLoader) return kFALSE;
1162 TString detStr = detectors;
1163 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1164 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1165 AliModule* det = (AliModule*) detArray->At(iDet);
1166 if (!det || !det->IsActive()) continue;
1167 if (IsSelected(det->GetName(), detStr)) {
1168 AliInfo(Form("creating summable digits for %s", det->GetName()));
1169 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1170 det->Hits2SDigits();
1171 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1172 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1176 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1177 AliError(Form("the following detectors were not found: %s",
1179 if (fStopOnError) return kFALSE;
1188 //_____________________________________________________________________________
1189 Bool_t AliSimulation::RunDigitization(const char* detectors,
1190 const char* excludeDetectors)
1192 // run the digitization and produce digits from sdigits
1194 AliCodeTimerAuto("",0)
1196 // initialize CDB storage, run number, set CDB lock
1198 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1201 delete AliRunLoader::Instance();
1206 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1207 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1208 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
1209 // manager->SetEmbeddingFlag(fEmbeddingFlag);
1210 manager->SetInputStream(0, fGAliceFileName.Data());
1211 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1212 const char* fileName = ((TObjString*)
1213 (fBkgrdFileNames->At(iStream-1)))->GetName();
1214 manager->SetInputStream(iStream, fileName);
1217 TString detStr = detectors;
1218 TString detExcl = excludeDetectors;
1219 manager->GetInputStream(0)->ImportgAlice();
1220 AliRunLoader* runLoader =
1221 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
1222 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1223 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1224 AliModule* det = (AliModule*) detArray->At(iDet);
1225 if (!det || !det->IsActive()) continue;
1226 if (IsSelected(det->GetName(), detStr) &&
1227 !IsSelected(det->GetName(), detExcl)) {
1228 AliDigitizer* digitizer = det->CreateDigitizer(manager);
1231 AliError(Form("no digitizer for %s", det->GetName()));
1232 if (fStopOnError) return kFALSE;
1234 digitizer->SetRegionOfInterest(fRegionOfInterest);
1239 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1240 AliError(Form("the following detectors were not found: %s",
1242 if (fStopOnError) return kFALSE;
1245 if (!manager->GetListOfTasks()->IsEmpty()) {
1246 AliInfo("executing digitization");
1255 //_____________________________________________________________________________
1256 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1258 // run the digitization and produce digits from hits
1260 AliCodeTimerAuto("",0)
1262 // initialize CDB storage, run number, set CDB lock
1264 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1267 AliRunLoader* runLoader = LoadRun("READ");
1268 if (!runLoader) return kFALSE;
1270 TString detStr = detectors;
1271 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1272 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1273 AliModule* det = (AliModule*) detArray->At(iDet);
1274 if (!det || !det->IsActive()) continue;
1275 if (IsSelected(det->GetName(), detStr)) {
1276 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1281 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1282 AliError(Form("the following detectors were not found: %s",
1284 if (fStopOnError) return kFALSE;
1290 //_____________________________________________________________________________
1291 Bool_t AliSimulation::WriteRawData(const char* detectors,
1292 const char* fileName,
1293 Bool_t deleteIntermediateFiles,
1296 // convert the digits to raw data
1297 // First DDL raw data files for the given detectors are created.
1298 // If a file name is given, the DDL files are then converted to a DATE file.
1299 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1301 // If the file name has the extension ".root", the DATE file is converted
1303 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1304 // 'selrawdata' flag can be used to enable writing of detectors raw data
1305 // accoring to the trigger cluster.
1307 AliCodeTimerAuto("",0)
1308 AliSysInfo::AddStamp("WriteRawData_Start");
1310 TString detStr = detectors;
1311 if (!WriteRawFiles(detStr.Data())) {
1312 if (fStopOnError) return kFALSE;
1314 AliSysInfo::AddStamp("WriteRawFiles");
1316 // run HLT simulation on simulated DDL raw files
1317 // and produce HLT ddl raw files to be included in date/root file
1318 // bugfix 2009-06-26: the decision whether to write HLT raw data
1319 // is taken in RunHLT. Here HLT always needs to be run in order to
1320 // create HLT digits, unless its switched off. This is due to the
1321 // special placement of the HLT between the generation of DDL files
1322 // and conversion to DATE/Root file.
1323 detStr.ReplaceAll("HLT", "");
1324 if (!fRunHLT.IsNull()) {
1326 if (fStopOnError) return kFALSE;
1329 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1331 TString dateFileName(fileName);
1332 if (!dateFileName.IsNull()) {
1333 Bool_t rootOutput = dateFileName.EndsWith(".root");
1334 if (rootOutput) dateFileName += ".date";
1335 TString selDateFileName;
1337 selDateFileName = "selected.";
1338 selDateFileName+= dateFileName;
1340 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1341 if (fStopOnError) return kFALSE;
1343 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1344 if (deleteIntermediateFiles) {
1345 AliRunLoader* runLoader = LoadRun("READ");
1346 if (runLoader) for (Int_t iEvent = 0;
1347 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1349 sprintf(command, "rm -r raw%d", iEvent);
1350 gSystem->Exec(command);
1355 if (!ConvertDateToRoot(dateFileName, fileName)) {
1356 if (fStopOnError) return kFALSE;
1358 AliSysInfo::AddStamp("ConvertDateToRoot");
1359 if (deleteIntermediateFiles) {
1360 gSystem->Unlink(dateFileName);
1363 TString selFileName = "selected.";
1364 selFileName += fileName;
1365 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1366 if (fStopOnError) return kFALSE;
1368 if (deleteIntermediateFiles) {
1369 gSystem->Unlink(selDateFileName);
1378 //_____________________________________________________________________________
1379 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1381 // convert the digits to raw data DDL files
1383 AliCodeTimerAuto("",0)
1385 AliRunLoader* runLoader = LoadRun("READ");
1386 if (!runLoader) return kFALSE;
1388 // write raw data to DDL files
1389 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1390 AliInfo(Form("processing event %d", iEvent));
1391 runLoader->GetEvent(iEvent);
1392 TString baseDir = gSystem->WorkingDirectory();
1394 sprintf(dirName, "raw%d", iEvent);
1395 gSystem->MakeDirectory(dirName);
1396 if (!gSystem->ChangeDirectory(dirName)) {
1397 AliError(Form("couldn't change to directory %s", dirName));
1398 if (fStopOnError) return kFALSE; else continue;
1401 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1404 TString detStr = detectors;
1405 if (IsSelected("HLT", detStr)) {
1406 // Do nothing. "HLT" will be removed from detStr and HLT raw
1407 // data files are generated in RunHLT.
1410 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1411 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1412 AliModule* det = (AliModule*) detArray->At(iDet);
1413 if (!det || !det->IsActive()) continue;
1414 if (IsSelected(det->GetName(), detStr)) {
1415 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1420 if (!WriteTriggerRawData())
1421 if (fStopOnError) return kFALSE;
1423 gSystem->ChangeDirectory(baseDir);
1424 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1425 AliError(Form("the following detectors were not found: %s",
1427 if (fStopOnError) return kFALSE;
1436 //_____________________________________________________________________________
1437 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1438 const char* selDateFileName)
1440 // convert raw data DDL files to a DATE file with the program "dateStream"
1441 // The second argument is not empty when the user decides to write
1442 // the detectors raw data according to the trigger cluster.
1444 AliCodeTimerAuto("",0)
1446 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1448 AliError("the program dateStream was not found");
1449 if (fStopOnError) return kFALSE;
1454 AliRunLoader* runLoader = LoadRun("READ");
1455 if (!runLoader) return kFALSE;
1457 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1458 Bool_t selrawdata = kFALSE;
1459 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1462 // Note the option -s. It is used in order to avoid
1463 // the generation of SOR/EOR events.
1464 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1465 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1466 FILE* pipe = gSystem->OpenPipe(command, "w");
1469 AliError(Form("Cannot execute command: %s",command));
1473 Int_t selEvents = 0;
1474 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1476 UInt_t detectorPattern = 0;
1477 runLoader->GetEvent(iEvent);
1478 if (!runLoader->LoadTrigger()) {
1479 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1480 detectorPattern = aCTP->GetClusterMask();
1481 // Check if the event was triggered by CTP
1483 if (aCTP->GetClassMask()) selEvents++;
1487 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1489 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1490 selrawdata = kFALSE;
1494 fprintf(pipe, "GDC DetectorPattern %u\n", detectorPattern);
1498 // loop over detectors and DDLs
1499 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1500 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1502 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1503 Int_t ldcID = Int_t(ldc + 0.0001);
1504 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1506 char rawFileName[256];
1507 sprintf(rawFileName, "raw%d/%s",
1508 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1510 // check existence and size of raw data file
1511 FILE* file = fopen(rawFileName, "rb");
1512 if (!file) continue;
1513 fseek(file, 0, SEEK_END);
1514 unsigned long size = ftell(file);
1516 if (!size) continue;
1518 if (ldcID != prevLDC) {
1519 fprintf(pipe, " LDC Id %d\n", ldcID);
1522 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1527 Int_t result = gSystem->ClosePipe(pipe);
1529 if (!(selrawdata && selEvents > 0)) {
1531 return (result == 0);
1534 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1536 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1537 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1538 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1540 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1542 // Get the trigger decision and cluster
1543 UInt_t detectorPattern = 0;
1545 runLoader->GetEvent(iEvent);
1546 if (!runLoader->LoadTrigger()) {
1547 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1548 if (aCTP->GetClassMask() == 0) continue;
1549 detectorPattern = aCTP->GetClusterMask();
1550 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1551 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1554 fprintf(pipe2, "GDC DetectorPattern %u\n", detectorPattern);
1558 // loop over detectors and DDLs
1559 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1560 // Write only raw data from detectors that
1561 // are contained in the trigger cluster(s)
1562 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1564 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1566 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1567 Int_t ldcID = Int_t(ldc + 0.0001);
1568 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1570 char rawFileName[256];
1571 sprintf(rawFileName, "raw%d/%s",
1572 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1574 // check existence and size of raw data file
1575 FILE* file = fopen(rawFileName, "rb");
1576 if (!file) continue;
1577 fseek(file, 0, SEEK_END);
1578 unsigned long size = ftell(file);
1580 if (!size) continue;
1582 if (ldcID != prevLDC) {
1583 fprintf(pipe2, " LDC Id %d\n", ldcID);
1586 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1591 Int_t result2 = gSystem->ClosePipe(pipe2);
1594 return ((result == 0) && (result2 == 0));
1597 //_____________________________________________________________________________
1598 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1599 const char* rootFileName)
1601 // convert a DATE file to a root file with the program "alimdc"
1604 const Int_t kDBSize = 2000000000;
1605 const Int_t kTagDBSize = 1000000000;
1606 const Bool_t kFilter = kFALSE;
1607 const Int_t kCompression = 1;
1609 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1611 AliError("the program alimdc was not found");
1612 if (fStopOnError) return kFALSE;
1617 AliInfo(Form("converting DATE file %s to root file %s",
1618 dateFileName, rootFileName));
1620 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1621 const char* tagDBFS = "/tmp/mdc1/tags";
1623 // User defined file system locations
1624 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1625 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1626 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1627 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1628 if (gSystem->Getenv("ALIMDC_TAGDB"))
1629 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1631 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1632 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1633 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1635 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1636 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1637 gSystem->Exec(Form("mkdir %s",tagDBFS));
1639 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1640 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1641 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1643 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1644 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1645 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1647 return (result == 0);
1651 //_____________________________________________________________________________
1652 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1654 // delete existing run loaders, open a new one and load gAlice
1656 delete AliRunLoader::Instance();
1657 AliRunLoader* runLoader =
1658 AliRunLoader::Open(fGAliceFileName.Data(),
1659 AliConfig::GetDefaultEventFolderName(), mode);
1661 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1664 runLoader->LoadgAlice();
1665 runLoader->LoadHeader();
1666 gAlice = runLoader->GetAliRun();
1668 AliError(Form("no gAlice object found in file %s",
1669 fGAliceFileName.Data()));
1675 //_____________________________________________________________________________
1676 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1678 // get or calculate the number of signal events per background event
1680 if (!fBkgrdFileNames) return 1;
1681 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1682 if (nBkgrdFiles == 0) return 1;
1684 // get the number of signal events
1686 AliRunLoader* runLoader =
1687 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1688 if (!runLoader) return 1;
1690 nEvents = runLoader->GetNumberOfEvents();
1695 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1696 // get the number of background events
1697 const char* fileName = ((TObjString*)
1698 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1699 AliRunLoader* runLoader =
1700 AliRunLoader::Open(fileName, "BKGRD");
1701 if (!runLoader) continue;
1702 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1705 // get or calculate the number of signal per background events
1706 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1707 if (nSignalPerBkgrd <= 0) {
1708 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1709 } else if (result && (result != nSignalPerBkgrd)) {
1710 AliInfo(Form("the number of signal events per background event "
1711 "will be changed from %d to %d for stream %d",
1712 nSignalPerBkgrd, result, iBkgrdFile+1));
1713 nSignalPerBkgrd = result;
1716 if (!result) result = nSignalPerBkgrd;
1717 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1718 AliWarning(Form("not enough background events (%d) for %d signal events "
1719 "using %d signal per background events for stream %d",
1720 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1727 //_____________________________________________________________________________
1728 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1730 // check whether detName is contained in detectors
1731 // if yes, it is removed from detectors
1733 // check if all detectors are selected
1734 if ((detectors.CompareTo("ALL") == 0) ||
1735 detectors.BeginsWith("ALL ") ||
1736 detectors.EndsWith(" ALL") ||
1737 detectors.Contains(" ALL ")) {
1742 // search for the given detector
1743 Bool_t result = kFALSE;
1744 if ((detectors.CompareTo(detName) == 0) ||
1745 detectors.BeginsWith(detName+" ") ||
1746 detectors.EndsWith(" "+detName) ||
1747 detectors.Contains(" "+detName+" ")) {
1748 detectors.ReplaceAll(detName, "");
1752 // clean up the detectors string
1753 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1754 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1755 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1760 //_____________________________________________________________________________
1761 Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName)
1764 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1765 // These can be used for embedding of MC tracks into RAW data using the standard
1766 // merging procedure.
1768 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1771 AliError("no gAlice object. Restart aliroot and try again.");
1774 if (gAlice->Modules()->GetEntries() > 0) {
1775 AliError("gAlice was already run. Restart aliroot and try again.");
1779 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1783 gROOT->LoadMacro(fConfigFileName.Data());
1784 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1786 if(AliCDBManager::Instance()->GetRun() >= 0) {
1787 SetRunNumber(AliCDBManager::Instance()->GetRun());
1789 AliWarning("Run number not initialized!!");
1792 AliRunLoader::Instance()->CdGAFile();
1794 AliPDG::AddParticlesToPdgDataBase();
1796 gAlice->GetMCApp()->Init();
1798 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1800 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1801 gAlice->InitLoaders();
1802 AliRunLoader::Instance()->MakeTree("E");
1803 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1804 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1805 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1807 // Save stuff at the beginning of the file to avoid file corruption
1808 AliRunLoader::Instance()->CdGAFile();
1813 //AliCDBManager* man = AliCDBManager::Instance();
1814 //man->SetRun(0); // Should this come from rawdata header ?
1818 // Get the runloader
1819 AliRunLoader* runLoader = AliRunLoader::Instance();
1821 // Open esd file if available
1822 TFile* esdFile = TFile::Open(esdFileName);
1824 AliESDEvent* esd = new AliESDEvent();
1825 esdFile->GetObject("esdTree", treeESD);
1826 if (treeESD) esd->ReadFromTree(treeESD);
1829 // Create the RawReader
1830 TString fileName(rawDirectory);
1831 AliRawReader* rawReader = 0x0;
1832 if (fileName.EndsWith("/")) {
1833 rawReader = new AliRawReaderFile(fileName);
1834 } else if (fileName.EndsWith(".root")) {
1835 rawReader = new AliRawReaderRoot(fileName);
1836 } else if (!fileName.IsNull()) {
1837 rawReader = new AliRawReaderDate(fileName);
1839 // if (!fEquipIdMap.IsNull() && fRawReader)
1840 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1842 // Get list of detectors
1843 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1846 AliHeader* header = runLoader->GetHeader();
1848 TString detStr = fMakeSDigits;
1852 if (!(rawReader->NextEvent())) break;
1855 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1856 AliModule* det = (AliModule*) detArray->At(iDet);
1857 if (!det || !det->IsActive()) continue;
1858 if (IsSelected(det->GetName(), detStr)) {
1859 AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1860 det->Raw2SDigits(rawReader);
1867 // If ESD information available obtain reconstructed vertex and store in header.
1869 treeESD->GetEvent(nev);
1870 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1871 Double_t position[3];
1872 esdVertex->GetXYZ(position);
1873 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1876 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1877 mcHeader->SetPrimaryVertex(mcV);
1878 header->Reset(0,nev);
1879 header->SetGenEventHeader(mcHeader);
1880 printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1885 runLoader->TreeE()->Fill();
1886 runLoader->SetNextEvent();
1892 runLoader->CdGAFile();
1893 runLoader->WriteHeader("OVERWRITE");
1894 runLoader->WriteRunLoader();
1899 //_____________________________________________________________________________
1900 void AliSimulation::FinishRun()
1903 // Called at the end of the run.
1908 AliDebug(1, "Finish Lego");
1909 AliRunLoader::Instance()->CdGAFile();
1913 // Clean detector information
1914 TIter next(gAlice->Modules());
1915 AliModule *detector;
1916 while((detector = dynamic_cast<AliModule*>(next()))) {
1917 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
1918 detector->FinishRun();
1921 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
1922 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
1924 // Write AliRun info and all detectors parameters
1925 AliRunLoader::Instance()->CdGAFile();
1926 gAlice->Write(0,TObject::kOverwrite);//write AliRun
1927 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
1929 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
1930 AliRunLoader::Instance()->Synchronize();
1933 //_____________________________________________________________________________
1934 Int_t AliSimulation::GetDetIndex(const char* detector)
1936 // return the detector index corresponding to detector
1938 for (index = 0; index < fgkNDetectors ; index++) {
1939 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
1945 //_____________________________________________________________________________
1946 Bool_t AliSimulation::CreateHLT()
1948 // Init the HLT simulation.
1949 // The function loads the library and creates the instance of AliHLTSimulation.
1950 // the main reason for the decoupled creation is to set the transient OCDB
1951 // objects before the OCDB is locked
1953 // load the library dynamically
1954 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
1956 // check for the library version
1957 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
1959 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
1962 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
1963 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
1966 // print compile info
1967 typedef void (*CompileInfo)( const char*& date, const char*& time);
1968 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
1970 const char* date="";
1971 const char* time="";
1972 (*fctInfo)(date, time);
1973 if (!date) date="unknown";
1974 if (!time) time="unknown";
1975 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
1977 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
1980 // create instance of the HLT simulation
1981 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
1982 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
1983 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
1987 TString specObjects;
1988 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
1989 if (specObjects.Length()>0) specObjects+=" ";
1990 specObjects+=fSpecCDBUri[i]->GetName();
1993 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
1994 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
1995 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2001 //_____________________________________________________________________________
2002 Bool_t AliSimulation::RunHLT()
2004 // Run the HLT simulation
2005 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2006 // Disabled if fRunHLT is empty, default vaule is "default".
2007 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2008 // The default simulation depends on the HLT component libraries and their
2009 // corresponding agents which define components and chains to run. See
2010 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2011 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2013 // The libraries to be loaded can be specified as an option.
2015 // AliSimulation sim;
2016 // sim.SetRunHLT("libAliHLTSample.so");
2018 // will only load <tt>libAliHLTSample.so</tt>
2020 // Other available options:
2021 // \li loglevel=<i>level</i> <br>
2022 // logging level for this processing
2024 // disable redirection of log messages to AliLog class
2025 // \li config=<i>macro</i>
2026 // configuration macro
2027 // \li chains=<i>configuration</i>
2028 // comma separated list of configurations to be run during simulation
2029 // \li rawfile=<i>file</i>
2030 // source for the RawReader to be created, the default is <i>./</i> if
2031 // raw data is simulated
2035 if (!fpHLT && !CreateHLT()) {
2038 AliHLTSimulation* pHLT=fpHLT;
2040 AliRunLoader* pRunLoader = LoadRun("READ");
2041 if (!pRunLoader) return kFALSE;
2043 // initialize CDB storage, run number, set CDB lock
2044 // thats for the case of running HLT simulation without all the other steps
2045 // multiple calls are handled by the function, so we can just call
2047 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2050 // init the HLT simulation
2052 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2053 TString detStr = fWriteRawData;
2054 if (!IsSelected("HLT", detStr)) {
2055 options+=" writerawfiles=";
2057 options+=" writerawfiles=HLT";
2060 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2061 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2062 // in order to get detector data. By default, RawReaderFile is used to read
2063 // the already simulated ddl files. Date and Root files from the raw data
2064 // are generated after the HLT simulation.
2065 options+=" rawfile=./";
2068 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2069 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2070 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2072 // run the HLT simulation
2073 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2074 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2075 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2079 // delete the instance
2080 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2081 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2082 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2086 return iResult>=0?kTRUE:kFALSE;
2089 //_____________________________________________________________________________
2090 Bool_t AliSimulation::RunQA()
2092 // run the QA on summable hits, digits or digits
2094 if(!gAlice) return kFALSE;
2095 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2097 TString detectorsw("") ;
2099 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2100 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2101 if ( detectorsw.IsNull() )
2106 //_____________________________________________________________________________
2107 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2109 // Allows to run QA for a selected set of detectors
2110 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2111 // all selected detectors run the same selected tasks
2113 if (!detAndAction.Contains(":")) {
2114 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2118 Int_t colon = detAndAction.Index(":") ;
2119 fQADetectors = detAndAction(0, colon) ;
2120 if (fQADetectors.Contains("ALL") ){
2121 TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2122 Int_t minus = fQADetectors.Last('-') ;
2123 TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2124 TString toRemove("") ;
2125 while (minus >= 0) {
2126 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
2127 toRemove = toRemove.Strip() ;
2128 toKeep.ReplaceAll(toRemove, "") ;
2129 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
2130 minus = fQADetectors.Last('-') ;
2132 fQADetectors = toKeep ;
2134 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2135 if (fQATasks.Contains("ALL") ) {
2136 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2138 fQATasks.ToUpper() ;
2140 if ( fQATasks.Contains("HIT") )
2141 tempo = Form("%d ", AliQAv1::kHITS) ;
2142 if ( fQATasks.Contains("SDIGIT") )
2143 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2144 if ( fQATasks.Contains("DIGIT") )
2145 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2147 if (fQATasks.IsNull()) {
2148 AliInfo("No QA requested\n") ;
2153 TString tempo(fQATasks) ;
2154 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2155 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2156 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2157 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2159 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2160 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2161 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2162 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2167 //_____________________________________________________________________________
2168 void AliSimulation::ProcessEnvironmentVars()
2170 // Extract run number and random generator seed from env variables
2172 AliInfo("Processing environment variables");
2174 // Random Number seed
2176 // first check that seed is not already set
2178 if (gSystem->Getenv("CONFIG_SEED")) {
2179 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2182 if (gSystem->Getenv("CONFIG_SEED")) {
2183 AliInfo(Form("Seed for random number generation already set (%d)"
2184 ": CONFIG_SEED variable ignored!", fSeed));
2188 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2192 // first check that run number is not already set
2194 if (gSystem->Getenv("DC_RUN")) {
2195 fRun = atoi(gSystem->Getenv("DC_RUN"));
2198 if (gSystem->Getenv("DC_RUN")) {
2199 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2203 AliInfo(Form("Run number = %d", fRun));
2206 //---------------------------------------------------------------------
2207 void AliSimulation::WriteGRPEntry()
2209 // Get the necessary information from galice (generator, trigger etc) and
2210 // write a GRP entry corresponding to the settings in the Config.C used
2211 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2214 AliInfo("Writing global run parameters entry into the OCDB");
2216 AliGRPObject* grpObj = new AliGRPObject();
2218 grpObj->SetRunType("PHYSICS");
2219 grpObj->SetTimeStart(0);
2221 grpObj->SetTimeStart(0);
2222 grpObj->SetTimeEnd(curtime.Convert());
2223 grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2225 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2227 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2230 gen->GetProjectile(projectile,a,z);
2232 gen->GetTarget(target,a,z);
2233 TString beamType = projectile + "-" + target;
2234 beamType.ReplaceAll(" ","");
2235 if (!beamType.CompareTo("-")) {
2236 grpObj->SetBeamType("UNKNOWN");
2239 grpObj->SetBeamType(beamType);
2240 // Heavy ion run, the event specie is set to kHighMult
2241 fEventSpecie = AliRecoParam::kHighMult;
2242 if ((strcmp(beamType,"p-p") == 0) ||
2243 (strcmp(beamType,"p-") == 0) ||
2244 (strcmp(beamType,"-p") == 0) ||
2245 (strcmp(beamType,"P-P") == 0) ||
2246 (strcmp(beamType,"P-") == 0) ||
2247 (strcmp(beamType,"-P") == 0)) {
2248 // Proton run, the event specie is set to kLowMult
2249 fEventSpecie = AliRecoParam::kLowMult;
2253 AliWarning("Unknown beam type and energy! Setting energy to 0");
2254 grpObj->SetBeamEnergy(0);
2255 grpObj->SetBeamType("UNKNOWN");
2258 UInt_t detectorPattern = 0;
2260 TObjArray *detArray = gAlice->Detectors();
2261 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2262 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2263 detectorPattern |= (1 << iDet);
2268 if (!fMakeTrigger.IsNull() || strcmp(gAlice->GetTriggerDescriptor(),""))
2269 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2272 if (!fRunHLT.IsNull())
2273 detectorPattern |= (1 << AliDAQ::kHLTId);
2275 grpObj->SetNumberOfDetectors((Char_t)nDets);
2276 grpObj->SetDetectorMask((Int_t)detectorPattern);
2277 grpObj->SetLHCPeriod("LHC08c");
2278 grpObj->SetLHCState("STABLE_BEAMS");
2279 grpObj->SetLHCLuminosity(0,(AliGRPObject::Stats)0);
2280 grpObj->SetBeamIntensity(0,(AliGRPObject::Stats)0);
2282 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2283 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2285 Float_t factorSol = field ? field->GetFactorSol() : 0;
2286 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2287 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2289 Float_t factorDip = field ? field->GetFactorDip() : 0;
2290 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2292 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2293 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2294 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2295 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2296 grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2297 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2299 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2301 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2303 // Now store the entry in OCDB
2304 AliCDBManager* man = AliCDBManager::Instance();
2306 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2307 AliCDBMetaData *metadata= new AliCDBMetaData();
2309 metadata->SetResponsible("alice-off@cern.ch");
2310 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2312 man->Put(grpObj,id,metadata);