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 **************************************************************************/
16 /* $Id: AliSimulation.cxx 64623 2013-10-21 13:38:58Z rgrosso $ */
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 // The trigger configuration is set by the method SetTriggerConfig(X) //
107 // X can take three kinds of values //
109 // - The exact string "none" - case insensitive. In this case, not trigger //
110 // information is generated from the digits. //
111 // - The empty string or "ocdb" - case insensitive. In this case the //
112 // trigger configuration is read from OCDB //
113 // - Some string - say "p-p" - in which case the configuration is read from //
114 // fixed files in $ALICE_ROOT/GRP/CTP/ - say $ALICE_ROOT/GRP/CTP/p-p.cfg //
116 // Default is to read from OCDB. //
118 ///////////////////////////////////////////////////////////////////////////////
121 #include <TGeoGlobalMagField.h>
122 #include <TGeoManager.h>
123 #include <TObjString.h>
126 #include <TVirtualMC.h>
127 #include <TVirtualMCApplication.h>
129 #include <TInterpreter.h>
131 #include "AliAlignObj.h"
132 #include "AliCDBEntry.h"
133 #include "AliCDBManager.h"
134 #include "AliGRPManager.h"
135 #include "AliCDBStorage.h"
136 #include "AliCTPRawData.h"
137 #include "AliCentralTrigger.h"
138 #include "AliCentralTrigger.h"
139 #include "AliCodeTimer.h"
141 #include "AliDigitizer.h"
142 #include "AliESDEvent.h"
143 #include "AliGRPObject.h"
144 #include "AliGenEventHeader.h"
145 #include "AliGenerator.h"
146 #include "AliGeomManager.h"
147 #include "AliHLTSimulation.h"
148 #include "AliHeader.h"
150 #include "AliLegoGenerator.h"
154 #include "AliModule.h"
156 #include "AliRawReaderDate.h"
157 #include "AliRawReaderFile.h"
158 #include "AliRawReaderRoot.h"
160 #include "AliDigitizationInput.h"
161 #include "AliRunLoader.h"
162 #include "AliStack.h"
163 #include "AliSimulation.h"
164 #include "AliSysInfo.h"
165 #include "AliVertexGenFile.h"
168 ClassImp(AliSimulation)
170 AliSimulation *AliSimulation::fgInstance = 0;
171 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD",
172 "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE","AD",
175 //_____________________________________________________________________________
176 AliSimulation::AliSimulation(const char* configFileName,
177 const char* name, const char* title) :
180 fRunGeneratorOnly(kFALSE),
181 fRunGeneration(kTRUE),
182 fRunSimulation(kTRUE),
183 fLoadAlignFromCDB(kTRUE),
184 fLoadAlObjsListOfDets("ALL"),
188 fMakeDigitsFromHits(""),
190 fRawDataFileName(""),
191 fDeleteIntermediateFiles(kFALSE),
192 fWriteSelRawData(kFALSE),
193 fStopOnError(kFALSE),
194 fUseMonitoring(kFALSE),
196 fConfigFileName(configFileName),
197 fGAliceFileName("galice.root"),
199 fBkgrdFileNames(NULL),
200 fAlignObjArray(NULL),
201 fUseBkgrdVertex(kTRUE),
202 fRegionOfInterest(kFALSE),
208 fInitCDBCalled(kFALSE),
209 fInitRunNumberCalled(kFALSE),
210 fSetRunNumberFromDataCalled(kFALSE),
211 fEmbeddingFlag(kFALSE),
214 fUseVertexFromCDB(0),
215 fUseMagFieldFromGRP(0),
216 fGRPWriteLocation(Form("local://%s", gSystem->pwd())),
217 fUseTimeStampFromCDB(0),
223 fEventSpecie(AliRecoParam::kDefault),
224 fWriteQAExpertData(kTRUE),
228 fWriteGRPEntry(kTRUE)
230 // create simulation object with default parameters
232 SetGAliceFile("galice.root");
235 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
236 qam->SetActiveDetectors(fQADetectors) ;
237 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
238 qam->SetTasks(fQATasks) ;
241 //_____________________________________________________________________________
242 AliSimulation::~AliSimulation()
246 fEventsPerFile.Delete();
247 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
248 // delete fAlignObjArray; fAlignObjArray=0;
250 if (fBkgrdFileNames) {
251 fBkgrdFileNames->Delete();
252 delete fBkgrdFileNames;
255 fSpecCDBUri.Delete();
256 if (fgInstance==this) fgInstance = 0;
258 AliQAManager::QAManager()->ShowQA() ;
259 AliQAManager::Destroy() ;
260 AliCodeTimer::Instance()->Print();
264 //_____________________________________________________________________________
265 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
267 // set the number of events for one run
272 //_____________________________________________________________________________
273 void AliSimulation::InitQA()
275 // activate a default CDB storage
276 // First check if we have any CDB storage set, because it is used
277 // to retrieve the calibration and alignment constants
279 if (fInitCDBCalled) return;
280 fInitCDBCalled = kTRUE;
282 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
283 qam->SetActiveDetectors(fQADetectors) ;
284 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
285 qam->SetTasks(fQATasks) ;
286 if (fWriteQAExpertData)
287 qam->SetWriteExpert() ;
289 if (qam->IsDefaultStorageSet()) {
290 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
291 AliWarning("Default QA reference storage has been already set !");
292 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
293 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
294 fQARefUri = qam->GetDefaultStorage()->GetURI();
296 if (fQARefUri.Length() > 0) {
297 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
298 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
299 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
301 fQARefUri="local://$ALICE_ROOT/QARef";
302 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
303 AliWarning("Default QA reference storage not yet set !!!!");
304 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
305 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
307 qam->SetDefaultStorage(fQARefUri);
311 //_____________________________________________________________________________
312 void AliSimulation::InitCDB()
314 // activate a default CDB storage
315 // First check if we have any CDB storage set, because it is used
316 // to retrieve the calibration and alignment constants
318 if (fInitCDBCalled) return;
319 fInitCDBCalled = kTRUE;
321 AliCDBManager* man = AliCDBManager::Instance();
322 if (man->IsDefaultStorageSet())
324 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
325 AliWarning("Default CDB storage has been already set !");
326 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
327 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
328 fCDBUri = man->GetDefaultStorage()->GetURI();
331 if (fCDBUri.Length() > 0)
333 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
334 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
335 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
337 fCDBUri="local://$ALICE_ROOT/OCDB";
338 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
339 AliWarning("Default CDB storage not yet set !!!!");
340 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
341 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
344 man->SetDefaultStorage(fCDBUri);
347 // Now activate the detector specific CDB storage locations
348 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
349 TObject* obj = fSpecCDBUri[i];
351 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
352 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
353 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
354 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
359 //_____________________________________________________________________________
360 void AliSimulation::InitRunNumber(){
361 // check run number. If not set, set it to 0 !!!!
363 if (fInitRunNumberCalled) return;
364 fInitRunNumberCalled = kTRUE;
366 AliCDBManager* man = AliCDBManager::Instance();
367 if (man->GetRun() >= 0)
369 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
370 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
374 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
377 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
386 //_____________________________________________________________________________
387 void AliSimulation::SetCDBLock() {
388 // Set CDB lock: from now on it is forbidden to reset the run number
389 // or the default storage or to activate any further storage!
391 ULong64_t key = AliCDBManager::Instance()->SetLock(1);
395 //_____________________________________________________________________________
396 void AliSimulation::SetDefaultStorage(const char* uri) {
397 // Store the desired default CDB storage location
398 // Activate it later within the Run() method
404 //_____________________________________________________________________________
405 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
406 // Store the desired default CDB storage location
407 // Activate it later within the Run() method
410 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
413 //_____________________________________________________________________________
414 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
415 // Store a detector-specific CDB storage location
416 // Activate it later within the Run() method
418 AliCDBPath aPath(calibType);
419 if(!aPath.IsValid()){
420 AliError(Form("Not a valid path: %s", calibType));
424 TObject* obj = fSpecCDBUri.FindObject(calibType);
425 if (obj) fSpecCDBUri.Remove(obj);
426 fSpecCDBUri.Add(new TNamed(calibType, uri));
430 //_____________________________________________________________________________
431 void AliSimulation::SetRunNumber(Int_t run)
434 // Activate it later within the Run() method
439 //_____________________________________________________________________________
440 void AliSimulation::SetSeed(Int_t seed)
443 // Activate it later within the Run() method
448 //_____________________________________________________________________________
449 Bool_t AliSimulation::SetRunNumberFromData()
451 // Set the CDB manager run number
452 // The run number is retrieved from gAlice
454 if (fSetRunNumberFromDataCalled) return kTRUE;
455 fSetRunNumberFromDataCalled = kTRUE;
457 AliCDBManager* man = AliCDBManager::Instance();
458 Int_t runData = -1, runCDB = -1;
460 AliRunLoader* runLoader = LoadRun("READ");
461 if (!runLoader) return kFALSE;
463 runData = runLoader->GetHeader()->GetRun();
467 runCDB = man->GetRun();
469 if (runCDB != runData) {
470 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
471 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
472 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
473 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
478 man->SetRun(runData);
481 if(man->GetRun() < 0) {
482 AliError("Run number not properly initalized!");
491 //_____________________________________________________________________________
492 void AliSimulation::SetConfigFile(const char* fileName)
494 // set the name of the config file
496 fConfigFileName = fileName;
499 //_____________________________________________________________________________
500 void AliSimulation::SetGAliceFile(const char* fileName)
502 // set the name of the galice file
503 // the path is converted to an absolute one if it is relative
505 fGAliceFileName = fileName;
506 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
507 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
509 fGAliceFileName = absFileName;
510 delete[] absFileName;
513 AliDebug(2, Form("galice file name set to %s", fileName));
516 //_____________________________________________________________________________
517 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
520 // set the number of events per file for the given detector and data type
521 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
523 TNamed* obj = new TNamed(detector, type);
524 obj->SetUniqueID(nEvents);
525 fEventsPerFile.Add(obj);
528 //_____________________________________________________________________________
529 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
531 // Read the alignment objects from CDB.
532 // Each detector is supposed to have the
533 // alignment objects in DET/Align/Data CDB path.
534 // All the detector objects are then collected,
535 // sorted by geometry level (starting from ALIC) and
536 // then applied to the TGeo geometry.
537 // Finally an overlaps check is performed.
539 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
540 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
544 // initialize CDB storage, run number, set CDB lock
546 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
549 Bool_t delRunLoader = kFALSE;
551 runLoader = LoadRun("READ");
552 if (!runLoader) return kFALSE;
553 delRunLoader = kTRUE;
556 // Export ideal geometry
557 if(!IsGeometryFromFile()) AliGeomManager::GetGeometry()->Export("geometry.root");
559 // Load alignment data from CDB and apply to geometry through AliGeomManager
560 if(fLoadAlignFromCDB){
562 TString detStr = fLoadAlObjsListOfDets;
563 TString loadAlObjsListOfDets = "";
565 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
566 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
567 AliModule* det = (AliModule*) detArray->At(iDet);
568 if (!det || !det->IsActive()) continue;
569 if (IsSelected(det->GetName(), detStr)) {
570 //add det to list of dets to be aligned from CDB
571 loadAlObjsListOfDets += det->GetName();
572 loadAlObjsListOfDets += " ";
574 } // end loop over detectors
575 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
576 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
578 // Check if the array with alignment objects was
579 // provided by the user. If yes, apply the objects
580 // to the present TGeo geometry
581 if (fAlignObjArray) {
582 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
583 AliError("The misalignment of one or more volumes failed!"
584 "Compare the list of simulated detectors and the list of detector alignment data!");
585 if (delRunLoader) delete runLoader;
591 // Update the internal geometry of modules (ITS needs it)
592 TString detStr = fLoadAlObjsListOfDets;
593 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
594 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
596 AliModule* det = (AliModule*) detArray->At(iDet);
597 if (!det || !det->IsActive()) continue;
598 if (IsSelected(det->GetName(), detStr)) {
599 det->UpdateInternalGeometry();
601 } // end loop over detectors
604 if (delRunLoader) delete runLoader;
609 //_____________________________________________________________________________
610 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
612 // add a file with background events for merging
614 TObjString* fileNameStr = new TObjString(fileName);
615 fileNameStr->SetUniqueID(nSignalPerBkgrd);
616 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
617 fBkgrdFileNames->Add(fileNameStr);
620 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
622 // add a file with background events for embeddin
623 MergeWith(fileName, nSignalPerBkgrd);
624 fEmbeddingFlag = kTRUE;
627 //_____________________________________________________________________________
628 Bool_t AliSimulation::Run(Int_t nEvents)
630 // run the generation, simulation and digitization
633 AliCodeTimerAuto("",0)
634 AliSysInfo::AddStamp("Start_Run");
636 // Load run number and seed from environmental vars
637 ProcessEnvironmentVars();
638 AliSysInfo::AddStamp("ProcessEnvironmentVars");
640 gRandom->SetSeed(fSeed);
642 if (nEvents > 0) fNEvents = nEvents;
644 // Run generator-only code on demand
645 if (fRunGeneratorOnly)
647 if(!RunGeneratorOnly())
649 if (fStopOnError) return kFALSE;
655 // create and setup the HLT instance
656 if (!fRunHLT.IsNull() && !CreateHLT()) {
657 if (fStopOnError) return kFALSE;
662 // generation and simulation -> hits
663 if (fRunGeneration) {
664 if (!RunSimulation()) if (fStopOnError) return kFALSE;
666 AliSysInfo::AddStamp("RunSimulation");
668 // initialize CDB storage from external environment
669 // (either CDB manager or AliSimulation setters),
670 // if not already done in RunSimulation()
672 AliSysInfo::AddStamp("InitCDB");
674 // Set run number in CDBManager from data
675 // From this point on the run number must be always loaded from data!
676 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
678 // Set CDB lock: from now on it is forbidden to reset the run number
679 // or the default storage or to activate any further storage!
682 // If RunSimulation was not called, load the geometry and misalign it
683 if (!AliGeomManager::GetGeometry()) {
684 // Initialize the geometry manager
685 AliGeomManager::LoadGeometry("geometry.root");
686 AliSysInfo::AddStamp("GetGeometry");
687 // // Check that the consistency of symbolic names for the activated subdetectors
688 // // in the geometry loaded by AliGeomManager
689 // AliRunLoader* runLoader = LoadRun("READ");
690 // if (!runLoader) return kFALSE;
692 // TString detsToBeChecked = "";
693 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
694 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
695 // AliModule* det = (AliModule*) detArray->At(iDet);
696 // if (!det || !det->IsActive()) continue;
697 // detsToBeChecked += det->GetName();
698 // detsToBeChecked += " ";
699 // } // end loop over detectors
700 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
701 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
702 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
704 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
706 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
708 AliSysInfo::AddStamp("MissalignGeometry");
711 // hits -> summable digits
712 AliSysInfo::AddStamp("Start_sdigitization");
713 if (!fMakeSDigits.IsNull()) {
714 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
717 AliSysInfo::AddStamp("Stop_sdigitization");
719 AliSysInfo::AddStamp("Start_digitization");
720 // summable digits -> digits
721 if (!fMakeDigits.IsNull()) {
722 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
723 if (fStopOnError) return kFALSE;
726 AliSysInfo::AddStamp("Stop_digitization");
731 if (!fMakeDigitsFromHits.IsNull()) {
732 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
733 AliWarning(Form("Merging and direct creation of digits from hits "
734 "was selected for some detectors. "
735 "No merging will be done for the following detectors: %s",
736 fMakeDigitsFromHits.Data()));
738 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
739 if (fStopOnError) return kFALSE;
743 AliSysInfo::AddStamp("Hits2Digits");
746 // digits -> trigger. Set trigger configuration to "none" - any
747 // case - to not generate the trigger information. Set the trigger
748 // configuration to some string X to read from file at
749 // $ALICE_ROOT/GRP/CTP/X. Set the trigger configuration to the
750 // empty string or "ocdb" - any case - to read from OCDB.
751 if (!fTriggerConfig.EqualTo("none",TString::kIgnoreCase) &&
752 !RunTrigger(fTriggerConfig,fMakeDigits)) {
753 if (fStopOnError) return kFALSE;
756 AliSysInfo::AddStamp("RunTrigger");
759 // digits -> raw data
760 if (!fWriteRawData.IsNull()) {
761 if (!WriteRawData(fWriteRawData, fRawDataFileName,
762 fDeleteIntermediateFiles,fWriteSelRawData)) {
763 if (fStopOnError) return kFALSE;
767 AliSysInfo::AddStamp("WriteRaw");
769 // run HLT simulation on simulated digit data if raw data is not
770 // simulated, otherwise its called as part of WriteRawData
771 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
773 if (fStopOnError) return kFALSE;
777 AliSysInfo::AddStamp("RunHLT");
781 Bool_t rv = RunQA() ;
787 AliSysInfo::AddStamp("RunQA");
791 TString snapshotFileOut("");
792 if(TString(gSystem->Getenv("OCDB_SNAPSHOT_CREATE")) == TString("kTRUE")){
793 AliInfo(" ******** Creating the snapshot! *********");
794 TString snapshotFile(gSystem->Getenv("OCDB_SNAPSHOT_FILENAME"));
795 if(!(snapshotFile.IsNull() || snapshotFile.IsWhitespace()))
796 snapshotFileOut = snapshotFile;
798 snapshotFileOut="OCDB.root";
799 AliCDBManager::Instance()->DumpToSnapshotFile(snapshotFileOut.Data(),kFALSE);
802 // Cleanup of CDB manager: cache and active storages!
803 AliCDBManager::Instance()->ClearCache();
808 //_______________________________________________________________________
809 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
810 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
811 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
814 // Generates lego plots of:
815 // - radiation length map phi vs theta
816 // - radiation length map phi vs eta
817 // - interaction length map
818 // - g/cm2 length map
820 // ntheta bins in theta, eta
821 // themin minimum angle in theta (degrees)
822 // themax maximum angle in theta (degrees)
824 // phimin minimum angle in phi (degrees)
825 // phimax maximum angle in phi (degrees)
826 // rmin minimum radius
827 // rmax maximum radius
830 // The number of events generated = ntheta*nphi
831 // run input parameters in macro setup (default="Config.C")
833 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
836 <img src="picts/AliRunLego1.gif">
841 <img src="picts/AliRunLego2.gif">
846 <img src="picts/AliRunLego3.gif">
851 // run the generation and simulation
853 AliCodeTimerAuto("",0)
855 // initialize CDB storage and run number from external environment
856 // (either CDB manager or AliSimulation setters)
862 AliError("no gAlice object. Restart aliroot and try again.");
865 if (gAlice->Modules()->GetEntries() > 0) {
866 AliError("gAlice was already run. Restart aliroot and try again.");
869 AliInfo(Form("initializing gAlice with config file %s",
870 fConfigFileName.Data()));
873 if (nev == -1) nev = nc1 * nc2;
875 // check if initialisation has been done
876 // If runloader has been initialized, set the number of events per file to nc1 * nc2
879 if (!gener) gener = new AliLegoGenerator();
881 // Configure Generator
883 gener->SetRadiusRange(rmin, rmax);
884 gener->SetZMax(zmax);
885 gener->SetCoor1Range(nc1, c1min, c1max);
886 gener->SetCoor2Range(nc2, c2min, c2max);
890 fLego = new AliLego("lego",gener);
892 //__________________________________________________________________________
896 // - cholm - Add this here for consitency
897 // If requested set the mag. field from the GRP entry.
898 // After this the field is loccked and cannot be changed by Config.C
899 if (fUseMagFieldFromGRP) {
903 AliInfo("Field is locked now. It cannot be changed in Config.C");
907 gROOT->LoadMacro(setup);
908 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
910 if(AliCDBManager::Instance()->GetRun() >= 0) {
911 SetRunNumber(AliCDBManager::Instance()->GetRun());
913 AliWarning("Run number not initialized!!");
916 AliRunLoader::Instance()->CdGAFile();
918 AliPDG::AddParticlesToPdgDataBase();
920 TVirtualMC::GetMC()->SetMagField(TGeoGlobalMagField::Instance()->GetField());
922 gAlice->GetMCApp()->Init();
925 //Must be here because some MCs (G4) adds detectors here and not in Config.C
926 gAlice->InitLoaders();
927 AliRunLoader::Instance()->MakeTree("E");
930 // Save stuff at the beginning of the file to avoid file corruption
931 AliRunLoader::Instance()->CdGAFile();
934 //Save current generator
935 AliGenerator *gen=gAlice->GetMCApp()->Generator();
936 gAlice->GetMCApp()->ResetGenerator(gener);
937 //Prepare MC for Lego Run
938 TVirtualMC::GetMC()->InitLego();
943 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
944 TVirtualMC::GetMC()->ProcessRun(nev);
946 // End of this run, close files
948 // Restore current generator
949 gAlice->GetMCApp()->ResetGenerator(gen);
950 // Delete Lego Object
956 //_____________________________________________________________________________
957 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
961 AliCodeTimerAuto("",0)
963 // initialize CDB storage from external environment
964 // (either CDB manager or AliSimulation setters),
965 // if not already done in RunSimulation()
968 // Set run number in CDBManager from data
969 // From this point on the run number must be always loaded from data!
970 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
972 // Set CDB lock: from now on it is forbidden to reset the run number
973 // or the default storage or to activate any further storage!
976 AliRunLoader* runLoader = LoadRun("READ");
977 if (!runLoader) return kFALSE;
978 TString trconfiguration = config;
980 if (trconfiguration.IsNull()) {
981 if(!fTriggerConfig.IsNull()) {
982 trconfiguration = fTriggerConfig;
985 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
988 runLoader->MakeTree( "GG" );
989 AliCentralTrigger* aCTP = runLoader->GetTrigger();
990 // Load Configuration
991 if (!aCTP->LoadConfiguration( trconfiguration ))
995 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
1007 //_____________________________________________________________________________
1008 Bool_t AliSimulation::WriteTriggerRawData()
1010 // Writes the CTP (trigger) DDL raw data
1011 // Details of the format are given in the
1012 // trigger TDR - pages 134 and 135.
1013 AliCTPRawData writer;
1015 writer.RawDataRun2();
1020 //_____________________________________________________________________________
1021 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
1023 // run the generation and simulation
1025 AliCodeTimerAuto("",0)
1027 // initialize CDB storage and run number from external environment
1028 // (either CDB manager or AliSimulation setters)
1029 AliSysInfo::AddStamp("RunSimulation_Begin");
1031 AliSysInfo::AddStamp("RunSimulation_InitCDB");
1035 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
1038 AliError("no gAlice object. Restart aliroot and try again.");
1041 if (gAlice->Modules()->GetEntries() > 0) {
1042 AliError("gAlice was already run. Restart aliroot and try again.");
1046 // Setup monitoring if requested
1047 gAlice->GetMCApp()->SetUseMonitoring(fUseMonitoring);
1049 AliInfo(Form("initializing gAlice with config file %s",
1050 fConfigFileName.Data()));
1053 // Initialize ALICE Simulation run
1058 // If requested set the mag. field from the GRP entry.
1059 // After this the field is loccked and cannot be changed by Config.C
1060 if (fUseMagFieldFromGRP) {
1062 grpM.ReadGRPEntry();
1064 AliInfo("Field is locked now. It cannot be changed in Config.C");
1068 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1069 gROOT->LoadMacro(fConfigFileName.Data());
1070 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1071 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1072 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1074 AliSysInfo::AddStamp("RunSimulation_Config");
1077 // If requested obtain the vertex position and vertex sigma_z from the CDB
1078 // This overwrites the settings from the Config.C
1079 if (fUseVertexFromCDB) {
1080 Double_t vtxPos[3] = {0., 0., 0.};
1081 Double_t vtxSig[3] = {0., 0., 0.};
1082 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1084 AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1086 if(vertex->GetXRes()>2.8) { // > pipe radius --> it's a dummy object, don't use it
1087 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1088 if (entry) vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1092 vertex->GetXYZ(vtxPos);
1093 vertex->GetSigmaXYZ(vtxSig);
1094 AliInfo("Overwriting Config.C vertex settings !");
1095 AliInfo(Form("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
1096 vtxPos[0], vtxPos[1], vtxPos[2], vtxSig[2]));
1098 AliGenerator *gen = gAlice->GetMCApp()->Generator();
1099 gen->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]); // vertex position
1100 gen->SetSigmaZ(vtxSig[2]);
1105 // If requested we take the SOR and EOR time-stamps from the GRP and use them
1106 // in order to generate the event time-stamps
1107 if (fUseTimeStampFromCDB) {
1109 grpM.ReadGRPEntry();
1110 const AliGRPObject *grpObj = grpM.GetGRPData();
1111 if (!grpObj || (grpObj->GetTimeEnd() <= grpObj->GetTimeStart())) {
1112 AliError("Missing GRP or bad SOR/EOR time-stamps! Switching off the time-stamp generation from GRP!");
1113 fTimeStart = fTimeEnd = 0;
1114 fUseTimeStampFromCDB = kFALSE;
1117 fTimeStart = grpObj->GetTimeStart();
1118 fTimeEnd = grpObj->GetTimeEnd();
1122 if(AliCDBManager::Instance()->GetRun() >= 0) {
1123 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
1124 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
1126 AliWarning("Run number not initialized!!");
1129 AliRunLoader::Instance()->CdGAFile();
1132 AliPDG::AddParticlesToPdgDataBase();
1134 TVirtualMC::GetMC()->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1135 AliSysInfo::AddStamp("RunSimulation_GetField");
1136 gAlice->GetMCApp()->Init();
1137 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1139 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1140 gAlice->InitLoaders();
1141 AliRunLoader::Instance()->MakeTree("E");
1142 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1143 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1144 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1146 // Save stuff at the beginning of the file to avoid file corruption
1147 AliRunLoader::Instance()->CdGAFile();
1149 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1150 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1151 //___________________________________________________________________________________________
1153 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1155 // Set run number in CDBManager
1156 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1158 AliRunLoader* runLoader = AliRunLoader::Instance();
1160 AliError(Form("gAlice has no run loader object. "
1161 "Check your config file: %s", fConfigFileName.Data()));
1164 SetGAliceFile(runLoader->GetFileName());
1166 // Misalign geometry
1167 #if ROOT_VERSION_CODE < 331527
1168 AliGeomManager::SetGeometry(gGeoManager);
1170 // Check that the consistency of symbolic names for the activated subdetectors
1171 // in the geometry loaded by AliGeomManager
1172 TString detsToBeChecked = "";
1173 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1174 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1175 AliModule* det = (AliModule*) detArray->At(iDet);
1176 if (!det || !det->IsActive()) continue;
1177 detsToBeChecked += det->GetName();
1178 detsToBeChecked += " ";
1179 } // end loop over detectors
1180 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1181 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1182 MisalignGeometry(runLoader);
1183 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1186 // AliRunLoader* runLoader = AliRunLoader::Instance();
1187 // if (!runLoader) {
1188 // AliError(Form("gAlice has no run loader object. "
1189 // "Check your config file: %s", fConfigFileName.Data()));
1192 // SetGAliceFile(runLoader->GetFileName());
1194 if (!gAlice->GetMCApp()->Generator()) {
1195 AliError(Form("gAlice has no generator object. "
1196 "Check your config file: %s", fConfigFileName.Data()));
1200 // Write GRP entry corresponding to the setting found in Cofig.C
1203 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1205 if (nEvents <= 0) nEvents = fNEvents;
1207 // get vertex from background file in case of merging
1208 if (fUseBkgrdVertex &&
1209 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1210 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1211 const char* fileName = ((TObjString*)
1212 (fBkgrdFileNames->At(0)))->GetName();
1213 AliInfo(Form("The vertex will be taken from the background "
1214 "file %s with nSignalPerBackground = %d",
1215 fileName, signalPerBkgrd));
1216 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1217 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1220 if (!fRunSimulation) {
1221 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1224 // set the number of events per file for given detectors and data types
1225 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1226 if (!fEventsPerFile[i]) continue;
1227 const char* detName = fEventsPerFile[i]->GetName();
1228 const char* typeName = fEventsPerFile[i]->GetTitle();
1229 TString loaderName(detName);
1230 loaderName += "Loader";
1231 AliLoader* loader = runLoader->GetLoader(loaderName);
1233 AliError(Form("RunSimulation no loader for %s found\n Number of events per file not set for %s %s",
1234 detName, typeName, detName));
1237 AliDataLoader* dataLoader =
1238 loader->GetDataLoader(typeName);
1240 AliError(Form("no data loader for %s found\n"
1241 "Number of events per file not set for %s %s",
1242 typeName, detName, typeName));
1245 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1246 AliDebug(1, Form("number of events per file set to %d for %s %s",
1247 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1250 AliInfo("running gAlice");
1251 AliSysInfo::AddStamp("Start_ProcessRun");
1253 // Create the Root Tree with one branch per detector
1254 //Hits moved to begin event -> now we are crating separate tree for each event
1255 TVirtualMC::GetMC()->ProcessRun(nEvents);
1257 // End of this run, close files
1258 if(nEvents>0) FinishRun();
1260 AliSysInfo::AddStamp("Stop_ProcessRun");
1266 //_____________________________________________________________________________
1267 Bool_t AliSimulation::RunGeneratorOnly()
1270 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1271 gROOT->LoadMacro(fConfigFileName.Data());
1272 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1273 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1274 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1277 // Setup the runloader and generator, check if everything is OK
1278 AliRunLoader* runLoader = AliRunLoader::Instance();
1279 AliGenerator* generator = gAlice->GetMCApp()->Generator();
1281 AliError(Form("gAlice has no run loader object. "
1282 "Check your config file: %s", fConfigFileName.Data()));
1286 AliError(Form("gAlice has no generator object. "
1287 "Check your config file: %s", fConfigFileName.Data()));
1291 runLoader->LoadKinematics("RECREATE");
1292 runLoader->MakeTree("E");
1294 // Create stack and header
1295 runLoader->MakeStack();
1296 AliStack* stack = runLoader->Stack();
1297 AliHeader* header = runLoader->GetHeader();
1299 // Intialize generator
1301 generator->SetStack(stack);
1303 // Run main generator loop
1305 for (Int_t iev=0; iev<fNEvents; iev++)
1308 header->Reset(0,iev);
1309 runLoader->SetEventNumber(iev);
1311 runLoader->MakeTree("K");
1314 generator->Generate();
1317 header->SetNprimary(stack->GetNprimary());
1318 header->SetNtrack(stack->GetNtrack());
1319 stack->FinishEvent();
1320 header->SetStack(stack);
1321 runLoader->TreeE()->Fill();
1322 runLoader->WriteKinematics("OVERWRITE");
1326 generator->FinishRun();
1328 runLoader->WriteHeader("OVERWRITE");
1335 //_____________________________________________________________________________
1336 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1338 // run the digitization and produce summable digits
1339 static Int_t eventNr=0;
1340 AliCodeTimerAuto("",0) ;
1342 // initialize CDB storage, run number, set CDB lock
1344 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1347 AliRunLoader* runLoader = LoadRun();
1348 if (!runLoader) return kFALSE;
1350 TString detStr = detectors;
1351 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1352 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1353 AliModule* det = (AliModule*) detArray->At(iDet);
1354 if (!det || !det->IsActive()) continue;
1355 if (IsSelected(det->GetName(), detStr)) {
1356 AliInfo(Form("creating summable digits for %s", det->GetName()));
1357 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1358 det->Hits2SDigits();
1359 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1360 AliSysInfo::AddStamp(Form("SDigit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1364 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1365 AliError(Form("the following detectors were not found: %s",
1367 if (fStopOnError) return kFALSE;
1376 //_____________________________________________________________________________
1377 Bool_t AliSimulation::RunDigitization(const char* detectors,
1378 const char* excludeDetectors)
1380 // run the digitization and produce digits from sdigits
1381 AliCodeTimerAuto("",0)
1383 // initialize CDB storage, run number, set CDB lock
1385 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1388 delete AliRunLoader::Instance();
1393 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1394 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1395 AliDigitizationInput digInp(nStreams, signalPerBkgrd);
1396 // digInp.SetEmbeddingFlag(fEmbeddingFlag);
1397 digInp.SetRegionOfInterest(fRegionOfInterest);
1398 digInp.SetInputStream(0, fGAliceFileName.Data());
1399 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1400 const char* fileName = ((TObjString*)(fBkgrdFileNames->At(iStream-1)))->GetName();
1401 digInp.SetInputStream(iStream, fileName);
1404 detArr.SetOwner(kTRUE);
1405 TString detStr = detectors;
1406 TString detExcl = excludeDetectors;
1407 if (!static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice()) {
1408 AliError("Error occured while getting gAlice from Input 0");
1411 AliRunLoader* runLoader = AliRunLoader::GetRunLoader(digInp.GetInputStream(0)->GetFolderName());
1412 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1413 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1414 AliModule* det = (AliModule*) detArray->At(iDet);
1415 if (!det || !det->IsActive()) continue;
1416 if (!IsSelected(det->GetName(), detStr) || IsSelected(det->GetName(), detExcl)) continue;
1417 AliDigitizer* digitizer = det->CreateDigitizer(&digInp);
1418 if (!digitizer || !digitizer->Init()) {
1419 AliError(Form("no digitizer for %s", det->GetName()));
1420 if (fStopOnError) return kFALSE;
1423 detArr.AddLast(digitizer);
1424 AliInfo(Form("Created digitizer from SDigits -> Digits for %s", det->GetName()));
1428 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1429 AliError(Form("the following detectors were not found: %s", detStr.Data()));
1430 if (fStopOnError) return kFALSE;
1433 Int_t ndigs = detArr.GetEntriesFast();
1434 Int_t eventsCreated = 0;
1435 AliRunLoader* outRl = digInp.GetOutRunLoader();
1436 while ((eventsCreated++ < fNEvents) || (fNEvents < 0)) {
1437 if (!digInp.ConnectInputTrees()) break;
1438 digInp.InitEvent(); //this must be after call of Connect Input tress.
1439 if (outRl) outRl->SetEventNumber(eventsCreated-1);
1440 static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice(); // use gAlice of the first input stream
1441 for (int id=0;id<ndigs;id++) {
1442 ((AliDigitizer*)detArr[id])->Digitize("");
1443 AliSysInfo::AddStamp(Form("Digit_%s_%d",detArr[id]->GetName(),eventsCreated), 0,2, eventsCreated);
1445 digInp.FinishEvent();
1447 digInp.FinishGlobal();
1452 //_____________________________________________________________________________
1453 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1455 // run the digitization and produce digits from hits
1457 AliCodeTimerAuto("",0)
1459 // initialize CDB storage, run number, set CDB lock
1461 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1464 AliRunLoader* runLoader = LoadRun("READ");
1465 if (!runLoader) return kFALSE;
1467 TString detStr = detectors;
1468 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1469 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1470 AliModule* det = (AliModule*) detArray->At(iDet);
1471 if (!det || !det->IsActive()) continue;
1472 if (IsSelected(det->GetName(), detStr)) {
1473 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1478 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1479 AliError(Form("the following detectors were not found: %s",
1481 if (fStopOnError) return kFALSE;
1487 //_____________________________________________________________________________
1488 Bool_t AliSimulation::WriteRawData(const char* detectors,
1489 const char* fileName,
1490 Bool_t deleteIntermediateFiles,
1493 // convert the digits to raw data
1494 // First DDL raw data files for the given detectors are created.
1495 // If a file name is given, the DDL files are then converted to a DATE file.
1496 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1498 // If the file name has the extension ".root", the DATE file is converted
1500 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1501 // 'selrawdata' flag can be used to enable writing of detectors raw data
1502 // accoring to the trigger cluster.
1504 AliCodeTimerAuto("",0)
1505 AliSysInfo::AddStamp("WriteRawData_Start");
1507 TString detStr = detectors;
1508 if (!WriteRawFiles(detStr.Data())) {
1509 if (fStopOnError) return kFALSE;
1511 AliSysInfo::AddStamp("WriteRawFiles");
1513 // run HLT simulation on simulated DDL raw files
1514 // and produce HLT ddl raw files to be included in date/root file
1515 // bugfix 2009-06-26: the decision whether to write HLT raw data
1516 // is taken in RunHLT. Here HLT always needs to be run in order to
1517 // create HLT digits, unless its switched off. This is due to the
1518 // special placement of the HLT between the generation of DDL files
1519 // and conversion to DATE/Root file.
1520 detStr.ReplaceAll("HLT", "");
1521 if (!fRunHLT.IsNull()) {
1523 if (fStopOnError) return kFALSE;
1526 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1528 TString dateFileName(fileName);
1529 if (!dateFileName.IsNull()) {
1530 Bool_t rootOutput = dateFileName.EndsWith(".root");
1531 if (rootOutput) dateFileName += ".date";
1532 TString selDateFileName;
1534 selDateFileName = "selected.";
1535 selDateFileName+= dateFileName;
1537 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1538 if (fStopOnError) return kFALSE;
1540 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1541 if (deleteIntermediateFiles) {
1542 AliRunLoader* runLoader = LoadRun("READ");
1543 if (runLoader) for (Int_t iEvent = 0;
1544 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1546 snprintf(command, 256, "rm -r raw%d", iEvent);
1547 gSystem->Exec(command);
1553 if (!ConvertDateToRoot(dateFileName, fileName)) {
1554 if (fStopOnError) return kFALSE;
1556 AliSysInfo::AddStamp("ConvertDateToRoot");
1557 if (deleteIntermediateFiles) {
1558 gSystem->Unlink(dateFileName);
1561 TString selFileName = "selected.";
1562 selFileName += fileName;
1563 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1564 if (fStopOnError) return kFALSE;
1566 if (deleteIntermediateFiles) {
1567 gSystem->Unlink(selDateFileName);
1576 //_____________________________________________________________________________
1577 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1579 // convert the digits to raw data DDL files
1581 AliCodeTimerAuto("",0)
1583 AliRunLoader* runLoader = LoadRun("READ");
1584 if (!runLoader) return kFALSE;
1586 // write raw data to DDL files
1587 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1588 AliInfo(Form("processing event %d", iEvent));
1589 runLoader->GetEvent(iEvent);
1590 TString baseDir = gSystem->WorkingDirectory();
1592 snprintf(dirName, 256, "raw%d", iEvent);
1593 gSystem->MakeDirectory(dirName);
1594 if (!gSystem->ChangeDirectory(dirName)) {
1595 AliError(Form("couldn't change to directory %s", dirName));
1596 if (fStopOnError) return kFALSE; else continue;
1599 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1602 TString detStr = detectors;
1603 if (IsSelected("HLT", detStr)) {
1604 // Do nothing. "HLT" will be removed from detStr and HLT raw
1605 // data files are generated in RunHLT.
1608 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1609 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1610 AliModule* det = (AliModule*) detArray->At(iDet);
1611 if (!det || !det->IsActive()) continue;
1612 if (IsSelected(det->GetName(), detStr)) {
1613 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1618 if (!WriteTriggerRawData())
1619 if (fStopOnError) return kFALSE;
1621 gSystem->ChangeDirectory(baseDir);
1622 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1623 AliError(Form("the following detectors were not found: %s",
1625 if (fStopOnError) return kFALSE;
1634 //_____________________________________________________________________________
1635 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1636 const char* selDateFileName)
1638 // convert raw data DDL files to a DATE file with the program "dateStream"
1639 // The second argument is not empty when the user decides to write
1640 // the detectors raw data according to the trigger cluster.
1642 AliCodeTimerAuto("",0)
1644 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1646 AliError("the program dateStream was not found");
1647 if (fStopOnError) return kFALSE;
1652 AliRunLoader* runLoader = LoadRun("READ");
1653 if (!runLoader) return kFALSE;
1655 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1656 Bool_t selrawdata = kFALSE;
1657 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1660 // Note the option -s. It is used in order to avoid
1661 // the generation of SOR/EOR events.
1662 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1663 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1664 FILE* pipe = gSystem->OpenPipe(command, "w");
1667 AliError(Form("Cannot execute command: %s",command));
1671 Int_t selEvents = 0;
1672 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1674 UInt_t detectorPattern = 0;
1675 runLoader->GetEvent(iEvent);
1676 if (!runLoader->LoadTrigger()) {
1677 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1678 detectorPattern = aCTP->GetClusterMask();
1679 // Check if the event was triggered by CTP
1681 if (aCTP->GetClassMask()) selEvents++;
1685 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1687 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1688 selrawdata = kFALSE;
1692 fprintf(pipe, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1696 // loop over detectors and DDLs
1697 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1698 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1700 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1701 Int_t ldcID = Int_t(ldc + 0.0001);
1702 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1704 char rawFileName[256];
1705 snprintf(rawFileName, 256, "raw%d/%s",
1706 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1708 // check existence and size of raw data file
1709 FILE* file = fopen(rawFileName, "rb");
1710 if (!file) continue;
1711 fseek(file, 0, SEEK_END);
1712 unsigned long size = ftell(file);
1714 if (!size) continue;
1716 if (ldcID != prevLDC) {
1717 fprintf(pipe, " LDC Id %d\n", ldcID);
1720 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1725 Int_t result = gSystem->ClosePipe(pipe);
1727 if (!(selrawdata && selEvents > 0)) {
1729 return (result == 0);
1732 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1734 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1735 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1736 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1738 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1740 // Get the trigger decision and cluster
1741 UInt_t detectorPattern = 0;
1743 runLoader->GetEvent(iEvent);
1744 if (!runLoader->LoadTrigger()) {
1745 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1746 if (aCTP->GetClassMask() == 0) continue;
1747 detectorPattern = aCTP->GetClusterMask();
1748 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1749 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1752 fprintf(pipe2, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1756 // loop over detectors and DDLs
1757 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1758 // Write only raw data from detectors that
1759 // are contained in the trigger cluster(s)
1760 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1762 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1764 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1765 Int_t ldcID = Int_t(ldc + 0.0001);
1766 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1768 char rawFileName[256];
1769 snprintf(rawFileName, 256, "raw%d/%s",
1770 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1772 // check existence and size of raw data file
1773 FILE* file = fopen(rawFileName, "rb");
1774 if (!file) continue;
1775 fseek(file, 0, SEEK_END);
1776 unsigned long size = ftell(file);
1778 if (!size) continue;
1780 if (ldcID != prevLDC) {
1781 fprintf(pipe2, " LDC Id %d\n", ldcID);
1784 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1789 Int_t result2 = gSystem->ClosePipe(pipe2);
1792 return ((result == 0) && (result2 == 0));
1795 //_____________________________________________________________________________
1796 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1797 const char* rootFileName)
1799 // convert a DATE file to a root file with the program "alimdc"
1802 const Int_t kDBSize = 2000000000;
1803 const Int_t kTagDBSize = 1000000000;
1804 const Bool_t kFilter = kFALSE;
1805 const Int_t kCompression = 1;
1807 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1809 AliError("the program alimdc was not found");
1810 if (fStopOnError) return kFALSE;
1815 AliInfo(Form("converting DATE file %s to root file %s",
1816 dateFileName, rootFileName));
1818 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1819 const char* tagDBFS = "/tmp/mdc1/tags";
1821 // User defined file system locations
1822 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1823 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1824 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1825 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1826 if (gSystem->Getenv("ALIMDC_TAGDB"))
1827 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1829 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1830 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1831 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1833 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1834 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1835 gSystem->Exec(Form("mkdir %s",tagDBFS));
1837 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1838 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1839 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1841 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1842 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1843 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1845 return (result == 0);
1849 //_____________________________________________________________________________
1850 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1852 // delete existing run loaders, open a new one and load gAlice
1854 delete AliRunLoader::Instance();
1855 AliRunLoader* runLoader =
1856 AliRunLoader::Open(fGAliceFileName.Data(),
1857 AliConfig::GetDefaultEventFolderName(), mode);
1859 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1862 runLoader->LoadgAlice();
1863 runLoader->LoadHeader();
1864 gAlice = runLoader->GetAliRun();
1866 AliError(Form("no gAlice object found in file %s",
1867 fGAliceFileName.Data()));
1873 //_____________________________________________________________________________
1874 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1876 // get or calculate the number of signal events per background event
1878 if (!fBkgrdFileNames) return 1;
1879 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1880 if (nBkgrdFiles == 0) return 1;
1882 // get the number of signal events
1884 AliRunLoader* runLoader =
1885 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1886 if (!runLoader) return 1;
1888 nEvents = runLoader->GetNumberOfEvents();
1893 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1894 // get the number of background events
1895 const char* fileName = ((TObjString*)
1896 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1897 AliRunLoader* runLoader =
1898 AliRunLoader::Open(fileName, "BKGRD");
1899 if (!runLoader) continue;
1900 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1903 // get or calculate the number of signal per background events
1904 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1905 if (nSignalPerBkgrd <= 0) {
1906 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1907 } else if (result && (result != nSignalPerBkgrd)) {
1908 AliInfo(Form("the number of signal events per background event "
1909 "will be changed from %d to %d for stream %d",
1910 nSignalPerBkgrd, result, iBkgrdFile+1));
1911 nSignalPerBkgrd = result;
1914 if (!result) result = nSignalPerBkgrd;
1915 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1916 AliWarning(Form("not enough background events (%d) for %d signal events "
1917 "using %d signal per background events for stream %d",
1918 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1925 //_____________________________________________________________________________
1926 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1928 // check whether detName is contained in detectors
1929 // if yes, it is removed from detectors
1931 // check if all detectors are selected
1932 if ((detectors.CompareTo("ALL") == 0) ||
1933 detectors.BeginsWith("ALL ") ||
1934 detectors.EndsWith(" ALL") ||
1935 detectors.Contains(" ALL ")) {
1940 // search for the given detector
1941 Bool_t result = kFALSE;
1942 if ((detectors.CompareTo(detName) == 0) ||
1943 detectors.BeginsWith(detName+" ") ||
1944 detectors.EndsWith(" "+detName) ||
1945 detectors.Contains(" "+detName+" ")) {
1946 detectors.ReplaceAll(detName, "");
1950 // clean up the detectors string
1951 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1952 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1953 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1958 //_____________________________________________________________________________
1959 Int_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName, Int_t N, Int_t nSkip)
1962 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1963 // These can be used for embedding of MC tracks into RAW data using the standard
1964 // merging procedure.
1966 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1969 AliError("no gAlice object. Restart aliroot and try again.");
1972 if (gAlice->Modules()->GetEntries() > 0) {
1973 AliError("gAlice was already run. Restart aliroot and try again.");
1977 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1981 gROOT->LoadMacro(fConfigFileName.Data());
1982 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1984 if(AliCDBManager::Instance()->GetRun() >= 0) {
1985 SetRunNumber(AliCDBManager::Instance()->GetRun());
1987 AliWarning("Run number not initialized!!");
1990 AliRunLoader::Instance()->CdGAFile();
1992 AliPDG::AddParticlesToPdgDataBase();
1994 TVirtualMC::GetMC()->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1996 gAlice->GetMCApp()->Init();
1998 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1999 gAlice->InitLoaders();
2000 AliRunLoader::Instance()->MakeTree("E");
2001 AliRunLoader::Instance()->LoadKinematics("RECREATE");
2002 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
2003 AliRunLoader::Instance()->LoadHits("all","RECREATE");
2006 // Save stuff at the beginning of the file to avoid file corruption
2007 AliRunLoader::Instance()->CdGAFile();
2012 //AliCDBManager* man = AliCDBManager::Instance();
2013 //man->SetRun(0); // Should this come from rawdata header ?
2017 // Get the runloader
2018 AliRunLoader* runLoader = AliRunLoader::Instance();
2020 // Open esd file if available
2023 AliESDEvent* esd = 0;
2024 if (esdFileName && (strlen(esdFileName)>0)) {
2025 esdFile = TFile::Open(esdFileName);
2027 esd = new AliESDEvent();
2028 esdFile->GetObject("esdTree", treeESD);
2030 esd->ReadFromTree(treeESD);
2032 AliInfo(Form("Asking to skip first %d ESDs events",nSkip));
2041 // Create the RawReader
2042 TString fileName(rawDirectory);
2043 AliRawReader* rawReader = AliRawReader::Create(fileName.Data());
2044 if (!rawReader) return (kFALSE);
2046 // if (!fEquipIdMap.IsNull() && fRawReader)
2047 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
2049 // Get list of detectors
2050 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
2053 AliHeader* header = runLoader->GetHeader();
2057 if (!(rawReader->NextEvent())) break;
2058 runLoader->SetEventNumber(nev);
2059 runLoader->GetHeader()->Reset(rawReader->GetRunNumber(),
2061 runLoader->GetEvent(nev);
2062 AliInfo(Form("We are at event %d",nev));
2065 TString detStr = fMakeSDigits;
2066 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
2067 AliModule* det = (AliModule*) detArray->At(iDet);
2068 if (!det || !det->IsActive()) continue;
2069 if (IsSelected(det->GetName(), detStr)) {
2070 AliInfo(Form("Calling Raw2SDigits for %s", det->GetName()));
2071 det->Raw2SDigits(rawReader);
2078 // If ESD information available obtain reconstructed vertex and store in header.
2080 AliInfo(Form("Selected event %d correspond to event %d in raw and to %d in esd",nev,rawReader->GetEventIndex(),nSkip+rawReader->GetEventIndex()));
2081 treeESD->GetEvent(nSkip+rawReader->GetEventIndex());
2082 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
2083 Double_t position[3];
2084 esdVertex->GetXYZ(position);
2085 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
2088 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
2089 mcHeader->SetPrimaryVertex(mcV);
2090 header->Reset(0,nev);
2091 header->SetGenEventHeader(mcHeader);
2092 AliInfo(Form("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]));
2096 runLoader->TreeE()->Fill();
2097 AliInfo(Form("Finished event %d",nev));
2106 runLoader->CdGAFile();
2107 runLoader->WriteHeader("OVERWRITE");
2108 runLoader->WriteRunLoader();
2113 //_____________________________________________________________________________
2114 void AliSimulation::FinishRun()
2117 // Called at the end of the run.
2122 AliDebug(1, "Finish Lego");
2123 AliRunLoader::Instance()->CdGAFile();
2127 // Clean detector information
2128 TIter next(gAlice->Modules());
2129 AliModule *detector;
2130 while((detector = dynamic_cast<AliModule*>(next()))) {
2131 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
2132 detector->FinishRun();
2135 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
2136 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
2138 // Write AliRun info and all detectors parameters
2139 AliRunLoader::Instance()->CdGAFile();
2140 gAlice->Write(0,TObject::kOverwrite);//write AliRun
2141 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
2143 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
2144 AliRunLoader::Instance()->Synchronize();
2147 //_____________________________________________________________________________
2148 Int_t AliSimulation::GetDetIndex(const char* detector)
2150 // return the detector index corresponding to detector
2152 for (index = 0; index < fgkNDetectors ; index++) {
2153 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2159 //_____________________________________________________________________________
2160 Bool_t AliSimulation::CreateHLT()
2162 // Init the HLT simulation.
2163 // The function loads the library and creates the instance of AliHLTSimulation.
2164 // the main reason for the decoupled creation is to set the transient OCDB
2165 // objects before the OCDB is locked
2167 // load the library dynamically
2168 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
2170 // check for the library version
2171 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
2173 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
2176 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
2177 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
2180 // print compile info
2181 typedef void (*CompileInfo)( const char*& date, const char*& time);
2182 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
2184 const char* date="";
2185 const char* time="";
2186 (*fctInfo)(date, time);
2187 if (!date) date="unknown";
2188 if (!time) time="unknown";
2189 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
2191 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
2194 // create instance of the HLT simulation
2195 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
2196 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
2197 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
2201 TString specObjects;
2202 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
2203 if (specObjects.Length()>0) specObjects+=" ";
2204 specObjects+=fSpecCDBUri[i]->GetName();
2207 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
2208 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
2209 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2215 //_____________________________________________________________________________
2216 Bool_t AliSimulation::RunHLT()
2218 // Run the HLT simulation
2219 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2220 // Disabled if fRunHLT is empty, default vaule is "default".
2221 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2222 // The default simulation depends on the HLT component libraries and their
2223 // corresponding agents which define components and chains to run. See
2224 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2225 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2227 // The libraries to be loaded can be specified as an option.
2229 // AliSimulation sim;
2230 // sim.SetRunHLT("libAliHLTSample.so");
2232 // will only load <tt>libAliHLTSample.so</tt>
2234 // Other available options:
2235 // \li loglevel=<i>level</i> <br>
2236 // logging level for this processing
2238 // disable redirection of log messages to AliLog class
2239 // \li config=<i>macro</i>
2240 // configuration macro
2241 // \li chains=<i>configuration</i>
2242 // comma separated list of configurations to be run during simulation
2243 // \li rawfile=<i>file</i>
2244 // source for the RawReader to be created, the default is <i>./</i> if
2245 // raw data is simulated
2249 if (!fpHLT && !CreateHLT()) {
2252 AliHLTSimulation* pHLT=fpHLT;
2254 AliRunLoader* pRunLoader = LoadRun("READ");
2255 if (!pRunLoader) return kFALSE;
2257 // initialize CDB storage, run number, set CDB lock
2258 // thats for the case of running HLT simulation without all the other steps
2259 // multiple calls are handled by the function, so we can just call
2261 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2264 // init the HLT simulation
2266 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2267 TString detStr = fWriteRawData;
2268 if (!IsSelected("HLT", detStr)) {
2269 options+=" writerawfiles=";
2271 options+=" writerawfiles=HLT";
2274 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2275 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2276 // in order to get detector data. By default, RawReaderFile is used to read
2277 // the already simulated ddl files. Date and Root files from the raw data
2278 // are generated after the HLT simulation.
2279 options+=" rawfile=./";
2282 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2283 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2284 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2286 // run the HLT simulation
2287 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2288 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2289 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2293 // delete the instance
2294 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2295 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2296 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2300 return iResult>=0?kTRUE:kFALSE;
2303 //_____________________________________________________________________________
2304 Bool_t AliSimulation::RunQA()
2306 // run the QA on summable hits, digits or digits
2308 //if(!gAlice) return kFALSE;
2309 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2311 TString detectorsw("") ;
2313 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2314 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2315 if ( detectorsw.IsNull() )
2320 //_____________________________________________________________________________
2321 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2323 // Allows to run QA for a selected set of detectors
2324 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2325 // all selected detectors run the same selected tasks
2327 if (!detAndAction.Contains(":")) {
2328 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2332 Int_t colon = detAndAction.Index(":") ;
2333 fQADetectors = detAndAction(0, colon) ;
2334 if (fQADetectors.Contains("ALL") ){
2335 TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2336 Int_t minus = fQADetectors.Last('-') ;
2337 TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2338 TString toRemove("") ;
2339 while (minus >= 0) {
2340 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
2341 toRemove = toRemove.Strip() ;
2342 toKeep.ReplaceAll(toRemove, "") ;
2343 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
2344 minus = fQADetectors.Last('-') ;
2346 fQADetectors = toKeep ;
2348 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2349 if (fQATasks.Contains("ALL") ) {
2350 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2352 fQATasks.ToUpper() ;
2354 if ( fQATasks.Contains("HIT") )
2355 tempo = Form("%d ", AliQAv1::kHITS) ;
2356 if ( fQATasks.Contains("SDIGIT") )
2357 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2358 if ( fQATasks.Contains("DIGIT") )
2359 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2361 if (fQATasks.IsNull()) {
2362 AliInfo("No QA requested\n") ;
2367 TString tempo(fQATasks) ;
2368 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2369 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2370 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2371 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2373 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2374 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2375 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2376 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2381 //_____________________________________________________________________________
2382 void AliSimulation::ProcessEnvironmentVars()
2384 // Extract run number and random generator seed from env variables
2386 AliInfo("Processing environment variables");
2388 // Random Number seed
2390 // first check that seed is not already set
2392 if (gSystem->Getenv("CONFIG_SEED")) {
2393 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2396 if (gSystem->Getenv("CONFIG_SEED")) {
2397 AliInfo(Form("Seed for random number generation already set (%d)"
2398 ": CONFIG_SEED variable ignored!", fSeed));
2402 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2406 // first check that run number is not already set
2408 if (gSystem->Getenv("DC_RUN")) {
2409 fRun = atoi(gSystem->Getenv("DC_RUN"));
2412 if (gSystem->Getenv("DC_RUN")) {
2413 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2417 AliInfo(Form("Run number = %d", fRun));
2420 //---------------------------------------------------------------------
2421 void AliSimulation::WriteGRPEntry()
2423 // Get the necessary information from galice (generator, trigger etc) and
2424 // write a GRP entry corresponding to the settings in the Config.C used
2425 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2428 AliInfo("Writing global run parameters entry into the OCDB");
2430 AliGRPObject* grpObj = new AliGRPObject();
2432 grpObj->SetRunType("PHYSICS");
2433 grpObj->SetTimeStart(fTimeStart);
2434 grpObj->SetTimeEnd(fTimeEnd);
2435 grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2437 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2443 gen->GetProjectile(projectile,a,z);
2445 gen->GetTarget(target,a,z);
2446 TString beamType = projectile + "-" + target;
2447 beamType.ReplaceAll(" ","");
2448 if (!beamType.CompareTo("-")) {
2449 grpObj->SetBeamType("UNKNOWN");
2450 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2453 grpObj->SetBeamType(beamType);
2455 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 * a / z);
2457 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 );
2459 // Heavy ion run, the event specie is set to kHighMult
2460 fEventSpecie = AliRecoParam::kHighMult;
2461 if ((strcmp(beamType,"p-p") == 0) ||
2462 (strcmp(beamType,"p-") == 0) ||
2463 (strcmp(beamType,"-p") == 0) ||
2464 (strcmp(beamType,"P-P") == 0) ||
2465 (strcmp(beamType,"P-") == 0) ||
2466 (strcmp(beamType,"-P") == 0)) {
2467 // Proton run, the event specie is set to kLowMult
2468 fEventSpecie = AliRecoParam::kLowMult;
2472 AliWarning("Unknown beam type and energy! Setting energy to 0");
2473 grpObj->SetBeamEnergy(0);
2474 grpObj->SetBeamType("UNKNOWN");
2477 UInt_t detectorPattern = 0;
2479 TObjArray *detArray = gAlice->Detectors();
2480 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2481 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2482 AliDebug(1, Form("Detector #%d found: %s", iDet, AliDAQ::OfflineModuleName(iDet)));
2483 detectorPattern |= (1 << iDet);
2488 if (!fTriggerConfig.IsNull())
2489 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2492 if (!fRunHLT.IsNull())
2493 detectorPattern |= (1 << AliDAQ::kHLTId);
2495 grpObj->SetNumberOfDetectors((Char_t)nDets);
2496 grpObj->SetDetectorMask((Int_t)detectorPattern);
2497 grpObj->SetLHCPeriod("LHC08c");
2498 grpObj->SetLHCState("STABLE_BEAMS");
2500 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2501 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2503 Float_t factorSol = field ? field->GetFactorSol() : 0;
2504 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2505 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2507 Float_t factorDip = field ? field->GetFactorDip() : 0;
2508 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2510 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2511 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2512 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2513 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2514 if (field) grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2515 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2517 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2519 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2521 // Now store the entry in OCDB
2522 AliCDBManager* man = AliCDBManager::Instance();
2524 man->SetLock(0, fKey);
2526 AliCDBStorage* sto = man->GetStorage(fGRPWriteLocation.Data());
2529 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2530 AliCDBMetaData *metadata= new AliCDBMetaData();
2532 metadata->SetResponsible("alice-off@cern.ch");
2533 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2535 sto->Put(grpObj,id,metadata);
2536 man->SetLock(1, fKey);
2539 //_____________________________________________________________________________
2540 time_t AliSimulation::GenerateTimeStamp() const
2542 // Generate event time-stamp according to
2543 // SOR/EOR time from GRP
2544 if (fUseTimeStampFromCDB)
2545 return fTimeStart + gRandom->Integer(fTimeEnd-fTimeStart);
2550 //_____________________________________________________________________________
2551 void AliSimulation::StoreUsedCDBMaps() const
2553 // write in galice.root maps with used CDB paths
2556 AliRunLoader* runLoader = LoadRun();
2558 AliError("Failed to open gAlice.root in write mode");
2562 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
2563 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
2565 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
2566 cdbMapCopy->SetOwner(1);
2567 // cdbMapCopy->SetName("cdbMap");
2568 TIter iter(cdbMap->GetTable());
2571 while((pair = dynamic_cast<TPair*> (iter.Next()))){
2572 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
2573 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
2574 if (keyStr && valStr)
2575 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
2578 TList *cdbListCopy = new TList();
2579 cdbListCopy->SetOwner(1);
2580 // cdbListCopy->SetName("cdbList");
2582 TIter iter2(cdbList);
2585 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
2586 cdbListCopy->Add(new TObjString(id->ToString().Data()));
2589 AliRunLoader::Instance()->CdGAFile();
2590 gDirectory->WriteObject(cdbMapCopy,"cdbMap","kSingleKey");
2591 gDirectory->WriteObject(cdbListCopy,"cdbList","kSingleKey");
2594 AliInfo(Form("Stored used OCDB entries as TMap %s and TList %s in %s",
2595 cdbMapCopy->GetName(),
2596 cdbListCopy->GetName(),
2597 fGAliceFileName.Data()));