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 "AliAlignObj.h"
119 #include "AliCDBEntry.h"
120 #include "AliCDBManager.h"
121 #include "AliCDBStorage.h"
122 #include "AliCTPRawData.h"
123 #include "AliCentralTrigger.h"
124 #include "AliCentralTrigger.h"
125 #include "AliCodeTimer.h"
127 #include "AliDigitizer.h"
128 #include "AliESDEvent.h"
129 #include "AliGRPObject.h"
130 #include "AliGenEventHeader.h"
131 #include "AliGenerator.h"
132 #include "AliGeomManager.h"
133 #include "AliHLTSimulation.h"
134 #include "AliHeader.h"
136 #include "AliLegoGenerator.h"
140 #include "AliModule.h"
142 #include "AliRawReaderDate.h"
143 #include "AliRawReaderFile.h"
144 #include "AliRawReaderRoot.h"
146 #include "AliRunDigitizer.h"
147 #include "AliRunLoader.h"
148 #include "AliSimulation.h"
149 #include "AliSysInfo.h"
150 #include "AliVertexGenFile.h"
152 ClassImp(AliSimulation)
154 AliSimulation *AliSimulation::fgInstance = 0;
155 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
157 //_____________________________________________________________________________
158 AliSimulation::AliSimulation(const char* configFileName,
159 const char* name, const char* title) :
162 fRunGeneration(kTRUE),
163 fRunSimulation(kTRUE),
164 fLoadAlignFromCDB(kTRUE),
165 fLoadAlObjsListOfDets("ALL"),
169 fMakeDigitsFromHits(""),
171 fRawDataFileName(""),
172 fDeleteIntermediateFiles(kFALSE),
173 fWriteSelRawData(kFALSE),
174 fStopOnError(kFALSE),
176 fConfigFileName(configFileName),
177 fGAliceFileName("galice.root"),
179 fBkgrdFileNames(NULL),
180 fAlignObjArray(NULL),
181 fUseBkgrdVertex(kTRUE),
182 fRegionOfInterest(kFALSE),
188 fInitCDBCalled(kFALSE),
189 fInitRunNumberCalled(kFALSE),
190 fSetRunNumberFromDataCalled(kFALSE),
191 fEmbeddingFlag(kFALSE),
196 fEventSpecie(AliRecoParam::kDefault),
197 fWriteQAExpertData(kTRUE),
200 fWriteGRPEntry(kTRUE)
202 // create simulation object with default parameters
204 SetGAliceFile("galice.root");
207 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
208 qam->SetActiveDetectors(fQADetectors) ;
209 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
210 qam->SetTasks(fQATasks) ;
213 //_____________________________________________________________________________
214 AliSimulation::~AliSimulation()
218 fEventsPerFile.Delete();
219 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
220 // delete fAlignObjArray; fAlignObjArray=0;
222 if (fBkgrdFileNames) {
223 fBkgrdFileNames->Delete();
224 delete fBkgrdFileNames;
227 fSpecCDBUri.Delete();
228 if (fgInstance==this) fgInstance = 0;
230 AliQAManager::QAManager()->ShowQA() ;
231 AliQAManager::Destroy() ;
232 AliCodeTimer::Instance()->Print();
236 //_____________________________________________________________________________
237 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
239 // set the number of events for one run
244 //_____________________________________________________________________________
245 void AliSimulation::InitQA()
247 // activate a default CDB storage
248 // First check if we have any CDB storage set, because it is used
249 // to retrieve the calibration and alignment constants
251 if (fInitCDBCalled) return;
252 fInitCDBCalled = kTRUE;
254 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
255 qam->SetActiveDetectors(fQADetectors) ;
256 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
257 qam->SetTasks(fQATasks) ;
258 if (fWriteQAExpertData)
259 qam->SetWriteExpert() ;
261 if (qam->IsDefaultStorageSet()) {
262 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
263 AliWarning("Default QA reference storage has been already set !");
264 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
265 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
266 fQARefUri = qam->GetDefaultStorage()->GetURI();
268 if (fQARefUri.Length() > 0) {
269 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
270 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
271 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
273 fQARefUri="local://$ALICE_ROOT/QARef";
274 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
275 AliWarning("Default QA reference storage not yet set !!!!");
276 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
277 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
279 qam->SetDefaultStorage(fQARefUri);
283 //_____________________________________________________________________________
284 void AliSimulation::InitCDB()
286 // activate a default CDB storage
287 // First check if we have any CDB storage set, because it is used
288 // to retrieve the calibration and alignment constants
290 if (fInitCDBCalled) return;
291 fInitCDBCalled = kTRUE;
293 AliCDBManager* man = AliCDBManager::Instance();
294 if (man->IsDefaultStorageSet())
296 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
297 AliWarning("Default CDB storage has been already set !");
298 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
299 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
300 fCDBUri = man->GetDefaultStorage()->GetURI();
303 if (fCDBUri.Length() > 0)
305 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
306 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
307 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
309 fCDBUri="local://$ALICE_ROOT/OCDB";
310 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
311 AliWarning("Default CDB storage not yet set !!!!");
312 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
313 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
316 man->SetDefaultStorage(fCDBUri);
319 // Now activate the detector specific CDB storage locations
320 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
321 TObject* obj = fSpecCDBUri[i];
323 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
324 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
325 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
326 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
331 //_____________________________________________________________________________
332 void AliSimulation::InitRunNumber(){
333 // check run number. If not set, set it to 0 !!!!
335 if (fInitRunNumberCalled) return;
336 fInitRunNumberCalled = kTRUE;
338 AliCDBManager* man = AliCDBManager::Instance();
339 if (man->GetRun() >= 0)
341 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
342 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
346 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
349 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
358 //_____________________________________________________________________________
359 void AliSimulation::SetCDBLock() {
360 // Set CDB lock: from now on it is forbidden to reset the run number
361 // or the default storage or to activate any further storage!
363 AliCDBManager::Instance()->SetLock(1);
366 //_____________________________________________________________________________
367 void AliSimulation::SetDefaultStorage(const char* uri) {
368 // Store the desired default CDB storage location
369 // Activate it later within the Run() method
375 //_____________________________________________________________________________
376 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
377 // Store the desired default CDB storage location
378 // Activate it later within the Run() method
381 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
384 //_____________________________________________________________________________
385 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
386 // Store a detector-specific CDB storage location
387 // Activate it later within the Run() method
389 AliCDBPath aPath(calibType);
390 if(!aPath.IsValid()){
391 AliError(Form("Not a valid path: %s", calibType));
395 TObject* obj = fSpecCDBUri.FindObject(calibType);
396 if (obj) fSpecCDBUri.Remove(obj);
397 fSpecCDBUri.Add(new TNamed(calibType, uri));
401 //_____________________________________________________________________________
402 void AliSimulation::SetRunNumber(Int_t run)
405 // Activate it later within the Run() method
410 //_____________________________________________________________________________
411 void AliSimulation::SetSeed(Int_t seed)
414 // Activate it later within the Run() method
419 //_____________________________________________________________________________
420 Bool_t AliSimulation::SetRunNumberFromData()
422 // Set the CDB manager run number
423 // The run number is retrieved from gAlice
425 if (fSetRunNumberFromDataCalled) return kTRUE;
426 fSetRunNumberFromDataCalled = kTRUE;
428 AliCDBManager* man = AliCDBManager::Instance();
429 Int_t runData = -1, runCDB = -1;
431 AliRunLoader* runLoader = LoadRun("READ");
432 if (!runLoader) return kFALSE;
434 runData = runLoader->GetHeader()->GetRun();
438 runCDB = man->GetRun();
440 if (runCDB != runData) {
441 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
442 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
443 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
444 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
449 man->SetRun(runData);
452 if(man->GetRun() < 0) {
453 AliError("Run number not properly initalized!");
462 //_____________________________________________________________________________
463 void AliSimulation::SetConfigFile(const char* fileName)
465 // set the name of the config file
467 fConfigFileName = fileName;
470 //_____________________________________________________________________________
471 void AliSimulation::SetGAliceFile(const char* fileName)
473 // set the name of the galice file
474 // the path is converted to an absolute one if it is relative
476 fGAliceFileName = fileName;
477 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
478 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
480 fGAliceFileName = absFileName;
481 delete[] absFileName;
484 AliDebug(2, Form("galice file name set to %s", fileName));
487 //_____________________________________________________________________________
488 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
491 // set the number of events per file for the given detector and data type
492 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
494 TNamed* obj = new TNamed(detector, type);
495 obj->SetUniqueID(nEvents);
496 fEventsPerFile.Add(obj);
499 //_____________________________________________________________________________
500 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
502 // Read the alignment objects from CDB.
503 // Each detector is supposed to have the
504 // alignment objects in DET/Align/Data CDB path.
505 // All the detector objects are then collected,
506 // sorted by geometry level (starting from ALIC) and
507 // then applied to the TGeo geometry.
508 // Finally an overlaps check is performed.
510 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
511 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
515 // initialize CDB storage, run number, set CDB lock
517 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
520 Bool_t delRunLoader = kFALSE;
522 runLoader = LoadRun("READ");
523 if (!runLoader) return kFALSE;
524 delRunLoader = kTRUE;
527 // Export ideal geometry
528 if(!gAlice->IsRootGeometry()) AliGeomManager::GetGeometry()->Export("geometry.root");
530 // Load alignment data from CDB and apply to geometry through AliGeomManager
531 if(fLoadAlignFromCDB){
533 TString detStr = fLoadAlObjsListOfDets;
534 TString loadAlObjsListOfDets = "";
536 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
537 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
538 AliModule* det = (AliModule*) detArray->At(iDet);
539 if (!det || !det->IsActive()) continue;
540 if (IsSelected(det->GetName(), detStr)) {
541 //add det to list of dets to be aligned from CDB
542 loadAlObjsListOfDets += det->GetName();
543 loadAlObjsListOfDets += " ";
545 } // end loop over detectors
546 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
547 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
549 // Check if the array with alignment objects was
550 // provided by the user. If yes, apply the objects
551 // to the present TGeo geometry
552 if (fAlignObjArray) {
553 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
554 AliError("The misalignment of one or more volumes failed!"
555 "Compare the list of simulated detectors and the list of detector alignment data!");
556 if (delRunLoader) delete runLoader;
562 // Update the internal geometry of modules (ITS needs it)
563 TString detStr = fLoadAlObjsListOfDets;
564 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
565 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
567 AliModule* det = (AliModule*) detArray->At(iDet);
568 if (!det || !det->IsActive()) continue;
569 if (IsSelected(det->GetName(), detStr)) {
570 det->UpdateInternalGeometry();
572 } // end loop over detectors
575 if (delRunLoader) delete runLoader;
580 //_____________________________________________________________________________
581 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
583 // add a file with background events for merging
585 TObjString* fileNameStr = new TObjString(fileName);
586 fileNameStr->SetUniqueID(nSignalPerBkgrd);
587 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
588 fBkgrdFileNames->Add(fileNameStr);
591 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
593 // add a file with background events for embeddin
594 MergeWith(fileName, nSignalPerBkgrd);
595 fEmbeddingFlag = kTRUE;
598 //_____________________________________________________________________________
599 Bool_t AliSimulation::Run(Int_t nEvents)
601 // run the generation, simulation and digitization
604 AliCodeTimerAuto("",0)
605 AliSysInfo::AddStamp("Start_Run");
607 // Load run number and seed from environmental vars
608 ProcessEnvironmentVars();
609 AliSysInfo::AddStamp("ProcessEnvironmentVars");
611 gRandom->SetSeed(fSeed);
613 if (nEvents > 0) fNEvents = nEvents;
615 // create and setup the HLT instance
616 if (!fRunHLT.IsNull() && !CreateHLT()) {
617 if (fStopOnError) return kFALSE;
622 // generation and simulation -> hits
623 if (fRunGeneration) {
624 if (!RunSimulation()) if (fStopOnError) return kFALSE;
626 AliSysInfo::AddStamp("RunSimulation");
628 // initialize CDB storage from external environment
629 // (either CDB manager or AliSimulation setters),
630 // if not already done in RunSimulation()
632 AliSysInfo::AddStamp("InitCDB");
634 // Set run number in CDBManager from data
635 // From this point on the run number must be always loaded from data!
636 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
638 // Set CDB lock: from now on it is forbidden to reset the run number
639 // or the default storage or to activate any further storage!
642 // If RunSimulation was not called, load the geometry and misalign it
643 if (!AliGeomManager::GetGeometry()) {
644 // Initialize the geometry manager
645 AliGeomManager::LoadGeometry("geometry.root");
646 AliSysInfo::AddStamp("GetGeometry");
649 // // Check that the consistency of symbolic names for the activated subdetectors
650 // // in the geometry loaded by AliGeomManager
651 // AliRunLoader* runLoader = LoadRun("READ");
652 // if (!runLoader) return kFALSE;
654 // TString detsToBeChecked = "";
655 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
656 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
657 // AliModule* det = (AliModule*) detArray->At(iDet);
658 // if (!det || !det->IsActive()) continue;
659 // detsToBeChecked += det->GetName();
660 // detsToBeChecked += " ";
661 // } // end loop over detectors
662 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
663 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
664 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
666 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
668 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
670 AliSysInfo::AddStamp("MissalignGeometry");
673 // hits -> summable digits
674 AliSysInfo::AddStamp("Start_sdigitization");
675 if (!fMakeSDigits.IsNull()) {
676 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
679 AliSysInfo::AddStamp("Stop_sdigitization");
681 AliSysInfo::AddStamp("Start_digitization");
682 // summable digits -> digits
683 if (!fMakeDigits.IsNull()) {
684 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
685 if (fStopOnError) return kFALSE;
688 AliSysInfo::AddStamp("Stop_digitization");
693 if (!fMakeDigitsFromHits.IsNull()) {
694 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
695 AliWarning(Form("Merging and direct creation of digits from hits "
696 "was selected for some detectors. "
697 "No merging will be done for the following detectors: %s",
698 fMakeDigitsFromHits.Data()));
700 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
701 if (fStopOnError) return kFALSE;
705 AliSysInfo::AddStamp("Hits2Digits");
709 if (!RunTrigger(fMakeTrigger,fMakeDigits)) {
710 if (fStopOnError) return kFALSE;
713 AliSysInfo::AddStamp("RunTrigger");
716 // digits -> raw data
717 if (!fWriteRawData.IsNull()) {
718 if (!WriteRawData(fWriteRawData, fRawDataFileName,
719 fDeleteIntermediateFiles,fWriteSelRawData)) {
720 if (fStopOnError) return kFALSE;
724 AliSysInfo::AddStamp("WriteRaw");
726 // run HLT simulation on simulated digit data if raw data is not
727 // simulated, otherwise its called as part of WriteRawData
728 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
730 if (fStopOnError) return kFALSE;
734 AliSysInfo::AddStamp("RunHLT");
738 Bool_t rv = RunQA() ;
744 AliSysInfo::AddStamp("RunQA");
746 // Cleanup of CDB manager: cache and active storages!
747 AliCDBManager::Instance()->ClearCache();
752 //_______________________________________________________________________
753 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
754 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
755 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
758 // Generates lego plots of:
759 // - radiation length map phi vs theta
760 // - radiation length map phi vs eta
761 // - interaction length map
762 // - g/cm2 length map
764 // ntheta bins in theta, eta
765 // themin minimum angle in theta (degrees)
766 // themax maximum angle in theta (degrees)
768 // phimin minimum angle in phi (degrees)
769 // phimax maximum angle in phi (degrees)
770 // rmin minimum radius
771 // rmax maximum radius
774 // The number of events generated = ntheta*nphi
775 // run input parameters in macro setup (default="Config.C")
777 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
780 <img src="picts/AliRunLego1.gif">
785 <img src="picts/AliRunLego2.gif">
790 <img src="picts/AliRunLego3.gif">
795 // run the generation and simulation
797 AliCodeTimerAuto("",0)
799 // initialize CDB storage and run number from external environment
800 // (either CDB manager or AliSimulation setters)
806 AliError("no gAlice object. Restart aliroot and try again.");
809 if (gAlice->Modules()->GetEntries() > 0) {
810 AliError("gAlice was already run. Restart aliroot and try again.");
814 AliInfo(Form("initializing gAlice with config file %s",
815 fConfigFileName.Data()));
818 if (nev == -1) nev = nc1 * nc2;
820 // check if initialisation has been done
821 // If runloader has been initialized, set the number of events per file to nc1 * nc2
824 if (!gener) gener = new AliLegoGenerator();
826 // Configure Generator
828 gener->SetRadiusRange(rmin, rmax);
829 gener->SetZMax(zmax);
830 gener->SetCoor1Range(nc1, c1min, c1max);
831 gener->SetCoor2Range(nc2, c2min, c2max);
835 fLego = new AliLego("lego",gener);
837 //__________________________________________________________________________
841 gROOT->LoadMacro(setup);
842 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
844 if(AliCDBManager::Instance()->GetRun() >= 0) {
845 SetRunNumber(AliCDBManager::Instance()->GetRun());
847 AliWarning("Run number not initialized!!");
850 AliRunLoader::Instance()->CdGAFile();
852 AliPDG::AddParticlesToPdgDataBase();
854 gAlice->GetMCApp()->Init();
856 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
859 //Must be here because some MCs (G4) adds detectors here and not in Config.C
860 gAlice->InitLoaders();
861 AliRunLoader::Instance()->MakeTree("E");
864 // Save stuff at the beginning of the file to avoid file corruption
865 AliRunLoader::Instance()->CdGAFile();
868 //Save current generator
869 AliGenerator *gen=gAlice->GetMCApp()->Generator();
870 gAlice->GetMCApp()->ResetGenerator(gener);
871 //Prepare MC for Lego Run
877 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
878 gMC->ProcessRun(nev);
880 // End of this run, close files
882 // Restore current generator
883 gAlice->GetMCApp()->ResetGenerator(gen);
884 // Delete Lego Object
890 //_____________________________________________________________________________
891 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
895 AliCodeTimerAuto("",0)
897 // initialize CDB storage from external environment
898 // (either CDB manager or AliSimulation setters),
899 // if not already done in RunSimulation()
902 // Set run number in CDBManager from data
903 // From this point on the run number must be always loaded from data!
904 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
906 // Set CDB lock: from now on it is forbidden to reset the run number
907 // or the default storage or to activate any further storage!
910 AliRunLoader* runLoader = LoadRun("READ");
911 if (!runLoader) return kFALSE;
912 TString trconfiguration = config;
914 if (trconfiguration.IsNull()) {
915 if (strcmp(gAlice->GetTriggerDescriptor(),"")) {
916 trconfiguration = gAlice->GetTriggerDescriptor();
919 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
922 runLoader->MakeTree( "GG" );
923 AliCentralTrigger* aCTP = runLoader->GetTrigger();
924 // Load Configuration
925 if (!aCTP->LoadConfiguration( trconfiguration ))
929 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
941 //_____________________________________________________________________________
942 Bool_t AliSimulation::WriteTriggerRawData()
944 // Writes the CTP (trigger) DDL raw data
945 // Details of the format are given in the
946 // trigger TDR - pages 134 and 135.
947 AliCTPRawData writer;
953 //_____________________________________________________________________________
954 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
956 // run the generation and simulation
958 AliCodeTimerAuto("",0)
960 // initialize CDB storage and run number from external environment
961 // (either CDB manager or AliSimulation setters)
962 AliSysInfo::AddStamp("RunSimulation_Begin");
964 AliSysInfo::AddStamp("RunSimulation_InitCDB");
967 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
970 AliError("no gAlice object. Restart aliroot and try again.");
973 if (gAlice->Modules()->GetEntries() > 0) {
974 AliError("gAlice was already run. Restart aliroot and try again.");
978 AliInfo(Form("initializing gAlice with config file %s",
979 fConfigFileName.Data()));
982 // Initialize ALICE Simulation run
987 gROOT->LoadMacro(fConfigFileName.Data());
988 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
989 AliSysInfo::AddStamp("RunSimulation_Config");
991 if(AliCDBManager::Instance()->GetRun() >= 0) {
992 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
993 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
995 AliWarning("Run number not initialized!!");
998 AliRunLoader::Instance()->CdGAFile();
1000 AliPDG::AddParticlesToPdgDataBase();
1002 gAlice->GetMCApp()->Init();
1003 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1005 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1006 AliSysInfo::AddStamp("RunSimulation_GetField");
1008 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1009 gAlice->InitLoaders();
1010 AliRunLoader::Instance()->MakeTree("E");
1011 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1012 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1013 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1015 // Save stuff at the beginning of the file to avoid file corruption
1016 AliRunLoader::Instance()->CdGAFile();
1018 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1019 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1020 //___________________________________________________________________________________________
1022 // Get the trigger descriptor string
1023 // Either from AliSimulation or from
1025 if (fMakeTrigger.IsNull()) {
1026 if (strcmp(gAlice->GetTriggerDescriptor(),""))
1027 fMakeTrigger = gAlice->GetTriggerDescriptor();
1030 gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
1031 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1033 // Set run number in CDBManager
1034 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1036 AliRunLoader* runLoader = AliRunLoader::Instance();
1038 AliError(Form("gAlice has no run loader object. "
1039 "Check your config file: %s", fConfigFileName.Data()));
1042 SetGAliceFile(runLoader->GetFileName());
1044 // Misalign geometry
1045 #if ROOT_VERSION_CODE < 331527
1046 AliGeomManager::SetGeometry(gGeoManager);
1048 // Check that the consistency of symbolic names for the activated subdetectors
1049 // in the geometry loaded by AliGeomManager
1050 TString detsToBeChecked = "";
1051 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1052 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1053 AliModule* det = (AliModule*) detArray->At(iDet);
1054 if (!det || !det->IsActive()) continue;
1055 detsToBeChecked += det->GetName();
1056 detsToBeChecked += " ";
1057 } // end loop over detectors
1058 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1059 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1060 MisalignGeometry(runLoader);
1061 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1064 // AliRunLoader* runLoader = AliRunLoader::Instance();
1065 // if (!runLoader) {
1066 // AliError(Form("gAlice has no run loader object. "
1067 // "Check your config file: %s", fConfigFileName.Data()));
1070 // SetGAliceFile(runLoader->GetFileName());
1072 if (!gAlice->GetMCApp()->Generator()) {
1073 AliError(Form("gAlice has no generator object. "
1074 "Check your config file: %s", fConfigFileName.Data()));
1078 // Write GRP entry corresponding to the setting found in Cofig.C
1081 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1083 if (nEvents <= 0) nEvents = fNEvents;
1085 // get vertex from background file in case of merging
1086 if (fUseBkgrdVertex &&
1087 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1088 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1089 const char* fileName = ((TObjString*)
1090 (fBkgrdFileNames->At(0)))->GetName();
1091 AliInfo(Form("The vertex will be taken from the background "
1092 "file %s with nSignalPerBackground = %d",
1093 fileName, signalPerBkgrd));
1094 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1095 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1098 if (!fRunSimulation) {
1099 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1102 // set the number of events per file for given detectors and data types
1103 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1104 if (!fEventsPerFile[i]) continue;
1105 const char* detName = fEventsPerFile[i]->GetName();
1106 const char* typeName = fEventsPerFile[i]->GetTitle();
1107 TString loaderName(detName);
1108 loaderName += "Loader";
1109 AliLoader* loader = runLoader->GetLoader(loaderName);
1111 AliError(Form("RunSimulation", "no loader for %s found\n"
1112 "Number of events per file not set for %s %s",
1113 detName, typeName, detName));
1116 AliDataLoader* dataLoader =
1117 loader->GetDataLoader(typeName);
1119 AliError(Form("no data loader for %s found\n"
1120 "Number of events per file not set for %s %s",
1121 typeName, detName, typeName));
1124 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1125 AliDebug(1, Form("number of events per file set to %d for %s %s",
1126 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1129 AliInfo("running gAlice");
1130 AliSysInfo::AddStamp("Start_ProcessRun");
1132 // Create the Root Tree with one branch per detector
1133 //Hits moved to begin event -> now we are crating separate tree for each event
1135 gMC->ProcessRun(nEvents);
1137 // End of this run, close files
1138 if(nEvents>0) FinishRun();
1140 AliSysInfo::AddStamp("Stop_ProcessRun");
1146 //_____________________________________________________________________________
1147 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1149 // run the digitization and produce summable digits
1150 static Int_t eventNr=0;
1151 AliCodeTimerAuto("",0) ;
1153 // initialize CDB storage, run number, set CDB lock
1155 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1158 AliRunLoader* runLoader = LoadRun();
1159 if (!runLoader) return kFALSE;
1161 TString detStr = detectors;
1162 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1163 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1164 AliModule* det = (AliModule*) detArray->At(iDet);
1165 if (!det || !det->IsActive()) continue;
1166 if (IsSelected(det->GetName(), detStr)) {
1167 AliInfo(Form("creating summable digits for %s", det->GetName()));
1168 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1169 det->Hits2SDigits();
1170 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1171 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1175 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1176 AliError(Form("the following detectors were not found: %s",
1178 if (fStopOnError) return kFALSE;
1187 //_____________________________________________________________________________
1188 Bool_t AliSimulation::RunDigitization(const char* detectors,
1189 const char* excludeDetectors)
1191 // run the digitization and produce digits from sdigits
1193 AliCodeTimerAuto("",0)
1195 // initialize CDB storage, run number, set CDB lock
1197 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1200 delete AliRunLoader::Instance();
1205 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1206 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1207 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
1208 // manager->SetEmbeddingFlag(fEmbeddingFlag);
1209 manager->SetInputStream(0, fGAliceFileName.Data());
1210 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1211 const char* fileName = ((TObjString*)
1212 (fBkgrdFileNames->At(iStream-1)))->GetName();
1213 manager->SetInputStream(iStream, fileName);
1216 TString detStr = detectors;
1217 TString detExcl = excludeDetectors;
1218 manager->GetInputStream(0)->ImportgAlice();
1219 AliRunLoader* runLoader =
1220 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
1221 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1222 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1223 AliModule* det = (AliModule*) detArray->At(iDet);
1224 if (!det || !det->IsActive()) continue;
1225 if (IsSelected(det->GetName(), detStr) &&
1226 !IsSelected(det->GetName(), detExcl)) {
1227 AliDigitizer* digitizer = det->CreateDigitizer(manager);
1230 AliError(Form("no digitizer for %s", det->GetName()));
1231 if (fStopOnError) return kFALSE;
1233 digitizer->SetRegionOfInterest(fRegionOfInterest);
1238 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1239 AliError(Form("the following detectors were not found: %s",
1241 if (fStopOnError) return kFALSE;
1244 if (!manager->GetListOfTasks()->IsEmpty()) {
1245 AliInfo("executing digitization");
1254 //_____________________________________________________________________________
1255 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1257 // run the digitization and produce digits from hits
1259 AliCodeTimerAuto("",0)
1261 // initialize CDB storage, run number, set CDB lock
1263 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1266 AliRunLoader* runLoader = LoadRun("READ");
1267 if (!runLoader) return kFALSE;
1269 TString detStr = detectors;
1270 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1271 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1272 AliModule* det = (AliModule*) detArray->At(iDet);
1273 if (!det || !det->IsActive()) continue;
1274 if (IsSelected(det->GetName(), detStr)) {
1275 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1280 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1281 AliError(Form("the following detectors were not found: %s",
1283 if (fStopOnError) return kFALSE;
1289 //_____________________________________________________________________________
1290 Bool_t AliSimulation::WriteRawData(const char* detectors,
1291 const char* fileName,
1292 Bool_t deleteIntermediateFiles,
1295 // convert the digits to raw data
1296 // First DDL raw data files for the given detectors are created.
1297 // If a file name is given, the DDL files are then converted to a DATE file.
1298 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1300 // If the file name has the extension ".root", the DATE file is converted
1302 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1303 // 'selrawdata' flag can be used to enable writing of detectors raw data
1304 // accoring to the trigger cluster.
1306 AliCodeTimerAuto("",0)
1307 AliSysInfo::AddStamp("WriteRawData_Start");
1309 TString detStr = detectors;
1310 if (!WriteRawFiles(detStr.Data())) {
1311 if (fStopOnError) return kFALSE;
1313 AliSysInfo::AddStamp("WriteRawFiles");
1315 // run HLT simulation on simulated DDL raw files
1316 // and produce HLT ddl raw files to be included in date/root file
1317 // bugfix 2009-06-26: the decision whether to write HLT raw data
1318 // is taken in RunHLT. Here HLT always needs to be run in order to
1319 // create HLT digits, unless its switched off. This is due to the
1320 // special placement of the HLT between the generation of DDL files
1321 // and conversion to DATE/Root file.
1322 detStr.ReplaceAll("HLT", "");
1323 if (!fRunHLT.IsNull()) {
1325 if (fStopOnError) return kFALSE;
1328 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1330 TString dateFileName(fileName);
1331 if (!dateFileName.IsNull()) {
1332 Bool_t rootOutput = dateFileName.EndsWith(".root");
1333 if (rootOutput) dateFileName += ".date";
1334 TString selDateFileName;
1336 selDateFileName = "selected.";
1337 selDateFileName+= dateFileName;
1339 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1340 if (fStopOnError) return kFALSE;
1342 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1343 if (deleteIntermediateFiles) {
1344 AliRunLoader* runLoader = LoadRun("READ");
1345 if (runLoader) for (Int_t iEvent = 0;
1346 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1348 sprintf(command, "rm -r raw%d", iEvent);
1349 gSystem->Exec(command);
1354 if (!ConvertDateToRoot(dateFileName, fileName)) {
1355 if (fStopOnError) return kFALSE;
1357 AliSysInfo::AddStamp("ConvertDateToRoot");
1358 if (deleteIntermediateFiles) {
1359 gSystem->Unlink(dateFileName);
1362 TString selFileName = "selected.";
1363 selFileName += fileName;
1364 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1365 if (fStopOnError) return kFALSE;
1367 if (deleteIntermediateFiles) {
1368 gSystem->Unlink(selDateFileName);
1377 //_____________________________________________________________________________
1378 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1380 // convert the digits to raw data DDL files
1382 AliCodeTimerAuto("",0)
1384 AliRunLoader* runLoader = LoadRun("READ");
1385 if (!runLoader) return kFALSE;
1387 // write raw data to DDL files
1388 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1389 AliInfo(Form("processing event %d", iEvent));
1390 runLoader->GetEvent(iEvent);
1391 TString baseDir = gSystem->WorkingDirectory();
1393 sprintf(dirName, "raw%d", iEvent);
1394 gSystem->MakeDirectory(dirName);
1395 if (!gSystem->ChangeDirectory(dirName)) {
1396 AliError(Form("couldn't change to directory %s", dirName));
1397 if (fStopOnError) return kFALSE; else continue;
1400 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1403 TString detStr = detectors;
1404 if (IsSelected("HLT", detStr)) {
1405 // Do nothing. "HLT" will be removed from detStr and HLT raw
1406 // data files are generated in RunHLT.
1409 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1410 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1411 AliModule* det = (AliModule*) detArray->At(iDet);
1412 if (!det || !det->IsActive()) continue;
1413 if (IsSelected(det->GetName(), detStr)) {
1414 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1419 if (!WriteTriggerRawData())
1420 if (fStopOnError) return kFALSE;
1422 gSystem->ChangeDirectory(baseDir);
1423 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1424 AliError(Form("the following detectors were not found: %s",
1426 if (fStopOnError) return kFALSE;
1435 //_____________________________________________________________________________
1436 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1437 const char* selDateFileName)
1439 // convert raw data DDL files to a DATE file with the program "dateStream"
1440 // The second argument is not empty when the user decides to write
1441 // the detectors raw data according to the trigger cluster.
1443 AliCodeTimerAuto("",0)
1445 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1447 AliError("the program dateStream was not found");
1448 if (fStopOnError) return kFALSE;
1453 AliRunLoader* runLoader = LoadRun("READ");
1454 if (!runLoader) return kFALSE;
1456 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1457 Bool_t selrawdata = kFALSE;
1458 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1461 // Note the option -s. It is used in order to avoid
1462 // the generation of SOR/EOR events.
1463 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1464 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1465 FILE* pipe = gSystem->OpenPipe(command, "w");
1468 AliError(Form("Cannot execute command: %s",command));
1472 Int_t selEvents = 0;
1473 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1475 UInt_t detectorPattern = 0;
1476 runLoader->GetEvent(iEvent);
1477 if (!runLoader->LoadTrigger()) {
1478 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1479 detectorPattern = aCTP->GetClusterMask();
1480 // Check if the event was triggered by CTP
1482 if (aCTP->GetClassMask()) selEvents++;
1486 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1488 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1489 selrawdata = kFALSE;
1493 fprintf(pipe, "GDC DetectorPattern %u\n", detectorPattern);
1497 // loop over detectors and DDLs
1498 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1499 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1501 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1502 Int_t ldcID = Int_t(ldc + 0.0001);
1503 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1505 char rawFileName[256];
1506 sprintf(rawFileName, "raw%d/%s",
1507 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1509 // check existence and size of raw data file
1510 FILE* file = fopen(rawFileName, "rb");
1511 if (!file) continue;
1512 fseek(file, 0, SEEK_END);
1513 unsigned long size = ftell(file);
1515 if (!size) continue;
1517 if (ldcID != prevLDC) {
1518 fprintf(pipe, " LDC Id %d\n", ldcID);
1521 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1526 Int_t result = gSystem->ClosePipe(pipe);
1528 if (!(selrawdata && selEvents > 0)) {
1530 return (result == 0);
1533 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1535 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1536 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1537 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1539 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1541 // Get the trigger decision and cluster
1542 UInt_t detectorPattern = 0;
1544 runLoader->GetEvent(iEvent);
1545 if (!runLoader->LoadTrigger()) {
1546 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1547 if (aCTP->GetClassMask() == 0) continue;
1548 detectorPattern = aCTP->GetClusterMask();
1549 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1550 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1553 fprintf(pipe2, "GDC DetectorPattern %u\n", detectorPattern);
1557 // loop over detectors and DDLs
1558 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1559 // Write only raw data from detectors that
1560 // are contained in the trigger cluster(s)
1561 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1563 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1565 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1566 Int_t ldcID = Int_t(ldc + 0.0001);
1567 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1569 char rawFileName[256];
1570 sprintf(rawFileName, "raw%d/%s",
1571 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1573 // check existence and size of raw data file
1574 FILE* file = fopen(rawFileName, "rb");
1575 if (!file) continue;
1576 fseek(file, 0, SEEK_END);
1577 unsigned long size = ftell(file);
1579 if (!size) continue;
1581 if (ldcID != prevLDC) {
1582 fprintf(pipe2, " LDC Id %d\n", ldcID);
1585 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1590 Int_t result2 = gSystem->ClosePipe(pipe2);
1593 return ((result == 0) && (result2 == 0));
1596 //_____________________________________________________________________________
1597 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1598 const char* rootFileName)
1600 // convert a DATE file to a root file with the program "alimdc"
1603 const Int_t kDBSize = 2000000000;
1604 const Int_t kTagDBSize = 1000000000;
1605 const Bool_t kFilter = kFALSE;
1606 const Int_t kCompression = 1;
1608 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1610 AliError("the program alimdc was not found");
1611 if (fStopOnError) return kFALSE;
1616 AliInfo(Form("converting DATE file %s to root file %s",
1617 dateFileName, rootFileName));
1619 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1620 const char* tagDBFS = "/tmp/mdc1/tags";
1622 // User defined file system locations
1623 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1624 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1625 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1626 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1627 if (gSystem->Getenv("ALIMDC_TAGDB"))
1628 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1630 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1631 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1632 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1634 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1635 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1636 gSystem->Exec(Form("mkdir %s",tagDBFS));
1638 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1639 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1640 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1642 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1643 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1644 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1646 return (result == 0);
1650 //_____________________________________________________________________________
1651 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1653 // delete existing run loaders, open a new one and load gAlice
1655 delete AliRunLoader::Instance();
1656 AliRunLoader* runLoader =
1657 AliRunLoader::Open(fGAliceFileName.Data(),
1658 AliConfig::GetDefaultEventFolderName(), mode);
1660 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1663 runLoader->LoadgAlice();
1664 runLoader->LoadHeader();
1665 gAlice = runLoader->GetAliRun();
1667 AliError(Form("no gAlice object found in file %s",
1668 fGAliceFileName.Data()));
1674 //_____________________________________________________________________________
1675 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1677 // get or calculate the number of signal events per background event
1679 if (!fBkgrdFileNames) return 1;
1680 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1681 if (nBkgrdFiles == 0) return 1;
1683 // get the number of signal events
1685 AliRunLoader* runLoader =
1686 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1687 if (!runLoader) return 1;
1689 nEvents = runLoader->GetNumberOfEvents();
1694 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1695 // get the number of background events
1696 const char* fileName = ((TObjString*)
1697 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1698 AliRunLoader* runLoader =
1699 AliRunLoader::Open(fileName, "BKGRD");
1700 if (!runLoader) continue;
1701 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1704 // get or calculate the number of signal per background events
1705 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1706 if (nSignalPerBkgrd <= 0) {
1707 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1708 } else if (result && (result != nSignalPerBkgrd)) {
1709 AliInfo(Form("the number of signal events per background event "
1710 "will be changed from %d to %d for stream %d",
1711 nSignalPerBkgrd, result, iBkgrdFile+1));
1712 nSignalPerBkgrd = result;
1715 if (!result) result = nSignalPerBkgrd;
1716 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1717 AliWarning(Form("not enough background events (%d) for %d signal events "
1718 "using %d signal per background events for stream %d",
1719 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1726 //_____________________________________________________________________________
1727 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1729 // check whether detName is contained in detectors
1730 // if yes, it is removed from detectors
1732 // check if all detectors are selected
1733 if ((detectors.CompareTo("ALL") == 0) ||
1734 detectors.BeginsWith("ALL ") ||
1735 detectors.EndsWith(" ALL") ||
1736 detectors.Contains(" ALL ")) {
1741 // search for the given detector
1742 Bool_t result = kFALSE;
1743 if ((detectors.CompareTo(detName) == 0) ||
1744 detectors.BeginsWith(detName+" ") ||
1745 detectors.EndsWith(" "+detName) ||
1746 detectors.Contains(" "+detName+" ")) {
1747 detectors.ReplaceAll(detName, "");
1751 // clean up the detectors string
1752 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1753 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1754 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1759 //_____________________________________________________________________________
1760 Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName)
1763 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1764 // These can be used for embedding of MC tracks into RAW data using the standard
1765 // merging procedure.
1767 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1770 AliError("no gAlice object. Restart aliroot and try again.");
1773 if (gAlice->Modules()->GetEntries() > 0) {
1774 AliError("gAlice was already run. Restart aliroot and try again.");
1778 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1782 gROOT->LoadMacro(fConfigFileName.Data());
1783 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1785 if(AliCDBManager::Instance()->GetRun() >= 0) {
1786 SetRunNumber(AliCDBManager::Instance()->GetRun());
1788 AliWarning("Run number not initialized!!");
1791 AliRunLoader::Instance()->CdGAFile();
1793 AliPDG::AddParticlesToPdgDataBase();
1795 gAlice->GetMCApp()->Init();
1797 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1799 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1800 gAlice->InitLoaders();
1801 AliRunLoader::Instance()->MakeTree("E");
1802 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1803 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1804 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1806 // Save stuff at the beginning of the file to avoid file corruption
1807 AliRunLoader::Instance()->CdGAFile();
1812 //AliCDBManager* man = AliCDBManager::Instance();
1813 //man->SetRun(0); // Should this come from rawdata header ?
1817 // Get the runloader
1818 AliRunLoader* runLoader = AliRunLoader::Instance();
1820 // Open esd file if available
1821 TFile* esdFile = TFile::Open(esdFileName);
1823 AliESDEvent* esd = new AliESDEvent();
1824 esdFile->GetObject("esdTree", treeESD);
1825 if (treeESD) esd->ReadFromTree(treeESD);
1828 // Create the RawReader
1829 TString fileName(rawDirectory);
1830 AliRawReader* rawReader = 0x0;
1831 if (fileName.EndsWith("/")) {
1832 rawReader = new AliRawReaderFile(fileName);
1833 } else if (fileName.EndsWith(".root")) {
1834 rawReader = new AliRawReaderRoot(fileName);
1835 } else if (!fileName.IsNull()) {
1836 rawReader = new AliRawReaderDate(fileName);
1838 // if (!fEquipIdMap.IsNull() && fRawReader)
1839 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1841 // Get list of detectors
1842 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1845 AliHeader* header = runLoader->GetHeader();
1847 TString detStr = fMakeSDigits;
1851 if (!(rawReader->NextEvent())) break;
1854 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1855 AliModule* det = (AliModule*) detArray->At(iDet);
1856 if (!det || !det->IsActive()) continue;
1857 if (IsSelected(det->GetName(), detStr)) {
1858 AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1859 det->Raw2SDigits(rawReader);
1866 // If ESD information available obtain reconstructed vertex and store in header.
1868 treeESD->GetEvent(nev);
1869 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1870 Double_t position[3];
1871 esdVertex->GetXYZ(position);
1872 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1875 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1876 mcHeader->SetPrimaryVertex(mcV);
1877 header->Reset(0,nev);
1878 header->SetGenEventHeader(mcHeader);
1879 printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1884 runLoader->TreeE()->Fill();
1885 runLoader->SetNextEvent();
1891 runLoader->CdGAFile();
1892 runLoader->WriteHeader("OVERWRITE");
1893 runLoader->WriteRunLoader();
1898 //_____________________________________________________________________________
1899 void AliSimulation::FinishRun()
1902 // Called at the end of the run.
1907 AliDebug(1, "Finish Lego");
1908 AliRunLoader::Instance()->CdGAFile();
1912 // Clean detector information
1913 TIter next(gAlice->Modules());
1914 AliModule *detector;
1915 while((detector = dynamic_cast<AliModule*>(next()))) {
1916 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
1917 detector->FinishRun();
1920 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
1921 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
1923 // Write AliRun info and all detectors parameters
1924 AliRunLoader::Instance()->CdGAFile();
1925 gAlice->Write(0,TObject::kOverwrite);//write AliRun
1926 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
1928 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
1929 AliRunLoader::Instance()->Synchronize();
1932 //_____________________________________________________________________________
1933 Int_t AliSimulation::GetDetIndex(const char* detector)
1935 // return the detector index corresponding to detector
1937 for (index = 0; index < fgkNDetectors ; index++) {
1938 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
1944 //_____________________________________________________________________________
1945 Bool_t AliSimulation::CreateHLT()
1947 // Init the HLT simulation.
1948 // The function loads the library and creates the instance of AliHLTSimulation.
1949 // the main reason for the decoupled creation is to set the transient OCDB
1950 // objects before the OCDB is locked
1952 // load the library dynamically
1953 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
1955 // check for the library version
1956 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
1958 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
1961 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
1962 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
1965 // print compile info
1966 typedef void (*CompileInfo)( const char*& date, const char*& time);
1967 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
1969 const char* date="";
1970 const char* time="";
1971 (*fctInfo)(date, time);
1972 if (!date) date="unknown";
1973 if (!time) time="unknown";
1974 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
1976 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
1979 // create instance of the HLT simulation
1980 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
1981 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
1982 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
1986 TString specObjects;
1987 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
1988 if (specObjects.Length()>0) specObjects+=" ";
1989 specObjects+=fSpecCDBUri[i]->GetName();
1992 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
1993 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
1994 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2000 //_____________________________________________________________________________
2001 Bool_t AliSimulation::RunHLT()
2003 // Run the HLT simulation
2004 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2005 // Disabled if fRunHLT is empty, default vaule is "default".
2006 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2007 // The default simulation depends on the HLT component libraries and their
2008 // corresponding agents which define components and chains to run. See
2009 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2010 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2012 // The libraries to be loaded can be specified as an option.
2014 // AliSimulation sim;
2015 // sim.SetRunHLT("libAliHLTSample.so");
2017 // will only load <tt>libAliHLTSample.so</tt>
2019 // Other available options:
2020 // \li loglevel=<i>level</i> <br>
2021 // logging level for this processing
2023 // disable redirection of log messages to AliLog class
2024 // \li config=<i>macro</i>
2025 // configuration macro
2026 // \li chains=<i>configuration</i>
2027 // comma separated list of configurations to be run during simulation
2028 // \li rawfile=<i>file</i>
2029 // source for the RawReader to be created, the default is <i>./</i> if
2030 // raw data is simulated
2034 if (!fpHLT && !CreateHLT()) {
2037 AliHLTSimulation* pHLT=fpHLT;
2039 AliRunLoader* pRunLoader = LoadRun("READ");
2040 if (!pRunLoader) return kFALSE;
2042 // initialize CDB storage, run number, set CDB lock
2043 // thats for the case of running HLT simulation without all the other steps
2044 // multiple calls are handled by the function, so we can just call
2046 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2049 // init the HLT simulation
2051 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2052 TString detStr = fWriteRawData;
2053 if (!IsSelected("HLT", detStr)) {
2054 options+=" writerawfiles=";
2056 options+=" writerawfiles=HLT";
2059 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2060 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2061 // in order to get detector data. By default, RawReaderFile is used to read
2062 // the already simulated ddl files. Date and Root files from the raw data
2063 // are generated after the HLT simulation.
2064 options+=" rawfile=./";
2067 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2068 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2069 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2071 // run the HLT simulation
2072 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2073 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2074 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2078 // delete the instance
2079 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2080 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2081 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2085 return iResult>=0?kTRUE:kFALSE;
2088 //_____________________________________________________________________________
2089 Bool_t AliSimulation::RunQA()
2091 // run the QA on summable hits, digits or digits
2093 if(!gAlice) return kFALSE;
2094 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2096 TString detectorsw("") ;
2098 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2099 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2100 if ( detectorsw.IsNull() )
2105 //_____________________________________________________________________________
2106 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2108 // Allows to run QA for a selected set of detectors
2109 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2110 // all selected detectors run the same selected tasks
2112 if (!detAndAction.Contains(":")) {
2113 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2117 Int_t colon = detAndAction.Index(":") ;
2118 fQADetectors = detAndAction(0, colon) ;
2119 if (fQADetectors.Contains("ALL") )
2120 fQADetectors = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2121 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2122 if (fQATasks.Contains("ALL") ) {
2123 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2125 fQATasks.ToUpper() ;
2127 if ( fQATasks.Contains("HIT") )
2128 tempo = Form("%d ", AliQAv1::kHITS) ;
2129 if ( fQATasks.Contains("SDIGIT") )
2130 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2131 if ( fQATasks.Contains("DIGIT") )
2132 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2134 if (fQATasks.IsNull()) {
2135 AliInfo("No QA requested\n") ;
2140 TString tempo(fQATasks) ;
2141 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2142 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2143 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2144 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2146 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2147 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2148 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2149 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2154 //_____________________________________________________________________________
2155 void AliSimulation::ProcessEnvironmentVars()
2157 // Extract run number and random generator seed from env variables
2159 AliInfo("Processing environment variables");
2161 // Random Number seed
2163 // first check that seed is not already set
2165 if (gSystem->Getenv("CONFIG_SEED")) {
2166 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2169 if (gSystem->Getenv("CONFIG_SEED")) {
2170 AliInfo(Form("Seed for random number generation already set (%d)"
2171 ": CONFIG_SEED variable ignored!", fSeed));
2175 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2179 // first check that run number is not already set
2181 if (gSystem->Getenv("DC_RUN")) {
2182 fRun = atoi(gSystem->Getenv("DC_RUN"));
2185 if (gSystem->Getenv("DC_RUN")) {
2186 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2190 AliInfo(Form("Run number = %d", fRun));
2193 //---------------------------------------------------------------------
2194 void AliSimulation::WriteGRPEntry()
2196 // Get the necessary information from galice (generator, trigger etc) and
2197 // write a GRP entry corresponding to the settings in the Config.C used
2198 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2201 AliInfo("Writing global run parameters entry into the OCDB");
2203 AliGRPObject* grpObj = new AliGRPObject();
2205 grpObj->SetRunType("PHYSICS");
2206 grpObj->SetTimeStart(0);
2207 grpObj->SetTimeEnd(9999);
2209 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2211 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/0.120);
2214 gen->GetProjectile(projectile,a,z);
2216 gen->GetTarget(target,a,z);
2217 TString beamType = projectile + "-" + target;
2218 beamType.ReplaceAll(" ","");
2219 if (!beamType.CompareTo("-")) {
2220 grpObj->SetBeamType("UNKNOWN");
2223 grpObj->SetBeamType(beamType);
2224 // Heavy ion run, the event specie is set to kHighMult
2225 fEventSpecie = AliRecoParam::kHighMult;
2226 if ((strcmp(beamType,"p-p") == 0) ||
2227 (strcmp(beamType,"p-") == 0) ||
2228 (strcmp(beamType,"-p") == 0) ||
2229 (strcmp(beamType,"P-P") == 0) ||
2230 (strcmp(beamType,"P-") == 0) ||
2231 (strcmp(beamType,"-P") == 0)) {
2232 // Proton run, the event specie is set to kLowMult
2233 fEventSpecie = AliRecoParam::kLowMult;
2237 AliWarning("Unknown beam type and energy! Setting energy to 0");
2238 grpObj->SetBeamEnergy(0);
2239 grpObj->SetBeamType("UNKNOWN");
2242 UInt_t detectorPattern = 0;
2244 TObjArray *detArray = gAlice->Detectors();
2245 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2246 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2247 detectorPattern |= (1 << iDet);
2252 if (!fMakeTrigger.IsNull() || strcmp(gAlice->GetTriggerDescriptor(),""))
2253 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2256 if (!fRunHLT.IsNull())
2257 detectorPattern |= (1 << AliDAQ::kHLTId);
2259 grpObj->SetNumberOfDetectors((Char_t)nDets);
2260 grpObj->SetDetectorMask((Int_t)detectorPattern);
2261 grpObj->SetLHCPeriod("LHC08c");
2262 grpObj->SetLHCState("STABLE_BEAMS");
2263 grpObj->SetLHCLuminosity(0,(AliGRPObject::Stats)0);
2264 grpObj->SetBeamIntensity(0,(AliGRPObject::Stats)0);
2266 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2267 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2269 Float_t factorSol = field ? field->GetFactorSol() : 0;
2270 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2271 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2273 Float_t factorDip = field ? field->GetFactorDip() : 0;
2274 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2276 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2277 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2278 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2279 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2280 grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2281 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2283 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2285 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2287 // Now store the entry in OCDB
2288 AliCDBManager* man = AliCDBManager::Instance();
2290 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2291 AliCDBMetaData *metadata= new AliCDBMetaData();
2293 metadata->SetResponsible("alice-off@cern.ch");
2294 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2296 man->Put(grpObj,id,metadata);