1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for running generation, simulation and digitization //
22 // Hits, sdigits and digits are created for all detectors by typing: //
24 // AliSimulation sim; //
27 // The Run method returns kTRUE in case of successful execution. //
28 // The number of events can be given as argument to the Run method or it //
31 // sim.SetNumberOfEvents(n); //
33 // The name of the configuration file can be passed as argument to the //
34 // AliSimulation constructor or can be specified by //
36 // sim.SetConfigFile("..."); //
38 // The generation of particles and the simulation of detector hits can be //
39 // switched on or off by //
41 // sim.SetRunGeneration(kTRUE); // generation of primary particles //
42 // sim.SetRunSimulation(kFALSE); // but no tracking //
44 // For which detectors sdigits and digits will be created, can be steered //
47 // sim.SetMakeSDigits("ALL"); // make sdigits for all detectors //
48 // sim.SetMakeDigits("ITS TPC"); // make digits only for ITS and TPC //
50 // The argument is a (case sensitive) string with the names of the //
51 // detectors separated by a space. An empty string ("") can be used to //
52 // disable the creation of sdigits or digits. The special string "ALL" //
53 // selects all available detectors. This is the default. //
55 // The creation of digits from hits instead of from sdigits can be selected //
58 // sim.SetMakeDigitsFromHits("TRD"); //
60 // The argument is again a string with the selected detectors. Be aware that //
61 // this feature is not available for all detectors and that merging is not //
62 // possible, when digits are created directly from hits. //
64 // Background events can be merged by calling //
66 // sim.MergeWith("background/galice.root", 2); //
68 // The first argument is the file name of the background galice file. The //
69 // second argument is the number of signal events per background event. //
70 // By default this number is calculated from the number of available //
71 // background events. MergeWith can be called several times to merge more //
72 // than two event streams. It is assumed that the sdigits were already //
73 // produced for the background events. //
75 // The output of raw data can be switched on by calling //
77 // sim.SetWriteRawData("MUON"); // write raw data for MUON //
79 // The default output format of the raw data are DDL files. They are //
80 // converted to a DATE file, if a file name is given as second argument. //
81 // For this conversion the program "dateStream" is required. If the file //
82 // name has the extension ".root", the DATE file is converted to a root //
83 // file. The program "alimdc" is used for this purpose. For the conversion //
84 // to DATE and root format the two conversion programs have to be installed. //
85 // Only the raw data in the final format is kept if the third argument is //
88 // The methods RunSimulation, RunSDigitization, RunDigitization, //
89 // RunHitsDigitization and WriteRawData can be used to run only parts of //
90 // the full simulation chain. The creation of raw data DDL files and their //
91 // conversion to the DATE or root format can be run directly by calling //
92 // the methods WriteRawFiles, ConvertRawFilesToDate and ConvertDateToRoot. //
94 // The default number of events per file, which is usually set in the //
95 // config file, can be changed for individual detectors and data types //
98 // sim.SetEventsPerFile("PHOS", "Reconstructed Points", 3); //
100 // The first argument is the detector, the second one the data type and the //
101 // last one the number of events per file. Valid data types are "Hits", //
102 // "Summable Digits", "Digits", "Reconstructed Points" and "Tracks". //
103 // The number of events per file has to be set before the simulation of //
104 // hits. Otherwise it has no effect. //
106 ///////////////////////////////////////////////////////////////////////////////
110 #include <TGeoGlobalMagField.h>
111 #include <TGeoManager.h>
112 #include <TObjString.h>
115 #include <TVirtualMC.h>
116 #include <TVirtualMCApplication.h>
118 #include <TInterpreter.h>
120 #include "AliAlignObj.h"
121 #include "AliCDBEntry.h"
122 #include "AliCDBManager.h"
123 #include "AliGRPManager.h"
124 #include "AliCDBStorage.h"
125 #include "AliCTPRawData.h"
126 #include "AliCentralTrigger.h"
127 #include "AliCentralTrigger.h"
128 #include "AliCodeTimer.h"
130 #include "AliDigitizer.h"
131 #include "AliESDEvent.h"
132 #include "AliGRPObject.h"
133 #include "AliGenEventHeader.h"
134 #include "AliGenerator.h"
135 #include "AliGeomManager.h"
136 #include "AliHLTSimulation.h"
137 #include "AliHeader.h"
139 #include "AliLegoGenerator.h"
143 #include "AliModule.h"
145 #include "AliRawReaderDate.h"
146 #include "AliRawReaderFile.h"
147 #include "AliRawReaderRoot.h"
149 #include "AliDigitizationInput.h"
150 #include "AliRunLoader.h"
151 #include "AliStack.h"
152 #include "AliSimulation.h"
153 #include "AliSysInfo.h"
154 #include "AliVertexGenFile.h"
157 ClassImp(AliSimulation)
159 AliSimulation *AliSimulation::fgInstance = 0;
160 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE"
161 // #ifdef MFT_UPGRADE
168 //_____________________________________________________________________________
169 AliSimulation::AliSimulation(const char* configFileName,
170 const char* name, const char* title) :
173 fRunGeneratorOnly(kFALSE),
174 fRunGeneration(kTRUE),
175 fRunSimulation(kTRUE),
176 fLoadAlignFromCDB(kTRUE),
177 fLoadAlObjsListOfDets("ALL"),
181 fMakeDigitsFromHits(""),
183 fRawDataFileName(""),
184 fDeleteIntermediateFiles(kFALSE),
185 fWriteSelRawData(kFALSE),
186 fStopOnError(kFALSE),
187 fUseMonitoring(kFALSE),
189 fConfigFileName(configFileName),
190 fGAliceFileName("galice.root"),
192 fBkgrdFileNames(NULL),
193 fAlignObjArray(NULL),
194 fUseBkgrdVertex(kTRUE),
195 fRegionOfInterest(kFALSE),
201 fInitCDBCalled(kFALSE),
202 fInitRunNumberCalled(kFALSE),
203 fSetRunNumberFromDataCalled(kFALSE),
204 fEmbeddingFlag(kFALSE),
207 fUseVertexFromCDB(0),
208 fUseMagFieldFromGRP(0),
209 fGRPWriteLocation(Form("local://%s", gSystem->pwd())),
210 fUseTimeStampFromCDB(0),
216 fEventSpecie(AliRecoParam::kDefault),
217 fWriteQAExpertData(kTRUE),
221 fWriteGRPEntry(kTRUE)
223 // create simulation object with default parameters
225 SetGAliceFile("galice.root");
228 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
229 qam->SetActiveDetectors(fQADetectors) ;
230 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
231 qam->SetTasks(fQATasks) ;
234 //_____________________________________________________________________________
235 AliSimulation::~AliSimulation()
239 fEventsPerFile.Delete();
240 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
241 // delete fAlignObjArray; fAlignObjArray=0;
243 if (fBkgrdFileNames) {
244 fBkgrdFileNames->Delete();
245 delete fBkgrdFileNames;
248 fSpecCDBUri.Delete();
249 if (fgInstance==this) fgInstance = 0;
251 AliQAManager::QAManager()->ShowQA() ;
252 AliQAManager::Destroy() ;
253 AliCodeTimer::Instance()->Print();
257 //_____________________________________________________________________________
258 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
260 // set the number of events for one run
265 //_____________________________________________________________________________
266 void AliSimulation::InitQA()
268 // activate a default CDB storage
269 // First check if we have any CDB storage set, because it is used
270 // to retrieve the calibration and alignment constants
272 if (fInitCDBCalled) return;
273 fInitCDBCalled = kTRUE;
275 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
276 qam->SetActiveDetectors(fQADetectors) ;
277 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
278 qam->SetTasks(fQATasks) ;
279 if (fWriteQAExpertData)
280 qam->SetWriteExpert() ;
282 if (qam->IsDefaultStorageSet()) {
283 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
284 AliWarning("Default QA reference storage has been already set !");
285 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
286 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
287 fQARefUri = qam->GetDefaultStorage()->GetURI();
289 if (fQARefUri.Length() > 0) {
290 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
291 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
292 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
294 fQARefUri="local://$ALICE_ROOT/QARef";
295 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
296 AliWarning("Default QA reference storage not yet set !!!!");
297 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
298 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
300 qam->SetDefaultStorage(fQARefUri);
304 //_____________________________________________________________________________
305 void AliSimulation::InitCDB()
307 // activate a default CDB storage
308 // First check if we have any CDB storage set, because it is used
309 // to retrieve the calibration and alignment constants
311 if (fInitCDBCalled) return;
312 fInitCDBCalled = kTRUE;
314 AliCDBManager* man = AliCDBManager::Instance();
315 if (man->IsDefaultStorageSet())
317 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
318 AliWarning("Default CDB storage has been already set !");
319 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
320 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
321 fCDBUri = man->GetDefaultStorage()->GetURI();
324 if (fCDBUri.Length() > 0)
326 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
327 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
328 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
330 fCDBUri="local://$ALICE_ROOT/OCDB";
331 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
332 AliWarning("Default CDB storage not yet set !!!!");
333 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
334 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
337 man->SetDefaultStorage(fCDBUri);
340 // Now activate the detector specific CDB storage locations
341 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
342 TObject* obj = fSpecCDBUri[i];
344 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
345 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
346 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
347 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
352 //_____________________________________________________________________________
353 void AliSimulation::InitRunNumber(){
354 // check run number. If not set, set it to 0 !!!!
356 if (fInitRunNumberCalled) return;
357 fInitRunNumberCalled = kTRUE;
359 AliCDBManager* man = AliCDBManager::Instance();
360 if (man->GetRun() >= 0)
362 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
363 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
367 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
370 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
379 //_____________________________________________________________________________
380 void AliSimulation::SetCDBLock() {
381 // Set CDB lock: from now on it is forbidden to reset the run number
382 // or the default storage or to activate any further storage!
384 ULong64_t key = AliCDBManager::Instance()->SetLock(1);
388 //_____________________________________________________________________________
389 void AliSimulation::SetDefaultStorage(const char* uri) {
390 // Store the desired default CDB storage location
391 // Activate it later within the Run() method
397 //_____________________________________________________________________________
398 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
399 // Store the desired default CDB storage location
400 // Activate it later within the Run() method
403 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
406 //_____________________________________________________________________________
407 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
408 // Store a detector-specific CDB storage location
409 // Activate it later within the Run() method
411 AliCDBPath aPath(calibType);
412 if(!aPath.IsValid()){
413 AliError(Form("Not a valid path: %s", calibType));
417 TObject* obj = fSpecCDBUri.FindObject(calibType);
418 if (obj) fSpecCDBUri.Remove(obj);
419 fSpecCDBUri.Add(new TNamed(calibType, uri));
423 //_____________________________________________________________________________
424 void AliSimulation::SetRunNumber(Int_t run)
427 // Activate it later within the Run() method
432 //_____________________________________________________________________________
433 void AliSimulation::SetSeed(Int_t seed)
436 // Activate it later within the Run() method
441 //_____________________________________________________________________________
442 Bool_t AliSimulation::SetRunNumberFromData()
444 // Set the CDB manager run number
445 // The run number is retrieved from gAlice
447 if (fSetRunNumberFromDataCalled) return kTRUE;
448 fSetRunNumberFromDataCalled = kTRUE;
450 AliCDBManager* man = AliCDBManager::Instance();
451 Int_t runData = -1, runCDB = -1;
453 AliRunLoader* runLoader = LoadRun("READ");
454 if (!runLoader) return kFALSE;
456 runData = runLoader->GetHeader()->GetRun();
460 runCDB = man->GetRun();
462 if (runCDB != runData) {
463 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
464 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
465 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
466 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
471 man->SetRun(runData);
474 if(man->GetRun() < 0) {
475 AliError("Run number not properly initalized!");
484 //_____________________________________________________________________________
485 void AliSimulation::SetConfigFile(const char* fileName)
487 // set the name of the config file
489 fConfigFileName = fileName;
492 //_____________________________________________________________________________
493 void AliSimulation::SetGAliceFile(const char* fileName)
495 // set the name of the galice file
496 // the path is converted to an absolute one if it is relative
498 fGAliceFileName = fileName;
499 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
500 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
502 fGAliceFileName = absFileName;
503 delete[] absFileName;
506 AliDebug(2, Form("galice file name set to %s", fileName));
509 //_____________________________________________________________________________
510 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
513 // set the number of events per file for the given detector and data type
514 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
516 TNamed* obj = new TNamed(detector, type);
517 obj->SetUniqueID(nEvents);
518 fEventsPerFile.Add(obj);
521 //_____________________________________________________________________________
522 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
524 // Read the alignment objects from CDB.
525 // Each detector is supposed to have the
526 // alignment objects in DET/Align/Data CDB path.
527 // All the detector objects are then collected,
528 // sorted by geometry level (starting from ALIC) and
529 // then applied to the TGeo geometry.
530 // Finally an overlaps check is performed.
532 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
533 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
537 // initialize CDB storage, run number, set CDB lock
539 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
542 Bool_t delRunLoader = kFALSE;
544 runLoader = LoadRun("READ");
545 if (!runLoader) return kFALSE;
546 delRunLoader = kTRUE;
549 // Export ideal geometry
550 if(!IsGeometryFromFile()) AliGeomManager::GetGeometry()->Export("geometry.root");
552 // Load alignment data from CDB and apply to geometry through AliGeomManager
553 if(fLoadAlignFromCDB){
555 TString detStr = fLoadAlObjsListOfDets;
556 TString loadAlObjsListOfDets = "";
558 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
559 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
560 AliModule* det = (AliModule*) detArray->At(iDet);
561 if (!det || !det->IsActive()) continue;
562 if (IsSelected(det->GetName(), detStr)) {
563 //add det to list of dets to be aligned from CDB
564 loadAlObjsListOfDets += det->GetName();
565 loadAlObjsListOfDets += " ";
567 } // end loop over detectors
568 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
569 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
571 // Check if the array with alignment objects was
572 // provided by the user. If yes, apply the objects
573 // to the present TGeo geometry
574 if (fAlignObjArray) {
575 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
576 AliError("The misalignment of one or more volumes failed!"
577 "Compare the list of simulated detectors and the list of detector alignment data!");
578 if (delRunLoader) delete runLoader;
584 // Update the internal geometry of modules (ITS needs it)
585 TString detStr = fLoadAlObjsListOfDets;
586 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
587 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
589 AliModule* det = (AliModule*) detArray->At(iDet);
590 if (!det || !det->IsActive()) continue;
591 if (IsSelected(det->GetName(), detStr)) {
592 det->UpdateInternalGeometry();
594 } // end loop over detectors
597 if (delRunLoader) delete runLoader;
602 //_____________________________________________________________________________
603 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
605 // add a file with background events for merging
607 TObjString* fileNameStr = new TObjString(fileName);
608 fileNameStr->SetUniqueID(nSignalPerBkgrd);
609 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
610 fBkgrdFileNames->Add(fileNameStr);
613 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
615 // add a file with background events for embeddin
616 MergeWith(fileName, nSignalPerBkgrd);
617 fEmbeddingFlag = kTRUE;
620 //_____________________________________________________________________________
621 Bool_t AliSimulation::Run(Int_t nEvents)
623 // run the generation, simulation and digitization
626 AliCodeTimerAuto("",0)
627 AliSysInfo::AddStamp("Start_Run");
629 // Load run number and seed from environmental vars
630 ProcessEnvironmentVars();
631 AliSysInfo::AddStamp("ProcessEnvironmentVars");
633 gRandom->SetSeed(fSeed);
635 if (nEvents > 0) fNEvents = nEvents;
637 // Run generator-only code on demand
638 if (fRunGeneratorOnly)
640 if(!RunGeneratorOnly())
642 if (fStopOnError) return kFALSE;
648 // create and setup the HLT instance
649 if (!fRunHLT.IsNull() && !CreateHLT()) {
650 if (fStopOnError) return kFALSE;
655 // generation and simulation -> hits
656 if (fRunGeneration) {
657 if (!RunSimulation()) if (fStopOnError) return kFALSE;
659 AliSysInfo::AddStamp("RunSimulation");
661 // initialize CDB storage from external environment
662 // (either CDB manager or AliSimulation setters),
663 // if not already done in RunSimulation()
665 AliSysInfo::AddStamp("InitCDB");
667 // Set run number in CDBManager from data
668 // From this point on the run number must be always loaded from data!
669 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
671 // Set CDB lock: from now on it is forbidden to reset the run number
672 // or the default storage or to activate any further storage!
675 // If RunSimulation was not called, load the geometry and misalign it
676 if (!AliGeomManager::GetGeometry()) {
677 // Initialize the geometry manager
678 AliGeomManager::LoadGeometry("geometry.root");
679 AliSysInfo::AddStamp("GetGeometry");
680 // // Check that the consistency of symbolic names for the activated subdetectors
681 // // in the geometry loaded by AliGeomManager
682 // AliRunLoader* runLoader = LoadRun("READ");
683 // if (!runLoader) return kFALSE;
685 // TString detsToBeChecked = "";
686 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
687 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
688 // AliModule* det = (AliModule*) detArray->At(iDet);
689 // if (!det || !det->IsActive()) continue;
690 // detsToBeChecked += det->GetName();
691 // detsToBeChecked += " ";
692 // } // end loop over detectors
693 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
694 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
695 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
697 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
699 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
701 AliSysInfo::AddStamp("MissalignGeometry");
704 // hits -> summable digits
705 AliSysInfo::AddStamp("Start_sdigitization");
706 if (!fMakeSDigits.IsNull()) {
707 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
710 AliSysInfo::AddStamp("Stop_sdigitization");
712 AliSysInfo::AddStamp("Start_digitization");
713 // summable digits -> digits
714 if (!fMakeDigits.IsNull()) {
715 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
716 if (fStopOnError) return kFALSE;
719 AliSysInfo::AddStamp("Stop_digitization");
724 if (!fMakeDigitsFromHits.IsNull()) {
725 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
726 AliWarning(Form("Merging and direct creation of digits from hits "
727 "was selected for some detectors. "
728 "No merging will be done for the following detectors: %s",
729 fMakeDigitsFromHits.Data()));
731 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
732 if (fStopOnError) return kFALSE;
736 AliSysInfo::AddStamp("Hits2Digits");
740 if (!fTriggerConfig.IsNull() && !RunTrigger(fTriggerConfig,fMakeDigits)) {
741 if (fStopOnError) return kFALSE;
744 AliSysInfo::AddStamp("RunTrigger");
747 // digits -> raw data
748 if (!fWriteRawData.IsNull()) {
749 if (!WriteRawData(fWriteRawData, fRawDataFileName,
750 fDeleteIntermediateFiles,fWriteSelRawData)) {
751 if (fStopOnError) return kFALSE;
755 AliSysInfo::AddStamp("WriteRaw");
757 // run HLT simulation on simulated digit data if raw data is not
758 // simulated, otherwise its called as part of WriteRawData
759 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
761 if (fStopOnError) return kFALSE;
765 AliSysInfo::AddStamp("RunHLT");
769 Bool_t rv = RunQA() ;
775 AliSysInfo::AddStamp("RunQA");
777 TString snapshotFileOut("");
778 if(TString(gSystem->Getenv("OCDB_SNAPSHOT_CREATE")) == TString("kTRUE")){
779 AliInfo(" ******** Creating the snapshot! *********");
780 TString snapshotFile(gSystem->Getenv("OCDB_SNAPSHOT_FILENAME"));
781 if(!(snapshotFile.IsNull() || snapshotFile.IsWhitespace()))
782 snapshotFileOut = snapshotFile;
784 snapshotFileOut="OCDB.root";
785 AliCDBManager::Instance()->DumpToSnapshotFile(snapshotFileOut.Data(),kFALSE);
788 // Cleanup of CDB manager: cache and active storages!
789 AliCDBManager::Instance()->ClearCache();
794 //_______________________________________________________________________
795 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
796 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
797 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
800 // Generates lego plots of:
801 // - radiation length map phi vs theta
802 // - radiation length map phi vs eta
803 // - interaction length map
804 // - g/cm2 length map
806 // ntheta bins in theta, eta
807 // themin minimum angle in theta (degrees)
808 // themax maximum angle in theta (degrees)
810 // phimin minimum angle in phi (degrees)
811 // phimax maximum angle in phi (degrees)
812 // rmin minimum radius
813 // rmax maximum radius
816 // The number of events generated = ntheta*nphi
817 // run input parameters in macro setup (default="Config.C")
819 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
822 <img src="picts/AliRunLego1.gif">
827 <img src="picts/AliRunLego2.gif">
832 <img src="picts/AliRunLego3.gif">
837 // run the generation and simulation
839 AliCodeTimerAuto("",0)
841 // initialize CDB storage and run number from external environment
842 // (either CDB manager or AliSimulation setters)
848 AliError("no gAlice object. Restart aliroot and try again.");
851 if (gAlice->Modules()->GetEntries() > 0) {
852 AliError("gAlice was already run. Restart aliroot and try again.");
856 AliInfo(Form("initializing gAlice with config file %s",
857 fConfigFileName.Data()));
860 if (nev == -1) nev = nc1 * nc2;
862 // check if initialisation has been done
863 // If runloader has been initialized, set the number of events per file to nc1 * nc2
866 if (!gener) gener = new AliLegoGenerator();
868 // Configure Generator
870 gener->SetRadiusRange(rmin, rmax);
871 gener->SetZMax(zmax);
872 gener->SetCoor1Range(nc1, c1min, c1max);
873 gener->SetCoor2Range(nc2, c2min, c2max);
877 fLego = new AliLego("lego",gener);
879 //__________________________________________________________________________
883 gROOT->LoadMacro(setup);
884 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
886 if(AliCDBManager::Instance()->GetRun() >= 0) {
887 SetRunNumber(AliCDBManager::Instance()->GetRun());
889 AliWarning("Run number not initialized!!");
892 AliRunLoader::Instance()->CdGAFile();
894 AliPDG::AddParticlesToPdgDataBase();
896 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
898 gAlice->GetMCApp()->Init();
901 //Must be here because some MCs (G4) adds detectors here and not in Config.C
902 gAlice->InitLoaders();
903 AliRunLoader::Instance()->MakeTree("E");
906 // Save stuff at the beginning of the file to avoid file corruption
907 AliRunLoader::Instance()->CdGAFile();
910 //Save current generator
911 AliGenerator *gen=gAlice->GetMCApp()->Generator();
912 gAlice->GetMCApp()->ResetGenerator(gener);
913 //Prepare MC for Lego Run
919 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
920 gMC->ProcessRun(nev);
922 // End of this run, close files
924 // Restore current generator
925 gAlice->GetMCApp()->ResetGenerator(gen);
926 // Delete Lego Object
932 //_____________________________________________________________________________
933 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
937 AliCodeTimerAuto("",0)
939 // initialize CDB storage from external environment
940 // (either CDB manager or AliSimulation setters),
941 // if not already done in RunSimulation()
944 // Set run number in CDBManager from data
945 // From this point on the run number must be always loaded from data!
946 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
948 // Set CDB lock: from now on it is forbidden to reset the run number
949 // or the default storage or to activate any further storage!
952 AliRunLoader* runLoader = LoadRun("READ");
953 if (!runLoader) return kFALSE;
954 TString trconfiguration = config;
956 if (trconfiguration.IsNull()) {
957 if(!fTriggerConfig.IsNull()) {
958 trconfiguration = fTriggerConfig;
961 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
964 runLoader->MakeTree( "GG" );
965 AliCentralTrigger* aCTP = runLoader->GetTrigger();
966 // Load Configuration
967 if (!aCTP->LoadConfiguration( trconfiguration ))
971 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
983 //_____________________________________________________________________________
984 Bool_t AliSimulation::WriteTriggerRawData()
986 // Writes the CTP (trigger) DDL raw data
987 // Details of the format are given in the
988 // trigger TDR - pages 134 and 135.
989 AliCTPRawData writer;
995 //_____________________________________________________________________________
996 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
998 // run the generation and simulation
1000 AliCodeTimerAuto("",0)
1002 // initialize CDB storage and run number from external environment
1003 // (either CDB manager or AliSimulation setters)
1004 AliSysInfo::AddStamp("RunSimulation_Begin");
1006 AliSysInfo::AddStamp("RunSimulation_InitCDB");
1010 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
1013 AliError("no gAlice object. Restart aliroot and try again.");
1016 if (gAlice->Modules()->GetEntries() > 0) {
1017 AliError("gAlice was already run. Restart aliroot and try again.");
1021 // Setup monitoring if requested
1022 gAlice->GetMCApp()->SetUseMonitoring(fUseMonitoring);
1024 AliInfo(Form("initializing gAlice with config file %s",
1025 fConfigFileName.Data()));
1028 // Initialize ALICE Simulation run
1033 // If requested set the mag. field from the GRP entry.
1034 // After this the field is loccked and cannot be changed by Config.C
1035 if (fUseMagFieldFromGRP) {
1037 grpM.ReadGRPEntry();
1039 AliInfo("Field is locked now. It cannot be changed in Config.C");
1043 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1044 gROOT->LoadMacro(fConfigFileName.Data());
1045 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1046 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1047 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1049 AliSysInfo::AddStamp("RunSimulation_Config");
1052 // If requested obtain the vertex position and vertex sigma_z from the CDB
1053 // This overwrites the settings from the Config.C
1054 if (fUseVertexFromCDB) {
1055 Double_t vtxPos[3] = {0., 0., 0.};
1056 Double_t vtxSig[3] = {0., 0., 0.};
1057 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1059 AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1061 if(vertex->GetXRes()>2.8) { // > pipe radius --> it's a dummy object, don't use it
1062 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1063 if (entry) vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
1067 vertex->GetXYZ(vtxPos);
1068 vertex->GetSigmaXYZ(vtxSig);
1069 AliInfo("Overwriting Config.C vertex settings !");
1070 AliInfo(Form("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
1071 vtxPos[0], vtxPos[1], vtxPos[2], vtxSig[2]));
1073 AliGenerator *gen = gAlice->GetMCApp()->Generator();
1074 gen->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]); // vertex position
1075 gen->SetSigmaZ(vtxSig[2]);
1080 // If requested we take the SOR and EOR time-stamps from the GRP and use them
1081 // in order to generate the event time-stamps
1082 if (fUseTimeStampFromCDB) {
1084 grpM.ReadGRPEntry();
1085 const AliGRPObject *grpObj = grpM.GetGRPData();
1086 if (!grpObj || (grpObj->GetTimeEnd() <= grpObj->GetTimeStart())) {
1087 AliError("Missing GRP or bad SOR/EOR time-stamps! Switching off the time-stamp generation from GRP!");
1088 fTimeStart = fTimeEnd = 0;
1089 fUseTimeStampFromCDB = kFALSE;
1092 fTimeStart = grpObj->GetTimeStart();
1093 fTimeEnd = grpObj->GetTimeEnd();
1097 if(AliCDBManager::Instance()->GetRun() >= 0) {
1098 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
1099 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
1101 AliWarning("Run number not initialized!!");
1104 AliRunLoader::Instance()->CdGAFile();
1107 AliPDG::AddParticlesToPdgDataBase();
1109 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1110 AliSysInfo::AddStamp("RunSimulation_GetField");
1111 gAlice->GetMCApp()->Init();
1112 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1114 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1115 gAlice->InitLoaders();
1116 AliRunLoader::Instance()->MakeTree("E");
1117 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1118 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1119 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1121 // Save stuff at the beginning of the file to avoid file corruption
1122 AliRunLoader::Instance()->CdGAFile();
1124 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1125 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1126 //___________________________________________________________________________________________
1128 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1130 // Set run number in CDBManager
1131 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1133 AliRunLoader* runLoader = AliRunLoader::Instance();
1135 AliError(Form("gAlice has no run loader object. "
1136 "Check your config file: %s", fConfigFileName.Data()));
1139 SetGAliceFile(runLoader->GetFileName());
1141 // Misalign geometry
1142 #if ROOT_VERSION_CODE < 331527
1143 AliGeomManager::SetGeometry(gGeoManager);
1145 // Check that the consistency of symbolic names for the activated subdetectors
1146 // in the geometry loaded by AliGeomManager
1147 TString detsToBeChecked = "";
1148 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1149 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1150 AliModule* det = (AliModule*) detArray->At(iDet);
1151 if (!det || !det->IsActive()) continue;
1152 detsToBeChecked += det->GetName();
1153 detsToBeChecked += " ";
1154 } // end loop over detectors
1155 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1156 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1157 MisalignGeometry(runLoader);
1158 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1161 // AliRunLoader* runLoader = AliRunLoader::Instance();
1162 // if (!runLoader) {
1163 // AliError(Form("gAlice has no run loader object. "
1164 // "Check your config file: %s", fConfigFileName.Data()));
1167 // SetGAliceFile(runLoader->GetFileName());
1169 if (!gAlice->GetMCApp()->Generator()) {
1170 AliError(Form("gAlice has no generator object. "
1171 "Check your config file: %s", fConfigFileName.Data()));
1175 // Write GRP entry corresponding to the setting found in Cofig.C
1178 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1180 if (nEvents <= 0) nEvents = fNEvents;
1182 // get vertex from background file in case of merging
1183 if (fUseBkgrdVertex &&
1184 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1185 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1186 const char* fileName = ((TObjString*)
1187 (fBkgrdFileNames->At(0)))->GetName();
1188 AliInfo(Form("The vertex will be taken from the background "
1189 "file %s with nSignalPerBackground = %d",
1190 fileName, signalPerBkgrd));
1191 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1192 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1195 if (!fRunSimulation) {
1196 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1199 // set the number of events per file for given detectors and data types
1200 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1201 if (!fEventsPerFile[i]) continue;
1202 const char* detName = fEventsPerFile[i]->GetName();
1203 const char* typeName = fEventsPerFile[i]->GetTitle();
1204 TString loaderName(detName);
1205 loaderName += "Loader";
1206 AliLoader* loader = runLoader->GetLoader(loaderName);
1208 AliError(Form("RunSimulation no loader for %s found\n Number of events per file not set for %s %s",
1209 detName, typeName, detName));
1212 AliDataLoader* dataLoader =
1213 loader->GetDataLoader(typeName);
1215 AliError(Form("no data loader for %s found\n"
1216 "Number of events per file not set for %s %s",
1217 typeName, detName, typeName));
1220 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1221 AliDebug(1, Form("number of events per file set to %d for %s %s",
1222 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1225 AliInfo("running gAlice");
1226 AliSysInfo::AddStamp("Start_ProcessRun");
1228 // Create the Root Tree with one branch per detector
1229 //Hits moved to begin event -> now we are crating separate tree for each event
1230 gMC->ProcessRun(nEvents);
1232 // End of this run, close files
1233 if(nEvents>0) FinishRun();
1235 AliSysInfo::AddStamp("Stop_ProcessRun");
1241 //_____________________________________________________________________________
1242 Bool_t AliSimulation::RunGeneratorOnly()
1245 TInterpreter::EErrorCode interpreterError=TInterpreter::kNoError;
1246 gROOT->LoadMacro(fConfigFileName.Data());
1247 Long_t interpreterResult=gInterpreter->ProcessLine(gAlice->GetConfigFunction(), &interpreterError);
1248 if (interpreterResult!=0 || interpreterError!=TInterpreter::kNoError) {
1249 AliFatal(Form("execution of config file \"%s\" failed with error %d", fConfigFileName.Data(), (int)interpreterError));
1252 // Setup the runloader and generator, check if everything is OK
1253 AliRunLoader* runLoader = AliRunLoader::Instance();
1254 AliGenerator* generator = gAlice->GetMCApp()->Generator();
1256 AliError(Form("gAlice has no run loader object. "
1257 "Check your config file: %s", fConfigFileName.Data()));
1261 AliError(Form("gAlice has no generator object. "
1262 "Check your config file: %s", fConfigFileName.Data()));
1266 runLoader->LoadKinematics("RECREATE");
1267 runLoader->MakeTree("E");
1269 // Create stack and header
1270 runLoader->MakeStack();
1271 AliStack* stack = runLoader->Stack();
1272 AliHeader* header = runLoader->GetHeader();
1274 // Intialize generator
1276 generator->SetStack(stack);
1278 // Run main generator loop
1280 for (Int_t iev=0; iev<fNEvents; iev++)
1283 header->Reset(0,iev);
1284 runLoader->SetEventNumber(iev);
1286 runLoader->MakeTree("K");
1289 generator->Generate();
1292 header->SetNprimary(stack->GetNprimary());
1293 header->SetNtrack(stack->GetNtrack());
1294 stack->FinishEvent();
1295 header->SetStack(stack);
1296 runLoader->TreeE()->Fill();
1297 runLoader->WriteKinematics("OVERWRITE");
1301 generator->FinishRun();
1303 runLoader->WriteHeader("OVERWRITE");
1310 //_____________________________________________________________________________
1311 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1313 // run the digitization and produce summable digits
1314 static Int_t eventNr=0;
1315 AliCodeTimerAuto("",0) ;
1317 // initialize CDB storage, run number, set CDB lock
1319 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1322 AliRunLoader* runLoader = LoadRun();
1323 if (!runLoader) return kFALSE;
1325 TString detStr = detectors;
1326 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1327 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1328 AliModule* det = (AliModule*) detArray->At(iDet);
1329 if (!det || !det->IsActive()) continue;
1330 if (IsSelected(det->GetName(), detStr)) {
1331 AliInfo(Form("creating summable digits for %s", det->GetName()));
1332 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1333 det->Hits2SDigits();
1334 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1335 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1339 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1340 AliError(Form("the following detectors were not found: %s",
1342 if (fStopOnError) return kFALSE;
1351 //_____________________________________________________________________________
1352 Bool_t AliSimulation::RunDigitization(const char* detectors,
1353 const char* excludeDetectors)
1355 // run the digitization and produce digits from sdigits
1357 AliCodeTimerAuto("",0)
1359 // initialize CDB storage, run number, set CDB lock
1361 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1364 delete AliRunLoader::Instance();
1369 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1370 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1371 AliDigitizationInput digInp(nStreams, signalPerBkgrd);
1372 // digInp.SetEmbeddingFlag(fEmbeddingFlag);
1373 digInp.SetRegionOfInterest(fRegionOfInterest);
1374 digInp.SetInputStream(0, fGAliceFileName.Data());
1375 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1376 const char* fileName = ((TObjString*)(fBkgrdFileNames->At(iStream-1)))->GetName();
1377 digInp.SetInputStream(iStream, fileName);
1380 detArr.SetOwner(kTRUE);
1381 TString detStr = detectors;
1382 TString detExcl = excludeDetectors;
1383 if (!static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice()) {
1384 AliError("Error occured while getting gAlice from Input 0");
1387 AliRunLoader* runLoader = AliRunLoader::GetRunLoader(digInp.GetInputStream(0)->GetFolderName());
1388 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1389 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1390 AliModule* det = (AliModule*) detArray->At(iDet);
1391 if (!det || !det->IsActive()) continue;
1392 if (!IsSelected(det->GetName(), detStr) || IsSelected(det->GetName(), detExcl)) continue;
1393 AliDigitizer* digitizer = det->CreateDigitizer(&digInp);
1394 if (!digitizer || !digitizer->Init()) {
1395 AliError(Form("no digitizer for %s", det->GetName()));
1396 if (fStopOnError) return kFALSE;
1399 detArr.AddLast(digitizer);
1400 AliInfo(Form("Created digitizer from SDigits -> Digits for %s", det->GetName()));
1403 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1404 AliError(Form("the following detectors were not found: %s", detStr.Data()));
1405 if (fStopOnError) return kFALSE;
1408 Int_t ndigs = detArr.GetEntriesFast();
1409 Int_t eventsCreated = 0;
1410 AliRunLoader* outRl = digInp.GetOutRunLoader();
1411 while ((eventsCreated++ < fNEvents) || (fNEvents < 0)) {
1412 if (!digInp.ConnectInputTrees()) break;
1413 digInp.InitEvent(); //this must be after call of Connect Input tress.
1414 if (outRl) outRl->SetEventNumber(eventsCreated-1);
1415 static_cast<AliStream*>(digInp.GetInputStream(0))->ImportgAlice(); // use gAlice of the first input stream
1416 for (int id=0;id<ndigs;id++) ((AliDigitizer*)detArr[id])->Digitize("");
1417 digInp.FinishEvent();
1419 digInp.FinishGlobal();
1424 //_____________________________________________________________________________
1425 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1427 // run the digitization and produce digits from hits
1429 AliCodeTimerAuto("",0)
1431 // initialize CDB storage, run number, set CDB lock
1433 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1436 AliRunLoader* runLoader = LoadRun("READ");
1437 if (!runLoader) return kFALSE;
1439 TString detStr = detectors;
1440 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1441 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1442 AliModule* det = (AliModule*) detArray->At(iDet);
1443 if (!det || !det->IsActive()) continue;
1444 if (IsSelected(det->GetName(), detStr)) {
1445 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1450 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1451 AliError(Form("the following detectors were not found: %s",
1453 if (fStopOnError) return kFALSE;
1459 //_____________________________________________________________________________
1460 Bool_t AliSimulation::WriteRawData(const char* detectors,
1461 const char* fileName,
1462 Bool_t deleteIntermediateFiles,
1465 // convert the digits to raw data
1466 // First DDL raw data files for the given detectors are created.
1467 // If a file name is given, the DDL files are then converted to a DATE file.
1468 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1470 // If the file name has the extension ".root", the DATE file is converted
1472 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1473 // 'selrawdata' flag can be used to enable writing of detectors raw data
1474 // accoring to the trigger cluster.
1476 AliCodeTimerAuto("",0)
1477 AliSysInfo::AddStamp("WriteRawData_Start");
1479 TString detStr = detectors;
1480 if (!WriteRawFiles(detStr.Data())) {
1481 if (fStopOnError) return kFALSE;
1483 AliSysInfo::AddStamp("WriteRawFiles");
1485 // run HLT simulation on simulated DDL raw files
1486 // and produce HLT ddl raw files to be included in date/root file
1487 // bugfix 2009-06-26: the decision whether to write HLT raw data
1488 // is taken in RunHLT. Here HLT always needs to be run in order to
1489 // create HLT digits, unless its switched off. This is due to the
1490 // special placement of the HLT between the generation of DDL files
1491 // and conversion to DATE/Root file.
1492 detStr.ReplaceAll("HLT", "");
1493 if (!fRunHLT.IsNull()) {
1495 if (fStopOnError) return kFALSE;
1498 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1500 TString dateFileName(fileName);
1501 if (!dateFileName.IsNull()) {
1502 Bool_t rootOutput = dateFileName.EndsWith(".root");
1503 if (rootOutput) dateFileName += ".date";
1504 TString selDateFileName;
1506 selDateFileName = "selected.";
1507 selDateFileName+= dateFileName;
1509 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1510 if (fStopOnError) return kFALSE;
1512 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1513 if (deleteIntermediateFiles) {
1514 AliRunLoader* runLoader = LoadRun("READ");
1515 if (runLoader) for (Int_t iEvent = 0;
1516 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1518 snprintf(command, 256, "rm -r raw%d", iEvent);
1519 gSystem->Exec(command);
1525 if (!ConvertDateToRoot(dateFileName, fileName)) {
1526 if (fStopOnError) return kFALSE;
1528 AliSysInfo::AddStamp("ConvertDateToRoot");
1529 if (deleteIntermediateFiles) {
1530 gSystem->Unlink(dateFileName);
1533 TString selFileName = "selected.";
1534 selFileName += fileName;
1535 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1536 if (fStopOnError) return kFALSE;
1538 if (deleteIntermediateFiles) {
1539 gSystem->Unlink(selDateFileName);
1548 //_____________________________________________________________________________
1549 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1551 // convert the digits to raw data DDL files
1553 AliCodeTimerAuto("",0)
1555 AliRunLoader* runLoader = LoadRun("READ");
1556 if (!runLoader) return kFALSE;
1558 // write raw data to DDL files
1559 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1560 AliInfo(Form("processing event %d", iEvent));
1561 runLoader->GetEvent(iEvent);
1562 TString baseDir = gSystem->WorkingDirectory();
1564 snprintf(dirName, 256, "raw%d", iEvent);
1565 gSystem->MakeDirectory(dirName);
1566 if (!gSystem->ChangeDirectory(dirName)) {
1567 AliError(Form("couldn't change to directory %s", dirName));
1568 if (fStopOnError) return kFALSE; else continue;
1571 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1574 TString detStr = detectors;
1575 if (IsSelected("HLT", detStr)) {
1576 // Do nothing. "HLT" will be removed from detStr and HLT raw
1577 // data files are generated in RunHLT.
1580 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1581 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1582 AliModule* det = (AliModule*) detArray->At(iDet);
1583 if (!det || !det->IsActive()) continue;
1584 if (IsSelected(det->GetName(), detStr)) {
1585 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1590 if (!WriteTriggerRawData())
1591 if (fStopOnError) return kFALSE;
1593 gSystem->ChangeDirectory(baseDir);
1594 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1595 AliError(Form("the following detectors were not found: %s",
1597 if (fStopOnError) return kFALSE;
1606 //_____________________________________________________________________________
1607 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1608 const char* selDateFileName)
1610 // convert raw data DDL files to a DATE file with the program "dateStream"
1611 // The second argument is not empty when the user decides to write
1612 // the detectors raw data according to the trigger cluster.
1614 AliCodeTimerAuto("",0)
1616 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1618 AliError("the program dateStream was not found");
1619 if (fStopOnError) return kFALSE;
1624 AliRunLoader* runLoader = LoadRun("READ");
1625 if (!runLoader) return kFALSE;
1627 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1628 Bool_t selrawdata = kFALSE;
1629 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1632 // Note the option -s. It is used in order to avoid
1633 // the generation of SOR/EOR events.
1634 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1635 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1636 FILE* pipe = gSystem->OpenPipe(command, "w");
1639 AliError(Form("Cannot execute command: %s",command));
1643 Int_t selEvents = 0;
1644 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1646 UInt_t detectorPattern = 0;
1647 runLoader->GetEvent(iEvent);
1648 if (!runLoader->LoadTrigger()) {
1649 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1650 detectorPattern = aCTP->GetClusterMask();
1651 // Check if the event was triggered by CTP
1653 if (aCTP->GetClassMask()) selEvents++;
1657 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1659 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1660 selrawdata = kFALSE;
1664 fprintf(pipe, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1668 // loop over detectors and DDLs
1669 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1670 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1672 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1673 Int_t ldcID = Int_t(ldc + 0.0001);
1674 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1676 char rawFileName[256];
1677 snprintf(rawFileName, 256, "raw%d/%s",
1678 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1680 // check existence and size of raw data file
1681 FILE* file = fopen(rawFileName, "rb");
1682 if (!file) continue;
1683 fseek(file, 0, SEEK_END);
1684 unsigned long size = ftell(file);
1686 if (!size) continue;
1688 if (ldcID != prevLDC) {
1689 fprintf(pipe, " LDC Id %d\n", ldcID);
1692 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1697 Int_t result = gSystem->ClosePipe(pipe);
1699 if (!(selrawdata && selEvents > 0)) {
1701 return (result == 0);
1704 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1706 snprintf(command, 256, "dateStream -c -s -D -o %s -# %d -C -run %d",
1707 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1708 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1710 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1712 // Get the trigger decision and cluster
1713 UInt_t detectorPattern = 0;
1715 runLoader->GetEvent(iEvent);
1716 if (!runLoader->LoadTrigger()) {
1717 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1718 if (aCTP->GetClassMask() == 0) continue;
1719 detectorPattern = aCTP->GetClusterMask();
1720 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1721 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1724 fprintf(pipe2, "GDC DetectorPattern %u Timestamp %ld\n", detectorPattern, runLoader->GetHeader()->GetTimeStamp());
1728 // loop over detectors and DDLs
1729 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1730 // Write only raw data from detectors that
1731 // are contained in the trigger cluster(s)
1732 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1734 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1736 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1737 Int_t ldcID = Int_t(ldc + 0.0001);
1738 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1740 char rawFileName[256];
1741 snprintf(rawFileName, 256, "raw%d/%s",
1742 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1744 // check existence and size of raw data file
1745 FILE* file = fopen(rawFileName, "rb");
1746 if (!file) continue;
1747 fseek(file, 0, SEEK_END);
1748 unsigned long size = ftell(file);
1750 if (!size) continue;
1752 if (ldcID != prevLDC) {
1753 fprintf(pipe2, " LDC Id %d\n", ldcID);
1756 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1761 Int_t result2 = gSystem->ClosePipe(pipe2);
1764 return ((result == 0) && (result2 == 0));
1767 //_____________________________________________________________________________
1768 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1769 const char* rootFileName)
1771 // convert a DATE file to a root file with the program "alimdc"
1774 const Int_t kDBSize = 2000000000;
1775 const Int_t kTagDBSize = 1000000000;
1776 const Bool_t kFilter = kFALSE;
1777 const Int_t kCompression = 1;
1779 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1781 AliError("the program alimdc was not found");
1782 if (fStopOnError) return kFALSE;
1787 AliInfo(Form("converting DATE file %s to root file %s",
1788 dateFileName, rootFileName));
1790 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1791 const char* tagDBFS = "/tmp/mdc1/tags";
1793 // User defined file system locations
1794 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1795 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1796 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1797 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1798 if (gSystem->Getenv("ALIMDC_TAGDB"))
1799 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1801 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1802 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1803 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1805 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1806 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1807 gSystem->Exec(Form("mkdir %s",tagDBFS));
1809 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1810 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1811 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1813 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1814 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1815 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1817 return (result == 0);
1821 //_____________________________________________________________________________
1822 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1824 // delete existing run loaders, open a new one and load gAlice
1826 delete AliRunLoader::Instance();
1827 AliRunLoader* runLoader =
1828 AliRunLoader::Open(fGAliceFileName.Data(),
1829 AliConfig::GetDefaultEventFolderName(), mode);
1831 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1834 runLoader->LoadgAlice();
1835 runLoader->LoadHeader();
1836 gAlice = runLoader->GetAliRun();
1838 AliError(Form("no gAlice object found in file %s",
1839 fGAliceFileName.Data()));
1845 //_____________________________________________________________________________
1846 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1848 // get or calculate the number of signal events per background event
1850 if (!fBkgrdFileNames) return 1;
1851 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1852 if (nBkgrdFiles == 0) return 1;
1854 // get the number of signal events
1856 AliRunLoader* runLoader =
1857 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1858 if (!runLoader) return 1;
1860 nEvents = runLoader->GetNumberOfEvents();
1865 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1866 // get the number of background events
1867 const char* fileName = ((TObjString*)
1868 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1869 AliRunLoader* runLoader =
1870 AliRunLoader::Open(fileName, "BKGRD");
1871 if (!runLoader) continue;
1872 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1875 // get or calculate the number of signal per background events
1876 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1877 if (nSignalPerBkgrd <= 0) {
1878 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1879 } else if (result && (result != nSignalPerBkgrd)) {
1880 AliInfo(Form("the number of signal events per background event "
1881 "will be changed from %d to %d for stream %d",
1882 nSignalPerBkgrd, result, iBkgrdFile+1));
1883 nSignalPerBkgrd = result;
1886 if (!result) result = nSignalPerBkgrd;
1887 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1888 AliWarning(Form("not enough background events (%d) for %d signal events "
1889 "using %d signal per background events for stream %d",
1890 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1897 //_____________________________________________________________________________
1898 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1900 // check whether detName is contained in detectors
1901 // if yes, it is removed from detectors
1903 // check if all detectors are selected
1904 if ((detectors.CompareTo("ALL") == 0) ||
1905 detectors.BeginsWith("ALL ") ||
1906 detectors.EndsWith(" ALL") ||
1907 detectors.Contains(" ALL ")) {
1912 // search for the given detector
1913 Bool_t result = kFALSE;
1914 if ((detectors.CompareTo(detName) == 0) ||
1915 detectors.BeginsWith(detName+" ") ||
1916 detectors.EndsWith(" "+detName) ||
1917 detectors.Contains(" "+detName+" ")) {
1918 detectors.ReplaceAll(detName, "");
1922 // clean up the detectors string
1923 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1924 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1925 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1930 //_____________________________________________________________________________
1931 Int_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName, Int_t N, Int_t nSkip)
1934 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1935 // These can be used for embedding of MC tracks into RAW data using the standard
1936 // merging procedure.
1938 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1941 AliError("no gAlice object. Restart aliroot and try again.");
1944 if (gAlice->Modules()->GetEntries() > 0) {
1945 AliError("gAlice was already run. Restart aliroot and try again.");
1949 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1953 gROOT->LoadMacro(fConfigFileName.Data());
1954 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1956 if(AliCDBManager::Instance()->GetRun() >= 0) {
1957 SetRunNumber(AliCDBManager::Instance()->GetRun());
1959 AliWarning("Run number not initialized!!");
1962 AliRunLoader::Instance()->CdGAFile();
1964 AliPDG::AddParticlesToPdgDataBase();
1966 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1968 gAlice->GetMCApp()->Init();
1970 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1971 gAlice->InitLoaders();
1972 AliRunLoader::Instance()->MakeTree("E");
1973 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1974 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1975 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1978 // Save stuff at the beginning of the file to avoid file corruption
1979 AliRunLoader::Instance()->CdGAFile();
1984 //AliCDBManager* man = AliCDBManager::Instance();
1985 //man->SetRun(0); // Should this come from rawdata header ?
1989 // Get the runloader
1990 AliRunLoader* runLoader = AliRunLoader::Instance();
1992 // Open esd file if available
1995 AliESDEvent* esd = 0;
1996 if (esdFileName && (strlen(esdFileName)>0)) {
1997 esdFile = TFile::Open(esdFileName);
1999 esd = new AliESDEvent();
2000 esdFile->GetObject("esdTree", treeESD);
2002 esd->ReadFromTree(treeESD);
2004 AliInfo(Form("Asking to skip first %d ESDs events",nSkip));
2013 // Create the RawReader
2014 TString fileName(rawDirectory);
2015 AliRawReader* rawReader = AliRawReader::Create(fileName.Data());
2016 if (!rawReader) return (kFALSE);
2018 // if (!fEquipIdMap.IsNull() && fRawReader)
2019 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
2021 // Get list of detectors
2022 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
2025 AliHeader* header = runLoader->GetHeader();
2029 if (!(rawReader->NextEvent())) break;
2030 runLoader->SetEventNumber(nev);
2031 runLoader->GetHeader()->Reset(rawReader->GetRunNumber(),
2033 runLoader->GetEvent(nev);
2034 AliInfo(Form("We are at event %d",nev));
2037 TString detStr = fMakeSDigits;
2038 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
2039 AliModule* det = (AliModule*) detArray->At(iDet);
2040 if (!det || !det->IsActive()) continue;
2041 if (IsSelected(det->GetName(), detStr)) {
2042 AliInfo(Form("Calling Raw2SDigits for %s", det->GetName()));
2043 det->Raw2SDigits(rawReader);
2050 // If ESD information available obtain reconstructed vertex and store in header.
2052 AliInfo(Form("Selected event %d correspond to event %d in raw and to %d in esd",nev,rawReader->GetEventIndex(),nSkip+rawReader->GetEventIndex()));
2053 treeESD->GetEvent(nSkip+rawReader->GetEventIndex());
2054 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
2055 Double_t position[3];
2056 esdVertex->GetXYZ(position);
2057 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
2060 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
2061 mcHeader->SetPrimaryVertex(mcV);
2062 header->Reset(0,nev);
2063 header->SetGenEventHeader(mcHeader);
2064 AliInfo(Form("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]));
2068 runLoader->TreeE()->Fill();
2069 AliInfo(Form("Finished event %d",nev));
2078 runLoader->CdGAFile();
2079 runLoader->WriteHeader("OVERWRITE");
2080 runLoader->WriteRunLoader();
2085 //_____________________________________________________________________________
2086 void AliSimulation::FinishRun()
2089 // Called at the end of the run.
2094 AliDebug(1, "Finish Lego");
2095 AliRunLoader::Instance()->CdGAFile();
2099 // Clean detector information
2100 TIter next(gAlice->Modules());
2101 AliModule *detector;
2102 while((detector = dynamic_cast<AliModule*>(next()))) {
2103 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
2104 detector->FinishRun();
2107 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
2108 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
2110 // Write AliRun info and all detectors parameters
2111 AliRunLoader::Instance()->CdGAFile();
2112 gAlice->Write(0,TObject::kOverwrite);//write AliRun
2113 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
2115 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
2116 AliRunLoader::Instance()->Synchronize();
2119 //_____________________________________________________________________________
2120 Int_t AliSimulation::GetDetIndex(const char* detector)
2122 // return the detector index corresponding to detector
2124 for (index = 0; index < fgkNDetectors ; index++) {
2125 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2131 //_____________________________________________________________________________
2132 Bool_t AliSimulation::CreateHLT()
2134 // Init the HLT simulation.
2135 // The function loads the library and creates the instance of AliHLTSimulation.
2136 // the main reason for the decoupled creation is to set the transient OCDB
2137 // objects before the OCDB is locked
2139 // load the library dynamically
2140 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
2142 // check for the library version
2143 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
2145 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
2148 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
2149 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
2152 // print compile info
2153 typedef void (*CompileInfo)( const char*& date, const char*& time);
2154 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
2156 const char* date="";
2157 const char* time="";
2158 (*fctInfo)(date, time);
2159 if (!date) date="unknown";
2160 if (!time) time="unknown";
2161 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
2163 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
2166 // create instance of the HLT simulation
2167 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
2168 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
2169 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
2173 TString specObjects;
2174 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
2175 if (specObjects.Length()>0) specObjects+=" ";
2176 specObjects+=fSpecCDBUri[i]->GetName();
2179 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
2180 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
2181 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2187 //_____________________________________________________________________________
2188 Bool_t AliSimulation::RunHLT()
2190 // Run the HLT simulation
2191 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2192 // Disabled if fRunHLT is empty, default vaule is "default".
2193 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2194 // The default simulation depends on the HLT component libraries and their
2195 // corresponding agents which define components and chains to run. See
2196 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2197 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2199 // The libraries to be loaded can be specified as an option.
2201 // AliSimulation sim;
2202 // sim.SetRunHLT("libAliHLTSample.so");
2204 // will only load <tt>libAliHLTSample.so</tt>
2206 // Other available options:
2207 // \li loglevel=<i>level</i> <br>
2208 // logging level for this processing
2210 // disable redirection of log messages to AliLog class
2211 // \li config=<i>macro</i>
2212 // configuration macro
2213 // \li chains=<i>configuration</i>
2214 // comma separated list of configurations to be run during simulation
2215 // \li rawfile=<i>file</i>
2216 // source for the RawReader to be created, the default is <i>./</i> if
2217 // raw data is simulated
2221 if (!fpHLT && !CreateHLT()) {
2224 AliHLTSimulation* pHLT=fpHLT;
2226 AliRunLoader* pRunLoader = LoadRun("READ");
2227 if (!pRunLoader) return kFALSE;
2229 // initialize CDB storage, run number, set CDB lock
2230 // thats for the case of running HLT simulation without all the other steps
2231 // multiple calls are handled by the function, so we can just call
2233 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2236 // init the HLT simulation
2238 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2239 TString detStr = fWriteRawData;
2240 if (!IsSelected("HLT", detStr)) {
2241 options+=" writerawfiles=";
2243 options+=" writerawfiles=HLT";
2246 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2247 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2248 // in order to get detector data. By default, RawReaderFile is used to read
2249 // the already simulated ddl files. Date and Root files from the raw data
2250 // are generated after the HLT simulation.
2251 options+=" rawfile=./";
2254 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2255 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2256 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2258 // run the HLT simulation
2259 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2260 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2261 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2265 // delete the instance
2266 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2267 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2268 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2272 return iResult>=0?kTRUE:kFALSE;
2275 //_____________________________________________________________________________
2276 Bool_t AliSimulation::RunQA()
2278 // run the QA on summable hits, digits or digits
2280 //if(!gAlice) return kFALSE;
2281 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2283 TString detectorsw("") ;
2285 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2286 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2287 if ( detectorsw.IsNull() )
2292 //_____________________________________________________________________________
2293 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2295 // Allows to run QA for a selected set of detectors
2296 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2297 // all selected detectors run the same selected tasks
2299 if (!detAndAction.Contains(":")) {
2300 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2304 Int_t colon = detAndAction.Index(":") ;
2305 fQADetectors = detAndAction(0, colon) ;
2306 if (fQADetectors.Contains("ALL") ){
2307 TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2308 Int_t minus = fQADetectors.Last('-') ;
2309 TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2310 TString toRemove("") ;
2311 while (minus >= 0) {
2312 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
2313 toRemove = toRemove.Strip() ;
2314 toKeep.ReplaceAll(toRemove, "") ;
2315 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
2316 minus = fQADetectors.Last('-') ;
2318 fQADetectors = toKeep ;
2320 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2321 if (fQATasks.Contains("ALL") ) {
2322 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2324 fQATasks.ToUpper() ;
2326 if ( fQATasks.Contains("HIT") )
2327 tempo = Form("%d ", AliQAv1::kHITS) ;
2328 if ( fQATasks.Contains("SDIGIT") )
2329 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2330 if ( fQATasks.Contains("DIGIT") )
2331 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2333 if (fQATasks.IsNull()) {
2334 AliInfo("No QA requested\n") ;
2339 TString tempo(fQATasks) ;
2340 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2341 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2342 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2343 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2345 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2346 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2347 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2348 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2353 //_____________________________________________________________________________
2354 void AliSimulation::ProcessEnvironmentVars()
2356 // Extract run number and random generator seed from env variables
2358 AliInfo("Processing environment variables");
2360 // Random Number seed
2362 // first check that seed is not already set
2364 if (gSystem->Getenv("CONFIG_SEED")) {
2365 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2368 if (gSystem->Getenv("CONFIG_SEED")) {
2369 AliInfo(Form("Seed for random number generation already set (%d)"
2370 ": CONFIG_SEED variable ignored!", fSeed));
2374 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2378 // first check that run number is not already set
2380 if (gSystem->Getenv("DC_RUN")) {
2381 fRun = atoi(gSystem->Getenv("DC_RUN"));
2384 if (gSystem->Getenv("DC_RUN")) {
2385 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2389 AliInfo(Form("Run number = %d", fRun));
2392 //---------------------------------------------------------------------
2393 void AliSimulation::WriteGRPEntry()
2395 // Get the necessary information from galice (generator, trigger etc) and
2396 // write a GRP entry corresponding to the settings in the Config.C used
2397 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2400 AliInfo("Writing global run parameters entry into the OCDB");
2402 AliGRPObject* grpObj = new AliGRPObject();
2404 grpObj->SetRunType("PHYSICS");
2405 grpObj->SetTimeStart(fTimeStart);
2406 grpObj->SetTimeEnd(fTimeEnd);
2407 grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2409 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2415 gen->GetProjectile(projectile,a,z);
2417 gen->GetTarget(target,a,z);
2418 TString beamType = projectile + "-" + target;
2419 beamType.ReplaceAll(" ","");
2420 if (!beamType.CompareTo("-")) {
2421 grpObj->SetBeamType("UNKNOWN");
2422 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2425 grpObj->SetBeamType(beamType);
2427 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 * a / z);
2429 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 );
2431 // Heavy ion run, the event specie is set to kHighMult
2432 fEventSpecie = AliRecoParam::kHighMult;
2433 if ((strcmp(beamType,"p-p") == 0) ||
2434 (strcmp(beamType,"p-") == 0) ||
2435 (strcmp(beamType,"-p") == 0) ||
2436 (strcmp(beamType,"P-P") == 0) ||
2437 (strcmp(beamType,"P-") == 0) ||
2438 (strcmp(beamType,"-P") == 0)) {
2439 // Proton run, the event specie is set to kLowMult
2440 fEventSpecie = AliRecoParam::kLowMult;
2444 AliWarning("Unknown beam type and energy! Setting energy to 0");
2445 grpObj->SetBeamEnergy(0);
2446 grpObj->SetBeamType("UNKNOWN");
2449 UInt_t detectorPattern = 0;
2451 TObjArray *detArray = gAlice->Detectors();
2452 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2453 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2454 AliDebug(1, Form("Detector #%d found: %s", iDet, AliDAQ::OfflineModuleName(iDet)));
2455 detectorPattern |= (1 << iDet);
2460 if (!fTriggerConfig.IsNull())
2461 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2464 if (!fRunHLT.IsNull())
2465 detectorPattern |= (1 << AliDAQ::kHLTId);
2467 grpObj->SetNumberOfDetectors((Char_t)nDets);
2468 grpObj->SetDetectorMask((Int_t)detectorPattern);
2469 grpObj->SetLHCPeriod("LHC08c");
2470 grpObj->SetLHCState("STABLE_BEAMS");
2472 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2473 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2475 Float_t factorSol = field ? field->GetFactorSol() : 0;
2476 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2477 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2479 Float_t factorDip = field ? field->GetFactorDip() : 0;
2480 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2482 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2483 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2484 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2485 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2486 if (field) grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2487 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2489 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2491 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2493 // Now store the entry in OCDB
2494 AliCDBManager* man = AliCDBManager::Instance();
2496 man->SetLock(0, fKey);
2498 AliCDBStorage* sto = man->GetStorage(fGRPWriteLocation.Data());
2501 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2502 AliCDBMetaData *metadata= new AliCDBMetaData();
2504 metadata->SetResponsible("alice-off@cern.ch");
2505 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2507 sto->Put(grpObj,id,metadata);
2508 man->SetLock(1, fKey);
2511 //_____________________________________________________________________________
2512 time_t AliSimulation::GenerateTimeStamp() const
2514 // Generate event time-stamp according to
2515 // SOR/EOR time from GRP
2516 if (fUseTimeStampFromCDB)
2517 return fTimeStart + gRandom->Integer(fTimeEnd-fTimeStart);