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>
119 #include "AliAlignObj.h"
120 #include "AliCDBEntry.h"
121 #include "AliCDBManager.h"
122 #include "AliCDBStorage.h"
123 #include "AliCTPRawData.h"
124 #include "AliCentralTrigger.h"
125 #include "AliCentralTrigger.h"
126 #include "AliCodeTimer.h"
128 #include "AliDigitizer.h"
129 #include "AliESDEvent.h"
130 #include "AliGRPObject.h"
131 #include "AliGenEventHeader.h"
132 #include "AliGenerator.h"
133 #include "AliGeomManager.h"
134 #include "AliHLTSimulation.h"
135 #include "AliHeader.h"
137 #include "AliLegoGenerator.h"
141 #include "AliModule.h"
143 #include "AliRawReaderDate.h"
144 #include "AliRawReaderFile.h"
145 #include "AliRawReaderRoot.h"
147 #include "AliRunDigitizer.h"
148 #include "AliRunLoader.h"
149 #include "AliSimulation.h"
150 #include "AliSysInfo.h"
151 #include "AliVertexGenFile.h"
153 ClassImp(AliSimulation)
155 AliSimulation *AliSimulation::fgInstance = 0;
156 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
158 //_____________________________________________________________________________
159 AliSimulation::AliSimulation(const char* configFileName,
160 const char* name, const char* title) :
163 fRunGeneration(kTRUE),
164 fRunSimulation(kTRUE),
165 fLoadAlignFromCDB(kTRUE),
166 fLoadAlObjsListOfDets("ALL"),
170 fMakeDigitsFromHits(""),
172 fRawDataFileName(""),
173 fDeleteIntermediateFiles(kFALSE),
174 fWriteSelRawData(kFALSE),
175 fStopOnError(kFALSE),
177 fConfigFileName(configFileName),
178 fGAliceFileName("galice.root"),
180 fBkgrdFileNames(NULL),
181 fAlignObjArray(NULL),
182 fUseBkgrdVertex(kTRUE),
183 fRegionOfInterest(kFALSE),
189 fInitCDBCalled(kFALSE),
190 fInitRunNumberCalled(kFALSE),
191 fSetRunNumberFromDataCalled(kFALSE),
192 fEmbeddingFlag(kFALSE),
197 fEventSpecie(AliRecoParam::kDefault),
198 fWriteQAExpertData(kTRUE),
202 fWriteGRPEntry(kTRUE)
204 // create simulation object with default parameters
206 SetGAliceFile("galice.root");
209 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
210 qam->SetActiveDetectors(fQADetectors) ;
211 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
212 qam->SetTasks(fQATasks) ;
215 //_____________________________________________________________________________
216 AliSimulation::~AliSimulation()
220 fEventsPerFile.Delete();
221 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
222 // delete fAlignObjArray; fAlignObjArray=0;
224 if (fBkgrdFileNames) {
225 fBkgrdFileNames->Delete();
226 delete fBkgrdFileNames;
229 fSpecCDBUri.Delete();
230 if (fgInstance==this) fgInstance = 0;
232 AliQAManager::QAManager()->ShowQA() ;
233 AliQAManager::Destroy() ;
234 AliCodeTimer::Instance()->Print();
238 //_____________________________________________________________________________
239 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
241 // set the number of events for one run
246 //_____________________________________________________________________________
247 void AliSimulation::InitQA()
249 // activate a default CDB storage
250 // First check if we have any CDB storage set, because it is used
251 // to retrieve the calibration and alignment constants
253 if (fInitCDBCalled) return;
254 fInitCDBCalled = kTRUE;
256 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kSIMMODE) ;
257 qam->SetActiveDetectors(fQADetectors) ;
258 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
259 qam->SetTasks(fQATasks) ;
260 if (fWriteQAExpertData)
261 qam->SetWriteExpert() ;
263 if (qam->IsDefaultStorageSet()) {
264 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
265 AliWarning("Default QA reference storage has been already set !");
266 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fQARefUri.Data()));
267 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
268 fQARefUri = qam->GetDefaultStorage()->GetURI();
270 if (fQARefUri.Length() > 0) {
271 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
272 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
273 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
275 fQARefUri="local://$ALICE_ROOT/QARef";
276 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
277 AliWarning("Default QA reference storage not yet set !!!!");
278 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
279 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
281 qam->SetDefaultStorage(fQARefUri);
285 //_____________________________________________________________________________
286 void AliSimulation::InitCDB()
288 // activate a default CDB storage
289 // First check if we have any CDB storage set, because it is used
290 // to retrieve the calibration and alignment constants
292 if (fInitCDBCalled) return;
293 fInitCDBCalled = kTRUE;
295 AliCDBManager* man = AliCDBManager::Instance();
296 if (man->IsDefaultStorageSet())
298 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
299 AliWarning("Default CDB storage has been already set !");
300 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
301 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
302 fCDBUri = man->GetDefaultStorage()->GetURI();
305 if (fCDBUri.Length() > 0)
307 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
308 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
309 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
311 fCDBUri="local://$ALICE_ROOT/OCDB";
312 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
313 AliWarning("Default CDB storage not yet set !!!!");
314 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
315 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
318 man->SetDefaultStorage(fCDBUri);
321 // Now activate the detector specific CDB storage locations
322 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
323 TObject* obj = fSpecCDBUri[i];
325 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
326 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
327 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
328 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
333 //_____________________________________________________________________________
334 void AliSimulation::InitRunNumber(){
335 // check run number. If not set, set it to 0 !!!!
337 if (fInitRunNumberCalled) return;
338 fInitRunNumberCalled = kTRUE;
340 AliCDBManager* man = AliCDBManager::Instance();
341 if (man->GetRun() >= 0)
343 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
344 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
348 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
351 AliWarning(Form("Run number not yet set !!!! Setting it now to: %d",
360 //_____________________________________________________________________________
361 void AliSimulation::SetCDBLock() {
362 // Set CDB lock: from now on it is forbidden to reset the run number
363 // or the default storage or to activate any further storage!
365 AliCDBManager::Instance()->SetLock(1);
368 //_____________________________________________________________________________
369 void AliSimulation::SetDefaultStorage(const char* uri) {
370 // Store the desired default CDB storage location
371 // Activate it later within the Run() method
377 //_____________________________________________________________________________
378 void AliSimulation::SetQARefDefaultStorage(const char* uri) {
379 // Store the desired default CDB storage location
380 // Activate it later within the Run() method
383 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
386 //_____________________________________________________________________________
387 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
388 // Store a detector-specific CDB storage location
389 // Activate it later within the Run() method
391 AliCDBPath aPath(calibType);
392 if(!aPath.IsValid()){
393 AliError(Form("Not a valid path: %s", calibType));
397 TObject* obj = fSpecCDBUri.FindObject(calibType);
398 if (obj) fSpecCDBUri.Remove(obj);
399 fSpecCDBUri.Add(new TNamed(calibType, uri));
403 //_____________________________________________________________________________
404 void AliSimulation::SetRunNumber(Int_t run)
407 // Activate it later within the Run() method
412 //_____________________________________________________________________________
413 void AliSimulation::SetSeed(Int_t seed)
416 // Activate it later within the Run() method
421 //_____________________________________________________________________________
422 Bool_t AliSimulation::SetRunNumberFromData()
424 // Set the CDB manager run number
425 // The run number is retrieved from gAlice
427 if (fSetRunNumberFromDataCalled) return kTRUE;
428 fSetRunNumberFromDataCalled = kTRUE;
430 AliCDBManager* man = AliCDBManager::Instance();
431 Int_t runData = -1, runCDB = -1;
433 AliRunLoader* runLoader = LoadRun("READ");
434 if (!runLoader) return kFALSE;
436 runData = runLoader->GetHeader()->GetRun();
440 runCDB = man->GetRun();
442 if (runCDB != runData) {
443 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
444 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
445 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
446 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
451 man->SetRun(runData);
454 if(man->GetRun() < 0) {
455 AliError("Run number not properly initalized!");
464 //_____________________________________________________________________________
465 void AliSimulation::SetConfigFile(const char* fileName)
467 // set the name of the config file
469 fConfigFileName = fileName;
472 //_____________________________________________________________________________
473 void AliSimulation::SetGAliceFile(const char* fileName)
475 // set the name of the galice file
476 // the path is converted to an absolute one if it is relative
478 fGAliceFileName = fileName;
479 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
480 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
482 fGAliceFileName = absFileName;
483 delete[] absFileName;
486 AliDebug(2, Form("galice file name set to %s", fileName));
489 //_____________________________________________________________________________
490 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
493 // set the number of events per file for the given detector and data type
494 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
496 TNamed* obj = new TNamed(detector, type);
497 obj->SetUniqueID(nEvents);
498 fEventsPerFile.Add(obj);
501 //_____________________________________________________________________________
502 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
504 // Read the alignment objects from CDB.
505 // Each detector is supposed to have the
506 // alignment objects in DET/Align/Data CDB path.
507 // All the detector objects are then collected,
508 // sorted by geometry level (starting from ALIC) and
509 // then applied to the TGeo geometry.
510 // Finally an overlaps check is performed.
512 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
513 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
517 // initialize CDB storage, run number, set CDB lock
519 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
522 Bool_t delRunLoader = kFALSE;
524 runLoader = LoadRun("READ");
525 if (!runLoader) return kFALSE;
526 delRunLoader = kTRUE;
529 // Export ideal geometry
530 if(!IsGeometryFromFile()) AliGeomManager::GetGeometry()->Export("geometry.root");
532 // Load alignment data from CDB and apply to geometry through AliGeomManager
533 if(fLoadAlignFromCDB){
535 TString detStr = fLoadAlObjsListOfDets;
536 TString loadAlObjsListOfDets = "";
538 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
539 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
540 AliModule* det = (AliModule*) detArray->At(iDet);
541 if (!det || !det->IsActive()) continue;
542 if (IsSelected(det->GetName(), detStr)) {
543 //add det to list of dets to be aligned from CDB
544 loadAlObjsListOfDets += det->GetName();
545 loadAlObjsListOfDets += " ";
547 } // end loop over detectors
548 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
549 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
551 // Check if the array with alignment objects was
552 // provided by the user. If yes, apply the objects
553 // to the present TGeo geometry
554 if (fAlignObjArray) {
555 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
556 AliError("The misalignment of one or more volumes failed!"
557 "Compare the list of simulated detectors and the list of detector alignment data!");
558 if (delRunLoader) delete runLoader;
564 // Update the internal geometry of modules (ITS needs it)
565 TString detStr = fLoadAlObjsListOfDets;
566 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
567 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
569 AliModule* det = (AliModule*) detArray->At(iDet);
570 if (!det || !det->IsActive()) continue;
571 if (IsSelected(det->GetName(), detStr)) {
572 det->UpdateInternalGeometry();
574 } // end loop over detectors
577 if (delRunLoader) delete runLoader;
582 //_____________________________________________________________________________
583 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
585 // add a file with background events for merging
587 TObjString* fileNameStr = new TObjString(fileName);
588 fileNameStr->SetUniqueID(nSignalPerBkgrd);
589 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
590 fBkgrdFileNames->Add(fileNameStr);
593 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
595 // add a file with background events for embeddin
596 MergeWith(fileName, nSignalPerBkgrd);
597 fEmbeddingFlag = kTRUE;
600 //_____________________________________________________________________________
601 Bool_t AliSimulation::Run(Int_t nEvents)
603 // run the generation, simulation and digitization
606 AliCodeTimerAuto("",0)
607 AliSysInfo::AddStamp("Start_Run");
609 // Load run number and seed from environmental vars
610 ProcessEnvironmentVars();
611 AliSysInfo::AddStamp("ProcessEnvironmentVars");
613 gRandom->SetSeed(fSeed);
615 if (nEvents > 0) fNEvents = nEvents;
617 // create and setup the HLT instance
618 if (!fRunHLT.IsNull() && !CreateHLT()) {
619 if (fStopOnError) return kFALSE;
624 // generation and simulation -> hits
625 if (fRunGeneration) {
626 if (!RunSimulation()) if (fStopOnError) return kFALSE;
628 AliSysInfo::AddStamp("RunSimulation");
630 // initialize CDB storage from external environment
631 // (either CDB manager or AliSimulation setters),
632 // if not already done in RunSimulation()
634 AliSysInfo::AddStamp("InitCDB");
636 // Set run number in CDBManager from data
637 // From this point on the run number must be always loaded from data!
638 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
640 // Set CDB lock: from now on it is forbidden to reset the run number
641 // or the default storage or to activate any further storage!
644 // If RunSimulation was not called, load the geometry and misalign it
645 if (!AliGeomManager::GetGeometry()) {
646 // Initialize the geometry manager
647 AliGeomManager::LoadGeometry("geometry.root");
648 AliSysInfo::AddStamp("GetGeometry");
651 // // Check that the consistency of symbolic names for the activated subdetectors
652 // // in the geometry loaded by AliGeomManager
653 // AliRunLoader* runLoader = LoadRun("READ");
654 // if (!runLoader) return kFALSE;
656 // TString detsToBeChecked = "";
657 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
658 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
659 // AliModule* det = (AliModule*) detArray->At(iDet);
660 // if (!det || !det->IsActive()) continue;
661 // detsToBeChecked += det->GetName();
662 // detsToBeChecked += " ";
663 // } // end loop over detectors
664 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
665 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
666 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
668 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
670 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
672 AliSysInfo::AddStamp("MissalignGeometry");
675 // hits -> summable digits
676 AliSysInfo::AddStamp("Start_sdigitization");
677 if (!fMakeSDigits.IsNull()) {
678 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
681 AliSysInfo::AddStamp("Stop_sdigitization");
683 AliSysInfo::AddStamp("Start_digitization");
684 // summable digits -> digits
685 if (!fMakeDigits.IsNull()) {
686 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
687 if (fStopOnError) return kFALSE;
690 AliSysInfo::AddStamp("Stop_digitization");
695 if (!fMakeDigitsFromHits.IsNull()) {
696 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
697 AliWarning(Form("Merging and direct creation of digits from hits "
698 "was selected for some detectors. "
699 "No merging will be done for the following detectors: %s",
700 fMakeDigitsFromHits.Data()));
702 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
703 if (fStopOnError) return kFALSE;
707 AliSysInfo::AddStamp("Hits2Digits");
711 if (!RunTrigger(fMakeTrigger,fMakeDigits)) {
712 if (fStopOnError) return kFALSE;
715 AliSysInfo::AddStamp("RunTrigger");
718 // digits -> raw data
719 if (!fWriteRawData.IsNull()) {
720 if (!WriteRawData(fWriteRawData, fRawDataFileName,
721 fDeleteIntermediateFiles,fWriteSelRawData)) {
722 if (fStopOnError) return kFALSE;
726 AliSysInfo::AddStamp("WriteRaw");
728 // run HLT simulation on simulated digit data if raw data is not
729 // simulated, otherwise its called as part of WriteRawData
730 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
732 if (fStopOnError) return kFALSE;
736 AliSysInfo::AddStamp("RunHLT");
740 Bool_t rv = RunQA() ;
746 AliSysInfo::AddStamp("RunQA");
748 // Cleanup of CDB manager: cache and active storages!
749 AliCDBManager::Instance()->ClearCache();
754 //_______________________________________________________________________
755 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
756 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
757 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
760 // Generates lego plots of:
761 // - radiation length map phi vs theta
762 // - radiation length map phi vs eta
763 // - interaction length map
764 // - g/cm2 length map
766 // ntheta bins in theta, eta
767 // themin minimum angle in theta (degrees)
768 // themax maximum angle in theta (degrees)
770 // phimin minimum angle in phi (degrees)
771 // phimax maximum angle in phi (degrees)
772 // rmin minimum radius
773 // rmax maximum radius
776 // The number of events generated = ntheta*nphi
777 // run input parameters in macro setup (default="Config.C")
779 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
782 <img src="picts/AliRunLego1.gif">
787 <img src="picts/AliRunLego2.gif">
792 <img src="picts/AliRunLego3.gif">
797 // run the generation and simulation
799 AliCodeTimerAuto("",0)
801 // initialize CDB storage and run number from external environment
802 // (either CDB manager or AliSimulation setters)
808 AliError("no gAlice object. Restart aliroot and try again.");
811 if (gAlice->Modules()->GetEntries() > 0) {
812 AliError("gAlice was already run. Restart aliroot and try again.");
816 AliInfo(Form("initializing gAlice with config file %s",
817 fConfigFileName.Data()));
820 if (nev == -1) nev = nc1 * nc2;
822 // check if initialisation has been done
823 // If runloader has been initialized, set the number of events per file to nc1 * nc2
826 if (!gener) gener = new AliLegoGenerator();
828 // Configure Generator
830 gener->SetRadiusRange(rmin, rmax);
831 gener->SetZMax(zmax);
832 gener->SetCoor1Range(nc1, c1min, c1max);
833 gener->SetCoor2Range(nc2, c2min, c2max);
837 fLego = new AliLego("lego",gener);
839 //__________________________________________________________________________
843 gROOT->LoadMacro(setup);
844 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
846 if(AliCDBManager::Instance()->GetRun() >= 0) {
847 SetRunNumber(AliCDBManager::Instance()->GetRun());
849 AliWarning("Run number not initialized!!");
852 AliRunLoader::Instance()->CdGAFile();
854 AliPDG::AddParticlesToPdgDataBase();
856 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
858 gAlice->GetMCApp()->Init();
861 //Must be here because some MCs (G4) adds detectors here and not in Config.C
862 gAlice->InitLoaders();
863 AliRunLoader::Instance()->MakeTree("E");
866 // Save stuff at the beginning of the file to avoid file corruption
867 AliRunLoader::Instance()->CdGAFile();
870 //Save current generator
871 AliGenerator *gen=gAlice->GetMCApp()->Generator();
872 gAlice->GetMCApp()->ResetGenerator(gener);
873 //Prepare MC for Lego Run
879 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
880 gMC->ProcessRun(nev);
882 // End of this run, close files
884 // Restore current generator
885 gAlice->GetMCApp()->ResetGenerator(gen);
886 // Delete Lego Object
892 //_____________________________________________________________________________
893 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
897 AliCodeTimerAuto("",0)
899 // initialize CDB storage from external environment
900 // (either CDB manager or AliSimulation setters),
901 // if not already done in RunSimulation()
904 // Set run number in CDBManager from data
905 // From this point on the run number must be always loaded from data!
906 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
908 // Set CDB lock: from now on it is forbidden to reset the run number
909 // or the default storage or to activate any further storage!
912 AliRunLoader* runLoader = LoadRun("READ");
913 if (!runLoader) return kFALSE;
914 TString trconfiguration = config;
916 if (trconfiguration.IsNull()) {
917 if (strcmp(gAlice->GetTriggerDescriptor(),"")) {
918 trconfiguration = gAlice->GetTriggerDescriptor();
921 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
924 runLoader->MakeTree( "GG" );
925 AliCentralTrigger* aCTP = runLoader->GetTrigger();
926 // Load Configuration
927 if (!aCTP->LoadConfiguration( trconfiguration ))
931 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
943 //_____________________________________________________________________________
944 Bool_t AliSimulation::WriteTriggerRawData()
946 // Writes the CTP (trigger) DDL raw data
947 // Details of the format are given in the
948 // trigger TDR - pages 134 and 135.
949 AliCTPRawData writer;
955 //_____________________________________________________________________________
956 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
958 // run the generation and simulation
960 AliCodeTimerAuto("",0)
962 // initialize CDB storage and run number from external environment
963 // (either CDB manager or AliSimulation setters)
964 AliSysInfo::AddStamp("RunSimulation_Begin");
966 AliSysInfo::AddStamp("RunSimulation_InitCDB");
969 AliSysInfo::AddStamp("RunSimulation_SetCDBLock");
972 AliError("no gAlice object. Restart aliroot and try again.");
975 if (gAlice->Modules()->GetEntries() > 0) {
976 AliError("gAlice was already run. Restart aliroot and try again.");
980 AliInfo(Form("initializing gAlice with config file %s",
981 fConfigFileName.Data()));
984 // Initialize ALICE Simulation run
989 gROOT->LoadMacro(fConfigFileName.Data());
990 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
991 AliSysInfo::AddStamp("RunSimulation_Config");
993 if(AliCDBManager::Instance()->GetRun() >= 0) {
994 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
995 AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
997 AliWarning("Run number not initialized!!");
1000 AliRunLoader::Instance()->CdGAFile();
1002 AliPDG::AddParticlesToPdgDataBase();
1004 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1005 AliSysInfo::AddStamp("RunSimulation_GetField");
1007 gAlice->GetMCApp()->Init();
1008 AliSysInfo::AddStamp("RunSimulation_InitMCApp");
1010 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1011 gAlice->InitLoaders();
1012 AliRunLoader::Instance()->MakeTree("E");
1013 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1014 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1015 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1017 // Save stuff at the beginning of the file to avoid file corruption
1018 AliRunLoader::Instance()->CdGAFile();
1020 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
1021 AliSysInfo::AddStamp("RunSimulation_InitLoaders");
1022 //___________________________________________________________________________________________
1024 // Get the trigger descriptor string
1025 // Either from AliSimulation or from
1027 if (fMakeTrigger.IsNull()) {
1028 if (strcmp(gAlice->GetTriggerDescriptor(),""))
1029 fMakeTrigger = gAlice->GetTriggerDescriptor();
1032 gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
1033 AliSysInfo::AddStamp("RunSimulation_TriggerDescriptor");
1035 // Set run number in CDBManager
1036 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
1038 AliRunLoader* runLoader = AliRunLoader::Instance();
1040 AliError(Form("gAlice has no run loader object. "
1041 "Check your config file: %s", fConfigFileName.Data()));
1044 SetGAliceFile(runLoader->GetFileName());
1046 // Misalign geometry
1047 #if ROOT_VERSION_CODE < 331527
1048 AliGeomManager::SetGeometry(gGeoManager);
1050 // Check that the consistency of symbolic names for the activated subdetectors
1051 // in the geometry loaded by AliGeomManager
1052 TString detsToBeChecked = "";
1053 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1054 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1055 AliModule* det = (AliModule*) detArray->At(iDet);
1056 if (!det || !det->IsActive()) continue;
1057 detsToBeChecked += det->GetName();
1058 detsToBeChecked += " ";
1059 } // end loop over detectors
1060 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
1061 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
1062 MisalignGeometry(runLoader);
1063 AliSysInfo::AddStamp("RunSimulation_MisalignGeometry");
1066 // AliRunLoader* runLoader = AliRunLoader::Instance();
1067 // if (!runLoader) {
1068 // AliError(Form("gAlice has no run loader object. "
1069 // "Check your config file: %s", fConfigFileName.Data()));
1072 // SetGAliceFile(runLoader->GetFileName());
1074 if (!gAlice->GetMCApp()->Generator()) {
1075 AliError(Form("gAlice has no generator object. "
1076 "Check your config file: %s", fConfigFileName.Data()));
1080 // Write GRP entry corresponding to the setting found in Cofig.C
1083 AliSysInfo::AddStamp("RunSimulation_WriteGRP");
1085 if (nEvents <= 0) nEvents = fNEvents;
1087 // get vertex from background file in case of merging
1088 if (fUseBkgrdVertex &&
1089 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1090 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1091 const char* fileName = ((TObjString*)
1092 (fBkgrdFileNames->At(0)))->GetName();
1093 AliInfo(Form("The vertex will be taken from the background "
1094 "file %s with nSignalPerBackground = %d",
1095 fileName, signalPerBkgrd));
1096 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1097 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1100 if (!fRunSimulation) {
1101 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1104 // set the number of events per file for given detectors and data types
1105 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1106 if (!fEventsPerFile[i]) continue;
1107 const char* detName = fEventsPerFile[i]->GetName();
1108 const char* typeName = fEventsPerFile[i]->GetTitle();
1109 TString loaderName(detName);
1110 loaderName += "Loader";
1111 AliLoader* loader = runLoader->GetLoader(loaderName);
1113 AliError(Form("RunSimulation", "no loader for %s found\n"
1114 "Number of events per file not set for %s %s",
1115 detName, typeName, detName));
1118 AliDataLoader* dataLoader =
1119 loader->GetDataLoader(typeName);
1121 AliError(Form("no data loader for %s found\n"
1122 "Number of events per file not set for %s %s",
1123 typeName, detName, typeName));
1126 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1127 AliDebug(1, Form("number of events per file set to %d for %s %s",
1128 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1131 AliInfo("running gAlice");
1132 AliSysInfo::AddStamp("Start_ProcessRun");
1134 // Create the Root Tree with one branch per detector
1135 //Hits moved to begin event -> now we are crating separate tree for each event
1137 gMC->ProcessRun(nEvents);
1139 // End of this run, close files
1140 if(nEvents>0) FinishRun();
1142 AliSysInfo::AddStamp("Stop_ProcessRun");
1148 //_____________________________________________________________________________
1149 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1151 // run the digitization and produce summable digits
1152 static Int_t eventNr=0;
1153 AliCodeTimerAuto("",0) ;
1155 // initialize CDB storage, run number, set CDB lock
1157 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1160 AliRunLoader* runLoader = LoadRun();
1161 if (!runLoader) return kFALSE;
1163 TString detStr = detectors;
1164 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1165 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1166 AliModule* det = (AliModule*) detArray->At(iDet);
1167 if (!det || !det->IsActive()) continue;
1168 if (IsSelected(det->GetName(), detStr)) {
1169 AliInfo(Form("creating summable digits for %s", det->GetName()));
1170 AliCodeTimerStart(Form("creating summable digits for %s", det->GetName()));
1171 det->Hits2SDigits();
1172 AliCodeTimerStop(Form("creating summable digits for %s", det->GetName()));
1173 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1177 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1178 AliError(Form("the following detectors were not found: %s",
1180 if (fStopOnError) return kFALSE;
1189 //_____________________________________________________________________________
1190 Bool_t AliSimulation::RunDigitization(const char* detectors,
1191 const char* excludeDetectors)
1193 // run the digitization and produce digits from sdigits
1195 AliCodeTimerAuto("",0)
1197 // initialize CDB storage, run number, set CDB lock
1199 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1202 delete AliRunLoader::Instance();
1207 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1208 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1209 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
1210 // manager->SetEmbeddingFlag(fEmbeddingFlag);
1211 manager->SetInputStream(0, fGAliceFileName.Data());
1212 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1213 const char* fileName = ((TObjString*)
1214 (fBkgrdFileNames->At(iStream-1)))->GetName();
1215 manager->SetInputStream(iStream, fileName);
1218 TString detStr = detectors;
1219 TString detExcl = excludeDetectors;
1220 manager->GetInputStream(0)->ImportgAlice();
1221 AliRunLoader* runLoader =
1222 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
1223 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1224 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1225 AliModule* det = (AliModule*) detArray->At(iDet);
1226 if (!det || !det->IsActive()) continue;
1227 if (IsSelected(det->GetName(), detStr) &&
1228 !IsSelected(det->GetName(), detExcl)) {
1229 AliDigitizer* digitizer = det->CreateDigitizer(manager);
1232 AliError(Form("no digitizer for %s", det->GetName()));
1233 if (fStopOnError) return kFALSE;
1235 digitizer->SetRegionOfInterest(fRegionOfInterest);
1240 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1241 AliError(Form("the following detectors were not found: %s",
1243 if (fStopOnError) return kFALSE;
1246 if (!manager->GetListOfTasks()->IsEmpty()) {
1247 AliInfo("executing digitization");
1256 //_____________________________________________________________________________
1257 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1259 // run the digitization and produce digits from hits
1261 AliCodeTimerAuto("",0)
1263 // initialize CDB storage, run number, set CDB lock
1265 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1268 AliRunLoader* runLoader = LoadRun("READ");
1269 if (!runLoader) return kFALSE;
1271 TString detStr = detectors;
1272 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1273 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1274 AliModule* det = (AliModule*) detArray->At(iDet);
1275 if (!det || !det->IsActive()) continue;
1276 if (IsSelected(det->GetName(), detStr)) {
1277 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1282 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1283 AliError(Form("the following detectors were not found: %s",
1285 if (fStopOnError) return kFALSE;
1291 //_____________________________________________________________________________
1292 Bool_t AliSimulation::WriteRawData(const char* detectors,
1293 const char* fileName,
1294 Bool_t deleteIntermediateFiles,
1297 // convert the digits to raw data
1298 // First DDL raw data files for the given detectors are created.
1299 // If a file name is given, the DDL files are then converted to a DATE file.
1300 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1302 // If the file name has the extension ".root", the DATE file is converted
1304 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1305 // 'selrawdata' flag can be used to enable writing of detectors raw data
1306 // accoring to the trigger cluster.
1308 AliCodeTimerAuto("",0)
1309 AliSysInfo::AddStamp("WriteRawData_Start");
1311 TString detStr = detectors;
1312 if (!WriteRawFiles(detStr.Data())) {
1313 if (fStopOnError) return kFALSE;
1315 AliSysInfo::AddStamp("WriteRawFiles");
1317 // run HLT simulation on simulated DDL raw files
1318 // and produce HLT ddl raw files to be included in date/root file
1319 // bugfix 2009-06-26: the decision whether to write HLT raw data
1320 // is taken in RunHLT. Here HLT always needs to be run in order to
1321 // create HLT digits, unless its switched off. This is due to the
1322 // special placement of the HLT between the generation of DDL files
1323 // and conversion to DATE/Root file.
1324 detStr.ReplaceAll("HLT", "");
1325 if (!fRunHLT.IsNull()) {
1327 if (fStopOnError) return kFALSE;
1330 AliSysInfo::AddStamp("WriteRawData_RunHLT");
1332 TString dateFileName(fileName);
1333 if (!dateFileName.IsNull()) {
1334 Bool_t rootOutput = dateFileName.EndsWith(".root");
1335 if (rootOutput) dateFileName += ".date";
1336 TString selDateFileName;
1338 selDateFileName = "selected.";
1339 selDateFileName+= dateFileName;
1341 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1342 if (fStopOnError) return kFALSE;
1344 AliSysInfo::AddStamp("ConvertRawFilesToDate");
1345 if (deleteIntermediateFiles) {
1346 AliRunLoader* runLoader = LoadRun("READ");
1347 if (runLoader) for (Int_t iEvent = 0;
1348 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1350 sprintf(command, "rm -r raw%d", iEvent);
1351 gSystem->Exec(command);
1356 if (!ConvertDateToRoot(dateFileName, fileName)) {
1357 if (fStopOnError) return kFALSE;
1359 AliSysInfo::AddStamp("ConvertDateToRoot");
1360 if (deleteIntermediateFiles) {
1361 gSystem->Unlink(dateFileName);
1364 TString selFileName = "selected.";
1365 selFileName += fileName;
1366 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1367 if (fStopOnError) return kFALSE;
1369 if (deleteIntermediateFiles) {
1370 gSystem->Unlink(selDateFileName);
1379 //_____________________________________________________________________________
1380 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1382 // convert the digits to raw data DDL files
1384 AliCodeTimerAuto("",0)
1386 AliRunLoader* runLoader = LoadRun("READ");
1387 if (!runLoader) return kFALSE;
1389 // write raw data to DDL files
1390 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1391 AliInfo(Form("processing event %d", iEvent));
1392 runLoader->GetEvent(iEvent);
1393 TString baseDir = gSystem->WorkingDirectory();
1395 sprintf(dirName, "raw%d", iEvent);
1396 gSystem->MakeDirectory(dirName);
1397 if (!gSystem->ChangeDirectory(dirName)) {
1398 AliError(Form("couldn't change to directory %s", dirName));
1399 if (fStopOnError) return kFALSE; else continue;
1402 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1405 TString detStr = detectors;
1406 if (IsSelected("HLT", detStr)) {
1407 // Do nothing. "HLT" will be removed from detStr and HLT raw
1408 // data files are generated in RunHLT.
1411 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1412 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1413 AliModule* det = (AliModule*) detArray->At(iDet);
1414 if (!det || !det->IsActive()) continue;
1415 if (IsSelected(det->GetName(), detStr)) {
1416 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1421 if (!WriteTriggerRawData())
1422 if (fStopOnError) return kFALSE;
1424 gSystem->ChangeDirectory(baseDir);
1425 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1426 AliError(Form("the following detectors were not found: %s",
1428 if (fStopOnError) return kFALSE;
1437 //_____________________________________________________________________________
1438 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1439 const char* selDateFileName)
1441 // convert raw data DDL files to a DATE file with the program "dateStream"
1442 // The second argument is not empty when the user decides to write
1443 // the detectors raw data according to the trigger cluster.
1445 AliCodeTimerAuto("",0)
1447 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1449 AliError("the program dateStream was not found");
1450 if (fStopOnError) return kFALSE;
1455 AliRunLoader* runLoader = LoadRun("READ");
1456 if (!runLoader) return kFALSE;
1458 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1459 Bool_t selrawdata = kFALSE;
1460 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1463 // Note the option -s. It is used in order to avoid
1464 // the generation of SOR/EOR events.
1465 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1466 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1467 FILE* pipe = gSystem->OpenPipe(command, "w");
1470 AliError(Form("Cannot execute command: %s",command));
1474 Int_t selEvents = 0;
1475 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1477 UInt_t detectorPattern = 0;
1478 runLoader->GetEvent(iEvent);
1479 if (!runLoader->LoadTrigger()) {
1480 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1481 detectorPattern = aCTP->GetClusterMask();
1482 // Check if the event was triggered by CTP
1484 if (aCTP->GetClassMask()) selEvents++;
1488 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1490 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1491 selrawdata = kFALSE;
1495 fprintf(pipe, "GDC DetectorPattern %u\n", detectorPattern);
1499 // loop over detectors and DDLs
1500 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1501 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1503 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1504 Int_t ldcID = Int_t(ldc + 0.0001);
1505 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1507 char rawFileName[256];
1508 sprintf(rawFileName, "raw%d/%s",
1509 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1511 // check existence and size of raw data file
1512 FILE* file = fopen(rawFileName, "rb");
1513 if (!file) continue;
1514 fseek(file, 0, SEEK_END);
1515 unsigned long size = ftell(file);
1517 if (!size) continue;
1519 if (ldcID != prevLDC) {
1520 fprintf(pipe, " LDC Id %d\n", ldcID);
1523 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1528 Int_t result = gSystem->ClosePipe(pipe);
1530 if (!(selrawdata && selEvents > 0)) {
1532 return (result == 0);
1535 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1537 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1538 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1539 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1541 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1543 // Get the trigger decision and cluster
1544 UInt_t detectorPattern = 0;
1546 runLoader->GetEvent(iEvent);
1547 if (!runLoader->LoadTrigger()) {
1548 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1549 if (aCTP->GetClassMask() == 0) continue;
1550 detectorPattern = aCTP->GetClusterMask();
1551 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1552 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1555 fprintf(pipe2, "GDC DetectorPattern %u\n", detectorPattern);
1559 // loop over detectors and DDLs
1560 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1561 // Write only raw data from detectors that
1562 // are contained in the trigger cluster(s)
1563 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1565 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1567 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1568 Int_t ldcID = Int_t(ldc + 0.0001);
1569 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1571 char rawFileName[256];
1572 sprintf(rawFileName, "raw%d/%s",
1573 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1575 // check existence and size of raw data file
1576 FILE* file = fopen(rawFileName, "rb");
1577 if (!file) continue;
1578 fseek(file, 0, SEEK_END);
1579 unsigned long size = ftell(file);
1581 if (!size) continue;
1583 if (ldcID != prevLDC) {
1584 fprintf(pipe2, " LDC Id %d\n", ldcID);
1587 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1592 Int_t result2 = gSystem->ClosePipe(pipe2);
1595 return ((result == 0) && (result2 == 0));
1598 //_____________________________________________________________________________
1599 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1600 const char* rootFileName)
1602 // convert a DATE file to a root file with the program "alimdc"
1605 const Int_t kDBSize = 2000000000;
1606 const Int_t kTagDBSize = 1000000000;
1607 const Bool_t kFilter = kFALSE;
1608 const Int_t kCompression = 1;
1610 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1612 AliError("the program alimdc was not found");
1613 if (fStopOnError) return kFALSE;
1618 AliInfo(Form("converting DATE file %s to root file %s",
1619 dateFileName, rootFileName));
1621 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1622 const char* tagDBFS = "/tmp/mdc1/tags";
1624 // User defined file system locations
1625 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1626 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1627 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1628 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1629 if (gSystem->Getenv("ALIMDC_TAGDB"))
1630 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1632 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1633 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1634 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1636 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1637 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1638 gSystem->Exec(Form("mkdir %s",tagDBFS));
1640 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1641 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1642 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1644 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1645 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1646 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1648 return (result == 0);
1652 //_____________________________________________________________________________
1653 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1655 // delete existing run loaders, open a new one and load gAlice
1657 delete AliRunLoader::Instance();
1658 AliRunLoader* runLoader =
1659 AliRunLoader::Open(fGAliceFileName.Data(),
1660 AliConfig::GetDefaultEventFolderName(), mode);
1662 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1665 runLoader->LoadgAlice();
1666 runLoader->LoadHeader();
1667 gAlice = runLoader->GetAliRun();
1669 AliError(Form("no gAlice object found in file %s",
1670 fGAliceFileName.Data()));
1676 //_____________________________________________________________________________
1677 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1679 // get or calculate the number of signal events per background event
1681 if (!fBkgrdFileNames) return 1;
1682 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1683 if (nBkgrdFiles == 0) return 1;
1685 // get the number of signal events
1687 AliRunLoader* runLoader =
1688 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1689 if (!runLoader) return 1;
1691 nEvents = runLoader->GetNumberOfEvents();
1696 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1697 // get the number of background events
1698 const char* fileName = ((TObjString*)
1699 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1700 AliRunLoader* runLoader =
1701 AliRunLoader::Open(fileName, "BKGRD");
1702 if (!runLoader) continue;
1703 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1706 // get or calculate the number of signal per background events
1707 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1708 if (nSignalPerBkgrd <= 0) {
1709 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1710 } else if (result && (result != nSignalPerBkgrd)) {
1711 AliInfo(Form("the number of signal events per background event "
1712 "will be changed from %d to %d for stream %d",
1713 nSignalPerBkgrd, result, iBkgrdFile+1));
1714 nSignalPerBkgrd = result;
1717 if (!result) result = nSignalPerBkgrd;
1718 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1719 AliWarning(Form("not enough background events (%d) for %d signal events "
1720 "using %d signal per background events for stream %d",
1721 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1728 //_____________________________________________________________________________
1729 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1731 // check whether detName is contained in detectors
1732 // if yes, it is removed from detectors
1734 // check if all detectors are selected
1735 if ((detectors.CompareTo("ALL") == 0) ||
1736 detectors.BeginsWith("ALL ") ||
1737 detectors.EndsWith(" ALL") ||
1738 detectors.Contains(" ALL ")) {
1743 // search for the given detector
1744 Bool_t result = kFALSE;
1745 if ((detectors.CompareTo(detName) == 0) ||
1746 detectors.BeginsWith(detName+" ") ||
1747 detectors.EndsWith(" "+detName) ||
1748 detectors.Contains(" "+detName+" ")) {
1749 detectors.ReplaceAll(detName, "");
1753 // clean up the detectors string
1754 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1755 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1756 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1761 //_____________________________________________________________________________
1762 Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName)
1765 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1766 // These can be used for embedding of MC tracks into RAW data using the standard
1767 // merging procedure.
1769 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1772 AliError("no gAlice object. Restart aliroot and try again.");
1775 if (gAlice->Modules()->GetEntries() > 0) {
1776 AliError("gAlice was already run. Restart aliroot and try again.");
1780 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1784 gROOT->LoadMacro(fConfigFileName.Data());
1785 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1787 if(AliCDBManager::Instance()->GetRun() >= 0) {
1788 SetRunNumber(AliCDBManager::Instance()->GetRun());
1790 AliWarning("Run number not initialized!!");
1793 AliRunLoader::Instance()->CdGAFile();
1795 AliPDG::AddParticlesToPdgDataBase();
1797 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1799 gAlice->GetMCApp()->Init();
1801 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1802 gAlice->InitLoaders();
1803 AliRunLoader::Instance()->MakeTree("E");
1804 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1805 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1806 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1808 // Save stuff at the beginning of the file to avoid file corruption
1809 AliRunLoader::Instance()->CdGAFile();
1814 //AliCDBManager* man = AliCDBManager::Instance();
1815 //man->SetRun(0); // Should this come from rawdata header ?
1819 // Get the runloader
1820 AliRunLoader* runLoader = AliRunLoader::Instance();
1822 // Open esd file if available
1823 TFile* esdFile = TFile::Open(esdFileName);
1825 AliESDEvent* esd = new AliESDEvent();
1826 esdFile->GetObject("esdTree", treeESD);
1827 if (treeESD) esd->ReadFromTree(treeESD);
1830 // Create the RawReader
1831 TString fileName(rawDirectory);
1832 AliRawReader* rawReader = 0x0;
1833 if (fileName.EndsWith("/")) {
1834 rawReader = new AliRawReaderFile(fileName);
1835 } else if (fileName.EndsWith(".root")) {
1836 rawReader = new AliRawReaderRoot(fileName);
1837 } else if (!fileName.IsNull()) {
1838 rawReader = new AliRawReaderDate(fileName);
1840 // if (!fEquipIdMap.IsNull() && fRawReader)
1841 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1843 // Get list of detectors
1844 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1847 AliHeader* header = runLoader->GetHeader();
1849 TString detStr = fMakeSDigits;
1853 if (!(rawReader->NextEvent())) break;
1856 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1857 AliModule* det = (AliModule*) detArray->At(iDet);
1858 if (!det || !det->IsActive()) continue;
1859 if (IsSelected(det->GetName(), detStr)) {
1860 AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1861 det->Raw2SDigits(rawReader);
1868 // If ESD information available obtain reconstructed vertex and store in header.
1870 treeESD->GetEvent(nev);
1871 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1872 Double_t position[3];
1873 esdVertex->GetXYZ(position);
1874 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1877 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1878 mcHeader->SetPrimaryVertex(mcV);
1879 header->Reset(0,nev);
1880 header->SetGenEventHeader(mcHeader);
1881 printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1886 runLoader->TreeE()->Fill();
1887 runLoader->SetNextEvent();
1893 runLoader->CdGAFile();
1894 runLoader->WriteHeader("OVERWRITE");
1895 runLoader->WriteRunLoader();
1900 //_____________________________________________________________________________
1901 void AliSimulation::FinishRun()
1904 // Called at the end of the run.
1909 AliDebug(1, "Finish Lego");
1910 AliRunLoader::Instance()->CdGAFile();
1914 // Clean detector information
1915 TIter next(gAlice->Modules());
1916 AliModule *detector;
1917 while((detector = dynamic_cast<AliModule*>(next()))) {
1918 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
1919 detector->FinishRun();
1922 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
1923 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
1925 // Write AliRun info and all detectors parameters
1926 AliRunLoader::Instance()->CdGAFile();
1927 gAlice->Write(0,TObject::kOverwrite);//write AliRun
1928 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
1930 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
1931 AliRunLoader::Instance()->Synchronize();
1934 //_____________________________________________________________________________
1935 Int_t AliSimulation::GetDetIndex(const char* detector)
1937 // return the detector index corresponding to detector
1939 for (index = 0; index < fgkNDetectors ; index++) {
1940 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
1946 //_____________________________________________________________________________
1947 Bool_t AliSimulation::CreateHLT()
1949 // Init the HLT simulation.
1950 // The function loads the library and creates the instance of AliHLTSimulation.
1951 // the main reason for the decoupled creation is to set the transient OCDB
1952 // objects before the OCDB is locked
1954 // load the library dynamically
1955 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
1957 // check for the library version
1958 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
1960 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
1963 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
1964 AliWarning(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
1967 // print compile info
1968 typedef void (*CompileInfo)( const char*& date, const char*& time);
1969 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
1971 const char* date="";
1972 const char* time="";
1973 (*fctInfo)(date, time);
1974 if (!date) date="unknown";
1975 if (!time) time="unknown";
1976 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
1978 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
1981 // create instance of the HLT simulation
1982 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
1983 if (fctCreate==NULL || (fpHLT=(fctCreate()))==NULL) {
1984 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
1988 TString specObjects;
1989 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
1990 if (specObjects.Length()>0) specObjects+=" ";
1991 specObjects+=fSpecCDBUri[i]->GetName();
1994 AliHLTSimulationSetup_t fctSetup=(AliHLTSimulationSetup_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_SETUP));
1995 if (fctSetup==NULL || fctSetup(fpHLT, this, specObjects.Data())<0) {
1996 AliWarning(Form("failed to setup HLT simulation (function %p)", fctSetup));
2002 //_____________________________________________________________________________
2003 Bool_t AliSimulation::RunHLT()
2005 // Run the HLT simulation
2006 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
2007 // Disabled if fRunHLT is empty, default vaule is "default".
2008 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
2009 // The default simulation depends on the HLT component libraries and their
2010 // corresponding agents which define components and chains to run. See
2011 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
2012 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
2014 // The libraries to be loaded can be specified as an option.
2016 // AliSimulation sim;
2017 // sim.SetRunHLT("libAliHLTSample.so");
2019 // will only load <tt>libAliHLTSample.so</tt>
2021 // Other available options:
2022 // \li loglevel=<i>level</i> <br>
2023 // logging level for this processing
2025 // disable redirection of log messages to AliLog class
2026 // \li config=<i>macro</i>
2027 // configuration macro
2028 // \li chains=<i>configuration</i>
2029 // comma separated list of configurations to be run during simulation
2030 // \li rawfile=<i>file</i>
2031 // source for the RawReader to be created, the default is <i>./</i> if
2032 // raw data is simulated
2036 if (!fpHLT && !CreateHLT()) {
2039 AliHLTSimulation* pHLT=fpHLT;
2041 AliRunLoader* pRunLoader = LoadRun("READ");
2042 if (!pRunLoader) return kFALSE;
2044 // initialize CDB storage, run number, set CDB lock
2045 // thats for the case of running HLT simulation without all the other steps
2046 // multiple calls are handled by the function, so we can just call
2048 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
2051 // init the HLT simulation
2053 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
2054 TString detStr = fWriteRawData;
2055 if (!IsSelected("HLT", detStr)) {
2056 options+=" writerawfiles=";
2058 options+=" writerawfiles=HLT";
2061 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
2062 // as a matter of fact, HLT will run reconstruction and needs the RawReader
2063 // in order to get detector data. By default, RawReaderFile is used to read
2064 // the already simulated ddl files. Date and Root files from the raw data
2065 // are generated after the HLT simulation.
2066 options+=" rawfile=./";
2069 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
2070 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
2071 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
2073 // run the HLT simulation
2074 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
2075 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
2076 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
2080 // delete the instance
2081 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
2082 if (fctDelete==NULL || fctDelete(pHLT)<0) {
2083 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
2087 return iResult>=0?kTRUE:kFALSE;
2090 //_____________________________________________________________________________
2091 Bool_t AliSimulation::RunQA()
2093 // run the QA on summable hits, digits or digits
2095 if(!gAlice) return kFALSE;
2096 AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
2098 TString detectorsw("") ;
2100 AliQAManager::QAManager()->SetEventSpecie(fEventSpecie) ;
2101 detectorsw = AliQAManager::QAManager()->Run(fQADetectors.Data()) ;
2102 if ( detectorsw.IsNull() )
2107 //_____________________________________________________________________________
2108 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2110 // Allows to run QA for a selected set of detectors
2111 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2112 // all selected detectors run the same selected tasks
2114 if (!detAndAction.Contains(":")) {
2115 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2119 Int_t colon = detAndAction.Index(":") ;
2120 fQADetectors = detAndAction(0, colon) ;
2121 if (fQADetectors.Contains("ALL") ){
2122 TString tmp = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2123 Int_t minus = fQADetectors.Last('-') ;
2124 TString toKeep = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2125 TString toRemove("") ;
2126 while (minus >= 0) {
2127 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
2128 toRemove = toRemove.Strip() ;
2129 toKeep.ReplaceAll(toRemove, "") ;
2130 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
2131 minus = fQADetectors.Last('-') ;
2133 fQADetectors = toKeep ;
2135 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2136 if (fQATasks.Contains("ALL") ) {
2137 fQATasks = Form("%d %d %d", AliQAv1::kHITS, AliQAv1::kSDIGITS, AliQAv1::kDIGITS) ;
2139 fQATasks.ToUpper() ;
2141 if ( fQATasks.Contains("HIT") )
2142 tempo = Form("%d ", AliQAv1::kHITS) ;
2143 if ( fQATasks.Contains("SDIGIT") )
2144 tempo += Form("%d ", AliQAv1::kSDIGITS) ;
2145 if ( fQATasks.Contains("DIGIT") )
2146 tempo += Form("%d ", AliQAv1::kDIGITS) ;
2148 if (fQATasks.IsNull()) {
2149 AliInfo("No QA requested\n") ;
2154 TString tempo(fQATasks) ;
2155 tempo.ReplaceAll(Form("%d", AliQAv1::kHITS), AliQAv1::GetTaskName(AliQAv1::kHITS)) ;
2156 tempo.ReplaceAll(Form("%d", AliQAv1::kSDIGITS), AliQAv1::GetTaskName(AliQAv1::kSDIGITS)) ;
2157 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITS), AliQAv1::GetTaskName(AliQAv1::kDIGITS)) ;
2158 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2160 AliQAManager::QAManager()->SetActiveDetectors(fQADetectors) ;
2161 AliQAManager::QAManager()->SetTasks(fQATasks) ;
2162 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++)
2163 AliQAManager::QAManager()->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
2168 //_____________________________________________________________________________
2169 void AliSimulation::ProcessEnvironmentVars()
2171 // Extract run number and random generator seed from env variables
2173 AliInfo("Processing environment variables");
2175 // Random Number seed
2177 // first check that seed is not already set
2179 if (gSystem->Getenv("CONFIG_SEED")) {
2180 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2183 if (gSystem->Getenv("CONFIG_SEED")) {
2184 AliInfo(Form("Seed for random number generation already set (%d)"
2185 ": CONFIG_SEED variable ignored!", fSeed));
2189 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2193 // first check that run number is not already set
2195 if (gSystem->Getenv("DC_RUN")) {
2196 fRun = atoi(gSystem->Getenv("DC_RUN"));
2199 if (gSystem->Getenv("DC_RUN")) {
2200 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2204 AliInfo(Form("Run number = %d", fRun));
2207 //---------------------------------------------------------------------
2208 void AliSimulation::WriteGRPEntry()
2210 // Get the necessary information from galice (generator, trigger etc) and
2211 // write a GRP entry corresponding to the settings in the Config.C used
2212 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2215 AliInfo("Writing global run parameters entry into the OCDB");
2217 AliGRPObject* grpObj = new AliGRPObject();
2219 grpObj->SetRunType("PHYSICS");
2220 grpObj->SetTimeStart(0);
2222 grpObj->SetTimeStart(0);
2223 grpObj->SetTimeEnd(curtime.Convert());
2224 grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
2226 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2228 grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
2231 gen->GetProjectile(projectile,a,z);
2233 gen->GetTarget(target,a,z);
2234 TString beamType = projectile + "-" + target;
2235 beamType.ReplaceAll(" ","");
2236 if (!beamType.CompareTo("-")) {
2237 grpObj->SetBeamType("UNKNOWN");
2240 grpObj->SetBeamType(beamType);
2241 // Heavy ion run, the event specie is set to kHighMult
2242 fEventSpecie = AliRecoParam::kHighMult;
2243 if ((strcmp(beamType,"p-p") == 0) ||
2244 (strcmp(beamType,"p-") == 0) ||
2245 (strcmp(beamType,"-p") == 0) ||
2246 (strcmp(beamType,"P-P") == 0) ||
2247 (strcmp(beamType,"P-") == 0) ||
2248 (strcmp(beamType,"-P") == 0)) {
2249 // Proton run, the event specie is set to kLowMult
2250 fEventSpecie = AliRecoParam::kLowMult;
2254 AliWarning("Unknown beam type and energy! Setting energy to 0");
2255 grpObj->SetBeamEnergy(0);
2256 grpObj->SetBeamType("UNKNOWN");
2259 UInt_t detectorPattern = 0;
2261 TObjArray *detArray = gAlice->Detectors();
2262 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2263 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2264 detectorPattern |= (1 << iDet);
2269 if (!fMakeTrigger.IsNull() || strcmp(gAlice->GetTriggerDescriptor(),""))
2270 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2273 if (!fRunHLT.IsNull())
2274 detectorPattern |= (1 << AliDAQ::kHLTId);
2276 grpObj->SetNumberOfDetectors((Char_t)nDets);
2277 grpObj->SetDetectorMask((Int_t)detectorPattern);
2278 grpObj->SetLHCPeriod("LHC08c");
2279 grpObj->SetLHCState("STABLE_BEAMS");
2281 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2282 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2284 Float_t factorSol = field ? field->GetFactorSol() : 0;
2285 Float_t currentSol = TMath::Abs(factorSol)>1E-6 ?
2286 TMath::Nint(TMath::Abs(solenoidField/factorSol))/5.*30000.*TMath::Abs(factorSol) : 0;
2288 Float_t factorDip = field ? field->GetFactorDip() : 0;
2289 Float_t currentDip = 6000.*TMath::Abs(factorDip);
2291 grpObj->SetL3Current(currentSol,(AliGRPObject::Stats)0);
2292 grpObj->SetDipoleCurrent(currentDip,(AliGRPObject::Stats)0);
2293 grpObj->SetL3Polarity(factorSol>0 ? 0:1);
2294 grpObj->SetDipolePolarity(factorDip>0 ? 0:1);
2295 grpObj->SetUniformBMap(field->IsUniform()); // for special MC with k5kGUniform map
2296 grpObj->SetPolarityConventionLHC(); // LHC convention +/+ current -> -/- field main components
2298 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2300 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2302 // Now store the entry in OCDB
2303 AliCDBManager* man = AliCDBManager::Instance();
2305 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2306 AliCDBMetaData *metadata= new AliCDBMetaData();
2308 metadata->SetResponsible("alice-off@cern.ch");
2309 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2311 man->Put(grpObj,id,metadata);