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 <TGeoManager.h>
111 #include <TObjString.h>
114 #include <TVirtualMC.h>
115 #include <TVirtualMCApplication.h>
117 #include "AliAlignObj.h"
118 #include "AliCDBEntry.h"
119 #include "AliCDBManager.h"
120 #include "AliCDBStorage.h"
121 #include "AliCTPRawData.h"
122 #include "AliCentralTrigger.h"
123 #include "AliCentralTrigger.h"
124 #include "AliCodeTimer.h"
126 #include "AliDigitizer.h"
128 #include "AliGRPObject.h"
129 #include "AliGenEventHeader.h"
130 #include "AliGenerator.h"
131 #include "AliGeomManager.h"
132 #include "AliHLTSimulation.h"
133 #include "AliHeader.h"
137 #include "AliModule.h"
139 #include "AliRawReaderDate.h"
140 #include "AliRawReaderFile.h"
141 #include "AliRawReaderRoot.h"
143 #include "AliRunDigitizer.h"
144 #include "AliRunLoader.h"
145 #include "AliSimulation.h"
146 #include "AliSysInfo.h"
147 #include "AliVertexGenFile.h"
148 #include "AliLegoGenerator.h"
151 ClassImp(AliSimulation)
153 AliSimulation *AliSimulation::fgInstance = 0;
154 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
156 //_____________________________________________________________________________
157 AliSimulation::AliSimulation(const char* configFileName,
158 const char* name, const char* title) :
161 fRunGeneration(kTRUE),
162 fRunSimulation(kTRUE),
163 fLoadAlignFromCDB(kTRUE),
164 fLoadAlObjsListOfDets("ALL"),
168 fMakeDigitsFromHits(""),
170 fRawDataFileName(""),
171 fDeleteIntermediateFiles(kFALSE),
172 fWriteSelRawData(kFALSE),
173 fStopOnError(kFALSE),
175 fConfigFileName(configFileName),
176 fGAliceFileName("galice.root"),
178 fBkgrdFileNames(NULL),
179 fAlignObjArray(NULL),
180 fUseBkgrdVertex(kTRUE),
181 fRegionOfInterest(kFALSE),
186 fInitCDBCalled(kFALSE),
187 fInitRunNumberCalled(kFALSE),
188 fSetRunNumberFromDataCalled(kFALSE),
189 fEmbeddingFlag(kFALSE),
195 fEventSpecie(AliRecoParam::kDefault),
197 fWriteGRPEntry(kTRUE)
199 // create simulation object with default parameters
201 SetGAliceFile("galice.root");
204 fQASteer = new AliQADataMakerSteer("sim") ;
205 fQASteer->SetActiveDetectors(fQADetectors) ;
206 fQATasks = Form("%d %d %d", AliQA::kHITS, AliQA::kSDIGITS, AliQA::kDIGITS) ;
207 fQASteer->SetTasks(fQATasks) ;
210 //_____________________________________________________________________________
211 AliSimulation::~AliSimulation()
215 fEventsPerFile.Delete();
216 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
217 // delete fAlignObjArray; fAlignObjArray=0;
219 if (fBkgrdFileNames) {
220 fBkgrdFileNames->Delete();
221 delete fBkgrdFileNames;
224 fSpecCDBUri.Delete();
225 if (fgInstance==this) fgInstance = 0;
229 AliCodeTimer::Instance()->Print();
233 //_____________________________________________________________________________
234 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
236 // set the number of events for one run
241 //_____________________________________________________________________________
242 void AliSimulation::InitCDB()
244 // activate a default CDB storage
245 // First check if we have any CDB storage set, because it is used
246 // to retrieve the calibration and alignment constants
248 if (fInitCDBCalled) return;
249 fInitCDBCalled = kTRUE;
251 AliCDBManager* man = AliCDBManager::Instance();
252 if (man->IsDefaultStorageSet())
254 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
255 AliWarning("Default CDB storage has been already set !");
256 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
257 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
258 fCDBUri = man->GetDefaultStorage()->GetURI();
261 if (fCDBUri.Length() > 0)
263 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
264 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
265 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
267 fCDBUri="local://$ALICE_ROOT";
268 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
269 AliWarning("Default CDB storage not yet set !!!!");
270 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
271 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
274 man->SetDefaultStorage(fCDBUri);
277 // Now activate the detector specific CDB storage locations
278 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
279 TObject* obj = fSpecCDBUri[i];
281 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
282 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
283 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
284 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
289 //_____________________________________________________________________________
290 void AliSimulation::InitRunNumber(){
291 // check run number. If not set, set it to 0 !!!!
293 if (fInitRunNumberCalled) return;
294 fInitRunNumberCalled = kTRUE;
296 AliCDBManager* man = AliCDBManager::Instance();
297 if (man->GetRun() >= 0)
299 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
300 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
304 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
305 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
306 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
309 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
310 AliWarning("Run number not yet set !!!!");
311 AliWarning(Form("Setting it now to: %d", fRun));
312 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
321 //_____________________________________________________________________________
322 void AliSimulation::SetCDBLock() {
323 // Set CDB lock: from now on it is forbidden to reset the run number
324 // or the default storage or to activate any further storage!
326 AliCDBManager::Instance()->SetLock(1);
329 //_____________________________________________________________________________
330 void AliSimulation::SetDefaultStorage(const char* uri) {
331 // Store the desired default CDB storage location
332 // Activate it later within the Run() method
338 //_____________________________________________________________________________
339 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
340 // Store a detector-specific CDB storage location
341 // Activate it later within the Run() method
343 AliCDBPath aPath(calibType);
344 if(!aPath.IsValid()){
345 AliError(Form("Not a valid path: %s", calibType));
349 TObject* obj = fSpecCDBUri.FindObject(calibType);
350 if (obj) fSpecCDBUri.Remove(obj);
351 fSpecCDBUri.Add(new TNamed(calibType, uri));
355 //_____________________________________________________________________________
356 void AliSimulation::SetRunNumber(Int_t run)
359 // Activate it later within the Run() method
364 //_____________________________________________________________________________
365 void AliSimulation::SetSeed(Int_t seed)
368 // Activate it later within the Run() method
373 //_____________________________________________________________________________
374 Bool_t AliSimulation::SetRunNumberFromData()
376 // Set the CDB manager run number
377 // The run number is retrieved from gAlice
379 if (fSetRunNumberFromDataCalled) return kTRUE;
380 fSetRunNumberFromDataCalled = kTRUE;
382 AliCDBManager* man = AliCDBManager::Instance();
383 Int_t runData = -1, runCDB = -1;
385 AliRunLoader* runLoader = LoadRun("READ");
386 if (!runLoader) return kFALSE;
388 runData = runLoader->GetHeader()->GetRun();
392 runCDB = man->GetRun();
394 if (runCDB != runData) {
395 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
396 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
397 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
398 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
403 man->SetRun(runData);
406 if(man->GetRun() < 0) {
407 AliError("Run number not properly initalized!");
416 //_____________________________________________________________________________
417 void AliSimulation::SetConfigFile(const char* fileName)
419 // set the name of the config file
421 fConfigFileName = fileName;
424 //_____________________________________________________________________________
425 void AliSimulation::SetGAliceFile(const char* fileName)
427 // set the name of the galice file
428 // the path is converted to an absolute one if it is relative
430 fGAliceFileName = fileName;
431 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
432 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
434 fGAliceFileName = absFileName;
435 delete[] absFileName;
438 AliDebug(2, Form("galice file name set to %s", fileName));
441 //_____________________________________________________________________________
442 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
445 // set the number of events per file for the given detector and data type
446 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
448 TNamed* obj = new TNamed(detector, type);
449 obj->SetUniqueID(nEvents);
450 fEventsPerFile.Add(obj);
453 //_____________________________________________________________________________
454 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
456 // Read the alignment objects from CDB.
457 // Each detector is supposed to have the
458 // alignment objects in DET/Align/Data CDB path.
459 // All the detector objects are then collected,
460 // sorted by geometry level (starting from ALIC) and
461 // then applied to the TGeo geometry.
462 // Finally an overlaps check is performed.
464 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
465 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
469 // initialize CDB storage, run number, set CDB lock
471 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
474 Bool_t delRunLoader = kFALSE;
476 runLoader = LoadRun("READ");
477 if (!runLoader) return kFALSE;
478 delRunLoader = kTRUE;
481 // Export ideal geometry
482 if(!gAlice->IsRootGeometry()) AliGeomManager::GetGeometry()->Export("geometry.root");
484 // Load alignment data from CDB and apply to geometry through AliGeomManager
485 if(fLoadAlignFromCDB){
487 TString detStr = fLoadAlObjsListOfDets;
488 TString loadAlObjsListOfDets = "";
490 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
491 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
492 AliModule* det = (AliModule*) detArray->At(iDet);
493 if (!det || !det->IsActive()) continue;
494 if (IsSelected(det->GetName(), detStr)) {
495 //add det to list of dets to be aligned from CDB
496 loadAlObjsListOfDets += det->GetName();
497 loadAlObjsListOfDets += " ";
499 } // end loop over detectors
500 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
501 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
503 // Check if the array with alignment objects was
504 // provided by the user. If yes, apply the objects
505 // to the present TGeo geometry
506 if (fAlignObjArray) {
507 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
508 AliError("The misalignment of one or more volumes failed!"
509 "Compare the list of simulated detectors and the list of detector alignment data!");
510 if (delRunLoader) delete runLoader;
516 // Update the internal geometry of modules (ITS needs it)
517 TString detStr = fLoadAlObjsListOfDets;
518 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
519 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
521 AliModule* det = (AliModule*) detArray->At(iDet);
522 if (!det || !det->IsActive()) continue;
523 if (IsSelected(det->GetName(), detStr)) {
524 det->UpdateInternalGeometry();
526 } // end loop over detectors
529 if (delRunLoader) delete runLoader;
534 //_____________________________________________________________________________
535 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
537 // add a file with background events for merging
539 TObjString* fileNameStr = new TObjString(fileName);
540 fileNameStr->SetUniqueID(nSignalPerBkgrd);
541 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
542 fBkgrdFileNames->Add(fileNameStr);
545 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
547 // add a file with background events for embeddin
548 MergeWith(fileName, nSignalPerBkgrd);
549 fEmbeddingFlag = kTRUE;
552 //_____________________________________________________________________________
553 Bool_t AliSimulation::Run(Int_t nEvents)
555 // run the generation, simulation and digitization
560 // Load run number and seed from environmental vars
561 ProcessEnvironmentVars();
563 gRandom->SetSeed(fSeed);
565 if (nEvents > 0) fNEvents = nEvents;
567 // generation and simulation -> hits
568 if (fRunGeneration) {
569 if (!RunSimulation()) if (fStopOnError) return kFALSE;
572 // initialize CDB storage from external environment
573 // (either CDB manager or AliSimulation setters),
574 // if not already done in RunSimulation()
577 // Set run number in CDBManager from data
578 // From this point on the run number must be always loaded from data!
579 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
581 // Set CDB lock: from now on it is forbidden to reset the run number
582 // or the default storage or to activate any further storage!
585 // If RunSimulation was not called, load the geometry and misalign it
586 if (!AliGeomManager::GetGeometry()) {
587 // Initialize the geometry manager
588 AliGeomManager::LoadGeometry("geometry.root");
590 // // Check that the consistency of symbolic names for the activated subdetectors
591 // // in the geometry loaded by AliGeomManager
592 // AliRunLoader* runLoader = LoadRun("READ");
593 // if (!runLoader) return kFALSE;
595 // TString detsToBeChecked = "";
596 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
597 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
598 // AliModule* det = (AliModule*) detArray->At(iDet);
599 // if (!det || !det->IsActive()) continue;
600 // detsToBeChecked += det->GetName();
601 // detsToBeChecked += " ";
602 // } // end loop over detectors
603 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
604 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
605 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
607 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
609 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
613 // hits -> summable digits
614 AliSysInfo::AddStamp("Start_sdigitization");
615 if (!fMakeSDigits.IsNull()) {
616 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
619 AliSysInfo::AddStamp("Stop_sdigitization");
621 AliSysInfo::AddStamp("Start_digitization");
622 // summable digits -> digits
623 if (!fMakeDigits.IsNull()) {
624 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
625 if (fStopOnError) return kFALSE;
628 AliSysInfo::AddStamp("Stop_digitization");
633 if (!fMakeDigitsFromHits.IsNull()) {
634 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
635 AliWarning(Form("Merging and direct creation of digits from hits "
636 "was selected for some detectors. "
637 "No merging will be done for the following detectors: %s",
638 fMakeDigitsFromHits.Data()));
640 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
641 if (fStopOnError) return kFALSE;
648 if (!RunTrigger(fMakeTrigger,fMakeDigits)) {
649 if (fStopOnError) return kFALSE;
654 // digits -> raw data
655 if (!fWriteRawData.IsNull()) {
656 if (!WriteRawData(fWriteRawData, fRawDataFileName,
657 fDeleteIntermediateFiles,fWriteSelRawData)) {
658 if (fStopOnError) return kFALSE;
663 // run HLT simulation on simulated digit data if raw data is not
664 // simulated, otherwise its called as part of WriteRawData
665 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
667 if (fStopOnError) return kFALSE;
673 Bool_t rv = RunQA() ;
679 // Cleanup of CDB manager: cache and active storages!
680 AliCDBManager::Instance()->ClearCache();
685 //_______________________________________________________________________
686 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
687 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
688 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
691 // Generates lego plots of:
692 // - radiation length map phi vs theta
693 // - radiation length map phi vs eta
694 // - interaction length map
695 // - g/cm2 length map
697 // ntheta bins in theta, eta
698 // themin minimum angle in theta (degrees)
699 // themax maximum angle in theta (degrees)
701 // phimin minimum angle in phi (degrees)
702 // phimax maximum angle in phi (degrees)
703 // rmin minimum radius
704 // rmax maximum radius
707 // The number of events generated = ntheta*nphi
708 // run input parameters in macro setup (default="Config.C")
710 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
713 <img src="picts/AliRunLego1.gif">
718 <img src="picts/AliRunLego2.gif">
723 <img src="picts/AliRunLego3.gif">
728 // run the generation and simulation
732 // initialize CDB storage and run number from external environment
733 // (either CDB manager or AliSimulation setters)
739 AliError("no gAlice object. Restart aliroot and try again.");
742 if (gAlice->Modules()->GetEntries() > 0) {
743 AliError("gAlice was already run. Restart aliroot and try again.");
747 AliInfo(Form("initializing gAlice with config file %s",
748 fConfigFileName.Data()));
751 if (nev == -1) nev = nc1 * nc2;
753 // check if initialisation has been done
754 // If runloader has been initialized, set the number of events per file to nc1 * nc2
757 if (!gener) gener = new AliLegoGenerator();
759 // Configure Generator
761 gener->SetRadiusRange(rmin, rmax);
762 gener->SetZMax(zmax);
763 gener->SetCoor1Range(nc1, c1min, c1max);
764 gener->SetCoor2Range(nc2, c2min, c2max);
768 fLego = new AliLego("lego",gener);
770 //__________________________________________________________________________
774 gROOT->LoadMacro(setup);
775 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
777 if(AliCDBManager::Instance()->GetRun() >= 0) {
778 SetRunNumber(AliCDBManager::Instance()->GetRun());
780 AliWarning("Run number not initialized!!");
783 AliRunLoader::GetRunLoader()->CdGAFile();
785 AliPDG::AddParticlesToPdgDataBase();
787 gAlice->GetMCApp()->Init();
789 //Must be here because some MCs (G4) adds detectors here and not in Config.C
790 gAlice->InitLoaders();
791 AliRunLoader::GetRunLoader()->MakeTree("E");
794 // Save stuff at the beginning of the file to avoid file corruption
795 AliRunLoader::GetRunLoader()->CdGAFile();
798 //Save current generator
799 AliGenerator *gen=gAlice->GetMCApp()->Generator();
800 gAlice->GetMCApp()->ResetGenerator(gener);
801 //Prepare MC for Lego Run
807 AliRunLoader::GetRunLoader()->SetNumberOfEventsPerFile(nev);
808 gMC->ProcessRun(nev);
810 // End of this run, close files
812 // Restore current generator
813 gAlice->GetMCApp()->ResetGenerator(gen);
814 // Delete Lego Object
820 //_____________________________________________________________________________
821 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
827 // initialize CDB storage from external environment
828 // (either CDB manager or AliSimulation setters),
829 // if not already done in RunSimulation()
832 // Set run number in CDBManager from data
833 // From this point on the run number must be always loaded from data!
834 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
836 // Set CDB lock: from now on it is forbidden to reset the run number
837 // or the default storage or to activate any further storage!
840 AliRunLoader* runLoader = LoadRun("READ");
841 if (!runLoader) return kFALSE;
842 TString trconfiguration = config;
844 if (trconfiguration.IsNull()) {
845 if (strcmp(gAlice->GetTriggerDescriptor(),"")) {
846 trconfiguration = gAlice->GetTriggerDescriptor();
849 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
852 runLoader->MakeTree( "GG" );
853 AliCentralTrigger* aCTP = runLoader->GetTrigger();
854 // Load Configuration
855 if (!aCTP->LoadConfiguration( trconfiguration ))
859 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
871 //_____________________________________________________________________________
872 Bool_t AliSimulation::WriteTriggerRawData()
874 // Writes the CTP (trigger) DDL raw data
875 // Details of the format are given in the
876 // trigger TDR - pages 134 and 135.
877 AliCTPRawData writer;
883 //_____________________________________________________________________________
884 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
886 // run the generation and simulation
890 // initialize CDB storage and run number from external environment
891 // (either CDB manager or AliSimulation setters)
897 AliError("no gAlice object. Restart aliroot and try again.");
900 if (gAlice->Modules()->GetEntries() > 0) {
901 AliError("gAlice was already run. Restart aliroot and try again.");
905 AliInfo(Form("initializing gAlice with config file %s",
906 fConfigFileName.Data()));
909 // Initialize ALICE Simulation run
914 gROOT->LoadMacro(fConfigFileName.Data());
915 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
917 if(AliCDBManager::Instance()->GetRun() >= 0) {
918 gAlice->SetRunNumber(AliCDBManager::Instance()->GetRun());
920 AliWarning("Run number not initialized!!");
923 AliRunLoader::GetRunLoader()->CdGAFile();
925 AliPDG::AddParticlesToPdgDataBase();
927 gAlice->GetMCApp()->Init();
929 //Must be here because some MCs (G4) adds detectors here and not in Config.C
930 gAlice->InitLoaders();
931 AliRunLoader::GetRunLoader()->MakeTree("E");
932 AliRunLoader::GetRunLoader()->LoadKinematics("RECREATE");
933 AliRunLoader::GetRunLoader()->LoadTrackRefs("RECREATE");
934 AliRunLoader::GetRunLoader()->LoadHits("all","RECREATE");
936 // Save stuff at the beginning of the file to avoid file corruption
937 AliRunLoader::GetRunLoader()->CdGAFile();
939 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
940 //___________________________________________________________________________________________
942 // Get the trigger descriptor string
943 // Either from AliSimulation or from
945 if (fMakeTrigger.IsNull()) {
946 if (strcmp(gAlice->GetTriggerDescriptor(),""))
947 fMakeTrigger = gAlice->GetTriggerDescriptor();
950 gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
952 // Set run number in CDBManager
953 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
955 AliRunLoader* runLoader = AliRunLoader::GetRunLoader();
957 AliError(Form("gAlice has no run loader object. "
958 "Check your config file: %s", fConfigFileName.Data()));
961 SetGAliceFile(runLoader->GetFileName());
964 #if ROOT_VERSION_CODE < 331527
965 AliGeomManager::SetGeometry(gGeoManager);
967 // Check that the consistency of symbolic names for the activated subdetectors
968 // in the geometry loaded by AliGeomManager
969 TString detsToBeChecked = "";
970 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
971 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
972 AliModule* det = (AliModule*) detArray->At(iDet);
973 if (!det || !det->IsActive()) continue;
974 detsToBeChecked += det->GetName();
975 detsToBeChecked += " ";
976 } // end loop over detectors
977 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
978 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
979 MisalignGeometry(runLoader);
982 // AliRunLoader* runLoader = AliRunLoader::GetRunLoader();
984 // AliError(Form("gAlice has no run loader object. "
985 // "Check your config file: %s", fConfigFileName.Data()));
988 // SetGAliceFile(runLoader->GetFileName());
990 if (!gAlice->GetMCApp()->Generator()) {
991 AliError(Form("gAlice has no generator object. "
992 "Check your config file: %s", fConfigFileName.Data()));
996 // Write GRP entry corresponding to the setting found in Cofig.C
1000 if (nEvents <= 0) nEvents = fNEvents;
1002 // get vertex from background file in case of merging
1003 if (fUseBkgrdVertex &&
1004 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1005 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1006 const char* fileName = ((TObjString*)
1007 (fBkgrdFileNames->At(0)))->GetName();
1008 AliInfo(Form("The vertex will be taken from the background "
1009 "file %s with nSignalPerBackground = %d",
1010 fileName, signalPerBkgrd));
1011 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1012 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1015 if (!fRunSimulation) {
1016 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1019 // set the number of events per file for given detectors and data types
1020 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1021 if (!fEventsPerFile[i]) continue;
1022 const char* detName = fEventsPerFile[i]->GetName();
1023 const char* typeName = fEventsPerFile[i]->GetTitle();
1024 TString loaderName(detName);
1025 loaderName += "Loader";
1026 AliLoader* loader = runLoader->GetLoader(loaderName);
1028 AliError(Form("RunSimulation", "no loader for %s found\n"
1029 "Number of events per file not set for %s %s",
1030 detName, typeName, detName));
1033 AliDataLoader* dataLoader =
1034 loader->GetDataLoader(typeName);
1036 AliError(Form("no data loader for %s found\n"
1037 "Number of events per file not set for %s %s",
1038 typeName, detName, typeName));
1041 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1042 AliDebug(1, Form("number of events per file set to %d for %s %s",
1043 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1046 AliInfo("running gAlice");
1047 AliSysInfo::AddStamp("Start_simulation");
1049 // Create the Root Tree with one branch per detector
1050 //Hits moved to begin event -> now we are crating separate tree for each event
1052 gMC->ProcessRun(nEvents);
1054 // End of this run, close files
1055 if(nEvents>0) FinishRun();
1057 AliSysInfo::AddStamp("Stop_simulation");
1063 //_____________________________________________________________________________
1064 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1066 // run the digitization and produce summable digits
1067 static Int_t eventNr=0;
1068 AliCodeTimerAuto("")
1070 // initialize CDB storage, run number, set CDB lock
1072 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1075 AliRunLoader* runLoader = LoadRun();
1076 if (!runLoader) return kFALSE;
1078 TString detStr = detectors;
1079 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1080 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1081 AliModule* det = (AliModule*) detArray->At(iDet);
1082 if (!det || !det->IsActive()) continue;
1083 if (IsSelected(det->GetName(), detStr)) {
1084 AliInfo(Form("creating summable digits for %s", det->GetName()));
1085 AliCodeTimerAuto(Form("creating summable digits for %s", det->GetName()));
1086 det->Hits2SDigits();
1087 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1091 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1092 AliError(Form("the following detectors were not found: %s",
1094 if (fStopOnError) return kFALSE;
1103 //_____________________________________________________________________________
1104 Bool_t AliSimulation::RunDigitization(const char* detectors,
1105 const char* excludeDetectors)
1107 // run the digitization and produce digits from sdigits
1109 AliCodeTimerAuto("")
1111 // initialize CDB storage, run number, set CDB lock
1113 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1116 delete AliRunLoader::GetRunLoader();
1121 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1122 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1123 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
1124 // manager->SetEmbeddingFlag(fEmbeddingFlag);
1125 manager->SetInputStream(0, fGAliceFileName.Data());
1126 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1127 const char* fileName = ((TObjString*)
1128 (fBkgrdFileNames->At(iStream-1)))->GetName();
1129 manager->SetInputStream(iStream, fileName);
1132 TString detStr = detectors;
1133 TString detExcl = excludeDetectors;
1134 manager->GetInputStream(0)->ImportgAlice();
1135 AliRunLoader* runLoader =
1136 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
1137 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1138 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1139 AliModule* det = (AliModule*) detArray->At(iDet);
1140 if (!det || !det->IsActive()) continue;
1141 if (IsSelected(det->GetName(), detStr) &&
1142 !IsSelected(det->GetName(), detExcl)) {
1143 AliDigitizer* digitizer = det->CreateDigitizer(manager);
1146 AliError(Form("no digitizer for %s", det->GetName()));
1147 if (fStopOnError) return kFALSE;
1149 digitizer->SetRegionOfInterest(fRegionOfInterest);
1154 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1155 AliError(Form("the following detectors were not found: %s",
1157 if (fStopOnError) return kFALSE;
1160 if (!manager->GetListOfTasks()->IsEmpty()) {
1161 AliInfo("executing digitization");
1170 //_____________________________________________________________________________
1171 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1173 // run the digitization and produce digits from hits
1175 AliCodeTimerAuto("")
1177 // initialize CDB storage, run number, set CDB lock
1179 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1182 AliRunLoader* runLoader = LoadRun("READ");
1183 if (!runLoader) return kFALSE;
1185 TString detStr = detectors;
1186 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1187 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1188 AliModule* det = (AliModule*) detArray->At(iDet);
1189 if (!det || !det->IsActive()) continue;
1190 if (IsSelected(det->GetName(), detStr)) {
1191 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1196 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1197 AliError(Form("the following detectors were not found: %s",
1199 if (fStopOnError) return kFALSE;
1205 //_____________________________________________________________________________
1206 Bool_t AliSimulation::WriteRawData(const char* detectors,
1207 const char* fileName,
1208 Bool_t deleteIntermediateFiles,
1211 // convert the digits to raw data
1212 // First DDL raw data files for the given detectors are created.
1213 // If a file name is given, the DDL files are then converted to a DATE file.
1214 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1216 // If the file name has the extension ".root", the DATE file is converted
1218 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1219 // 'selrawdata' flag can be used to enable writing of detectors raw data
1220 // accoring to the trigger cluster.
1222 AliCodeTimerAuto("")
1224 TString detStr = detectors;
1225 if (!WriteRawFiles(detStr.Data())) {
1226 if (fStopOnError) return kFALSE;
1229 // run HLT simulation on simulated DDL raw files
1230 // and produce HLT ddl raw files to be included in date/root file
1231 if (IsSelected("HLT", detStr) && !fRunHLT.IsNull()) {
1233 if (fStopOnError) return kFALSE;
1237 TString dateFileName(fileName);
1238 if (!dateFileName.IsNull()) {
1239 Bool_t rootOutput = dateFileName.EndsWith(".root");
1240 if (rootOutput) dateFileName += ".date";
1241 TString selDateFileName;
1243 selDateFileName = "selected.";
1244 selDateFileName+= dateFileName;
1246 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1247 if (fStopOnError) return kFALSE;
1249 if (deleteIntermediateFiles) {
1250 AliRunLoader* runLoader = LoadRun("READ");
1251 if (runLoader) for (Int_t iEvent = 0;
1252 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1254 sprintf(command, "rm -r raw%d", iEvent);
1255 gSystem->Exec(command);
1260 if (!ConvertDateToRoot(dateFileName, fileName)) {
1261 if (fStopOnError) return kFALSE;
1263 if (deleteIntermediateFiles) {
1264 gSystem->Unlink(dateFileName);
1267 TString selFileName = "selected.";
1268 selFileName += fileName;
1269 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1270 if (fStopOnError) return kFALSE;
1272 if (deleteIntermediateFiles) {
1273 gSystem->Unlink(selDateFileName);
1282 //_____________________________________________________________________________
1283 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1285 // convert the digits to raw data DDL files
1287 AliCodeTimerAuto("")
1289 AliRunLoader* runLoader = LoadRun("READ");
1290 if (!runLoader) return kFALSE;
1292 // write raw data to DDL files
1293 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1294 AliInfo(Form("processing event %d", iEvent));
1295 runLoader->GetEvent(iEvent);
1296 TString baseDir = gSystem->WorkingDirectory();
1298 sprintf(dirName, "raw%d", iEvent);
1299 gSystem->MakeDirectory(dirName);
1300 if (!gSystem->ChangeDirectory(dirName)) {
1301 AliError(Form("couldn't change to directory %s", dirName));
1302 if (fStopOnError) return kFALSE; else continue;
1305 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1308 TString detStr = detectors;
1309 if (IsSelected("HLT", detStr)) {
1310 // Do nothing. "HLT" will be removed from detStr and HLT raw
1311 // data files are generated in RunHLT.
1314 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1315 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1316 AliModule* det = (AliModule*) detArray->At(iDet);
1317 if (!det || !det->IsActive()) continue;
1318 if (IsSelected(det->GetName(), detStr)) {
1319 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1324 if (!WriteTriggerRawData())
1325 if (fStopOnError) return kFALSE;
1327 gSystem->ChangeDirectory(baseDir);
1328 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1329 AliError(Form("the following detectors were not found: %s",
1331 if (fStopOnError) return kFALSE;
1340 //_____________________________________________________________________________
1341 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1342 const char* selDateFileName)
1344 // convert raw data DDL files to a DATE file with the program "dateStream"
1345 // The second argument is not empty when the user decides to write
1346 // the detectors raw data according to the trigger cluster.
1348 AliCodeTimerAuto("")
1350 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1352 AliError("the program dateStream was not found");
1353 if (fStopOnError) return kFALSE;
1358 AliRunLoader* runLoader = LoadRun("READ");
1359 if (!runLoader) return kFALSE;
1361 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1362 Bool_t selrawdata = kFALSE;
1363 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1366 // Note the option -s. It is used in order to avoid
1367 // the generation of SOR/EOR events.
1368 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1369 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1370 FILE* pipe = gSystem->OpenPipe(command, "w");
1373 AliError(Form("Cannot execute command: %s",command));
1377 Int_t selEvents = 0;
1378 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1380 UInt_t detectorPattern = 0;
1381 runLoader->GetEvent(iEvent);
1382 if (!runLoader->LoadTrigger()) {
1383 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1384 detectorPattern = aCTP->GetClusterMask();
1385 // Check if the event was triggered by CTP
1387 if (aCTP->GetClassMask()) selEvents++;
1391 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1393 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1394 selrawdata = kFALSE;
1398 fprintf(pipe, "GDC DetectorPattern %u\n", detectorPattern);
1402 // loop over detectors and DDLs
1403 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1404 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1406 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1407 Int_t ldcID = Int_t(ldc + 0.0001);
1408 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1410 char rawFileName[256];
1411 sprintf(rawFileName, "raw%d/%s",
1412 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1414 // check existence and size of raw data file
1415 FILE* file = fopen(rawFileName, "rb");
1416 if (!file) continue;
1417 fseek(file, 0, SEEK_END);
1418 unsigned long size = ftell(file);
1420 if (!size) continue;
1422 if (ldcID != prevLDC) {
1423 fprintf(pipe, " LDC Id %d\n", ldcID);
1426 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1431 Int_t result = gSystem->ClosePipe(pipe);
1433 if (!(selrawdata && selEvents > 0)) {
1435 return (result == 0);
1438 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1440 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1441 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1442 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1444 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1446 // Get the trigger decision and cluster
1447 UInt_t detectorPattern = 0;
1449 runLoader->GetEvent(iEvent);
1450 if (!runLoader->LoadTrigger()) {
1451 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1452 if (aCTP->GetClassMask() == 0) continue;
1453 detectorPattern = aCTP->GetClusterMask();
1454 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1455 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1458 fprintf(pipe2, "GDC DetectorPattern %u\n", detectorPattern);
1462 // loop over detectors and DDLs
1463 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1464 // Write only raw data from detectors that
1465 // are contained in the trigger cluster(s)
1466 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1468 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1470 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1471 Int_t ldcID = Int_t(ldc + 0.0001);
1472 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1474 char rawFileName[256];
1475 sprintf(rawFileName, "raw%d/%s",
1476 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1478 // check existence and size of raw data file
1479 FILE* file = fopen(rawFileName, "rb");
1480 if (!file) continue;
1481 fseek(file, 0, SEEK_END);
1482 unsigned long size = ftell(file);
1484 if (!size) continue;
1486 if (ldcID != prevLDC) {
1487 fprintf(pipe2, " LDC Id %d\n", ldcID);
1490 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1495 Int_t result2 = gSystem->ClosePipe(pipe2);
1498 return ((result == 0) && (result2 == 0));
1501 //_____________________________________________________________________________
1502 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1503 const char* rootFileName)
1505 // convert a DATE file to a root file with the program "alimdc"
1508 const Int_t kDBSize = 2000000000;
1509 const Int_t kTagDBSize = 1000000000;
1510 const Bool_t kFilter = kFALSE;
1511 const Int_t kCompression = 1;
1513 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1515 AliError("the program alimdc was not found");
1516 if (fStopOnError) return kFALSE;
1521 AliInfo(Form("converting DATE file %s to root file %s",
1522 dateFileName, rootFileName));
1524 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1525 const char* tagDBFS = "/tmp/mdc1/tags";
1527 // User defined file system locations
1528 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1529 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1530 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1531 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1532 if (gSystem->Getenv("ALIMDC_TAGDB"))
1533 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1535 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1536 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1537 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1539 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1540 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1541 gSystem->Exec(Form("mkdir %s",tagDBFS));
1543 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1544 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1545 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1547 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1548 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1549 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1551 return (result == 0);
1555 //_____________________________________________________________________________
1556 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1558 // delete existing run loaders, open a new one and load gAlice
1560 delete AliRunLoader::GetRunLoader();
1561 AliRunLoader* runLoader =
1562 AliRunLoader::Open(fGAliceFileName.Data(),
1563 AliConfig::GetDefaultEventFolderName(), mode);
1565 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1568 runLoader->LoadgAlice();
1569 runLoader->LoadHeader();
1570 gAlice = runLoader->GetAliRun();
1572 AliError(Form("no gAlice object found in file %s",
1573 fGAliceFileName.Data()));
1579 //_____________________________________________________________________________
1580 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1582 // get or calculate the number of signal events per background event
1584 if (!fBkgrdFileNames) return 1;
1585 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1586 if (nBkgrdFiles == 0) return 1;
1588 // get the number of signal events
1590 AliRunLoader* runLoader =
1591 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1592 if (!runLoader) return 1;
1594 nEvents = runLoader->GetNumberOfEvents();
1599 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1600 // get the number of background events
1601 const char* fileName = ((TObjString*)
1602 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1603 AliRunLoader* runLoader =
1604 AliRunLoader::Open(fileName, "BKGRD");
1605 if (!runLoader) continue;
1606 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1609 // get or calculate the number of signal per background events
1610 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1611 if (nSignalPerBkgrd <= 0) {
1612 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1613 } else if (result && (result != nSignalPerBkgrd)) {
1614 AliInfo(Form("the number of signal events per background event "
1615 "will be changed from %d to %d for stream %d",
1616 nSignalPerBkgrd, result, iBkgrdFile+1));
1617 nSignalPerBkgrd = result;
1620 if (!result) result = nSignalPerBkgrd;
1621 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1622 AliWarning(Form("not enough background events (%d) for %d signal events "
1623 "using %d signal per background events for stream %d",
1624 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1631 //_____________________________________________________________________________
1632 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1634 // check whether detName is contained in detectors
1635 // if yes, it is removed from detectors
1637 // check if all detectors are selected
1638 if ((detectors.CompareTo("ALL") == 0) ||
1639 detectors.BeginsWith("ALL ") ||
1640 detectors.EndsWith(" ALL") ||
1641 detectors.Contains(" ALL ")) {
1646 // search for the given detector
1647 Bool_t result = kFALSE;
1648 if ((detectors.CompareTo(detName) == 0) ||
1649 detectors.BeginsWith(detName+" ") ||
1650 detectors.EndsWith(" "+detName) ||
1651 detectors.Contains(" "+detName+" ")) {
1652 detectors.ReplaceAll(detName, "");
1656 // clean up the detectors string
1657 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1658 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1659 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1664 //_____________________________________________________________________________
1665 Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName)
1668 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1669 // These can be used for embedding of MC tracks into RAW data using the standard
1670 // merging procedure.
1672 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1675 AliError("no gAlice object. Restart aliroot and try again.");
1678 if (gAlice->Modules()->GetEntries() > 0) {
1679 AliError("gAlice was already run. Restart aliroot and try again.");
1683 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1687 gROOT->LoadMacro(fConfigFileName.Data());
1688 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1690 if(AliCDBManager::Instance()->GetRun() >= 0) {
1691 SetRunNumber(AliCDBManager::Instance()->GetRun());
1693 AliWarning("Run number not initialized!!");
1696 AliRunLoader::GetRunLoader()->CdGAFile();
1698 AliPDG::AddParticlesToPdgDataBase();
1700 gAlice->GetMCApp()->Init();
1702 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1703 gAlice->InitLoaders();
1704 AliRunLoader::GetRunLoader()->MakeTree("E");
1705 AliRunLoader::GetRunLoader()->LoadKinematics("RECREATE");
1706 AliRunLoader::GetRunLoader()->LoadTrackRefs("RECREATE");
1707 AliRunLoader::GetRunLoader()->LoadHits("all","RECREATE");
1709 // Save stuff at the beginning of the file to avoid file corruption
1710 AliRunLoader::GetRunLoader()->CdGAFile();
1715 //AliCDBManager* man = AliCDBManager::Instance();
1716 //man->SetRun(0); // Should this come from rawdata header ?
1720 // Get the runloader
1721 AliRunLoader* runLoader = AliRunLoader::GetRunLoader();
1723 // Open esd file if available
1724 TFile* esdFile = TFile::Open(esdFileName);
1725 Bool_t esdOK = (esdFile != 0);
1726 AliESD* esd = new AliESD;
1729 treeESD = (TTree*) esdFile->Get("esdTree");
1731 AliWarning("No ESD tree found");
1734 treeESD->SetBranchAddress("ESD", &esd);
1738 // Create the RawReader
1739 TString fileName(rawDirectory);
1740 AliRawReader* rawReader = 0x0;
1741 if (fileName.EndsWith("/")) {
1742 rawReader = new AliRawReaderFile(fileName);
1743 } else if (fileName.EndsWith(".root")) {
1744 rawReader = new AliRawReaderRoot(fileName);
1745 } else if (!fileName.IsNull()) {
1746 rawReader = new AliRawReaderDate(fileName);
1748 // if (!fEquipIdMap.IsNull() && fRawReader)
1749 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1751 // Get list of detectors
1752 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1755 AliHeader* header = runLoader->GetHeader();
1757 TString detStr = fMakeSDigits;
1761 if (!(rawReader->NextEvent())) break;
1764 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1765 AliModule* det = (AliModule*) detArray->At(iDet);
1766 if (!det || !det->IsActive()) continue;
1767 if (IsSelected(det->GetName(), detStr)) {
1768 AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1769 det->Raw2SDigits(rawReader);
1776 // If ESD information available obtain reconstructed vertex and store in header.
1778 treeESD->GetEvent(nev);
1779 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1780 Double_t position[3];
1781 esdVertex->GetXYZ(position);
1782 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1785 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1786 mcHeader->SetPrimaryVertex(mcV);
1787 header->Reset(0,nev);
1788 header->SetGenEventHeader(mcHeader);
1789 printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1794 runLoader->TreeE()->Fill();
1795 runLoader->SetNextEvent();
1801 runLoader->CdGAFile();
1802 runLoader->WriteHeader("OVERWRITE");
1803 runLoader->WriteRunLoader();
1808 //_____________________________________________________________________________
1809 void AliSimulation::FinishRun()
1812 // Called at the end of the run.
1817 AliDebug(1, "Finish Lego");
1818 AliRunLoader::GetRunLoader()->CdGAFile();
1822 // Clean detector information
1823 TIter next(gAlice->Modules());
1824 AliModule *detector;
1825 while((detector = dynamic_cast<AliModule*>(next()))) {
1826 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
1827 detector->FinishRun();
1830 AliDebug(1, "AliRunLoader::GetRunLoader()->WriteHeader(OVERWRITE)");
1831 AliRunLoader::GetRunLoader()->WriteHeader("OVERWRITE");
1833 // Write AliRun info and all detectors parameters
1834 AliRunLoader::GetRunLoader()->CdGAFile();
1835 gAlice->Write(0,TObject::kOverwrite);//write AliRun
1836 AliRunLoader::GetRunLoader()->Write(0,TObject::kOverwrite);//write RunLoader itself
1838 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
1839 AliRunLoader::GetRunLoader()->Synchronize();
1842 //_____________________________________________________________________________
1843 Int_t AliSimulation::GetDetIndex(const char* detector)
1845 // return the detector index corresponding to detector
1847 for (index = 0; index < fgkNDetectors ; index++) {
1848 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
1854 //_____________________________________________________________________________
1855 Bool_t AliSimulation::RunHLT()
1857 // Run the HLT simulation
1858 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
1859 // Disabled if fRunHLT is empty, default vaule is "default".
1860 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
1861 // The default simulation depends on the HLT component libraries and their
1862 // corresponding agents which define components and chains to run. See
1863 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
1864 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
1866 // The libraries to be loaded can be specified as an option.
1868 // AliSimulation sim;
1869 // sim.SetRunHLT("libAliHLTSample.so");
1871 // will only load <tt>libAliHLTSample.so</tt>
1873 // Other available options:
1874 // \li loglevel=<i>level</i> <br>
1875 // logging level for this processing
1877 // disable redirection of log messages to AliLog class
1878 // \li config=<i>macro</i>
1879 // configuration macro
1880 // \li chains=<i>configuration</i>
1881 // comma separated list of configurations to be run during simulation
1882 // \li rawfile=<i>file</i>
1883 // source for the RawReader to be created, the default is <i>./</i> if
1884 // raw data is simulated
1887 AliRunLoader* pRunLoader = LoadRun("READ");
1888 if (!pRunLoader) return kFALSE;
1890 // initialize CDB storage, run number, set CDB lock
1892 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1895 // load the library dynamically
1896 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
1898 // check for the library version
1899 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
1901 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
1904 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
1905 AliError(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
1909 // print compile info
1910 typedef void (*CompileInfo)( char*& date, char*& time);
1911 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
1915 (*fctInfo)(date,time);
1916 if (!date) {date=new Char_t[8]; strcpy(date,"unknown");}
1917 if (!time) {time=new Char_t[8]; strcpy(time,"unknown");}
1918 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
1922 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
1925 // create instance of the HLT simulation
1926 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
1927 AliHLTSimulation* pHLT=NULL;
1928 if (fctCreate==NULL || (pHLT=(fctCreate()))==NULL) {
1929 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
1933 // init the HLT simulation
1935 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
1936 TString detStr = fWriteRawData;
1937 if (!IsSelected("HLT", detStr)) {
1938 options+=" writerawfiles=";
1940 options+=" writerawfiles=HLT";
1943 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
1944 // as a matter of fact, HLT will run reconstruction and needs the RawReader
1945 // in order to get detector data. By default, RawReaderFile is used to read
1946 // the already simulated ddl files. Date and Root files from the raw data
1947 // are generated after the HLT simulation.
1948 options+=" rawfile=./";
1951 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
1952 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
1953 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
1955 // run the HLT simulation
1956 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
1957 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
1958 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
1962 // delete the instance
1963 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
1964 if (fctDelete==NULL || fctDelete(pHLT)<0) {
1965 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
1969 return iResult>=0?kTRUE:kFALSE;
1972 //_____________________________________________________________________________
1973 Bool_t AliSimulation::RunQA()
1975 // run the QA on summable hits, digits or digits
1977 if(!gAlice) return kFALSE;
1978 fQASteer->SetRunLoader(AliRunLoader::GetRunLoader()) ;
1980 TString detectorsw("") ;
1982 fQASteer->SetEventSpecie(fEventSpecie) ;
1983 detectorsw = fQASteer->Run(fQADetectors.Data()) ;
1984 if ( detectorsw.IsNull() )
1989 //_____________________________________________________________________________
1990 Bool_t AliSimulation::SetRunQA(TString detAndAction)
1992 // Allows to run QA for a selected set of detectors
1993 // and a selected set of tasks among HITS, SDIGITS and DIGITS
1994 // all selected detectors run the same selected tasks
1996 if (!detAndAction.Contains(":")) {
1997 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2001 Int_t colon = detAndAction.Index(":") ;
2002 fQADetectors = detAndAction(0, colon) ;
2003 if (fQADetectors.Contains("ALL") )
2004 fQADetectors = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2005 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2006 if (fQATasks.Contains("ALL") ) {
2007 fQATasks = Form("%d %d %d", AliQA::kHITS, AliQA::kSDIGITS, AliQA::kDIGITS) ;
2009 fQATasks.ToUpper() ;
2011 if ( fQATasks.Contains("HIT") )
2012 tempo = Form("%d ", AliQA::kHITS) ;
2013 if ( fQATasks.Contains("SDIGIT") )
2014 tempo += Form("%d ", AliQA::kSDIGITS) ;
2015 if ( fQATasks.Contains("DIGIT") )
2016 tempo += Form("%d ", AliQA::kDIGITS) ;
2018 if (fQATasks.IsNull()) {
2019 AliInfo("No QA requested\n") ;
2024 TString tempo(fQATasks) ;
2025 tempo.ReplaceAll(Form("%d", AliQA::kHITS), AliQA::GetTaskName(AliQA::kHITS)) ;
2026 tempo.ReplaceAll(Form("%d", AliQA::kSDIGITS), AliQA::GetTaskName(AliQA::kSDIGITS)) ;
2027 tempo.ReplaceAll(Form("%d", AliQA::kDIGITS), AliQA::GetTaskName(AliQA::kDIGITS)) ;
2028 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2030 fQASteer->SetActiveDetectors(fQADetectors) ;
2031 fQASteer->SetTasks(fQATasks) ;
2032 for (Int_t det = 0 ; det < AliQA::kNDET ; det++)
2033 fQASteer->SetWriteExpert(AliQA::DETECTORINDEX_t(det)) ;
2038 //_____________________________________________________________________________
2039 void AliSimulation::ProcessEnvironmentVars()
2041 // Extract run number and random generator seed from env variables
2043 AliInfo("Processing environment variables");
2045 // Random Number seed
2047 // first check that seed is not already set
2049 if (gSystem->Getenv("CONFIG_SEED")) {
2050 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2053 if (gSystem->Getenv("CONFIG_SEED")) {
2054 AliInfo(Form("Seed for random number generation already set (%d)"
2055 ": CONFIG_SEED variable ignored!", fSeed));
2059 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2063 // first check that run number is not already set
2065 if (gSystem->Getenv("DC_RUN")) {
2066 fRun = atoi(gSystem->Getenv("DC_RUN"));
2069 if (gSystem->Getenv("DC_RUN")) {
2070 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2074 AliInfo(Form("Run number = %d", fRun));
2077 //---------------------------------------------------------------------
2078 void AliSimulation::WriteGRPEntry()
2080 // Get the necessary information from galice (generator, trigger etc) and
2081 // write a GRP entry corresponding to the settings in the Config.C used
2082 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2085 AliInfo("Writing global run parameters entry into the OCDB");
2087 AliGRPObject* grpObj = new AliGRPObject();
2089 grpObj->SetRunType("PHYSICS");
2090 grpObj->SetTimeStart(0);
2091 grpObj->SetTimeEnd(9999);
2093 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2095 grpObj->SetBeamEnergy(gen->GetEnergyCMS());
2098 gen->GetProjectile(projectile,a,z);
2100 gen->GetTarget(target,a,z);
2101 TString beamType = projectile + "-" + target;
2102 beamType.ReplaceAll(" ","");
2103 if (!beamType.CompareTo("-")) {
2104 grpObj->SetBeamType("UNKNOWN");
2107 grpObj->SetBeamType(beamType);
2108 // Heavy ion run, the event specie is set to kHighMult
2109 fEventSpecie = AliRecoParam::kHighMult;
2110 if ((strcmp(beamType,"p-p") == 0) ||
2111 (strcmp(beamType,"p-") == 0) ||
2112 (strcmp(beamType,"-p") == 0) ||
2113 (strcmp(beamType,"P-P") == 0) ||
2114 (strcmp(beamType,"P-") == 0) ||
2115 (strcmp(beamType,"-P") == 0)) {
2116 // Proton run, the event specie is set to kLowMult
2117 fEventSpecie = AliRecoParam::kLowMult;
2121 AliWarning("Unknown beam type and energy! Setting energy to 0");
2122 grpObj->SetBeamEnergy(0);
2123 grpObj->SetBeamType("UNKNOWN");
2126 UInt_t detectorPattern = 0;
2128 TObjArray *detArray = gAlice->Detectors();
2129 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2130 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2131 detectorPattern |= (1 << iDet);
2136 if (!fMakeTrigger.IsNull() || strcmp(gAlice->GetTriggerDescriptor(),""))
2137 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2140 if (!fRunHLT.IsNull())
2141 detectorPattern |= (1 << AliDAQ::kHLTId);
2143 grpObj->SetNumberOfDetectors((Char_t)nDets);
2144 grpObj->SetDetectorMask((Int_t)detectorPattern);
2145 grpObj->SetLHCPeriod("LHC08c");
2146 grpObj->SetLHCState("STABLE_BEAMS");
2147 grpObj->SetLHCLuminosity(0,(AliGRPObject::Stats)0);
2148 grpObj->SetBeamIntensity(0,(AliGRPObject::Stats)0);
2150 AliMagF *field = gAlice->Field();
2151 Float_t solenoidField = TMath::Abs(field->SolenoidField());
2152 Float_t factor = field->Factor();
2153 Float_t l3current = TMath::Abs(factor)*solenoidField*30000./5.;
2154 grpObj->SetL3Current(l3current,(AliGRPObject::Stats)0);
2157 grpObj->SetL3Polarity(0);
2158 grpObj->SetDipolePolarity(0);
2161 grpObj->SetL3Polarity(1);
2162 grpObj->SetDipolePolarity(1);
2165 if (TMath::Abs(factor) != 0)
2166 grpObj->SetDipoleCurrent(6000,(AliGRPObject::Stats)0);
2168 grpObj->SetDipoleCurrent(0,(AliGRPObject::Stats)0);
2170 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2172 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2174 // Now store the entry in OCDB
2175 AliCDBManager* man = AliCDBManager::Instance();
2177 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2178 AliCDBMetaData *metadata= new AliCDBMetaData();
2180 metadata->SetResponsible("alice-off@cern.ch");
2181 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2183 man->Put(grpObj,id,metadata);