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 63204 2013-06-26 13:33:28Z 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 ///////////////////////////////////////////////////////////////////////////////
109 #include <TGeoGlobalMagField.h>
110 #include <TGeoManager.h>
111 #include <TObjString.h>
114 #include <TVirtualMC.h>
115 #include <TVirtualMCApplication.h>
117 #include <TInterpreter.h>
119 #include "AliAlignObj.h"
120 #include "AliCDBEntry.h"
121 #include "AliCDBManager.h"
122 #include "AliGRPManager.h"
123 #include "AliCDBStorage.h"
124 #include "AliCTPRawData.h"
125 #include "AliCentralTrigger.h"
126 #include "AliCentralTrigger.h"
127 #include "AliCodeTimer.h"
129 #include "AliDigitizer.h"
130 #include "AliESDEvent.h"
131 #include "AliGRPObject.h"
132 #include "AliGenEventHeader.h"
133 #include "AliGenerator.h"
134 #include "AliGeomManager.h"
135 #include "AliHLTSimulation.h"
136 #include "AliHeader.h"
138 #include "AliLegoGenerator.h"
142 #include "AliModule.h"
144 #include "AliRawReaderDate.h"
145 #include "AliRawReaderFile.h"
146 #include "AliRawReaderRoot.h"
148 #include "AliDigitizationInput.h"
149 #include "AliRunLoader.h"
150 #include "AliStack.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","AD",
160 // #ifdef MFT_UPGRADE
167 //_____________________________________________________________________________
168 AliSimulation::AliSimulation(const char* configFileName,
169 const char* name, const char* title) :
172 fRunGeneratorOnly(kFALSE),
173 fRunGeneration(kTRUE),
174 fRunSimulation(kTRUE),
175 fLoadAlignFromCDB(kTRUE),
176 fLoadAlObjsListOfDets("ALL"),
180 fMakeDigitsFromHits(""),
182 fRawDataFileName(""),
183 fDeleteIntermediateFiles(kFALSE),
184 fWriteSelRawData(kFALSE),
185 fStopOnError(kFALSE),
186 fUseMonitoring(kFALSE),
188 fConfigFileName(configFileName),
189 fGAliceFileName("galice.root"),
191 fBkgrdFileNames(NULL),
192 fAlignObjArray(NULL),
193 fUseBkgrdVertex(kTRUE),
194 fRegionOfInterest(kFALSE),
200 fInitCDBCalled(kFALSE),
201 fInitRunNumberCalled(kFALSE),
202 fSetRunNumberFromDataCalled(kFALSE),
203 fEmbeddingFlag(kFALSE),
206 fUseVertexFromCDB(0),
207 fUseMagFieldFromGRP(0),
208 fGRPWriteLocation(Form("local://%s", gSystem->pwd())),
209 fUseTimeStampFromCDB(0),
215 fEventSpecie(AliRecoParam::kDefault),
216 fWriteQAExpertData(kTRUE),
220 fWriteGRPEntry(kTRUE)
222 // create simulation object with default parameters
224 SetGAliceFile("galice.root");
227 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
228 qam->SetActiveDetectors(fQADetectors) ;
229 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
230 qam->SetTasks(fQATasks) ;
233 //_____________________________________________________________________________
234 AliSimulation::~AliSimulation()
238 fEventsPerFile.Delete();
239 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
240 // delete fAlignObjArray; fAlignObjArray=0;
242 if (fBkgrdFileNames) {
243 fBkgrdFileNames->Delete();
244 delete fBkgrdFileNames;
247 fSpecCDBUri.Delete();
248 if (fgInstance==this) fgInstance = 0;
250 AliQAManager::QAManager()->ShowQA() ;
251 AliQAManager::Destroy() ;
252 AliCodeTimer::Instance()->Print();
256 //_____________________________________________________________________________
257 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
259 // set the number of events for one run
264 //_____________________________________________________________________________
265 void AliSimulation::InitQA()
267 // activate a default CDB storage
268 // First check if we have any CDB storage set, because it is used
269 // to retrieve the calibration and alignment constants
271 if (fInitCDBCalled) return;
272 fInitCDBCalled = kTRUE;
274 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
275 qam->SetActiveDetectors(fQADetectors) ;
276 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
277 qam->SetTasks(fQATasks) ;
278 if (fWriteQAExpertData)
279 qam->SetWriteExpert() ;
281 if (qam->IsDefaultStorageSet()) {
282 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
283 AliWarning("Default QA reference storage has been already set !");
284 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
285 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
286 fQARefUri = qam->GetDefaultStorage()->GetURI();
288 if (fQARefUri.Length() > 0) {
289 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
290 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
291 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
293 fQARefUri="local://$ALICE_ROOT/QARef";
294 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
295 AliWarning("Default QA reference storage not yet set !!!!");
296 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
297 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
299 qam->SetDefaultStorage(fQARefUri);
303 //_____________________________________________________________________________
304 void AliSimulation::InitCDB()
306 // activate a default CDB storage
307 // First check if we have any CDB storage set, because it is used
308 // to retrieve the calibration and alignment constants
310 if (fInitCDBCalled) return;
311 fInitCDBCalled = kTRUE;
313 AliCDBManager* man = AliCDBManager::Instance();
314 if (man->IsDefaultStorageSet())
316 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
317 AliWarning("Default CDB storage has been already set !");
318 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
319 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
320 fCDBUri = man->GetDefaultStorage()->GetURI();
323 if (fCDBUri.Length() > 0)
325 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
326 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
327 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
329 fCDBUri="local://$ALICE_ROOT/OCDB";
330 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
331 AliWarning("Default CDB storage not yet set !!!!");
332 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
333 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
336 man->SetDefaultStorage(fCDBUri);
339 // Now activate the detector specific CDB storage locations
340 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
341 TObject* obj = fSpecCDBUri[i];
343 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
344 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
345 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
346 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
351 //_____________________________________________________________________________
352 void AliSimulation::InitRunNumber(){
353 // check run number. If not set, set it to 0 !!!!
355 if (fInitRunNumberCalled) return;
356 fInitRunNumberCalled = kTRUE;
358 AliCDBManager* man = AliCDBManager::Instance();
359 if (man->GetRun() >= 0)
361 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
362 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
366 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
369 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
378 //_____________________________________________________________________________
379 void AliSimulation::SetCDBLock() {
380 // Set CDB lock: from now on it is forbidden to reset the run number
381 // or the default storage or to activate any further storage!
383 ULong64_t key = AliCDBManager::Instance()->SetLock(1);
387 //_____________________________________________________________________________
388 void AliSimulation::SetDefaultStorage(const char* uri) {
389 // Store the desired default CDB storage location
390 // Activate it later within the Run() method
396 //_____________________________________________________________________________
397 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
398 // Store the desired default CDB storage location
399 // Activate it later within the Run() method
402 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
405 //_____________________________________________________________________________
406 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
407 // Store a detector-specific CDB storage location
408 // Activate it later within the Run() method
410 AliCDBPath aPath(calibType);
411 if(!aPath.IsValid()){
412 AliError(Form("Not a valid path: %s", calibType));
416 TObject* obj = fSpecCDBUri.FindObject(calibType);
417 if (obj) fSpecCDBUri.Remove(obj);
418 fSpecCDBUri.Add(new TNamed(calibType, uri));
422 //_____________________________________________________________________________
423 void AliSimulation::SetRunNumber(Int_t run)
426 // Activate it later within the Run() method
431 //_____________________________________________________________________________
432 void AliSimulation::SetSeed(Int_t seed)
435 // Activate it later within the Run() method
440 //_____________________________________________________________________________
441 Bool_t AliSimulation::SetRunNumberFromData()
443 // Set the CDB manager run number
444 // The run number is retrieved from gAlice
446 if (fSetRunNumberFromDataCalled) return kTRUE;
447 fSetRunNumberFromDataCalled = kTRUE;
449 AliCDBManager* man = AliCDBManager::Instance();
450 Int_t runData = -1, runCDB = -1;
452 AliRunLoader* runLoader = LoadRun("READ");
453 if (!runLoader) return kFALSE;
455 runData = runLoader->GetHeader()->GetRun();
459 runCDB = man->GetRun();
461 if (runCDB != runData) {
462 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
463 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
464 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
465 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
470 man->SetRun(runData);
473 if(man->GetRun() < 0) {
474 AliError("Run number not properly initalized!");
483 //_____________________________________________________________________________
484 void AliSimulation::SetConfigFile(const char* fileName)
486 // set the name of the config file
488 fConfigFileName = fileName;
491 //_____________________________________________________________________________
492 void AliSimulation::SetGAliceFile(const char* fileName)
494 // set the name of the galice file
495 // the path is converted to an absolute one if it is relative
497 fGAliceFileName = fileName;
498 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
499 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
501 fGAliceFileName = absFileName;
502 delete[] absFileName;
505 AliDebug(2, Form("galice file name set to %s", fileName));
508 //_____________________________________________________________________________
509 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
512 // set the number of events per file for the given detector and data type
513 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
515 TNamed* obj = new TNamed(detector, type);
516 obj->SetUniqueID(nEvents);
517 fEventsPerFile.Add(obj);
520 //_____________________________________________________________________________
521 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
523 // Read the alignment objects from CDB.
524 // Each detector is supposed to have the
525 // alignment objects in DET/Align/Data CDB path.
526 // All the detector objects are then collected,
527 // sorted by geometry level (starting from ALIC) and
528 // then applied to the TGeo geometry.
529 // Finally an overlaps check is performed.
531 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
532 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
536 // initialize CDB storage, run number, set CDB lock
538 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
541 Bool_t delRunLoader = kFALSE;
543 runLoader = LoadRun("READ");
544 if (!runLoader) return kFALSE;
545 delRunLoader = kTRUE;
548 // Export ideal geometry
549 if(!IsGeometryFromFile()) AliGeomManager::GetGeometry()->Export("geometry.root");
551 // Load alignment data from CDB and apply to geometry through AliGeomManager
552 if(fLoadAlignFromCDB){
554 TString detStr = fLoadAlObjsListOfDets;
555 TString loadAlObjsListOfDets = "";
557 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
558 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
559 AliModule* det = (AliModule*) detArray->At(iDet);
560 if (!det || !det->IsActive()) continue;
561 if (IsSelected(det->GetName(), detStr)) {
562 //add det to list of dets to be aligned from CDB
563 loadAlObjsListOfDets += det->GetName();
564 loadAlObjsListOfDets += " ";
566 } // end loop over detectors
567 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
568 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
570 // Check if the array with alignment objects was
571 // provided by the user. If yes, apply the objects
572 // to the present TGeo geometry
573 if (fAlignObjArray) {
574 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
575 AliError("The misalignment of one or more volumes failed!"
576 "Compare the list of simulated detectors and the list of detector alignment data!");
577 if (delRunLoader) delete runLoader;
583 // Update the internal geometry of modules (ITS needs it)
584 TString detStr = fLoadAlObjsListOfDets;
585 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
586 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
588 AliModule* det = (AliModule*) detArray->At(iDet);
589 if (!det || !det->IsActive()) continue;
590 if (IsSelected(det->GetName(), detStr)) {
591 det->UpdateInternalGeometry();
593 } // end loop over detectors
596 if (delRunLoader) delete runLoader;
601 //_____________________________________________________________________________
602 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
604 // add a file with background events for merging
606 TObjString* fileNameStr = new TObjString(fileName);
607 fileNameStr->SetUniqueID(nSignalPerBkgrd);
608 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
609 fBkgrdFileNames->Add(fileNameStr);
612 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
614 // add a file with background events for embeddin
615 MergeWith(fileName, nSignalPerBkgrd);
616 fEmbeddingFlag = kTRUE;
619 //_____________________________________________________________________________
620 Bool_t AliSimulation::Run(Int_t nEvents)
622 // run the generation, simulation and digitization
625 AliCodeTimerAuto("",0)
626 AliSysInfo::AddStamp("Start_Run");
628 // Load run number and seed from environmental vars
629 ProcessEnvironmentVars();
630 AliSysInfo::AddStamp("ProcessEnvironmentVars");
632 gRandom->SetSeed(fSeed);
634 if (nEvents > 0) fNEvents = nEvents;
636 // Run generator-only code on demand
637 if (fRunGeneratorOnly)
639 if(!RunGeneratorOnly())
641 if (fStopOnError) return kFALSE;
647 // create and setup the HLT instance
648 if (!fRunHLT.IsNull() && !CreateHLT()) {
649 if (fStopOnError) return kFALSE;
654 // generation and simulation -> hits
655 if (fRunGeneration) {
656 if (!RunSimulation()) if (fStopOnError) return kFALSE;
658 AliSysInfo::AddStamp("RunSimulation");
660 // initialize CDB storage from external environment
661 // (either CDB manager or AliSimulation setters),
662 // if not already done in RunSimulation()
664 AliSysInfo::AddStamp("InitCDB");
666 // Set run number in CDBManager from data
667 // From this point on the run number must be always loaded from data!
668 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
670 // Set CDB lock: from now on it is forbidden to reset the run number
671 // or the default storage or to activate any further storage!
674 // If RunSimulation was not called, load the geometry and misalign it
675 if (!AliGeomManager::GetGeometry()) {
676 // Initialize the geometry manager
677 AliGeomManager::LoadGeometry("geometry.root");
678 AliSysInfo::AddStamp("GetGeometry");
679 // // Check that the consistency of symbolic names for the activated subdetectors
680 // // in the geometry loaded by AliGeomManager
681 // AliRunLoader* runLoader = LoadRun("READ");
682 // if (!runLoader) return kFALSE;
684 // TString detsToBeChecked = "";
685 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
686 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
687 // AliModule* det = (AliModule*) detArray->At(iDet);
688 // if (!det || !det->IsActive()) continue;
689 // detsToBeChecked += det->GetName();
690 // detsToBeChecked += " ";
691 // } // end loop over detectors
692 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
693 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
694 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
696 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
698 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
700 AliSysInfo::AddStamp("MissalignGeometry");
703 // hits -> summable digits
704 AliSysInfo::AddStamp("Start_sdigitization");
705 if (!fMakeSDigits.IsNull()) {
706 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
709 AliSysInfo::AddStamp("Stop_sdigitization");
711 AliSysInfo::AddStamp("Start_digitization");
712 // summable digits -> digits
713 if (!fMakeDigits.IsNull()) {
714 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
715 if (fStopOnError) return kFALSE;
718 AliSysInfo::AddStamp("Stop_digitization");
723 if (!fMakeDigitsFromHits.IsNull()) {
724 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
725 AliWarning(Form("Merging and direct creation of digits from hits "
726 "was selected for some detectors. "
727 "No merging will be done for the following detectors: %s",
728 fMakeDigitsFromHits.Data()));
730 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
731 if (fStopOnError) return kFALSE;
735 AliSysInfo::AddStamp("Hits2Digits");
739 if (!fTriggerConfig.IsNull() && !RunTrigger(fTriggerConfig,fMakeDigits)) {
740 if (fStopOnError) return kFALSE;
743 AliSysInfo::AddStamp("RunTrigger");
746 // digits -> raw data
747 if (!fWriteRawData.IsNull()) {
748 if (!WriteRawData(fWriteRawData, fRawDataFileName,
749 fDeleteIntermediateFiles,fWriteSelRawData)) {
750 if (fStopOnError) return kFALSE;
754 AliSysInfo::AddStamp("WriteRaw");
756 // run HLT simulation on simulated digit data if raw data is not
757 // simulated, otherwise its called as part of WriteRawData
758 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
760 if (fStopOnError) return kFALSE;
764 AliSysInfo::AddStamp("RunHLT");
768 Bool_t rv = RunQA() ;
774 AliSysInfo::AddStamp("RunQA");
776 TString snapshotFileOut("");
777 if(TString(gSystem->Getenv("OCDB_SNAPSHOT_CREATE")) == TString("kTRUE")){
778 AliInfo(" ******** Creating the snapshot! *********");
779 TString snapshotFile(gSystem->Getenv("OCDB_SNAPSHOT_FILENAME"));
780 if(!(snapshotFile.IsNull() || snapshotFile.IsWhitespace()))
781 snapshotFileOut = snapshotFile;
783 snapshotFileOut="OCDB.root";
784 AliCDBManager::Instance()->DumpToSnapshotFile(snapshotFileOut.Data(),kFALSE);
787 // Cleanup of CDB manager: cache and active storages!
788 AliCDBManager::Instance()->ClearCache();
793 //_______________________________________________________________________
794 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
795 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
796 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
799 // Generates lego plots of:
800 // - radiation length map phi vs theta
801 // - radiation length map phi vs eta
802 // - interaction length map
803 // - g/cm2 length map
805 // ntheta bins in theta, eta
806 // themin minimum angle in theta (degrees)
807 // themax maximum angle in theta (degrees)
809 // phimin minimum angle in phi (degrees)
810 // phimax maximum angle in phi (degrees)
811 // rmin minimum radius
812 // rmax maximum radius
815 // The number of events generated = ntheta*nphi
816 // run input parameters in macro setup (default="Config.C")
818 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
821 <img src="picts/AliRunLego1.gif">
826 <img src="picts/AliRunLego2.gif">
831 <img src="picts/AliRunLego3.gif">
836 // run the generation and simulation
838 AliCodeTimerAuto("",0)
840 // initialize CDB storage and run number from external environment
841 // (either CDB manager or AliSimulation setters)
847 AliError("no gAlice object. Restart aliroot and try again.");
850 if (gAlice->Modules()->GetEntries() > 0) {
851 AliError("gAlice was already run. Restart aliroot and try again.");
855 AliInfo(Form("initializing gAlice with config file %s",
856 fConfigFileName.Data()));
859 if (nev == -1) nev = nc1 * nc2;
861 // check if initialisation has been done
862 // If runloader has been initialized, set the number of events per file to nc1 * nc2
865 if (!gener) gener = new AliLegoGenerator();
867 // Configure Generator
869 gener->SetRadiusRange(rmin, rmax);
870 gener->SetZMax(zmax);
871 gener->SetCoor1Range(nc1, c1min, c1max);
872 gener->SetCoor2Range(nc2, c2min, c2max);
876 fLego = new AliLego("lego",gener);
878 //__________________________________________________________________________
882 gROOT->LoadMacro(setup);
883 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
885 if(AliCDBManager::Instance()->GetRun() >= 0) {
886 SetRunNumber(AliCDBManager::Instance()->GetRun());
888 AliWarning("Run number not initialized!!");
891 AliRunLoader::Instance()->CdGAFile();
893 AliPDG::AddParticlesToPdgDataBase();
895 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
897 gAlice->GetMCApp()->Init();
900 //Must be here because some MCs (G4) adds detectors here and not in Config.C
901 gAlice->InitLoaders();
902 AliRunLoader::Instance()->MakeTree("E");
905 // Save stuff at the beginning of the file to avoid file corruption
906 AliRunLoader::Instance()->CdGAFile();
909 //Save current generator
910 AliGenerator *gen=gAlice->GetMCApp()->Generator();
911 gAlice->GetMCApp()->ResetGenerator(gener);
912 //Prepare MC for Lego Run
918 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
919 gMC->ProcessRun(nev);
921 // End of this run, close files
923 // Restore current generator
924 gAlice->GetMCApp()->ResetGenerator(gen);
925 // Delete Lego Object
931 //_____________________________________________________________________________
932 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
936 AliCodeTimerAuto("",0)
938 // initialize CDB storage from external environment
939 // (either CDB manager or AliSimulation setters),
940 // if not already done in RunSimulation()
943 // Set run number in CDBManager from data
944 // From this point on the run number must be always loaded from data!
945 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
947 // Set CDB lock: from now on it is forbidden to reset the run number
948 // or the default storage or to activate any further storage!
951 AliRunLoader* runLoader = LoadRun("READ");
952 if (!runLoader) return kFALSE;
953 TString trconfiguration = config;
955 if (trconfiguration.IsNull()) {
956 if(!fTriggerConfig.IsNull()) {
957 trconfiguration = fTriggerConfig;
960 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
963 runLoader->MakeTree( "GG" );
964 AliCentralTrigger* aCTP = runLoader->GetTrigger();
965 // Load Configuration
966 if (!aCTP->LoadConfiguration( trconfiguration ))
970 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
982 //_____________________________________________________________________________
983 Bool_t AliSimulation::WriteTriggerRawData()
985 // Writes the CTP (trigger) DDL raw data
986 // Details of the format are given in the
987 // trigger TDR - pages 134 and 135.
988 AliCTPRawData writer;
994 //_____________________________________________________________________________
995 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
997 // run the generation and simulation
999 AliCodeTimerAuto("",0)
1001 // initialize CDB storage and run number from external environment
1002 // (either CDB manager or AliSimulation setters)
1003 AliSysInfo::AddStamp("RunSimulation_Begin");
1005 AliSysInfo::AddStamp("RunSimulation_InitCDB");
1009 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
1012 AliError("no gAlice object. Restart aliroot and try again.");
1015 if (gAlice->Modules()->GetEntries() > 0) {
1016 AliError("gAlice was already run. Restart aliroot and try again.");
1020 // Setup monitoring if requested
1021 gAlice->GetMCApp()->SetUseMonitoring(fUseMonitoring);
1023 AliInfo(Form("initializing gAlice with config file %s",
1024 fConfigFileName.Data()));
1027 // Initialize ALICE Simulation run
1032 // If requested set the mag. field from the GRP entry.
1033 // After this the field is loccked and cannot be changed by Config.C
1034 if (fUseMagFieldFromGRP) {
1036 grpM.ReadGRPEntry();
1038 AliInfo("Field is locked now. It cannot be changed in Config.C");
1042 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1043 gROOT->LoadMacro(fConfigFileName.Data());
1044 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1045 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1046 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1048 AliSysInfo::AddStamp("RunSimulation_Config");
1051 // If requested obtain the vertex position and vertex sigma_z from the CDB
1052 // This overwrites the settings from the Config.C
1053 if (fUseVertexFromCDB) {
1054 Double_t vtxPos[3] = {0., 0., 0.};
1055 Double_t vtxSig[3] = {0., 0., 0.};
1056 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1058 AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1060 if(vertex->GetXRes()>2.8) { // > pipe radius --> it's a dummy object, don't use it
1061 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1062 if (entry) vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1066 vertex->GetXYZ(vtxPos);
1067 vertex->GetSigmaXYZ(vtxSig);
1068 AliInfo("Overwriting Config.C vertex settings !");
1069 AliInfo(Form("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
1070 vtxPos[0], vtxPos[1], vtxPos[2], vtxSig[2]));
1072 AliGenerator *gen = gAlice->GetMCApp()->Generator();
1073 gen->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]); // vertex position
1074 gen->SetSigmaZ(vtxSig[2]);
1079 // If requested we take the SOR and EOR time-stamps from the GRP and use them
1080 // in order to generate the event time-stamps
1081 if (fUseTimeStampFromCDB) {
1083 grpM.ReadGRPEntry();
1084 const AliGRPObject *grpObj = grpM.GetGRPData();
1085 if (!grpObj || (grpObj->GetTimeEnd() <= grpObj->GetTimeStart())) {
1086 AliError("Missing GRP or bad SOR/EOR time-stamps! Switching off the time-stamp generation from GRP!");
1087 fTimeStart = fTimeEnd = 0;
1088 fUseTimeStampFromCDB = kFALSE;
1091 fTimeStart = grpObj->GetTimeStart();
1092 fTimeEnd = grpObj->GetTimeEnd();
1096 if(AliCDBManager::Instance()->GetRun() >= 0) {
1097 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
1098 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
1100 AliWarning("Run number not initialized!!");
1103 AliRunLoader::Instance()->CdGAFile();
1106 AliPDG::AddParticlesToPdgDataBase();
1108 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1109 AliSysInfo::AddStamp("RunSimulation_GetField");
1110 gAlice->GetMCApp()->Init();
1111 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1113 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1114 gAlice->InitLoaders();
1115 AliRunLoader::Instance()->MakeTree("E");
1116 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1117 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1118 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1120 // Save stuff at the beginning of the file to avoid file corruption
1121 AliRunLoader::Instance()->CdGAFile();
1123 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1124 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1125 //___________________________________________________________________________________________
1127 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1129 // Set run number in CDBManager
1130 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1132 AliRunLoader* runLoader = AliRunLoader::Instance();
1134 AliError(Form("gAlice has no run loader object. "
1135 "Check your config file: %s", fConfigFileName.Data()));
1138 SetGAliceFile(runLoader->GetFileName());
1140 // Misalign geometry
1141 #if ROOT_VERSION_CODE < 331527
1142 AliGeomManager::SetGeometry(gGeoManager);
1144 // Check that the consistency of symbolic names for the activated subdetectors
1145 // in the geometry loaded by AliGeomManager
1146 TString detsToBeChecked = "";
1147 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1148 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1149 AliModule* det = (AliModule*) detArray->At(iDet);
1150 if (!det || !det->IsActive()) continue;
1151 detsToBeChecked += det->GetName();
1152 detsToBeChecked += " ";
1153 } // end loop over detectors
1154 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1155 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1156 MisalignGeometry(runLoader);
1157 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1160 // AliRunLoader* runLoader = AliRunLoader::Instance();
1161 // if (!runLoader) {
1162 // AliError(Form("gAlice has no run loader object. "
1163 // "Check your config file: %s", fConfigFileName.Data()));
1166 // SetGAliceFile(runLoader->GetFileName());
1168 if (!gAlice->GetMCApp()->Generator()) {
1169 AliError(Form("gAlice has no generator object. "
1170 "Check your config file: %s", fConfigFileName.Data()));
1174 // Write GRP entry corresponding to the setting found in Cofig.C
1177 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1179 if (nEvents <= 0) nEvents = fNEvents;
1181 // get vertex from background file in case of merging
1182 if (fUseBkgrdVertex &&
1183 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1184 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1185 const char* fileName = ((TObjString*)
1186 (fBkgrdFileNames->At(0)))->GetName();
1187 AliInfo(Form("The vertex will be taken from the background "
1188 "file %s with nSignalPerBackground = %d",
1189 fileName, signalPerBkgrd));
1190 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1191 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1194 if (!fRunSimulation) {
1195 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1198 // set the number of events per file for given detectors and data types
1199 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1200 if (!fEventsPerFile[i]) continue;
1201 const char* detName = fEventsPerFile[i]->GetName();
1202 const char* typeName = fEventsPerFile[i]->GetTitle();
1203 TString loaderName(detName);
1204 loaderName += "Loader";
1205 AliLoader* loader = runLoader->GetLoader(loaderName);
1207 AliError(Form("RunSimulation no loader for %s found\n Number of events per file not set for %s %s",
1208 detName, typeName, detName));
1211 AliDataLoader* dataLoader =
1212 loader->GetDataLoader(typeName);
1214 AliError(Form("no data loader for %s found\n"
1215 "Number of events per file not set for %s %s",
1216 typeName, detName, typeName));
1219 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1220 AliDebug(1, Form("number of events per file set to %d for %s %s",
1221 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1224 AliInfo("running gAlice");
1225 AliSysInfo::AddStamp("Start_ProcessRun");
1227 // Create the Root Tree with one branch per detector
1228 //Hits moved to begin event -> now we are crating separate tree for each event
1229 gMC->ProcessRun(nEvents);
1231 // End of this run, close files
1232 if(nEvents>0) FinishRun();
1234 AliSysInfo::AddStamp("Stop_ProcessRun");
1240 //_____________________________________________________________________________
1241 Bool_t AliSimulation::RunGeneratorOnly()
1244 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1245 gROOT->LoadMacro(fConfigFileName.Data());
1246 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1247 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1248 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1251 // Setup the runloader and generator, check if everything is OK
1252 AliRunLoader* runLoader = AliRunLoader::Instance();
1253 AliGenerator* generator = gAlice->GetMCApp()->Generator();
1255 AliError(Form("gAlice has no run loader object. "
1256 "Check your config file: %s", fConfigFileName.Data()));
1260 AliError(Form("gAlice has no generator object. "
1261 "Check your config file: %s", fConfigFileName.Data()));
1265 runLoader->LoadKinematics("RECREATE");
1266 runLoader->MakeTree("E");
1268 // Create stack and header
1269 runLoader->MakeStack();
1270 AliStack* stack = runLoader->Stack();
1271 AliHeader* header = runLoader->GetHeader();
1273 // Intialize generator
1275 generator->SetStack(stack);
1277 // Run main generator loop
1279 for (Int_t iev=0; iev<fNEvents; iev++)
1282 header->Reset(0,iev);
1283 runLoader->SetEventNumber(iev);
1285 runLoader->MakeTree("K");
1288 generator->Generate();
1291 header->SetNprimary(stack->GetNprimary());
1292 header->SetNtrack(stack->GetNtrack());
1293 stack->FinishEvent();
1294 header->SetStack(stack);
1295 runLoader->TreeE()->Fill();
1296 runLoader->WriteKinematics("OVERWRITE");
1300 generator->FinishRun();
1302 runLoader->WriteHeader("OVERWRITE");
1309 //_____________________________________________________________________________
1310 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1312 // run the digitization and produce summable digits
1313 static Int_t eventNr=0;
1314 AliCodeTimerAuto("",0) ;
1316 // initialize CDB storage, run number, set CDB lock
1318 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1321 AliRunLoader* runLoader = LoadRun();
1322 if (!runLoader) return kFALSE;
1324 TString detStr = detectors;
1325 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1326 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1327 AliModule* det = (AliModule*) detArray->At(iDet);
1328 if (!det || !det->IsActive()) continue;
1329 if (IsSelected(det->GetName(), detStr)) {
1330 AliInfo(Form("creating summable digits for %s", det->GetName()));
1331 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1332 det->Hits2SDigits();
1333 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1334 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1338 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1339 AliError(Form("the following detectors were not found: %s",
1341 if (fStopOnError) return kFALSE;
1350 //_____________________________________________________________________________
1351 Bool_t AliSimulation::RunDigitization(const char* detectors,
1352 const char* excludeDetectors)
1354 // run the digitization and produce digits from sdigits
1356 AliCodeTimerAuto("",0)
1358 // initialize CDB storage, run number, set CDB lock
1360 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1363 delete AliRunLoader::Instance();
1368 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1369 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1370 AliDigitizationInput digInp(nStreams, signalPerBkgrd);
1371 // digInp.SetEmbeddingFlag(fEmbeddingFlag);
1372 digInp.SetRegionOfInterest(fRegionOfInterest);
1373 digInp.SetInputStream(0, fGAliceFileName.Data());
1374 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1375 const char* fileName = ((TObjString*)(fBkgrdFileNames->At(iStream-1)))->GetName();
1376 digInp.SetInputStream(iStream, fileName);
1379 detArr.SetOwner(kTRUE);
1380 TString detStr = detectors;
1381 TString detExcl = excludeDetectors;
1382 if (!static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice()) {
1383 AliError("Error occured while getting gAlice from Input 0");
1386 AliRunLoader* runLoader = AliRunLoader::GetRunLoader(digInp.GetInputStream(0)->GetFolderName());
1387 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1388 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1389 AliModule* det = (AliModule*) detArray->At(iDet);
1390 if (!det || !det->IsActive()) continue;
1391 if (!IsSelected(det->GetName(), detStr) || IsSelected(det->GetName(), detExcl)) continue;
1392 AliDigitizer* digitizer = det->CreateDigitizer(&digInp);
1393 if (!digitizer || !digitizer->Init()) {
1394 AliError(Form("no digitizer for %s", det->GetName()));
1395 if (fStopOnError) return kFALSE;
1398 detArr.AddLast(digitizer);
1399 AliInfo(Form("Created digitizer from SDigits -> Digits for %s", det->GetName()));
1402 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1403 AliError(Form("the following detectors were not found: %s", detStr.Data()));
1404 if (fStopOnError) return kFALSE;
1407 Int_t ndigs = detArr.GetEntriesFast();
1408 Int_t eventsCreated = 0;
1409 AliRunLoader* outRl = digInp.GetOutRunLoader();
1410 while ((eventsCreated++ < fNEvents) || (fNEvents < 0)) {
1411 if (!digInp.ConnectInputTrees()) break;
1412 digInp.InitEvent(); //this must be after call of Connect Input tress.
1413 if (outRl) outRl->SetEventNumber(eventsCreated-1);
1414 static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice(); // use gAlice of the first input stream
1415 for (int id=0;id<ndigs;id++) ((AliDigitizer*)detArr[id])->Digitize("");
1416 digInp.FinishEvent();
1418 digInp.FinishGlobal();
1423 //_____________________________________________________________________________
1424 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1426 // run the digitization and produce digits from hits
1428 AliCodeTimerAuto("",0)
1430 // initialize CDB storage, run number, set CDB lock
1432 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1435 AliRunLoader* runLoader = LoadRun("READ");
1436 if (!runLoader) return kFALSE;
1438 TString detStr = detectors;
1439 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1440 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1441 AliModule* det = (AliModule*) detArray->At(iDet);
1442 if (!det || !det->IsActive()) continue;
1443 if (IsSelected(det->GetName(), detStr)) {
1444 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1449 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1450 AliError(Form("the following detectors were not found: %s",
1452 if (fStopOnError) return kFALSE;
1458 //_____________________________________________________________________________
1459 Bool_t AliSimulation::WriteRawData(const char* detectors,
1460 const char* fileName,
1461 Bool_t deleteIntermediateFiles,
1464 // convert the digits to raw data
1465 // First DDL raw data files for the given detectors are created.
1466 // If a file name is given, the DDL files are then converted to a DATE file.
1467 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1469 // If the file name has the extension ".root", the DATE file is converted
1471 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1472 // 'selrawdata' flag can be used to enable writing of detectors raw data
1473 // accoring to the trigger cluster.
1475 AliCodeTimerAuto("",0)
1476 AliSysInfo::AddStamp("WriteRawData_Start");
1478 TString detStr = detectors;
1479 if (!WriteRawFiles(detStr.Data())) {
1480 if (fStopOnError) return kFALSE;
1482 AliSysInfo::AddStamp("WriteRawFiles");
1484 // run HLT simulation on simulated DDL raw files
1485 // and produce HLT ddl raw files to be included in date/root file
1486 // bugfix 2009-06-26: the decision whether to write HLT raw data
1487 // is taken in RunHLT. Here HLT always needs to be run in order to
1488 // create HLT digits, unless its switched off. This is due to the
1489 // special placement of the HLT between the generation of DDL files
1490 // and conversion to DATE/Root file.
1491 detStr.ReplaceAll("HLT", "");
1492 if (!fRunHLT.IsNull()) {
1494 if (fStopOnError) return kFALSE;
1497 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1499 TString dateFileName(fileName);
1500 if (!dateFileName.IsNull()) {
1501 Bool_t rootOutput = dateFileName.EndsWith(".root");
1502 if (rootOutput) dateFileName += ".date";
1503 TString selDateFileName;
1505 selDateFileName = "selected.";
1506 selDateFileName+= dateFileName;
1508 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1509 if (fStopOnError) return kFALSE;
1511 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1512 if (deleteIntermediateFiles) {
1513 AliRunLoader* runLoader = LoadRun("READ");
1514 if (runLoader) for (Int_t iEvent = 0;
1515 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1517 snprintf(command, 256, "rm -r raw%d", iEvent);
1518 gSystem->Exec(command);
1524 if (!ConvertDateToRoot(dateFileName, fileName)) {
1525 if (fStopOnError) return kFALSE;
1527 AliSysInfo::AddStamp("ConvertDateToRoot");
1528 if (deleteIntermediateFiles) {
1529 gSystem->Unlink(dateFileName);
1532 TString selFileName = "selected.";
1533 selFileName += fileName;
1534 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1535 if (fStopOnError) return kFALSE;
1537 if (deleteIntermediateFiles) {
1538 gSystem->Unlink(selDateFileName);
1547 //_____________________________________________________________________________
1548 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1550 // convert the digits to raw data DDL files
1552 AliCodeTimerAuto("",0)
1554 AliRunLoader* runLoader = LoadRun("READ");
1555 if (!runLoader) return kFALSE;
1557 // write raw data to DDL files
1558 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1559 AliInfo(Form("processing event %d", iEvent));
1560 runLoader->GetEvent(iEvent);
1561 TString baseDir = gSystem->WorkingDirectory();
1563 snprintf(dirName, 256, "raw%d", iEvent);
1564 gSystem->MakeDirectory(dirName);
1565 if (!gSystem->ChangeDirectory(dirName)) {
1566 AliError(Form("couldn't change to directory %s", dirName));
1567 if (fStopOnError) return kFALSE; else continue;
1570 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1573 TString detStr = detectors;
1574 if (IsSelected("HLT", detStr)) {
1575 // Do nothing. "HLT" will be removed from detStr and HLT raw
1576 // data files are generated in RunHLT.
1579 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1580 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1581 AliModule* det = (AliModule*) detArray->At(iDet);
1582 if (!det || !det->IsActive()) continue;
1583 if (IsSelected(det->GetName(), detStr)) {
1584 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1589 if (!WriteTriggerRawData())
1590 if (fStopOnError) return kFALSE;
1592 gSystem->ChangeDirectory(baseDir);
1593 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1594 AliError(Form("the following detectors were not found: %s",
1596 if (fStopOnError) return kFALSE;
1605 //_____________________________________________________________________________
1606 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1607 const char* selDateFileName)
1609 // convert raw data DDL files to a DATE file with the program "dateStream"
1610 // The second argument is not empty when the user decides to write
1611 // the detectors raw data according to the trigger cluster.
1613 AliCodeTimerAuto("",0)
1615 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1617 AliError("the program dateStream was not found");
1618 if (fStopOnError) return kFALSE;
1623 AliRunLoader* runLoader = LoadRun("READ");
1624 if (!runLoader) return kFALSE;
1626 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1627 Bool_t selrawdata = kFALSE;
1628 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1631 // Note the option -s. It is used in order to avoid
1632 // the generation of SOR/EOR events.
1633 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1634 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1635 FILE* pipe = gSystem->OpenPipe(command, "w");
1638 AliError(Form("Cannot execute command: %s",command));
1642 Int_t selEvents = 0;
1643 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1645 UInt_t detectorPattern = 0;
1646 runLoader->GetEvent(iEvent);
1647 if (!runLoader->LoadTrigger()) {
1648 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1649 detectorPattern = aCTP->GetClusterMask();
1650 // Check if the event was triggered by CTP
1652 if (aCTP->GetClassMask()) selEvents++;
1656 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1658 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1659 selrawdata = kFALSE;
1663 fprintf(pipe, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1667 // loop over detectors and DDLs
1668 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1669 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1671 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1672 Int_t ldcID = Int_t(ldc + 0.0001);
1673 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1675 char rawFileName[256];
1676 snprintf(rawFileName, 256, "raw%d/%s",
1677 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1679 // check existence and size of raw data file
1680 FILE* file = fopen(rawFileName, "rb");
1681 if (!file) continue;
1682 fseek(file, 0, SEEK_END);
1683 unsigned long size = ftell(file);
1685 if (!size) continue;
1687 if (ldcID != prevLDC) {
1688 fprintf(pipe, " LDC Id %d\n", ldcID);
1691 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1696 Int_t result = gSystem->ClosePipe(pipe);
1698 if (!(selrawdata && selEvents > 0)) {
1700 return (result == 0);
1703 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1705 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1706 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1707 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1709 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1711 // Get the trigger decision and cluster
1712 UInt_t detectorPattern = 0;
1714 runLoader->GetEvent(iEvent);
1715 if (!runLoader->LoadTrigger()) {
1716 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1717 if (aCTP->GetClassMask() == 0) continue;
1718 detectorPattern = aCTP->GetClusterMask();
1719 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1720 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1723 fprintf(pipe2, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1727 // loop over detectors and DDLs
1728 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1729 // Write only raw data from detectors that
1730 // are contained in the trigger cluster(s)
1731 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1733 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1735 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1736 Int_t ldcID = Int_t(ldc + 0.0001);
1737 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1739 char rawFileName[256];
1740 snprintf(rawFileName, 256, "raw%d/%s",
1741 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1743 // check existence and size of raw data file
1744 FILE* file = fopen(rawFileName, "rb");
1745 if (!file) continue;
1746 fseek(file, 0, SEEK_END);
1747 unsigned long size = ftell(file);
1749 if (!size) continue;
1751 if (ldcID != prevLDC) {
1752 fprintf(pipe2, " LDC Id %d\n", ldcID);
1755 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1760 Int_t result2 = gSystem->ClosePipe(pipe2);
1763 return ((result == 0) && (result2 == 0));
1766 //_____________________________________________________________________________
1767 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1768 const char* rootFileName)
1770 // convert a DATE file to a root file with the program "alimdc"
1773 const Int_t kDBSize = 2000000000;
1774 const Int_t kTagDBSize = 1000000000;
1775 const Bool_t kFilter = kFALSE;
1776 const Int_t kCompression = 1;
1778 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1780 AliError("the program alimdc was not found");
1781 if (fStopOnError) return kFALSE;
1786 AliInfo(Form("converting DATE file %s to root file %s",
1787 dateFileName, rootFileName));
1789 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1790 const char* tagDBFS = "/tmp/mdc1/tags";
1792 // User defined file system locations
1793 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1794 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1795 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1796 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1797 if (gSystem->Getenv("ALIMDC_TAGDB"))
1798 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1800 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1801 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1802 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1804 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1805 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1806 gSystem->Exec(Form("mkdir %s",tagDBFS));
1808 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1809 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1810 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1812 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1813 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1814 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1816 return (result == 0);
1820 //_____________________________________________________________________________
1821 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1823 // delete existing run loaders, open a new one and load gAlice
1825 delete AliRunLoader::Instance();
1826 AliRunLoader* runLoader =
1827 AliRunLoader::Open(fGAliceFileName.Data(),
1828 AliConfig::GetDefaultEventFolderName(), mode);
1830 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1833 runLoader->LoadgAlice();
1834 runLoader->LoadHeader();
1835 gAlice = runLoader->GetAliRun();
1837 AliError(Form("no gAlice object found in file %s",
1838 fGAliceFileName.Data()));
1844 //_____________________________________________________________________________
1845 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1847 // get or calculate the number of signal events per background event
1849 if (!fBkgrdFileNames) return 1;
1850 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1851 if (nBkgrdFiles == 0) return 1;
1853 // get the number of signal events
1855 AliRunLoader* runLoader =
1856 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1857 if (!runLoader) return 1;
1859 nEvents = runLoader->GetNumberOfEvents();
1864 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1865 // get the number of background events
1866 const char* fileName = ((TObjString*)
1867 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1868 AliRunLoader* runLoader =
1869 AliRunLoader::Open(fileName, "BKGRD");
1870 if (!runLoader) continue;
1871 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1874 // get or calculate the number of signal per background events
1875 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1876 if (nSignalPerBkgrd <= 0) {
1877 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1878 } else if (result && (result != nSignalPerBkgrd)) {
1879 AliInfo(Form("the number of signal events per background event "
1880 "will be changed from %d to %d for stream %d",
1881 nSignalPerBkgrd, result, iBkgrdFile+1));
1882 nSignalPerBkgrd = result;
1885 if (!result) result = nSignalPerBkgrd;
1886 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1887 AliWarning(Form("not enough background events (%d) for %d signal events "
1888 "using %d signal per background events for stream %d",
1889 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1896 //_____________________________________________________________________________
1897 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1899 // check whether detName is contained in detectors
1900 // if yes, it is removed from detectors
1902 // check if all detectors are selected
1903 if ((detectors.CompareTo("ALL") == 0) ||
1904 detectors.BeginsWith("ALL ") ||
1905 detectors.EndsWith(" ALL") ||
1906 detectors.Contains(" ALL ")) {
1911 // search for the given detector
1912 Bool_t result = kFALSE;
1913 if ((detectors.CompareTo(detName) == 0) ||
1914 detectors.BeginsWith(detName+" ") ||
1915 detectors.EndsWith(" "+detName) ||
1916 detectors.Contains(" "+detName+" ")) {
1917 detectors.ReplaceAll(detName, "");
1921 // clean up the detectors string
1922 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1923 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1924 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1929 //_____________________________________________________________________________
1930 Int_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName, Int_t N, Int_t nSkip)
1933 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1934 // These can be used for embedding of MC tracks into RAW data using the standard
1935 // merging procedure.
1937 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1940 AliError("no gAlice object. Restart aliroot and try again.");
1943 if (gAlice->Modules()->GetEntries() > 0) {
1944 AliError("gAlice was already run. Restart aliroot and try again.");
1948 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1952 gROOT->LoadMacro(fConfigFileName.Data());
1953 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1955 if(AliCDBManager::Instance()->GetRun() >= 0) {
1956 SetRunNumber(AliCDBManager::Instance()->GetRun());
1958 AliWarning("Run number not initialized!!");
1961 AliRunLoader::Instance()->CdGAFile();
1963 AliPDG::AddParticlesToPdgDataBase();
1965 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1967 gAlice->GetMCApp()->Init();
1969 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1970 gAlice->InitLoaders();
1971 AliRunLoader::Instance()->MakeTree("E");
1972 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1973 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1974 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1977 // Save stuff at the beginning of the file to avoid file corruption
1978 AliRunLoader::Instance()->CdGAFile();
1983 //AliCDBManager* man = AliCDBManager::Instance();
1984 //man->SetRun(0); // Should this come from rawdata header ?
1988 // Get the runloader
1989 AliRunLoader* runLoader = AliRunLoader::Instance();
1991 // Open esd file if available
1994 AliESDEvent* esd = 0;
1995 if (esdFileName && (strlen(esdFileName)>0)) {
1996 esdFile = TFile::Open(esdFileName);
1998 esd = new AliESDEvent();
1999 esdFile->GetObject("esdTree", treeESD);
2001 esd->ReadFromTree(treeESD);
2003 AliInfo(Form("Asking to skip first %d ESDs events",nSkip));
2012 // Create the RawReader
2013 TString fileName(rawDirectory);
2014 AliRawReader* rawReader = AliRawReader::Create(fileName.Data());
2015 if (!rawReader) return (kFALSE);
2017 // if (!fEquipIdMap.IsNull() && fRawReader)
2018 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
2020 // Get list of detectors
2021 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
2024 AliHeader* header = runLoader->GetHeader();
2028 if (!(rawReader->NextEvent())) break;
2029 runLoader->SetEventNumber(nev);
2030 runLoader->GetHeader()->Reset(rawReader->GetRunNumber(),
2032 runLoader->GetEvent(nev);
2033 AliInfo(Form("We are at event %d",nev));
2036 TString detStr = fMakeSDigits;
2037 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
2038 AliModule* det = (AliModule*) detArray->At(iDet);
2039 if (!det || !det->IsActive()) continue;
2040 if (IsSelected(det->GetName(), detStr)) {
2041 AliInfo(Form("Calling Raw2SDigits for %s", det->GetName()));
2042 det->Raw2SDigits(rawReader);
2049 // If ESD information available obtain reconstructed vertex and store in header.
2051 AliInfo(Form("Selected event %d correspond to event %d in raw and to %d in esd",nev,rawReader->GetEventIndex(),nSkip+rawReader->GetEventIndex()));
2052 treeESD->GetEvent(nSkip+rawReader->GetEventIndex());
2053 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
2054 Double_t position[3];
2055 esdVertex->GetXYZ(position);
2056 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
2059 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
2060 mcHeader->SetPrimaryVertex(mcV);
2061 header->Reset(0,nev);
2062 header->SetGenEventHeader(mcHeader);
2063 AliInfo(Form("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]));
2067 runLoader->TreeE()->Fill();
2068 AliInfo(Form("Finished event %d",nev));
2077 runLoader->CdGAFile();
2078 runLoader->WriteHeader("OVERWRITE");
2079 runLoader->WriteRunLoader();
2084 //_____________________________________________________________________________
2085 void AliSimulation::FinishRun()
2088 // Called at the end of the run.
2093 AliDebug(1, "Finish Lego");
2094 AliRunLoader::Instance()->CdGAFile();
2098 // Clean detector information
2099 TIter next(gAlice->Modules());
2100 AliModule *detector;
2101 while((detector = dynamic_cast<AliModule*>(next()))) {
2102 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
2103 detector->FinishRun();
2106 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
2107 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
2109 // Write AliRun info and all detectors parameters
2110 AliRunLoader::Instance()->CdGAFile();
2111 gAlice->Write(0,TObject::kOverwrite);//write AliRun
2112 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
2114 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
2115 AliRunLoader::Instance()->Synchronize();
2118 //_____________________________________________________________________________
2119 Int_t AliSimulation::GetDetIndex(const char* detector)
2121 // return the detector index corresponding to detector
2123 for (index = 0; index < fgkNDetectors ; index++) {
2124 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2130 //_____________________________________________________________________________
2131 Bool_t AliSimulation::CreateHLT()
2133 // Init the HLT simulation.
2134 // The function loads the library and creates the instance of AliHLTSimulation.
2135 // the main reason for the decoupled creation is to set the transient OCDB
2136 // objects before the OCDB is locked
2138 // load the library dynamically
2139 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
2141 // check for the library version
2142 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
2144 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
2147 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
2148 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
2151 // print compile info
2152 typedef void (*CompileInfo)( const char*& date, const char*& time);
2153 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
2155 const char* date="";
2156 const char* time="";
2157 (*fctInfo)(date, time);
2158 if (!date) date="unknown";
2159 if (!time) time="unknown";
2160 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
2162 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
2165 // create instance of the HLT simulation
2166 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
2167 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
2168 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
2172 TString specObjects;
2173 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
2174 if (specObjects.Length()>0) specObjects+=" ";
2175 specObjects+=fSpecCDBUri[i]->GetName();
2178 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
2179 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
2180 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2186 //_____________________________________________________________________________
2187 Bool_t AliSimulation::RunHLT()
2189 // Run the HLT simulation
2190 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2191 // Disabled if fRunHLT is empty, default vaule is "default".
2192 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2193 // The default simulation depends on the HLT component libraries and their
2194 // corresponding agents which define components and chains to run. See
2195 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2196 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2198 // The libraries to be loaded can be specified as an option.
2200 // AliSimulation sim;
2201 // sim.SetRunHLT("libAliHLTSample.so");
2203 // will only load <tt>libAliHLTSample.so</tt>
2205 // Other available options:
2206 // \li loglevel=<i>level</i> <br>
2207 // logging level for this processing
2209 // disable redirection of log messages to AliLog class
2210 // \li config=<i>macro</i>
2211 // configuration macro
2212 // \li chains=<i>configuration</i>
2213 // comma separated list of configurations to be run during simulation
2214 // \li rawfile=<i>file</i>
2215 // source for the RawReader to be created, the default is <i>./</i> if
2216 // raw data is simulated
2220 if (!fpHLT && !CreateHLT()) {
2223 AliHLTSimulation* pHLT=fpHLT;
2225 AliRunLoader* pRunLoader = LoadRun("READ");
2226 if (!pRunLoader) return kFALSE;
2228 // initialize CDB storage, run number, set CDB lock
2229 // thats for the case of running HLT simulation without all the other steps
2230 // multiple calls are handled by the function, so we can just call
2232 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2235 // init the HLT simulation
2237 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2238 TString detStr = fWriteRawData;
2239 if (!IsSelected("HLT", detStr)) {
2240 options+=" writerawfiles=";
2242 options+=" writerawfiles=HLT";
2245 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2246 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2247 // in order to get detector data. By default, RawReaderFile is used to read
2248 // the already simulated ddl files. Date and Root files from the raw data
2249 // are generated after the HLT simulation.
2250 options+=" rawfile=./";
2253 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2254 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2255 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2257 // run the HLT simulation
2258 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2259 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2260 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2264 // delete the instance
2265 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2266 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2267 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2271 return iResult>=0?kTRUE:kFALSE;
2274 //_____________________________________________________________________________
2275 Bool_t AliSimulation::RunQA()
2277 // run the QA on summable hits, digits or digits
2279 //if(!gAlice) return kFALSE;
2280 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2282 TString detectorsw("") ;
2284 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2285 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2286 if ( detectorsw.IsNull() )
2291 //_____________________________________________________________________________
2292 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2294 // Allows to run QA for a selected set of detectors
2295 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2296 // all selected detectors run the same selected tasks
2298 if (!detAndAction.Contains(":")) {
2299 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2303 Int_t colon = detAndAction.Index(":") ;
2304 fQADetectors = detAndAction(0, colon) ;
2305 if (fQADetectors.Contains("ALL") ){
2306 TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2307 Int_t minus = fQADetectors.Last('-') ;
2308 TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2309 TString toRemove("") ;
2310 while (minus >= 0) {
2311 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
2312 toRemove = toRemove.Strip() ;
2313 toKeep.ReplaceAll(toRemove, "") ;
2314 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
2315 minus = fQADetectors.Last('-') ;
2317 fQADetectors = toKeep ;
2319 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2320 if (fQATasks.Contains("ALL") ) {
2321 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2323 fQATasks.ToUpper() ;
2325 if ( fQATasks.Contains("HIT") )
2326 tempo = Form("%d ", AliQAv1::kHITS) ;
2327 if ( fQATasks.Contains("SDIGIT") )
2328 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2329 if ( fQATasks.Contains("DIGIT") )
2330 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2332 if (fQATasks.IsNull()) {
2333 AliInfo("No QA requested\n") ;
2338 TString tempo(fQATasks) ;
2339 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2340 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2341 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2342 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2344 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2345 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2346 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2347 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2352 //_____________________________________________________________________________
2353 void AliSimulation::ProcessEnvironmentVars()
2355 // Extract run number and random generator seed from env variables
2357 AliInfo("Processing environment variables");
2359 // Random Number seed
2361 // first check that seed is not already set
2363 if (gSystem->Getenv("CONFIG_SEED")) {
2364 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2367 if (gSystem->Getenv("CONFIG_SEED")) {
2368 AliInfo(Form("Seed for random number generation already set (%d)"
2369 ": CONFIG_SEED variable ignored!", fSeed));
2373 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2377 // first check that run number is not already set
2379 if (gSystem->Getenv("DC_RUN")) {
2380 fRun = atoi(gSystem->Getenv("DC_RUN"));
2383 if (gSystem->Getenv("DC_RUN")) {
2384 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2388 AliInfo(Form("Run number = %d", fRun));
2391 //---------------------------------------------------------------------
2392 void AliSimulation::WriteGRPEntry()
2394 // Get the necessary information from galice (generator, trigger etc) and
2395 // write a GRP entry corresponding to the settings in the Config.C used
2396 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2399 AliInfo("Writing global run parameters entry into the OCDB");
2401 AliGRPObject* grpObj = new AliGRPObject();
2403 grpObj->SetRunType("PHYSICS");
2404 grpObj->SetTimeStart(fTimeStart);
2405 grpObj->SetTimeEnd(fTimeEnd);
2406 grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2408 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2414 gen->GetProjectile(projectile,a,z);
2416 gen->GetTarget(target,a,z);
2417 TString beamType = projectile + "-" + target;
2418 beamType.ReplaceAll(" ","");
2419 if (!beamType.CompareTo("-")) {
2420 grpObj->SetBeamType("UNKNOWN");
2421 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2424 grpObj->SetBeamType(beamType);
2426 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 * a / z);
2428 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 );
2430 // Heavy ion run, the event specie is set to kHighMult
2431 fEventSpecie = AliRecoParam::kHighMult;
2432 if ((strcmp(beamType,"p-p") == 0) ||
2433 (strcmp(beamType,"p-") == 0) ||
2434 (strcmp(beamType,"-p") == 0) ||
2435 (strcmp(beamType,"P-P") == 0) ||
2436 (strcmp(beamType,"P-") == 0) ||
2437 (strcmp(beamType,"-P") == 0)) {
2438 // Proton run, the event specie is set to kLowMult
2439 fEventSpecie = AliRecoParam::kLowMult;
2443 AliWarning("Unknown beam type and energy! Setting energy to 0");
2444 grpObj->SetBeamEnergy(0);
2445 grpObj->SetBeamType("UNKNOWN");
2448 UInt_t detectorPattern = 0;
2450 TObjArray *detArray = gAlice->Detectors();
2451 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2452 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2453 AliDebug(1, Form("Detector #%d found: %s", iDet, AliDAQ::OfflineModuleName(iDet)));
2454 detectorPattern |= (1 << iDet);
2459 if (!fTriggerConfig.IsNull())
2460 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2463 if (!fRunHLT.IsNull())
2464 detectorPattern |= (1 << AliDAQ::kHLTId);
2466 grpObj->SetNumberOfDetectors((Char_t)nDets);
2467 grpObj->SetDetectorMask((Int_t)detectorPattern);
2468 grpObj->SetLHCPeriod("LHC08c");
2469 grpObj->SetLHCState("STABLE_BEAMS");
2471 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2472 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2474 Float_t factorSol = field ? field->GetFactorSol() : 0;
2475 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2476 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2478 Float_t factorDip = field ? field->GetFactorDip() : 0;
2479 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2481 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2482 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2483 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2484 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2485 if (field) grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2486 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2488 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2490 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2492 // Now store the entry in OCDB
2493 AliCDBManager* man = AliCDBManager::Instance();
2495 man->SetLock(0, fKey);
2497 AliCDBStorage* sto = man->GetStorage(fGRPWriteLocation.Data());
2500 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2501 AliCDBMetaData *metadata= new AliCDBMetaData();
2503 metadata->SetResponsible("alice-off@cern.ch");
2504 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2506 sto->Put(grpObj,id,metadata);
2507 man->SetLock(1, fKey);
2510 //_____________________________________________________________________________
2511 time_t AliSimulation::GenerateTimeStamp() const
2513 // Generate event time-stamp according to
2514 // SOR/EOR time from GRP
2515 if (fUseTimeStampFromCDB)
2516 return fTimeStart + gRandom->Integer(fTimeEnd-fTimeStart);