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>
118 #include <TInterpreter.h>
120 #include "AliAlignObj.h"
121 #include "AliCDBEntry.h"
122 #include "AliCDBManager.h"
123 #include "AliGRPManager.h"
124 #include "AliCDBStorage.h"
125 #include "AliCTPRawData.h"
126 #include "AliCentralTrigger.h"
127 #include "AliCentralTrigger.h"
128 #include "AliCodeTimer.h"
130 #include "AliDigitizer.h"
131 #include "AliESDEvent.h"
132 #include "AliGRPObject.h"
133 #include "AliGenEventHeader.h"
134 #include "AliGenerator.h"
135 #include "AliGeomManager.h"
136 #include "AliHLTSimulation.h"
137 #include "AliHeader.h"
139 #include "AliLegoGenerator.h"
143 #include "AliModule.h"
145 #include "AliRawReaderDate.h"
146 #include "AliRawReaderFile.h"
147 #include "AliRawReaderRoot.h"
149 #include "AliDigitizationInput.h"
150 #include "AliRunLoader.h"
151 #include "AliSimulation.h"
152 #include "AliSysInfo.h"
153 #include "AliVertexGenFile.h"
156 ClassImp(AliSimulation)
158 AliSimulation *AliSimulation::fgInstance = 0;
159 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE"
160 // #ifdef MFT_UPGRADE
167 //_____________________________________________________________________________
168 AliSimulation::AliSimulation(const char* configFileName,
169 const char* name, const char* title) :
172 fRunGeneration(kTRUE),
173 fRunSimulation(kTRUE),
174 fLoadAlignFromCDB(kTRUE),
175 fLoadAlObjsListOfDets("ALL"),
179 fMakeDigitsFromHits(""),
181 fRawDataFileName(""),
182 fDeleteIntermediateFiles(kFALSE),
183 fWriteSelRawData(kFALSE),
184 fStopOnError(kFALSE),
186 fConfigFileName(configFileName),
187 fGAliceFileName("galice.root"),
189 fBkgrdFileNames(NULL),
190 fAlignObjArray(NULL),
191 fUseBkgrdVertex(kTRUE),
192 fRegionOfInterest(kFALSE),
198 fInitCDBCalled(kFALSE),
199 fInitRunNumberCalled(kFALSE),
200 fSetRunNumberFromDataCalled(kFALSE),
201 fEmbeddingFlag(kFALSE),
204 fUseVertexFromCDB(0),
205 fUseMagFieldFromGRP(0),
206 fGRPWriteLocation(Form("local://%s", gSystem->pwd())),
207 fUseTimeStampFromCDB(0),
213 fEventSpecie(AliRecoParam::kDefault),
214 fWriteQAExpertData(kTRUE),
218 fWriteGRPEntry(kTRUE)
220 // create simulation object with default parameters
222 SetGAliceFile("galice.root");
225 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
226 qam->SetActiveDetectors(fQADetectors) ;
227 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
228 qam->SetTasks(fQATasks) ;
231 //_____________________________________________________________________________
232 AliSimulation::~AliSimulation()
236 fEventsPerFile.Delete();
237 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
238 // delete fAlignObjArray; fAlignObjArray=0;
240 if (fBkgrdFileNames) {
241 fBkgrdFileNames->Delete();
242 delete fBkgrdFileNames;
245 fSpecCDBUri.Delete();
246 if (fgInstance==this) fgInstance = 0;
248 AliQAManager::QAManager()->ShowQA() ;
249 AliQAManager::Destroy() ;
250 AliCodeTimer::Instance()->Print();
254 //_____________________________________________________________________________
255 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
257 // set the number of events for one run
262 //_____________________________________________________________________________
263 void AliSimulation::InitQA()
265 // activate a default CDB storage
266 // First check if we have any CDB storage set, because it is used
267 // to retrieve the calibration and alignment constants
269 if (fInitCDBCalled) return;
270 fInitCDBCalled = kTRUE;
272 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
273 qam->SetActiveDetectors(fQADetectors) ;
274 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
275 qam->SetTasks(fQATasks) ;
276 if (fWriteQAExpertData)
277 qam->SetWriteExpert() ;
279 if (qam->IsDefaultStorageSet()) {
280 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
281 AliWarning("Default QA reference storage has been already set !");
282 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
283 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
284 fQARefUri = qam->GetDefaultStorage()->GetURI();
286 if (fQARefUri.Length() > 0) {
287 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
288 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
289 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
291 fQARefUri="local://$ALICE_ROOT/QARef";
292 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
293 AliWarning("Default QA reference storage not yet set !!!!");
294 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
295 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
297 qam->SetDefaultStorage(fQARefUri);
301 //_____________________________________________________________________________
302 void AliSimulation::InitCDB()
304 // activate a default CDB storage
305 // First check if we have any CDB storage set, because it is used
306 // to retrieve the calibration and alignment constants
308 if (fInitCDBCalled) return;
309 fInitCDBCalled = kTRUE;
311 AliCDBManager* man = AliCDBManager::Instance();
312 if (man->IsDefaultStorageSet())
314 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
315 AliWarning("Default CDB storage has been already set !");
316 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
317 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
318 fCDBUri = man->GetDefaultStorage()->GetURI();
321 if (fCDBUri.Length() > 0)
323 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
324 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
325 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
327 fCDBUri="local://$ALICE_ROOT/OCDB";
328 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
329 AliWarning("Default CDB storage not yet set !!!!");
330 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
331 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
334 man->SetDefaultStorage(fCDBUri);
337 // Now activate the detector specific CDB storage locations
338 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
339 TObject* obj = fSpecCDBUri[i];
341 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
342 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
343 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
344 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
349 //_____________________________________________________________________________
350 void AliSimulation::InitRunNumber(){
351 // check run number. If not set, set it to 0 !!!!
353 if (fInitRunNumberCalled) return;
354 fInitRunNumberCalled = kTRUE;
356 AliCDBManager* man = AliCDBManager::Instance();
357 if (man->GetRun() >= 0)
359 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
360 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
364 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
367 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
376 //_____________________________________________________________________________
377 void AliSimulation::SetCDBLock() {
378 // Set CDB lock: from now on it is forbidden to reset the run number
379 // or the default storage or to activate any further storage!
381 ULong_t key = AliCDBManager::Instance()->SetLock(1);
385 //_____________________________________________________________________________
386 void AliSimulation::SetDefaultStorage(const char* uri) {
387 // Store the desired default CDB storage location
388 // Activate it later within the Run() method
394 //_____________________________________________________________________________
395 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
396 // Store the desired default CDB storage location
397 // Activate it later within the Run() method
400 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
403 //_____________________________________________________________________________
404 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
405 // Store a detector-specific CDB storage location
406 // Activate it later within the Run() method
408 AliCDBPath aPath(calibType);
409 if(!aPath.IsValid()){
410 AliError(Form("Not a valid path: %s", calibType));
414 TObject* obj = fSpecCDBUri.FindObject(calibType);
415 if (obj) fSpecCDBUri.Remove(obj);
416 fSpecCDBUri.Add(new TNamed(calibType, uri));
420 //_____________________________________________________________________________
421 void AliSimulation::SetRunNumber(Int_t run)
424 // Activate it later within the Run() method
429 //_____________________________________________________________________________
430 void AliSimulation::SetSeed(Int_t seed)
433 // Activate it later within the Run() method
438 //_____________________________________________________________________________
439 Bool_t AliSimulation::SetRunNumberFromData()
441 // Set the CDB manager run number
442 // The run number is retrieved from gAlice
444 if (fSetRunNumberFromDataCalled) return kTRUE;
445 fSetRunNumberFromDataCalled = kTRUE;
447 AliCDBManager* man = AliCDBManager::Instance();
448 Int_t runData = -1, runCDB = -1;
450 AliRunLoader* runLoader = LoadRun("READ");
451 if (!runLoader) return kFALSE;
453 runData = runLoader->GetHeader()->GetRun();
457 runCDB = man->GetRun();
459 if (runCDB != runData) {
460 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
461 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
462 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
463 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
468 man->SetRun(runData);
471 if(man->GetRun() < 0) {
472 AliError("Run number not properly initalized!");
481 //_____________________________________________________________________________
482 void AliSimulation::SetConfigFile(const char* fileName)
484 // set the name of the config file
486 fConfigFileName = fileName;
489 //_____________________________________________________________________________
490 void AliSimulation::SetGAliceFile(const char* fileName)
492 // set the name of the galice file
493 // the path is converted to an absolute one if it is relative
495 fGAliceFileName = fileName;
496 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
497 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
499 fGAliceFileName = absFileName;
500 delete[] absFileName;
503 AliDebug(2, Form("galice file name set to %s", fileName));
506 //_____________________________________________________________________________
507 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
510 // set the number of events per file for the given detector and data type
511 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
513 TNamed* obj = new TNamed(detector, type);
514 obj->SetUniqueID(nEvents);
515 fEventsPerFile.Add(obj);
518 //_____________________________________________________________________________
519 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
521 // Read the alignment objects from CDB.
522 // Each detector is supposed to have the
523 // alignment objects in DET/Align/Data CDB path.
524 // All the detector objects are then collected,
525 // sorted by geometry level (starting from ALIC) and
526 // then applied to the TGeo geometry.
527 // Finally an overlaps check is performed.
529 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
530 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
534 // initialize CDB storage, run number, set CDB lock
536 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
539 Bool_t delRunLoader = kFALSE;
541 runLoader = LoadRun("READ");
542 if (!runLoader) return kFALSE;
543 delRunLoader = kTRUE;
546 // Export ideal geometry
547 if(!IsGeometryFromFile()) AliGeomManager::GetGeometry()->Export("geometry.root");
549 // Load alignment data from CDB and apply to geometry through AliGeomManager
550 if(fLoadAlignFromCDB){
552 TString detStr = fLoadAlObjsListOfDets;
553 TString loadAlObjsListOfDets = "";
555 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
556 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
557 AliModule* det = (AliModule*) detArray->At(iDet);
558 if (!det || !det->IsActive()) continue;
559 if (IsSelected(det->GetName(), detStr)) {
560 //add det to list of dets to be aligned from CDB
561 loadAlObjsListOfDets += det->GetName();
562 loadAlObjsListOfDets += " ";
564 } // end loop over detectors
565 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
566 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
568 // Check if the array with alignment objects was
569 // provided by the user. If yes, apply the objects
570 // to the present TGeo geometry
571 if (fAlignObjArray) {
572 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
573 AliError("The misalignment of one or more volumes failed!"
574 "Compare the list of simulated detectors and the list of detector alignment data!");
575 if (delRunLoader) delete runLoader;
581 // Update the internal geometry of modules (ITS needs it)
582 TString detStr = fLoadAlObjsListOfDets;
583 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
584 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
586 AliModule* det = (AliModule*) detArray->At(iDet);
587 if (!det || !det->IsActive()) continue;
588 if (IsSelected(det->GetName(), detStr)) {
589 det->UpdateInternalGeometry();
591 } // end loop over detectors
594 if (delRunLoader) delete runLoader;
599 //_____________________________________________________________________________
600 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
602 // add a file with background events for merging
604 TObjString* fileNameStr = new TObjString(fileName);
605 fileNameStr->SetUniqueID(nSignalPerBkgrd);
606 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
607 fBkgrdFileNames->Add(fileNameStr);
610 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
612 // add a file with background events for embeddin
613 MergeWith(fileName, nSignalPerBkgrd);
614 fEmbeddingFlag = kTRUE;
617 //_____________________________________________________________________________
618 Bool_t AliSimulation::Run(Int_t nEvents)
620 // run the generation, simulation and digitization
623 AliCodeTimerAuto("",0)
624 AliSysInfo::AddStamp("Start_Run");
626 // Load run number and seed from environmental vars
627 ProcessEnvironmentVars();
628 AliSysInfo::AddStamp("ProcessEnvironmentVars");
630 gRandom->SetSeed(fSeed);
632 if (nEvents > 0) fNEvents = nEvents;
634 // create and setup the HLT instance
635 if (!fRunHLT.IsNull() && !CreateHLT()) {
636 if (fStopOnError) return kFALSE;
641 // generation and simulation -> hits
642 if (fRunGeneration) {
643 if (!RunSimulation()) if (fStopOnError) return kFALSE;
645 AliSysInfo::AddStamp("RunSimulation");
647 // initialize CDB storage from external environment
648 // (either CDB manager or AliSimulation setters),
649 // if not already done in RunSimulation()
651 AliSysInfo::AddStamp("InitCDB");
653 // Set run number in CDBManager from data
654 // From this point on the run number must be always loaded from data!
655 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
657 // Set CDB lock: from now on it is forbidden to reset the run number
658 // or the default storage or to activate any further storage!
661 // If RunSimulation was not called, load the geometry and misalign it
662 if (!AliGeomManager::GetGeometry()) {
663 // Initialize the geometry manager
664 AliGeomManager::LoadGeometry("geometry.root");
665 AliSysInfo::AddStamp("GetGeometry");
666 // // Check that the consistency of symbolic names for the activated subdetectors
667 // // in the geometry loaded by AliGeomManager
668 // AliRunLoader* runLoader = LoadRun("READ");
669 // if (!runLoader) return kFALSE;
671 // TString detsToBeChecked = "";
672 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
673 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
674 // AliModule* det = (AliModule*) detArray->At(iDet);
675 // if (!det || !det->IsActive()) continue;
676 // detsToBeChecked += det->GetName();
677 // detsToBeChecked += " ";
678 // } // end loop over detectors
679 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
680 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
681 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
683 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
685 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
687 AliSysInfo::AddStamp("MissalignGeometry");
690 // hits -> summable digits
691 AliSysInfo::AddStamp("Start_sdigitization");
692 if (!fMakeSDigits.IsNull()) {
693 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
696 AliSysInfo::AddStamp("Stop_sdigitization");
698 AliSysInfo::AddStamp("Start_digitization");
699 // summable digits -> digits
700 if (!fMakeDigits.IsNull()) {
701 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
702 if (fStopOnError) return kFALSE;
705 AliSysInfo::AddStamp("Stop_digitization");
710 if (!fMakeDigitsFromHits.IsNull()) {
711 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
712 AliWarning(Form("Merging and direct creation of digits from hits "
713 "was selected for some detectors. "
714 "No merging will be done for the following detectors: %s",
715 fMakeDigitsFromHits.Data()));
717 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
718 if (fStopOnError) return kFALSE;
722 AliSysInfo::AddStamp("Hits2Digits");
726 if (!fTriggerConfig.IsNull() && !RunTrigger(fTriggerConfig,fMakeDigits)) {
727 if (fStopOnError) return kFALSE;
730 AliSysInfo::AddStamp("RunTrigger");
733 // digits -> raw data
734 if (!fWriteRawData.IsNull()) {
735 if (!WriteRawData(fWriteRawData, fRawDataFileName,
736 fDeleteIntermediateFiles,fWriteSelRawData)) {
737 if (fStopOnError) return kFALSE;
741 AliSysInfo::AddStamp("WriteRaw");
743 // run HLT simulation on simulated digit data if raw data is not
744 // simulated, otherwise its called as part of WriteRawData
745 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
747 if (fStopOnError) return kFALSE;
751 AliSysInfo::AddStamp("RunHLT");
755 Bool_t rv = RunQA() ;
761 AliSysInfo::AddStamp("RunQA");
763 TString snapshotFileOut("");
764 if(TString(gSystem->Getenv("OCDB_SNAPSHOT_CREATE")) == TString("kTRUE")){
765 AliInfo(" ******** Creating the snapshot! *********");
766 TString snapshotFile(gSystem->Getenv("OCDB_SNAPSHOT_FILENAME"));
767 if(!(snapshotFile.IsNull() || snapshotFile.IsWhitespace()))
768 snapshotFileOut = snapshotFile;
770 snapshotFileOut="OCDB.root";
771 AliCDBManager::Instance()->DumpToSnapshotFile(snapshotFileOut.Data(),kFALSE);
774 // Cleanup of CDB manager: cache and active storages!
775 AliCDBManager::Instance()->ClearCache();
780 //_______________________________________________________________________
781 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
782 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
783 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
786 // Generates lego plots of:
787 // - radiation length map phi vs theta
788 // - radiation length map phi vs eta
789 // - interaction length map
790 // - g/cm2 length map
792 // ntheta bins in theta, eta
793 // themin minimum angle in theta (degrees)
794 // themax maximum angle in theta (degrees)
796 // phimin minimum angle in phi (degrees)
797 // phimax maximum angle in phi (degrees)
798 // rmin minimum radius
799 // rmax maximum radius
802 // The number of events generated = ntheta*nphi
803 // run input parameters in macro setup (default="Config.C")
805 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
808 <img src="picts/AliRunLego1.gif">
813 <img src="picts/AliRunLego2.gif">
818 <img src="picts/AliRunLego3.gif">
823 // run the generation and simulation
825 AliCodeTimerAuto("",0)
827 // initialize CDB storage and run number from external environment
828 // (either CDB manager or AliSimulation setters)
834 AliError("no gAlice object. Restart aliroot and try again.");
837 if (gAlice->Modules()->GetEntries() > 0) {
838 AliError("gAlice was already run. Restart aliroot and try again.");
842 AliInfo(Form("initializing gAlice with config file %s",
843 fConfigFileName.Data()));
846 if (nev == -1) nev = nc1 * nc2;
848 // check if initialisation has been done
849 // If runloader has been initialized, set the number of events per file to nc1 * nc2
852 if (!gener) gener = new AliLegoGenerator();
854 // Configure Generator
856 gener->SetRadiusRange(rmin, rmax);
857 gener->SetZMax(zmax);
858 gener->SetCoor1Range(nc1, c1min, c1max);
859 gener->SetCoor2Range(nc2, c2min, c2max);
863 fLego = new AliLego("lego",gener);
865 //__________________________________________________________________________
869 gROOT->LoadMacro(setup);
870 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
872 if(AliCDBManager::Instance()->GetRun() >= 0) {
873 SetRunNumber(AliCDBManager::Instance()->GetRun());
875 AliWarning("Run number not initialized!!");
878 AliRunLoader::Instance()->CdGAFile();
880 AliPDG::AddParticlesToPdgDataBase();
882 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
884 gAlice->GetMCApp()->Init();
887 //Must be here because some MCs (G4) adds detectors here and not in Config.C
888 gAlice->InitLoaders();
889 AliRunLoader::Instance()->MakeTree("E");
892 // Save stuff at the beginning of the file to avoid file corruption
893 AliRunLoader::Instance()->CdGAFile();
896 //Save current generator
897 AliGenerator *gen=gAlice->GetMCApp()->Generator();
898 gAlice->GetMCApp()->ResetGenerator(gener);
899 //Prepare MC for Lego Run
905 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
906 gMC->ProcessRun(nev);
908 // End of this run, close files
910 // Restore current generator
911 gAlice->GetMCApp()->ResetGenerator(gen);
912 // Delete Lego Object
918 //_____________________________________________________________________________
919 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
923 AliCodeTimerAuto("",0)
925 // initialize CDB storage from external environment
926 // (either CDB manager or AliSimulation setters),
927 // if not already done in RunSimulation()
930 // Set run number in CDBManager from data
931 // From this point on the run number must be always loaded from data!
932 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
934 // Set CDB lock: from now on it is forbidden to reset the run number
935 // or the default storage or to activate any further storage!
938 AliRunLoader* runLoader = LoadRun("READ");
939 if (!runLoader) return kFALSE;
940 TString trconfiguration = config;
942 if (trconfiguration.IsNull()) {
943 if(!fTriggerConfig.IsNull()) {
944 trconfiguration = fTriggerConfig;
947 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
950 runLoader->MakeTree( "GG" );
951 AliCentralTrigger* aCTP = runLoader->GetTrigger();
952 // Load Configuration
953 if (!aCTP->LoadConfiguration( trconfiguration ))
957 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
969 //_____________________________________________________________________________
970 Bool_t AliSimulation::WriteTriggerRawData()
972 // Writes the CTP (trigger) DDL raw data
973 // Details of the format are given in the
974 // trigger TDR - pages 134 and 135.
975 AliCTPRawData writer;
981 //_____________________________________________________________________________
982 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
984 // run the generation and simulation
986 AliCodeTimerAuto("",0)
988 // initialize CDB storage and run number from external environment
989 // (either CDB manager or AliSimulation setters)
990 AliSysInfo::AddStamp("RunSimulation_Begin");
992 AliSysInfo::AddStamp("RunSimulation_InitCDB");
996 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
999 AliError("no gAlice object. Restart aliroot and try again.");
1002 if (gAlice->Modules()->GetEntries() > 0) {
1003 AliError("gAlice was already run. Restart aliroot and try again.");
1007 AliInfo(Form("initializing gAlice with config file %s",
1008 fConfigFileName.Data()));
1011 // Initialize ALICE Simulation run
1016 // If requested set the mag. field from the GRP entry.
1017 // After this the field is loccked and cannot be changed by Config.C
1018 if (fUseMagFieldFromGRP) {
1020 grpM.ReadGRPEntry();
1022 AliInfo("Field is locked now. It cannot be changed in Config.C");
1026 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1027 gROOT->LoadMacro(fConfigFileName.Data());
1028 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1029 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1030 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1032 AliSysInfo::AddStamp("RunSimulation_Config");
1035 // If requested obtain the vertex position and vertex sigma_z from the CDB
1036 // This overwrites the settings from the Config.C
1037 if (fUseVertexFromCDB) {
1038 Double_t vtxPos[3] = {0., 0., 0.};
1039 Double_t vtxSig[3] = {0., 0., 0.};
1040 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1042 AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1044 if(vertex->GetXRes()>2.8) { // > pipe radius --> it's a dummy object, don't use it
1045 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1046 if (entry) vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1050 vertex->GetXYZ(vtxPos);
1051 vertex->GetSigmaXYZ(vtxSig);
1052 AliInfo("Overwriting Config.C vertex settings !");
1053 AliInfo(Form("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
1054 vtxPos[0], vtxPos[1], vtxPos[2], vtxSig[2]));
1056 AliGenerator *gen = gAlice->GetMCApp()->Generator();
1057 gen->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]); // vertex position
1058 gen->SetSigmaZ(vtxSig[2]);
1063 // If requested we take the SOR and EOR time-stamps from the GRP and use them
1064 // in order to generate the event time-stamps
1065 if (fUseTimeStampFromCDB) {
1067 grpM.ReadGRPEntry();
1068 const AliGRPObject *grpObj = grpM.GetGRPData();
1069 if (!grpObj || (grpObj->GetTimeEnd() <= grpObj->GetTimeStart())) {
1070 AliError("Missing GRP or bad SOR/EOR time-stamps! Switching off the time-stamp generation from GRP!");
1071 fTimeStart = fTimeEnd = 0;
1072 fUseTimeStampFromCDB = kFALSE;
1075 fTimeStart = grpObj->GetTimeStart();
1076 fTimeEnd = grpObj->GetTimeEnd();
1080 if(AliCDBManager::Instance()->GetRun() >= 0) {
1081 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
1082 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
1084 AliWarning("Run number not initialized!!");
1087 AliRunLoader::Instance()->CdGAFile();
1090 AliPDG::AddParticlesToPdgDataBase();
1092 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1093 AliSysInfo::AddStamp("RunSimulation_GetField");
1095 gAlice->GetMCApp()->Init();
1096 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1098 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1099 gAlice->InitLoaders();
1100 AliRunLoader::Instance()->MakeTree("E");
1101 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1102 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1103 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1105 // Save stuff at the beginning of the file to avoid file corruption
1106 AliRunLoader::Instance()->CdGAFile();
1108 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1109 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1110 //___________________________________________________________________________________________
1112 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1114 // Set run number in CDBManager
1115 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1117 AliRunLoader* runLoader = AliRunLoader::Instance();
1119 AliError(Form("gAlice has no run loader object. "
1120 "Check your config file: %s", fConfigFileName.Data()));
1123 SetGAliceFile(runLoader->GetFileName());
1125 // Misalign geometry
1126 #if ROOT_VERSION_CODE < 331527
1127 AliGeomManager::SetGeometry(gGeoManager);
1129 // Check that the consistency of symbolic names for the activated subdetectors
1130 // in the geometry loaded by AliGeomManager
1131 TString detsToBeChecked = "";
1132 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1133 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1134 AliModule* det = (AliModule*) detArray->At(iDet);
1135 if (!det || !det->IsActive()) continue;
1136 detsToBeChecked += det->GetName();
1137 detsToBeChecked += " ";
1138 } // end loop over detectors
1139 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1140 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1141 MisalignGeometry(runLoader);
1142 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1145 // AliRunLoader* runLoader = AliRunLoader::Instance();
1146 // if (!runLoader) {
1147 // AliError(Form("gAlice has no run loader object. "
1148 // "Check your config file: %s", fConfigFileName.Data()));
1151 // SetGAliceFile(runLoader->GetFileName());
1153 if (!gAlice->GetMCApp()->Generator()) {
1154 AliError(Form("gAlice has no generator object. "
1155 "Check your config file: %s", fConfigFileName.Data()));
1159 // Write GRP entry corresponding to the setting found in Cofig.C
1162 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1164 if (nEvents <= 0) nEvents = fNEvents;
1166 // get vertex from background file in case of merging
1167 if (fUseBkgrdVertex &&
1168 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1169 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1170 const char* fileName = ((TObjString*)
1171 (fBkgrdFileNames->At(0)))->GetName();
1172 AliInfo(Form("The vertex will be taken from the background "
1173 "file %s with nSignalPerBackground = %d",
1174 fileName, signalPerBkgrd));
1175 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1176 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1179 if (!fRunSimulation) {
1180 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1183 // set the number of events per file for given detectors and data types
1184 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1185 if (!fEventsPerFile[i]) continue;
1186 const char* detName = fEventsPerFile[i]->GetName();
1187 const char* typeName = fEventsPerFile[i]->GetTitle();
1188 TString loaderName(detName);
1189 loaderName += "Loader";
1190 AliLoader* loader = runLoader->GetLoader(loaderName);
1192 AliError(Form("RunSimulation no loader for %s found\n Number of events per file not set for %s %s",
1193 detName, typeName, detName));
1196 AliDataLoader* dataLoader =
1197 loader->GetDataLoader(typeName);
1199 AliError(Form("no data loader for %s found\n"
1200 "Number of events per file not set for %s %s",
1201 typeName, detName, typeName));
1204 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1205 AliDebug(1, Form("number of events per file set to %d for %s %s",
1206 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1209 AliInfo("running gAlice");
1210 AliSysInfo::AddStamp("Start_ProcessRun");
1212 // Create the Root Tree with one branch per detector
1213 //Hits moved to begin event -> now we are crating separate tree for each event
1215 gMC->ProcessRun(nEvents);
1217 // End of this run, close files
1218 if(nEvents>0) FinishRun();
1220 AliSysInfo::AddStamp("Stop_ProcessRun");
1226 //_____________________________________________________________________________
1227 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1229 // run the digitization and produce summable digits
1230 static Int_t eventNr=0;
1231 AliCodeTimerAuto("",0) ;
1233 // initialize CDB storage, run number, set CDB lock
1235 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1238 AliRunLoader* runLoader = LoadRun();
1239 if (!runLoader) return kFALSE;
1241 TString detStr = detectors;
1242 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1243 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1244 AliModule* det = (AliModule*) detArray->At(iDet);
1245 if (!det || !det->IsActive()) continue;
1246 if (IsSelected(det->GetName(), detStr)) {
1247 AliInfo(Form("creating summable digits for %s", det->GetName()));
1248 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1249 det->Hits2SDigits();
1250 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1251 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1255 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1256 AliError(Form("the following detectors were not found: %s",
1258 if (fStopOnError) return kFALSE;
1267 //_____________________________________________________________________________
1268 Bool_t AliSimulation::RunDigitization(const char* detectors,
1269 const char* excludeDetectors)
1271 // run the digitization and produce digits from sdigits
1273 AliCodeTimerAuto("",0)
1275 // initialize CDB storage, run number, set CDB lock
1277 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1280 delete AliRunLoader::Instance();
1285 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1286 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1287 AliDigitizationInput digInp(nStreams, signalPerBkgrd);
1288 // digInp.SetEmbeddingFlag(fEmbeddingFlag);
1289 digInp.SetRegionOfInterest(fRegionOfInterest);
1290 digInp.SetInputStream(0, fGAliceFileName.Data());
1291 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1292 const char* fileName = ((TObjString*)(fBkgrdFileNames->At(iStream-1)))->GetName();
1293 digInp.SetInputStream(iStream, fileName);
1296 detArr.SetOwner(kTRUE);
1297 TString detStr = detectors;
1298 TString detExcl = excludeDetectors;
1299 if (!static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice()) {
1300 AliError("Error occured while getting gAlice from Input 0");
1303 AliRunLoader* runLoader = AliRunLoader::GetRunLoader(digInp.GetInputStream(0)->GetFolderName());
1304 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1305 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1306 AliModule* det = (AliModule*) detArray->At(iDet);
1307 if (!det || !det->IsActive()) continue;
1308 if (!IsSelected(det->GetName(), detStr) || IsSelected(det->GetName(), detExcl)) continue;
1309 AliDigitizer* digitizer = det->CreateDigitizer(&digInp);
1310 if (!digitizer || !digitizer->Init()) {
1311 AliError(Form("no digitizer for %s", det->GetName()));
1312 if (fStopOnError) return kFALSE;
1315 detArr.AddLast(digitizer);
1316 AliInfo(Form("Created digitizer from SDigits -> Digits for %s", det->GetName()));
1319 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1320 AliError(Form("the following detectors were not found: %s", detStr.Data()));
1321 if (fStopOnError) return kFALSE;
1324 Int_t ndigs = detArr.GetEntriesFast();
1325 Int_t eventsCreated = 0;
1326 AliRunLoader* outRl = digInp.GetOutRunLoader();
1327 while ((eventsCreated++ < fNEvents) || (fNEvents < 0)) {
1328 if (!digInp.ConnectInputTrees()) break;
1329 digInp.InitEvent(); //this must be after call of Connect Input tress.
1330 if (outRl) outRl->SetEventNumber(eventsCreated-1);
1331 static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice(); // use gAlice of the first input stream
1332 for (int id=0;id<ndigs;id++) ((AliDigitizer*)detArr[id])->Digitize("");
1333 digInp.FinishEvent();
1335 digInp.FinishGlobal();
1340 //_____________________________________________________________________________
1341 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1343 // run the digitization and produce digits from hits
1345 AliCodeTimerAuto("",0)
1347 // initialize CDB storage, run number, set CDB lock
1349 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1352 AliRunLoader* runLoader = LoadRun("READ");
1353 if (!runLoader) return kFALSE;
1355 TString detStr = detectors;
1356 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1357 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1358 AliModule* det = (AliModule*) detArray->At(iDet);
1359 if (!det || !det->IsActive()) continue;
1360 if (IsSelected(det->GetName(), detStr)) {
1361 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1366 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1367 AliError(Form("the following detectors were not found: %s",
1369 if (fStopOnError) return kFALSE;
1375 //_____________________________________________________________________________
1376 Bool_t AliSimulation::WriteRawData(const char* detectors,
1377 const char* fileName,
1378 Bool_t deleteIntermediateFiles,
1381 // convert the digits to raw data
1382 // First DDL raw data files for the given detectors are created.
1383 // If a file name is given, the DDL files are then converted to a DATE file.
1384 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1386 // If the file name has the extension ".root", the DATE file is converted
1388 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1389 // 'selrawdata' flag can be used to enable writing of detectors raw data
1390 // accoring to the trigger cluster.
1392 AliCodeTimerAuto("",0)
1393 AliSysInfo::AddStamp("WriteRawData_Start");
1395 TString detStr = detectors;
1396 if (!WriteRawFiles(detStr.Data())) {
1397 if (fStopOnError) return kFALSE;
1399 AliSysInfo::AddStamp("WriteRawFiles");
1401 // run HLT simulation on simulated DDL raw files
1402 // and produce HLT ddl raw files to be included in date/root file
1403 // bugfix 2009-06-26: the decision whether to write HLT raw data
1404 // is taken in RunHLT. Here HLT always needs to be run in order to
1405 // create HLT digits, unless its switched off. This is due to the
1406 // special placement of the HLT between the generation of DDL files
1407 // and conversion to DATE/Root file.
1408 detStr.ReplaceAll("HLT", "");
1409 if (!fRunHLT.IsNull()) {
1411 if (fStopOnError) return kFALSE;
1414 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1416 TString dateFileName(fileName);
1417 if (!dateFileName.IsNull()) {
1418 Bool_t rootOutput = dateFileName.EndsWith(".root");
1419 if (rootOutput) dateFileName += ".date";
1420 TString selDateFileName;
1422 selDateFileName = "selected.";
1423 selDateFileName+= dateFileName;
1425 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1426 if (fStopOnError) return kFALSE;
1428 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1429 if (deleteIntermediateFiles) {
1430 AliRunLoader* runLoader = LoadRun("READ");
1431 if (runLoader) for (Int_t iEvent = 0;
1432 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1434 snprintf(command, 256, "rm -r raw%d", iEvent);
1435 gSystem->Exec(command);
1441 if (!ConvertDateToRoot(dateFileName, fileName)) {
1442 if (fStopOnError) return kFALSE;
1444 AliSysInfo::AddStamp("ConvertDateToRoot");
1445 if (deleteIntermediateFiles) {
1446 gSystem->Unlink(dateFileName);
1449 TString selFileName = "selected.";
1450 selFileName += fileName;
1451 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1452 if (fStopOnError) return kFALSE;
1454 if (deleteIntermediateFiles) {
1455 gSystem->Unlink(selDateFileName);
1464 //_____________________________________________________________________________
1465 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1467 // convert the digits to raw data DDL files
1469 AliCodeTimerAuto("",0)
1471 AliRunLoader* runLoader = LoadRun("READ");
1472 if (!runLoader) return kFALSE;
1474 // write raw data to DDL files
1475 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1476 AliInfo(Form("processing event %d", iEvent));
1477 runLoader->GetEvent(iEvent);
1478 TString baseDir = gSystem->WorkingDirectory();
1480 snprintf(dirName, 256, "raw%d", iEvent);
1481 gSystem->MakeDirectory(dirName);
1482 if (!gSystem->ChangeDirectory(dirName)) {
1483 AliError(Form("couldn't change to directory %s", dirName));
1484 if (fStopOnError) return kFALSE; else continue;
1487 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1490 TString detStr = detectors;
1491 if (IsSelected("HLT", detStr)) {
1492 // Do nothing. "HLT" will be removed from detStr and HLT raw
1493 // data files are generated in RunHLT.
1496 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1497 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1498 AliModule* det = (AliModule*) detArray->At(iDet);
1499 if (!det || !det->IsActive()) continue;
1500 if (IsSelected(det->GetName(), detStr)) {
1501 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1506 if (!WriteTriggerRawData())
1507 if (fStopOnError) return kFALSE;
1509 gSystem->ChangeDirectory(baseDir);
1510 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1511 AliError(Form("the following detectors were not found: %s",
1513 if (fStopOnError) return kFALSE;
1522 //_____________________________________________________________________________
1523 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1524 const char* selDateFileName)
1526 // convert raw data DDL files to a DATE file with the program "dateStream"
1527 // The second argument is not empty when the user decides to write
1528 // the detectors raw data according to the trigger cluster.
1530 AliCodeTimerAuto("",0)
1532 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1534 AliError("the program dateStream was not found");
1535 if (fStopOnError) return kFALSE;
1540 AliRunLoader* runLoader = LoadRun("READ");
1541 if (!runLoader) return kFALSE;
1543 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1544 Bool_t selrawdata = kFALSE;
1545 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1548 // Note the option -s. It is used in order to avoid
1549 // the generation of SOR/EOR events.
1550 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1551 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1552 FILE* pipe = gSystem->OpenPipe(command, "w");
1555 AliError(Form("Cannot execute command: %s",command));
1559 Int_t selEvents = 0;
1560 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1562 UInt_t detectorPattern = 0;
1563 runLoader->GetEvent(iEvent);
1564 if (!runLoader->LoadTrigger()) {
1565 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1566 detectorPattern = aCTP->GetClusterMask();
1567 // Check if the event was triggered by CTP
1569 if (aCTP->GetClassMask()) selEvents++;
1573 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1575 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1576 selrawdata = kFALSE;
1580 fprintf(pipe, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1584 // loop over detectors and DDLs
1585 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1586 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1588 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1589 Int_t ldcID = Int_t(ldc + 0.0001);
1590 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1592 char rawFileName[256];
1593 snprintf(rawFileName, 256, "raw%d/%s",
1594 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1596 // check existence and size of raw data file
1597 FILE* file = fopen(rawFileName, "rb");
1598 if (!file) continue;
1599 fseek(file, 0, SEEK_END);
1600 unsigned long size = ftell(file);
1602 if (!size) continue;
1604 if (ldcID != prevLDC) {
1605 fprintf(pipe, " LDC Id %d\n", ldcID);
1608 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1613 Int_t result = gSystem->ClosePipe(pipe);
1615 if (!(selrawdata && selEvents > 0)) {
1617 return (result == 0);
1620 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1622 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1623 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1624 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1626 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1628 // Get the trigger decision and cluster
1629 UInt_t detectorPattern = 0;
1631 runLoader->GetEvent(iEvent);
1632 if (!runLoader->LoadTrigger()) {
1633 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1634 if (aCTP->GetClassMask() == 0) continue;
1635 detectorPattern = aCTP->GetClusterMask();
1636 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1637 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1640 fprintf(pipe2, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1644 // loop over detectors and DDLs
1645 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1646 // Write only raw data from detectors that
1647 // are contained in the trigger cluster(s)
1648 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1650 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1652 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1653 Int_t ldcID = Int_t(ldc + 0.0001);
1654 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1656 char rawFileName[256];
1657 snprintf(rawFileName, 256, "raw%d/%s",
1658 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1660 // check existence and size of raw data file
1661 FILE* file = fopen(rawFileName, "rb");
1662 if (!file) continue;
1663 fseek(file, 0, SEEK_END);
1664 unsigned long size = ftell(file);
1666 if (!size) continue;
1668 if (ldcID != prevLDC) {
1669 fprintf(pipe2, " LDC Id %d\n", ldcID);
1672 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1677 Int_t result2 = gSystem->ClosePipe(pipe2);
1680 return ((result == 0) && (result2 == 0));
1683 //_____________________________________________________________________________
1684 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1685 const char* rootFileName)
1687 // convert a DATE file to a root file with the program "alimdc"
1690 const Int_t kDBSize = 2000000000;
1691 const Int_t kTagDBSize = 1000000000;
1692 const Bool_t kFilter = kFALSE;
1693 const Int_t kCompression = 1;
1695 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1697 AliError("the program alimdc was not found");
1698 if (fStopOnError) return kFALSE;
1703 AliInfo(Form("converting DATE file %s to root file %s",
1704 dateFileName, rootFileName));
1706 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1707 const char* tagDBFS = "/tmp/mdc1/tags";
1709 // User defined file system locations
1710 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1711 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1712 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1713 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1714 if (gSystem->Getenv("ALIMDC_TAGDB"))
1715 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1717 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1718 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1719 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1721 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1722 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1723 gSystem->Exec(Form("mkdir %s",tagDBFS));
1725 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1726 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1727 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1729 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1730 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1731 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1733 return (result == 0);
1737 //_____________________________________________________________________________
1738 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1740 // delete existing run loaders, open a new one and load gAlice
1742 delete AliRunLoader::Instance();
1743 AliRunLoader* runLoader =
1744 AliRunLoader::Open(fGAliceFileName.Data(),
1745 AliConfig::GetDefaultEventFolderName(), mode);
1747 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1750 runLoader->LoadgAlice();
1751 runLoader->LoadHeader();
1752 gAlice = runLoader->GetAliRun();
1754 AliError(Form("no gAlice object found in file %s",
1755 fGAliceFileName.Data()));
1761 //_____________________________________________________________________________
1762 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1764 // get or calculate the number of signal events per background event
1766 if (!fBkgrdFileNames) return 1;
1767 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1768 if (nBkgrdFiles == 0) return 1;
1770 // get the number of signal events
1772 AliRunLoader* runLoader =
1773 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1774 if (!runLoader) return 1;
1776 nEvents = runLoader->GetNumberOfEvents();
1781 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1782 // get the number of background events
1783 const char* fileName = ((TObjString*)
1784 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1785 AliRunLoader* runLoader =
1786 AliRunLoader::Open(fileName, "BKGRD");
1787 if (!runLoader) continue;
1788 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1791 // get or calculate the number of signal per background events
1792 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1793 if (nSignalPerBkgrd <= 0) {
1794 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1795 } else if (result && (result != nSignalPerBkgrd)) {
1796 AliInfo(Form("the number of signal events per background event "
1797 "will be changed from %d to %d for stream %d",
1798 nSignalPerBkgrd, result, iBkgrdFile+1));
1799 nSignalPerBkgrd = result;
1802 if (!result) result = nSignalPerBkgrd;
1803 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1804 AliWarning(Form("not enough background events (%d) for %d signal events "
1805 "using %d signal per background events for stream %d",
1806 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1813 //_____________________________________________________________________________
1814 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1816 // check whether detName is contained in detectors
1817 // if yes, it is removed from detectors
1819 // check if all detectors are selected
1820 if ((detectors.CompareTo("ALL") == 0) ||
1821 detectors.BeginsWith("ALL ") ||
1822 detectors.EndsWith(" ALL") ||
1823 detectors.Contains(" ALL ")) {
1828 // search for the given detector
1829 Bool_t result = kFALSE;
1830 if ((detectors.CompareTo(detName) == 0) ||
1831 detectors.BeginsWith(detName+" ") ||
1832 detectors.EndsWith(" "+detName) ||
1833 detectors.Contains(" "+detName+" ")) {
1834 detectors.ReplaceAll(detName, "");
1838 // clean up the detectors string
1839 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1840 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1841 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1846 //_____________________________________________________________________________
1847 Int_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName, Int_t N, Int_t nSkip)
1850 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1851 // These can be used for embedding of MC tracks into RAW data using the standard
1852 // merging procedure.
1854 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1857 AliError("no gAlice object. Restart aliroot and try again.");
1860 if (gAlice->Modules()->GetEntries() > 0) {
1861 AliError("gAlice was already run. Restart aliroot and try again.");
1865 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1869 gROOT->LoadMacro(fConfigFileName.Data());
1870 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1872 if(AliCDBManager::Instance()->GetRun() >= 0) {
1873 SetRunNumber(AliCDBManager::Instance()->GetRun());
1875 AliWarning("Run number not initialized!!");
1878 AliRunLoader::Instance()->CdGAFile();
1880 AliPDG::AddParticlesToPdgDataBase();
1882 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1884 gAlice->GetMCApp()->Init();
1886 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1887 gAlice->InitLoaders();
1888 AliRunLoader::Instance()->MakeTree("E");
1889 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1890 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1891 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1894 // Save stuff at the beginning of the file to avoid file corruption
1895 AliRunLoader::Instance()->CdGAFile();
1900 //AliCDBManager* man = AliCDBManager::Instance();
1901 //man->SetRun(0); // Should this come from rawdata header ?
1905 // Get the runloader
1906 AliRunLoader* runLoader = AliRunLoader::Instance();
1908 // Open esd file if available
1911 AliESDEvent* esd = 0;
1912 if (esdFileName && (strlen(esdFileName)>0)) {
1913 esdFile = TFile::Open(esdFileName);
1915 esd = new AliESDEvent();
1916 esdFile->GetObject("esdTree", treeESD);
1918 esd->ReadFromTree(treeESD);
1920 AliInfo(Form("Asking to skip first %d ESDs events",nSkip));
1929 // Create the RawReader
1930 TString fileName(rawDirectory);
1931 AliRawReader* rawReader = AliRawReader::Create(fileName.Data());
1932 if (!rawReader) return (kFALSE);
1934 // if (!fEquipIdMap.IsNull() && fRawReader)
1935 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1937 // Get list of detectors
1938 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1941 AliHeader* header = runLoader->GetHeader();
1945 if (!(rawReader->NextEvent())) break;
1946 runLoader->SetEventNumber(nev);
1947 runLoader->GetHeader()->Reset(rawReader->GetRunNumber(),
1949 runLoader->GetEvent(nev);
1950 AliInfo(Form("We are at event %d",nev));
1953 TString detStr = fMakeSDigits;
1954 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1955 AliModule* det = (AliModule*) detArray->At(iDet);
1956 if (!det || !det->IsActive()) continue;
1957 if (IsSelected(det->GetName(), detStr)) {
1958 AliInfo(Form("Calling Raw2SDigits for %s", det->GetName()));
1959 det->Raw2SDigits(rawReader);
1966 // If ESD information available obtain reconstructed vertex and store in header.
1968 AliInfo(Form("Selected event %d correspond to event %d in raw and to %d in esd",nev,rawReader->GetEventIndex(),nSkip+rawReader->GetEventIndex()));
1969 treeESD->GetEvent(nSkip+rawReader->GetEventIndex());
1970 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1971 Double_t position[3];
1972 esdVertex->GetXYZ(position);
1973 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1976 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1977 mcHeader->SetPrimaryVertex(mcV);
1978 header->Reset(0,nev);
1979 header->SetGenEventHeader(mcHeader);
1980 AliInfo(Form("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]));
1984 runLoader->TreeE()->Fill();
1985 AliInfo(Form("Finished event %d",nev));
1994 runLoader->CdGAFile();
1995 runLoader->WriteHeader("OVERWRITE");
1996 runLoader->WriteRunLoader();
2001 //_____________________________________________________________________________
2002 void AliSimulation::FinishRun()
2005 // Called at the end of the run.
2010 AliDebug(1, "Finish Lego");
2011 AliRunLoader::Instance()->CdGAFile();
2015 // Clean detector information
2016 TIter next(gAlice->Modules());
2017 AliModule *detector;
2018 while((detector = dynamic_cast<AliModule*>(next()))) {
2019 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
2020 detector->FinishRun();
2023 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
2024 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
2026 // Write AliRun info and all detectors parameters
2027 AliRunLoader::Instance()->CdGAFile();
2028 gAlice->Write(0,TObject::kOverwrite);//write AliRun
2029 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
2031 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
2032 AliRunLoader::Instance()->Synchronize();
2035 //_____________________________________________________________________________
2036 Int_t AliSimulation::GetDetIndex(const char* detector)
2038 // return the detector index corresponding to detector
2040 for (index = 0; index < fgkNDetectors ; index++) {
2041 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2047 //_____________________________________________________________________________
2048 Bool_t AliSimulation::CreateHLT()
2050 // Init the HLT simulation.
2051 // The function loads the library and creates the instance of AliHLTSimulation.
2052 // the main reason for the decoupled creation is to set the transient OCDB
2053 // objects before the OCDB is locked
2055 // load the library dynamically
2056 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
2058 // check for the library version
2059 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
2061 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
2064 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
2065 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
2068 // print compile info
2069 typedef void (*CompileInfo)( const char*& date, const char*& time);
2070 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
2072 const char* date="";
2073 const char* time="";
2074 (*fctInfo)(date, time);
2075 if (!date) date="unknown";
2076 if (!time) time="unknown";
2077 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
2079 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
2082 // create instance of the HLT simulation
2083 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
2084 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
2085 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
2089 TString specObjects;
2090 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
2091 if (specObjects.Length()>0) specObjects+=" ";
2092 specObjects+=fSpecCDBUri[i]->GetName();
2095 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
2096 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
2097 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2103 //_____________________________________________________________________________
2104 Bool_t AliSimulation::RunHLT()
2106 // Run the HLT simulation
2107 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2108 // Disabled if fRunHLT is empty, default vaule is "default".
2109 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2110 // The default simulation depends on the HLT component libraries and their
2111 // corresponding agents which define components and chains to run. See
2112 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2113 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2115 // The libraries to be loaded can be specified as an option.
2117 // AliSimulation sim;
2118 // sim.SetRunHLT("libAliHLTSample.so");
2120 // will only load <tt>libAliHLTSample.so</tt>
2122 // Other available options:
2123 // \li loglevel=<i>level</i> <br>
2124 // logging level for this processing
2126 // disable redirection of log messages to AliLog class
2127 // \li config=<i>macro</i>
2128 // configuration macro
2129 // \li chains=<i>configuration</i>
2130 // comma separated list of configurations to be run during simulation
2131 // \li rawfile=<i>file</i>
2132 // source for the RawReader to be created, the default is <i>./</i> if
2133 // raw data is simulated
2137 if (!fpHLT && !CreateHLT()) {
2140 AliHLTSimulation* pHLT=fpHLT;
2142 AliRunLoader* pRunLoader = LoadRun("READ");
2143 if (!pRunLoader) return kFALSE;
2145 // initialize CDB storage, run number, set CDB lock
2146 // thats for the case of running HLT simulation without all the other steps
2147 // multiple calls are handled by the function, so we can just call
2149 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2152 // init the HLT simulation
2154 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2155 TString detStr = fWriteRawData;
2156 if (!IsSelected("HLT", detStr)) {
2157 options+=" writerawfiles=";
2159 options+=" writerawfiles=HLT";
2162 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2163 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2164 // in order to get detector data. By default, RawReaderFile is used to read
2165 // the already simulated ddl files. Date and Root files from the raw data
2166 // are generated after the HLT simulation.
2167 options+=" rawfile=./";
2170 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2171 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2172 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2174 // run the HLT simulation
2175 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2176 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2177 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2181 // delete the instance
2182 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2183 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2184 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2188 return iResult>=0?kTRUE:kFALSE;
2191 //_____________________________________________________________________________
2192 Bool_t AliSimulation::RunQA()
2194 // run the QA on summable hits, digits or digits
2196 //if(!gAlice) return kFALSE;
2197 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2199 TString detectorsw("") ;
2201 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2202 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2203 if ( detectorsw.IsNull() )
2208 //_____________________________________________________________________________
2209 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2211 // Allows to run QA for a selected set of detectors
2212 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2213 // all selected detectors run the same selected tasks
2215 if (!detAndAction.Contains(":")) {
2216 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2220 Int_t colon = detAndAction.Index(":") ;
2221 fQADetectors = detAndAction(0, colon) ;
2222 if (fQADetectors.Contains("ALL") ){
2223 TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2224 Int_t minus = fQADetectors.Last('-') ;
2225 TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2226 TString toRemove("") ;
2227 while (minus >= 0) {
2228 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
2229 toRemove = toRemove.Strip() ;
2230 toKeep.ReplaceAll(toRemove, "") ;
2231 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
2232 minus = fQADetectors.Last('-') ;
2234 fQADetectors = toKeep ;
2236 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2237 if (fQATasks.Contains("ALL") ) {
2238 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2240 fQATasks.ToUpper() ;
2242 if ( fQATasks.Contains("HIT") )
2243 tempo = Form("%d ", AliQAv1::kHITS) ;
2244 if ( fQATasks.Contains("SDIGIT") )
2245 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2246 if ( fQATasks.Contains("DIGIT") )
2247 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2249 if (fQATasks.IsNull()) {
2250 AliInfo("No QA requested\n") ;
2255 TString tempo(fQATasks) ;
2256 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2257 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2258 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2259 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2261 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2262 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2263 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2264 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2269 //_____________________________________________________________________________
2270 void AliSimulation::ProcessEnvironmentVars()
2272 // Extract run number and random generator seed from env variables
2274 AliInfo("Processing environment variables");
2276 // Random Number seed
2278 // first check that seed is not already set
2280 if (gSystem->Getenv("CONFIG_SEED")) {
2281 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2284 if (gSystem->Getenv("CONFIG_SEED")) {
2285 AliInfo(Form("Seed for random number generation already set (%d)"
2286 ": CONFIG_SEED variable ignored!", fSeed));
2290 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2294 // first check that run number is not already set
2296 if (gSystem->Getenv("DC_RUN")) {
2297 fRun = atoi(gSystem->Getenv("DC_RUN"));
2300 if (gSystem->Getenv("DC_RUN")) {
2301 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2305 AliInfo(Form("Run number = %d", fRun));
2308 //---------------------------------------------------------------------
2309 void AliSimulation::WriteGRPEntry()
2311 // Get the necessary information from galice (generator, trigger etc) and
2312 // write a GRP entry corresponding to the settings in the Config.C used
2313 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2316 AliInfo("Writing global run parameters entry into the OCDB");
2318 AliGRPObject* grpObj = new AliGRPObject();
2320 grpObj->SetRunType("PHYSICS");
2321 grpObj->SetTimeStart(fTimeStart);
2322 grpObj->SetTimeEnd(fTimeEnd);
2323 grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2325 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2331 gen->GetProjectile(projectile,a,z);
2333 gen->GetTarget(target,a,z);
2334 TString beamType = projectile + "-" + target;
2335 beamType.ReplaceAll(" ","");
2336 if (!beamType.CompareTo("-")) {
2337 grpObj->SetBeamType("UNKNOWN");
2338 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2341 grpObj->SetBeamType(beamType);
2343 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 * a / z);
2345 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 );
2347 // Heavy ion run, the event specie is set to kHighMult
2348 fEventSpecie = AliRecoParam::kHighMult;
2349 if ((strcmp(beamType,"p-p") == 0) ||
2350 (strcmp(beamType,"p-") == 0) ||
2351 (strcmp(beamType,"-p") == 0) ||
2352 (strcmp(beamType,"P-P") == 0) ||
2353 (strcmp(beamType,"P-") == 0) ||
2354 (strcmp(beamType,"-P") == 0)) {
2355 // Proton run, the event specie is set to kLowMult
2356 fEventSpecie = AliRecoParam::kLowMult;
2360 AliWarning("Unknown beam type and energy! Setting energy to 0");
2361 grpObj->SetBeamEnergy(0);
2362 grpObj->SetBeamType("UNKNOWN");
2365 UInt_t detectorPattern = 0;
2367 TObjArray *detArray = gAlice->Detectors();
2368 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2369 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2370 AliDebug(1, Form("Detector #%d found: %s", iDet, AliDAQ::OfflineModuleName(iDet)));
2371 detectorPattern |= (1 << iDet);
2376 if (!fTriggerConfig.IsNull())
2377 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2380 if (!fRunHLT.IsNull())
2381 detectorPattern |= (1 << AliDAQ::kHLTId);
2383 grpObj->SetNumberOfDetectors((Char_t)nDets);
2384 grpObj->SetDetectorMask((Int_t)detectorPattern);
2385 grpObj->SetLHCPeriod("LHC08c");
2386 grpObj->SetLHCState("STABLE_BEAMS");
2388 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2389 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2391 Float_t factorSol = field ? field->GetFactorSol() : 0;
2392 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2393 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2395 Float_t factorDip = field ? field->GetFactorDip() : 0;
2396 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2398 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2399 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2400 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2401 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2402 if (field) grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2403 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2405 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2407 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2409 // Now store the entry in OCDB
2410 AliCDBManager* man = AliCDBManager::Instance();
2412 man->SetLock(0, fKey);
2414 AliCDBStorage* sto = man->GetStorage(fGRPWriteLocation.Data());
2417 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2418 AliCDBMetaData *metadata= new AliCDBMetaData();
2420 metadata->SetResponsible("alice-off@cern.ch");
2421 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2423 sto->Put(grpObj,id,metadata);
2424 man->SetLock(1, fKey);
2427 //_____________________________________________________________________________
2428 time_t AliSimulation::GenerateTimeStamp() const
2430 // Generate event time-stamp according to
2431 // SOR/EOR time from GRP
2432 if (fUseTimeStampFromCDB)
2433 return fTimeStart + gRandom->Integer(fTimeEnd-fTimeStart);