1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for running generation, simulation and digitization //
22 // Hits, sdigits and digits are created for all detectors by typing: //
24 // AliSimulation sim; //
27 // The Run method returns kTRUE in case of successful execution. //
28 // The number of events can be given as argument to the Run method or it //
31 // sim.SetNumberOfEvents(n); //
33 // The name of the configuration file can be passed as argument to the //
34 // AliSimulation constructor or can be specified by //
36 // sim.SetConfigFile("..."); //
38 // The generation of particles and the simulation of detector hits can be //
39 // switched on or off by //
41 // sim.SetRunGeneration(kTRUE); // generation of primary particles //
42 // sim.SetRunSimulation(kFALSE); // but no tracking //
44 // For which detectors sdigits and digits will be created, can be steered //
47 // sim.SetMakeSDigits("ALL"); // make sdigits for all detectors //
48 // sim.SetMakeDigits("ITS TPC"); // make digits only for ITS and TPC //
50 // The argument is a (case sensitive) string with the names of the //
51 // detectors separated by a space. An empty string ("") can be used to //
52 // disable the creation of sdigits or digits. The special string "ALL" //
53 // selects all available detectors. This is the default. //
55 // The creation of digits from hits instead of from sdigits can be selected //
58 // sim.SetMakeDigitsFromHits("TRD"); //
60 // The argument is again a string with the selected detectors. Be aware that //
61 // this feature is not available for all detectors and that merging is not //
62 // possible, when digits are created directly from hits. //
64 // Background events can be merged by calling //
66 // sim.MergeWith("background/galice.root", 2); //
68 // The first argument is the file name of the background galice file. The //
69 // second argument is the number of signal events per background event. //
70 // By default this number is calculated from the number of available //
71 // background events. MergeWith can be called several times to merge more //
72 // than two event streams. It is assumed that the sdigits were already //
73 // produced for the background events. //
75 // The output of raw data can be switched on by calling //
77 // sim.SetWriteRawData("MUON"); // write raw data for MUON //
79 // The default output format of the raw data are DDL files. They are //
80 // converted to a DATE file, if a file name is given as second argument. //
81 // For this conversion the program "dateStream" is required. If the file //
82 // name has the extension ".root", the DATE file is converted to a root //
83 // file. The program "alimdc" is used for this purpose. For the conversion //
84 // to DATE and root format the two conversion programs have to be installed. //
85 // Only the raw data in the final format is kept if the third argument is //
88 // The methods RunSimulation, RunSDigitization, RunDigitization, //
89 // RunHitsDigitization and WriteRawData can be used to run only parts of //
90 // the full simulation chain. The creation of raw data DDL files and their //
91 // conversion to the DATE or root format can be run directly by calling //
92 // the methods WriteRawFiles, ConvertRawFilesToDate and ConvertDateToRoot. //
94 // The default number of events per file, which is usually set in the //
95 // config file, can be changed for individual detectors and data types //
98 // sim.SetEventsPerFile("PHOS", "Reconstructed Points", 3); //
100 // The first argument is the detector, the second one the data type and the //
101 // last one the number of events per file. Valid data types are "Hits", //
102 // "Summable Digits", "Digits", "Reconstructed Points" and "Tracks". //
103 // The number of events per file has to be set before the simulation of //
104 // hits. Otherwise it has no effect. //
106 ///////////////////////////////////////////////////////////////////////////////
110 #include <TGeoGlobalMagField.h>
111 #include <TGeoManager.h>
112 #include <TObjString.h>
115 #include <TVirtualMC.h>
116 #include <TVirtualMCApplication.h>
118 #include "AliAlignObj.h"
119 #include "AliCDBEntry.h"
120 #include "AliCDBManager.h"
121 #include "AliCDBStorage.h"
122 #include "AliCTPRawData.h"
123 #include "AliCentralTrigger.h"
124 #include "AliCentralTrigger.h"
125 #include "AliCodeTimer.h"
127 #include "AliDigitizer.h"
129 #include "AliGRPObject.h"
130 #include "AliGenEventHeader.h"
131 #include "AliGenerator.h"
132 #include "AliGeomManager.h"
133 #include "AliHLTSimulation.h"
134 #include "AliHeader.h"
136 #include "AliLegoGenerator.h"
140 #include "AliModule.h"
142 #include "AliRawReaderDate.h"
143 #include "AliRawReaderFile.h"
144 #include "AliRawReaderRoot.h"
146 #include "AliRunDigitizer.h"
147 #include "AliRunLoader.h"
148 #include "AliSimulation.h"
149 #include "AliSysInfo.h"
150 #include "AliVertexGenFile.h"
152 ClassImp(AliSimulation)
154 AliSimulation *AliSimulation::fgInstance = 0;
155 const char* AliSimulation::fgkDetectorName[AliSimulation::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
157 //_____________________________________________________________________________
158 AliSimulation::AliSimulation(const char* configFileName,
159 const char* name, const char* title) :
162 fRunGeneration(kTRUE),
163 fRunSimulation(kTRUE),
164 fLoadAlignFromCDB(kTRUE),
165 fLoadAlObjsListOfDets("ALL"),
169 fMakeDigitsFromHits(""),
171 fRawDataFileName(""),
172 fDeleteIntermediateFiles(kFALSE),
173 fWriteSelRawData(kFALSE),
174 fStopOnError(kFALSE),
176 fConfigFileName(configFileName),
177 fGAliceFileName("galice.root"),
179 fBkgrdFileNames(NULL),
180 fAlignObjArray(NULL),
181 fUseBkgrdVertex(kTRUE),
182 fRegionOfInterest(kFALSE),
187 fInitCDBCalled(kFALSE),
188 fInitRunNumberCalled(kFALSE),
189 fSetRunNumberFromDataCalled(kFALSE),
190 fEmbeddingFlag(kFALSE),
196 fEventSpecie(AliRecoParam::kDefault),
198 fWriteGRPEntry(kTRUE)
200 // create simulation object with default parameters
202 SetGAliceFile("galice.root");
205 fQAManager = AliQAManager::QAManager("sim") ;
206 fQAManager->SetActiveDetectors(fQADetectors) ;
207 fQATasks = Form("%d %d %d", AliQA::kHITS, AliQA::kSDIGITS, AliQA::kDIGITS) ;
208 fQAManager->SetTasks(fQATasks) ;
211 //_____________________________________________________________________________
212 AliSimulation::~AliSimulation()
216 fEventsPerFile.Delete();
217 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
218 // delete fAlignObjArray; fAlignObjArray=0;
220 if (fBkgrdFileNames) {
221 fBkgrdFileNames->Delete();
222 delete fBkgrdFileNames;
225 fSpecCDBUri.Delete();
226 if (fgInstance==this) fgInstance = 0;
230 AliCodeTimer::Instance()->Print();
234 //_____________________________________________________________________________
235 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
237 // set the number of events for one run
242 //_____________________________________________________________________________
243 void AliSimulation::InitCDB()
245 // activate a default CDB storage
246 // First check if we have any CDB storage set, because it is used
247 // to retrieve the calibration and alignment constants
249 if (fInitCDBCalled) return;
250 fInitCDBCalled = kTRUE;
252 AliCDBManager* man = AliCDBManager::Instance();
253 if (man->IsDefaultStorageSet())
255 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
256 AliWarning("Default CDB storage has been already set !");
257 AliWarning(Form("Ignoring the default storage declared in AliSimulation: %s",fCDBUri.Data()));
258 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
259 fCDBUri = man->GetDefaultStorage()->GetURI();
262 if (fCDBUri.Length() > 0)
264 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
265 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
266 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
268 fCDBUri="local://$ALICE_ROOT/OCDB";
269 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
270 AliWarning("Default CDB storage not yet set !!!!");
271 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
272 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
275 man->SetDefaultStorage(fCDBUri);
278 // Now activate the detector specific CDB storage locations
279 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
280 TObject* obj = fSpecCDBUri[i];
282 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
283 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
284 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
285 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
290 //_____________________________________________________________________________
291 void AliSimulation::InitRunNumber(){
292 // check run number. If not set, set it to 0 !!!!
294 if (fInitRunNumberCalled) return;
295 fInitRunNumberCalled = kTRUE;
297 AliCDBManager* man = AliCDBManager::Instance();
298 if (man->GetRun() >= 0)
300 AliFatal(Form("Run number cannot be set in AliCDBManager before start of simulation: "
301 "Use external variable DC_RUN or AliSimulation::SetRun()!"));
305 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
306 AliDebug(2, Form("Setting CDB run number to: %d",fRun));
307 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
310 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
311 AliWarning("Run number not yet set !!!!");
312 AliWarning(Form("Setting it now to: %d", fRun));
313 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
322 //_____________________________________________________________________________
323 void AliSimulation::SetCDBLock() {
324 // Set CDB lock: from now on it is forbidden to reset the run number
325 // or the default storage or to activate any further storage!
327 AliCDBManager::Instance()->SetLock(1);
330 //_____________________________________________________________________________
331 void AliSimulation::SetDefaultStorage(const char* uri) {
332 // Store the desired default CDB storage location
333 // Activate it later within the Run() method
339 //_____________________________________________________________________________
340 void AliSimulation::SetSpecificStorage(const char* calibType, const char* uri) {
341 // Store a detector-specific CDB storage location
342 // Activate it later within the Run() method
344 AliCDBPath aPath(calibType);
345 if(!aPath.IsValid()){
346 AliError(Form("Not a valid path: %s", calibType));
350 TObject* obj = fSpecCDBUri.FindObject(calibType);
351 if (obj) fSpecCDBUri.Remove(obj);
352 fSpecCDBUri.Add(new TNamed(calibType, uri));
356 //_____________________________________________________________________________
357 void AliSimulation::SetRunNumber(Int_t run)
360 // Activate it later within the Run() method
365 //_____________________________________________________________________________
366 void AliSimulation::SetSeed(Int_t seed)
369 // Activate it later within the Run() method
374 //_____________________________________________________________________________
375 Bool_t AliSimulation::SetRunNumberFromData()
377 // Set the CDB manager run number
378 // The run number is retrieved from gAlice
380 if (fSetRunNumberFromDataCalled) return kTRUE;
381 fSetRunNumberFromDataCalled = kTRUE;
383 AliCDBManager* man = AliCDBManager::Instance();
384 Int_t runData = -1, runCDB = -1;
386 AliRunLoader* runLoader = LoadRun("READ");
387 if (!runLoader) return kFALSE;
389 runData = runLoader->GetHeader()->GetRun();
393 runCDB = man->GetRun();
395 if (runCDB != runData) {
396 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
397 AliWarning(Form("A run number was previously set in AliCDBManager: %d !", runCDB));
398 AliWarning(Form("It will be replaced with the run number got from run header: %d !", runData));
399 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
404 man->SetRun(runData);
407 if(man->GetRun() < 0) {
408 AliError("Run number not properly initalized!");
417 //_____________________________________________________________________________
418 void AliSimulation::SetConfigFile(const char* fileName)
420 // set the name of the config file
422 fConfigFileName = fileName;
425 //_____________________________________________________________________________
426 void AliSimulation::SetGAliceFile(const char* fileName)
428 // set the name of the galice file
429 // the path is converted to an absolute one if it is relative
431 fGAliceFileName = fileName;
432 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
433 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
435 fGAliceFileName = absFileName;
436 delete[] absFileName;
439 AliDebug(2, Form("galice file name set to %s", fileName));
442 //_____________________________________________________________________________
443 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
446 // set the number of events per file for the given detector and data type
447 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
449 TNamed* obj = new TNamed(detector, type);
450 obj->SetUniqueID(nEvents);
451 fEventsPerFile.Add(obj);
454 //_____________________________________________________________________________
455 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
457 // Read the alignment objects from CDB.
458 // Each detector is supposed to have the
459 // alignment objects in DET/Align/Data CDB path.
460 // All the detector objects are then collected,
461 // sorted by geometry level (starting from ALIC) and
462 // then applied to the TGeo geometry.
463 // Finally an overlaps check is performed.
465 if (!AliGeomManager::GetGeometry() || !AliGeomManager::GetGeometry()->IsClosed()) {
466 AliError("Can't apply the misalignment! Geometry is not loaded or it is still opened!");
470 // initialize CDB storage, run number, set CDB lock
472 // if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
475 Bool_t delRunLoader = kFALSE;
477 runLoader = LoadRun("READ");
478 if (!runLoader) return kFALSE;
479 delRunLoader = kTRUE;
482 // Export ideal geometry
483 if(!gAlice->IsRootGeometry()) AliGeomManager::GetGeometry()->Export("geometry.root");
485 // Load alignment data from CDB and apply to geometry through AliGeomManager
486 if(fLoadAlignFromCDB){
488 TString detStr = fLoadAlObjsListOfDets;
489 TString loadAlObjsListOfDets = "";
491 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
492 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
493 AliModule* det = (AliModule*) detArray->At(iDet);
494 if (!det || !det->IsActive()) continue;
495 if (IsSelected(det->GetName(), detStr)) {
496 //add det to list of dets to be aligned from CDB
497 loadAlObjsListOfDets += det->GetName();
498 loadAlObjsListOfDets += " ";
500 } // end loop over detectors
501 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
502 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
504 // Check if the array with alignment objects was
505 // provided by the user. If yes, apply the objects
506 // to the present TGeo geometry
507 if (fAlignObjArray) {
508 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
509 AliError("The misalignment of one or more volumes failed!"
510 "Compare the list of simulated detectors and the list of detector alignment data!");
511 if (delRunLoader) delete runLoader;
517 // Update the internal geometry of modules (ITS needs it)
518 TString detStr = fLoadAlObjsListOfDets;
519 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
520 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
522 AliModule* det = (AliModule*) detArray->At(iDet);
523 if (!det || !det->IsActive()) continue;
524 if (IsSelected(det->GetName(), detStr)) {
525 det->UpdateInternalGeometry();
527 } // end loop over detectors
530 if (delRunLoader) delete runLoader;
535 //_____________________________________________________________________________
536 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
538 // add a file with background events for merging
540 TObjString* fileNameStr = new TObjString(fileName);
541 fileNameStr->SetUniqueID(nSignalPerBkgrd);
542 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
543 fBkgrdFileNames->Add(fileNameStr);
546 void AliSimulation::EmbedInto(const char* fileName, Int_t nSignalPerBkgrd)
548 // add a file with background events for embeddin
549 MergeWith(fileName, nSignalPerBkgrd);
550 fEmbeddingFlag = kTRUE;
553 //_____________________________________________________________________________
554 Bool_t AliSimulation::Run(Int_t nEvents)
556 // run the generation, simulation and digitization
561 // Load run number and seed from environmental vars
562 ProcessEnvironmentVars();
564 gRandom->SetSeed(fSeed);
566 if (nEvents > 0) fNEvents = nEvents;
568 // generation and simulation -> hits
569 if (fRunGeneration) {
570 if (!RunSimulation()) if (fStopOnError) return kFALSE;
573 // initialize CDB storage from external environment
574 // (either CDB manager or AliSimulation setters),
575 // if not already done in RunSimulation()
578 // Set run number in CDBManager from data
579 // From this point on the run number must be always loaded from data!
580 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
582 // Set CDB lock: from now on it is forbidden to reset the run number
583 // or the default storage or to activate any further storage!
586 // If RunSimulation was not called, load the geometry and misalign it
587 if (!AliGeomManager::GetGeometry()) {
588 // Initialize the geometry manager
589 AliGeomManager::LoadGeometry("geometry.root");
591 // // Check that the consistency of symbolic names for the activated subdetectors
592 // // in the geometry loaded by AliGeomManager
593 // AliRunLoader* runLoader = LoadRun("READ");
594 // if (!runLoader) return kFALSE;
596 // TString detsToBeChecked = "";
597 // TObjArray* detArray = runLoader->GetAliRun()->Detectors();
598 // for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
599 // AliModule* det = (AliModule*) detArray->At(iDet);
600 // if (!det || !det->IsActive()) continue;
601 // detsToBeChecked += det->GetName();
602 // detsToBeChecked += " ";
603 // } // end loop over detectors
604 // if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
605 if(!AliGeomManager::CheckSymNamesLUT("ALL"))
606 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
608 if (!AliGeomManager::GetGeometry()) if (fStopOnError) return kFALSE;
610 if(!MisalignGeometry()) if (fStopOnError) return kFALSE;
614 // hits -> summable digits
615 AliSysInfo::AddStamp("Start_sdigitization");
616 if (!fMakeSDigits.IsNull()) {
617 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
620 AliSysInfo::AddStamp("Stop_sdigitization");
622 AliSysInfo::AddStamp("Start_digitization");
623 // summable digits -> digits
624 if (!fMakeDigits.IsNull()) {
625 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
626 if (fStopOnError) return kFALSE;
629 AliSysInfo::AddStamp("Stop_digitization");
634 if (!fMakeDigitsFromHits.IsNull()) {
635 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
636 AliWarning(Form("Merging and direct creation of digits from hits "
637 "was selected for some detectors. "
638 "No merging will be done for the following detectors: %s",
639 fMakeDigitsFromHits.Data()));
641 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
642 if (fStopOnError) return kFALSE;
649 if (!RunTrigger(fMakeTrigger,fMakeDigits)) {
650 if (fStopOnError) return kFALSE;
655 // digits -> raw data
656 if (!fWriteRawData.IsNull()) {
657 if (!WriteRawData(fWriteRawData, fRawDataFileName,
658 fDeleteIntermediateFiles,fWriteSelRawData)) {
659 if (fStopOnError) return kFALSE;
664 // run HLT simulation on simulated digit data if raw data is not
665 // simulated, otherwise its called as part of WriteRawData
666 if (!fRunHLT.IsNull() && fWriteRawData.IsNull()) {
668 if (fStopOnError) return kFALSE;
674 Bool_t rv = RunQA() ;
680 // Cleanup of CDB manager: cache and active storages!
681 AliCDBManager::Instance()->ClearCache();
686 //_______________________________________________________________________
687 Bool_t AliSimulation::RunLego(const char *setup, Int_t nc1, Float_t c1min,
688 Float_t c1max,Int_t nc2,Float_t c2min,Float_t c2max,
689 Float_t rmin,Float_t rmax,Float_t zmax, AliLegoGenerator* gener, Int_t nev)
692 // Generates lego plots of:
693 // - radiation length map phi vs theta
694 // - radiation length map phi vs eta
695 // - interaction length map
696 // - g/cm2 length map
698 // ntheta bins in theta, eta
699 // themin minimum angle in theta (degrees)
700 // themax maximum angle in theta (degrees)
702 // phimin minimum angle in phi (degrees)
703 // phimax maximum angle in phi (degrees)
704 // rmin minimum radius
705 // rmax maximum radius
708 // The number of events generated = ntheta*nphi
709 // run input parameters in macro setup (default="Config.C")
711 // Use macro "lego.C" to visualize the 3 lego plots in spherical coordinates
714 <img src="picts/AliRunLego1.gif">
719 <img src="picts/AliRunLego2.gif">
724 <img src="picts/AliRunLego3.gif">
729 // run the generation and simulation
733 // initialize CDB storage and run number from external environment
734 // (either CDB manager or AliSimulation setters)
740 AliError("no gAlice object. Restart aliroot and try again.");
743 if (gAlice->Modules()->GetEntries() > 0) {
744 AliError("gAlice was already run. Restart aliroot and try again.");
748 AliInfo(Form("initializing gAlice with config file %s",
749 fConfigFileName.Data()));
752 if (nev == -1) nev = nc1 * nc2;
754 // check if initialisation has been done
755 // If runloader has been initialized, set the number of events per file to nc1 * nc2
758 if (!gener) gener = new AliLegoGenerator();
760 // Configure Generator
762 gener->SetRadiusRange(rmin, rmax);
763 gener->SetZMax(zmax);
764 gener->SetCoor1Range(nc1, c1min, c1max);
765 gener->SetCoor2Range(nc2, c2min, c2max);
769 fLego = new AliLego("lego",gener);
771 //__________________________________________________________________________
775 gROOT->LoadMacro(setup);
776 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
778 if(AliCDBManager::Instance()->GetRun() >= 0) {
779 SetRunNumber(AliCDBManager::Instance()->GetRun());
781 AliWarning("Run number not initialized!!");
784 AliRunLoader::Instance()->CdGAFile();
786 AliPDG::AddParticlesToPdgDataBase();
788 gAlice->GetMCApp()->Init();
790 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
793 //Must be here because some MCs (G4) adds detectors here and not in Config.C
794 gAlice->InitLoaders();
795 AliRunLoader::Instance()->MakeTree("E");
798 // Save stuff at the beginning of the file to avoid file corruption
799 AliRunLoader::Instance()->CdGAFile();
802 //Save current generator
803 AliGenerator *gen=gAlice->GetMCApp()->Generator();
804 gAlice->GetMCApp()->ResetGenerator(gener);
805 //Prepare MC for Lego Run
811 AliRunLoader::Instance()->SetNumberOfEventsPerFile(nev);
812 gMC->ProcessRun(nev);
814 // End of this run, close files
816 // Restore current generator
817 gAlice->GetMCApp()->ResetGenerator(gen);
818 // Delete Lego Object
824 //_____________________________________________________________________________
825 Bool_t AliSimulation::RunTrigger(const char* config, const char* detectors)
831 // initialize CDB storage from external environment
832 // (either CDB manager or AliSimulation setters),
833 // if not already done in RunSimulation()
836 // Set run number in CDBManager from data
837 // From this point on the run number must be always loaded from data!
838 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
840 // Set CDB lock: from now on it is forbidden to reset the run number
841 // or the default storage or to activate any further storage!
844 AliRunLoader* runLoader = LoadRun("READ");
845 if (!runLoader) return kFALSE;
846 TString trconfiguration = config;
848 if (trconfiguration.IsNull()) {
849 if (strcmp(gAlice->GetTriggerDescriptor(),"")) {
850 trconfiguration = gAlice->GetTriggerDescriptor();
853 AliWarning("No trigger descriptor is specified. Loading the one that is in the CDB.");
856 runLoader->MakeTree( "GG" );
857 AliCentralTrigger* aCTP = runLoader->GetTrigger();
858 // Load Configuration
859 if (!aCTP->LoadConfiguration( trconfiguration ))
863 if( !aCTP->RunTrigger( runLoader , detectors ) ) {
875 //_____________________________________________________________________________
876 Bool_t AliSimulation::WriteTriggerRawData()
878 // Writes the CTP (trigger) DDL raw data
879 // Details of the format are given in the
880 // trigger TDR - pages 134 and 135.
881 AliCTPRawData writer;
887 //_____________________________________________________________________________
888 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
890 // run the generation and simulation
894 // initialize CDB storage and run number from external environment
895 // (either CDB manager or AliSimulation setters)
901 AliError("no gAlice object. Restart aliroot and try again.");
904 if (gAlice->Modules()->GetEntries() > 0) {
905 AliError("gAlice was already run. Restart aliroot and try again.");
909 AliInfo(Form("initializing gAlice with config file %s",
910 fConfigFileName.Data()));
913 // Initialize ALICE Simulation run
918 gROOT->LoadMacro(fConfigFileName.Data());
919 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
921 if(AliCDBManager::Instance()->GetRun() >= 0) {
922 AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
924 AliWarning("Run number not initialized!!");
927 AliRunLoader::Instance()->CdGAFile();
929 AliPDG::AddParticlesToPdgDataBase();
931 gAlice->GetMCApp()->Init();
933 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
935 //Must be here because some MCs (G4) adds detectors here and not in Config.C
936 gAlice->InitLoaders();
937 AliRunLoader::Instance()->MakeTree("E");
938 AliRunLoader::Instance()->LoadKinematics("RECREATE");
939 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
940 AliRunLoader::Instance()->LoadHits("all","RECREATE");
942 // Save stuff at the beginning of the file to avoid file corruption
943 AliRunLoader::Instance()->CdGAFile();
945 gAlice->SetEventNrInRun(-1); //important - we start Begin event from increasing current number in run
946 //___________________________________________________________________________________________
948 // Get the trigger descriptor string
949 // Either from AliSimulation or from
951 if (fMakeTrigger.IsNull()) {
952 if (strcmp(gAlice->GetTriggerDescriptor(),""))
953 fMakeTrigger = gAlice->GetTriggerDescriptor();
956 gAlice->SetTriggerDescriptor(fMakeTrigger.Data());
958 // Set run number in CDBManager
959 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
961 AliRunLoader* runLoader = AliRunLoader::Instance();
963 AliError(Form("gAlice has no run loader object. "
964 "Check your config file: %s", fConfigFileName.Data()));
967 SetGAliceFile(runLoader->GetFileName());
970 #if ROOT_VERSION_CODE < 331527
971 AliGeomManager::SetGeometry(gGeoManager);
973 // Check that the consistency of symbolic names for the activated subdetectors
974 // in the geometry loaded by AliGeomManager
975 TString detsToBeChecked = "";
976 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
977 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
978 AliModule* det = (AliModule*) detArray->At(iDet);
979 if (!det || !det->IsActive()) continue;
980 detsToBeChecked += det->GetName();
981 detsToBeChecked += " ";
982 } // end loop over detectors
983 if(!AliGeomManager::CheckSymNamesLUT(detsToBeChecked.Data()))
984 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
985 MisalignGeometry(runLoader);
988 // AliRunLoader* runLoader = AliRunLoader::Instance();
990 // AliError(Form("gAlice has no run loader object. "
991 // "Check your config file: %s", fConfigFileName.Data()));
994 // SetGAliceFile(runLoader->GetFileName());
996 if (!gAlice->GetMCApp()->Generator()) {
997 AliError(Form("gAlice has no generator object. "
998 "Check your config file: %s", fConfigFileName.Data()));
1002 // Write GRP entry corresponding to the setting found in Cofig.C
1006 if (nEvents <= 0) nEvents = fNEvents;
1008 // get vertex from background file in case of merging
1009 if (fUseBkgrdVertex &&
1010 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
1011 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
1012 const char* fileName = ((TObjString*)
1013 (fBkgrdFileNames->At(0)))->GetName();
1014 AliInfo(Form("The vertex will be taken from the background "
1015 "file %s with nSignalPerBackground = %d",
1016 fileName, signalPerBkgrd));
1017 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
1018 gAlice->GetMCApp()->Generator()->SetVertexGenerator(vtxGen);
1021 if (!fRunSimulation) {
1022 gAlice->GetMCApp()->Generator()->SetTrackingFlag(0);
1025 // set the number of events per file for given detectors and data types
1026 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
1027 if (!fEventsPerFile[i]) continue;
1028 const char* detName = fEventsPerFile[i]->GetName();
1029 const char* typeName = fEventsPerFile[i]->GetTitle();
1030 TString loaderName(detName);
1031 loaderName += "Loader";
1032 AliLoader* loader = runLoader->GetLoader(loaderName);
1034 AliError(Form("RunSimulation", "no loader for %s found\n"
1035 "Number of events per file not set for %s %s",
1036 detName, typeName, detName));
1039 AliDataLoader* dataLoader =
1040 loader->GetDataLoader(typeName);
1042 AliError(Form("no data loader for %s found\n"
1043 "Number of events per file not set for %s %s",
1044 typeName, detName, typeName));
1047 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
1048 AliDebug(1, Form("number of events per file set to %d for %s %s",
1049 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
1052 AliInfo("running gAlice");
1053 AliSysInfo::AddStamp("Start_simulation");
1055 // Create the Root Tree with one branch per detector
1056 //Hits moved to begin event -> now we are crating separate tree for each event
1058 gMC->ProcessRun(nEvents);
1060 // End of this run, close files
1061 if(nEvents>0) FinishRun();
1063 AliSysInfo::AddStamp("Stop_simulation");
1069 //_____________________________________________________________________________
1070 Bool_t AliSimulation::RunSDigitization(const char* detectors)
1072 // run the digitization and produce summable digits
1073 static Int_t eventNr=0;
1074 AliCodeTimerAuto("")
1076 // initialize CDB storage, run number, set CDB lock
1078 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1081 AliRunLoader* runLoader = LoadRun();
1082 if (!runLoader) return kFALSE;
1084 TString detStr = detectors;
1085 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1086 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1087 AliModule* det = (AliModule*) detArray->At(iDet);
1088 if (!det || !det->IsActive()) continue;
1089 if (IsSelected(det->GetName(), detStr)) {
1090 AliInfo(Form("creating summable digits for %s", det->GetName()));
1091 AliCodeTimerAuto(Form("creating summable digits for %s", det->GetName()));
1092 det->Hits2SDigits();
1093 AliSysInfo::AddStamp(Form("Digit_%s_%d",det->GetName(),eventNr), 0,1, eventNr);
1097 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1098 AliError(Form("the following detectors were not found: %s",
1100 if (fStopOnError) return kFALSE;
1109 //_____________________________________________________________________________
1110 Bool_t AliSimulation::RunDigitization(const char* detectors,
1111 const char* excludeDetectors)
1113 // run the digitization and produce digits from sdigits
1115 AliCodeTimerAuto("")
1117 // initialize CDB storage, run number, set CDB lock
1119 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1122 delete AliRunLoader::Instance();
1127 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
1128 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
1129 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
1130 // manager->SetEmbeddingFlag(fEmbeddingFlag);
1131 manager->SetInputStream(0, fGAliceFileName.Data());
1132 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
1133 const char* fileName = ((TObjString*)
1134 (fBkgrdFileNames->At(iStream-1)))->GetName();
1135 manager->SetInputStream(iStream, fileName);
1138 TString detStr = detectors;
1139 TString detExcl = excludeDetectors;
1140 manager->GetInputStream(0)->ImportgAlice();
1141 AliRunLoader* runLoader =
1142 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
1143 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1144 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1145 AliModule* det = (AliModule*) detArray->At(iDet);
1146 if (!det || !det->IsActive()) continue;
1147 if (IsSelected(det->GetName(), detStr) &&
1148 !IsSelected(det->GetName(), detExcl)) {
1149 AliDigitizer* digitizer = det->CreateDigitizer(manager);
1152 AliError(Form("no digitizer for %s", det->GetName()));
1153 if (fStopOnError) return kFALSE;
1155 digitizer->SetRegionOfInterest(fRegionOfInterest);
1160 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1161 AliError(Form("the following detectors were not found: %s",
1163 if (fStopOnError) return kFALSE;
1166 if (!manager->GetListOfTasks()->IsEmpty()) {
1167 AliInfo("executing digitization");
1176 //_____________________________________________________________________________
1177 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
1179 // run the digitization and produce digits from hits
1181 AliCodeTimerAuto("")
1183 // initialize CDB storage, run number, set CDB lock
1185 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1188 AliRunLoader* runLoader = LoadRun("READ");
1189 if (!runLoader) return kFALSE;
1191 TString detStr = detectors;
1192 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1193 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1194 AliModule* det = (AliModule*) detArray->At(iDet);
1195 if (!det || !det->IsActive()) continue;
1196 if (IsSelected(det->GetName(), detStr)) {
1197 AliInfo(Form("creating digits from hits for %s", det->GetName()));
1202 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1203 AliError(Form("the following detectors were not found: %s",
1205 if (fStopOnError) return kFALSE;
1211 //_____________________________________________________________________________
1212 Bool_t AliSimulation::WriteRawData(const char* detectors,
1213 const char* fileName,
1214 Bool_t deleteIntermediateFiles,
1217 // convert the digits to raw data
1218 // First DDL raw data files for the given detectors are created.
1219 // If a file name is given, the DDL files are then converted to a DATE file.
1220 // If deleteIntermediateFiles is true, the DDL raw files are deleted
1222 // If the file name has the extension ".root", the DATE file is converted
1224 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
1225 // 'selrawdata' flag can be used to enable writing of detectors raw data
1226 // accoring to the trigger cluster.
1228 AliCodeTimerAuto("")
1230 TString detStr = detectors;
1231 if (!WriteRawFiles(detStr.Data())) {
1232 if (fStopOnError) return kFALSE;
1235 // run HLT simulation on simulated DDL raw files
1236 // and produce HLT ddl raw files to be included in date/root file
1237 if (IsSelected("HLT", detStr) && !fRunHLT.IsNull()) {
1239 if (fStopOnError) return kFALSE;
1243 TString dateFileName(fileName);
1244 if (!dateFileName.IsNull()) {
1245 Bool_t rootOutput = dateFileName.EndsWith(".root");
1246 if (rootOutput) dateFileName += ".date";
1247 TString selDateFileName;
1249 selDateFileName = "selected.";
1250 selDateFileName+= dateFileName;
1252 if (!ConvertRawFilesToDate(dateFileName,selDateFileName)) {
1253 if (fStopOnError) return kFALSE;
1255 if (deleteIntermediateFiles) {
1256 AliRunLoader* runLoader = LoadRun("READ");
1257 if (runLoader) for (Int_t iEvent = 0;
1258 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1260 sprintf(command, "rm -r raw%d", iEvent);
1261 gSystem->Exec(command);
1266 if (!ConvertDateToRoot(dateFileName, fileName)) {
1267 if (fStopOnError) return kFALSE;
1269 if (deleteIntermediateFiles) {
1270 gSystem->Unlink(dateFileName);
1273 TString selFileName = "selected.";
1274 selFileName += fileName;
1275 if (!ConvertDateToRoot(selDateFileName, selFileName)) {
1276 if (fStopOnError) return kFALSE;
1278 if (deleteIntermediateFiles) {
1279 gSystem->Unlink(selDateFileName);
1288 //_____________________________________________________________________________
1289 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1291 // convert the digits to raw data DDL files
1293 AliCodeTimerAuto("")
1295 AliRunLoader* runLoader = LoadRun("READ");
1296 if (!runLoader) return kFALSE;
1298 // write raw data to DDL files
1299 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1300 AliInfo(Form("processing event %d", iEvent));
1301 runLoader->GetEvent(iEvent);
1302 TString baseDir = gSystem->WorkingDirectory();
1304 sprintf(dirName, "raw%d", iEvent);
1305 gSystem->MakeDirectory(dirName);
1306 if (!gSystem->ChangeDirectory(dirName)) {
1307 AliError(Form("couldn't change to directory %s", dirName));
1308 if (fStopOnError) return kFALSE; else continue;
1311 ofstream runNbFile(Form("run%u",runLoader->GetHeader()->GetRun()));
1314 TString detStr = detectors;
1315 if (IsSelected("HLT", detStr)) {
1316 // Do nothing. "HLT" will be removed from detStr and HLT raw
1317 // data files are generated in RunHLT.
1320 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1321 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1322 AliModule* det = (AliModule*) detArray->At(iDet);
1323 if (!det || !det->IsActive()) continue;
1324 if (IsSelected(det->GetName(), detStr)) {
1325 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1330 if (!WriteTriggerRawData())
1331 if (fStopOnError) return kFALSE;
1333 gSystem->ChangeDirectory(baseDir);
1334 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1335 AliError(Form("the following detectors were not found: %s",
1337 if (fStopOnError) return kFALSE;
1346 //_____________________________________________________________________________
1347 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName,
1348 const char* selDateFileName)
1350 // convert raw data DDL files to a DATE file with the program "dateStream"
1351 // The second argument is not empty when the user decides to write
1352 // the detectors raw data according to the trigger cluster.
1354 AliCodeTimerAuto("")
1356 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1358 AliError("the program dateStream was not found");
1359 if (fStopOnError) return kFALSE;
1364 AliRunLoader* runLoader = LoadRun("READ");
1365 if (!runLoader) return kFALSE;
1367 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1368 Bool_t selrawdata = kFALSE;
1369 if (strcmp(selDateFileName,"") != 0) selrawdata = kTRUE;
1372 // Note the option -s. It is used in order to avoid
1373 // the generation of SOR/EOR events.
1374 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1375 dateFileName, runLoader->GetNumberOfEvents(),runLoader->GetHeader()->GetRun());
1376 FILE* pipe = gSystem->OpenPipe(command, "w");
1379 AliError(Form("Cannot execute command: %s",command));
1383 Int_t selEvents = 0;
1384 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1386 UInt_t detectorPattern = 0;
1387 runLoader->GetEvent(iEvent);
1388 if (!runLoader->LoadTrigger()) {
1389 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1390 detectorPattern = aCTP->GetClusterMask();
1391 // Check if the event was triggered by CTP
1393 if (aCTP->GetClassMask()) selEvents++;
1397 AliWarning("No trigger can be loaded! Some fields in the event header will be empty !");
1399 AliWarning("No trigger can be loaded! Writing of selected raw data is abandoned !");
1400 selrawdata = kFALSE;
1404 fprintf(pipe, "GDC DetectorPattern %u\n", detectorPattern);
1408 // loop over detectors and DDLs
1409 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1410 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1412 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1413 Int_t ldcID = Int_t(ldc + 0.0001);
1414 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1416 char rawFileName[256];
1417 sprintf(rawFileName, "raw%d/%s",
1418 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1420 // check existence and size of raw data file
1421 FILE* file = fopen(rawFileName, "rb");
1422 if (!file) continue;
1423 fseek(file, 0, SEEK_END);
1424 unsigned long size = ftell(file);
1426 if (!size) continue;
1428 if (ldcID != prevLDC) {
1429 fprintf(pipe, " LDC Id %d\n", ldcID);
1432 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1437 Int_t result = gSystem->ClosePipe(pipe);
1439 if (!(selrawdata && selEvents > 0)) {
1441 return (result == 0);
1444 AliInfo(Form("converting selected by trigger cluster raw data DDL files to DATE file %s", selDateFileName));
1446 sprintf(command, "dateStream -c -s -D -o %s -# %d -C -run %d",
1447 selDateFileName,selEvents,runLoader->GetHeader()->GetRun());
1448 FILE* pipe2 = gSystem->OpenPipe(command, "w");
1450 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1452 // Get the trigger decision and cluster
1453 UInt_t detectorPattern = 0;
1455 runLoader->GetEvent(iEvent);
1456 if (!runLoader->LoadTrigger()) {
1457 AliCentralTrigger *aCTP = runLoader->GetTrigger();
1458 if (aCTP->GetClassMask() == 0) continue;
1459 detectorPattern = aCTP->GetClusterMask();
1460 detClust = AliDAQ::ListOfTriggeredDetectors(detectorPattern);
1461 AliInfo(Form("List of detectors to be read out: %s",detClust.Data()));
1464 fprintf(pipe2, "GDC DetectorPattern %u\n", detectorPattern);
1468 // loop over detectors and DDLs
1469 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors; iDet++) {
1470 // Write only raw data from detectors that
1471 // are contained in the trigger cluster(s)
1472 if (!IsSelected(AliDAQ::DetectorName(iDet),detClust)) continue;
1474 for (Int_t iDDL = 0; iDDL < AliDAQ::NumberOfDdls(iDet); iDDL++) {
1476 Int_t ddlID = AliDAQ::DdlID(iDet,iDDL);
1477 Int_t ldcID = Int_t(ldc + 0.0001);
1478 ldc += AliDAQ::NumberOfLdcs(iDet) / AliDAQ::NumberOfDdls(iDet);
1480 char rawFileName[256];
1481 sprintf(rawFileName, "raw%d/%s",
1482 iEvent, AliDAQ::DdlFileName(iDet,iDDL));
1484 // check existence and size of raw data file
1485 FILE* file = fopen(rawFileName, "rb");
1486 if (!file) continue;
1487 fseek(file, 0, SEEK_END);
1488 unsigned long size = ftell(file);
1490 if (!size) continue;
1492 if (ldcID != prevLDC) {
1493 fprintf(pipe2, " LDC Id %d\n", ldcID);
1496 fprintf(pipe2, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1501 Int_t result2 = gSystem->ClosePipe(pipe2);
1504 return ((result == 0) && (result2 == 0));
1507 //_____________________________________________________________________________
1508 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1509 const char* rootFileName)
1511 // convert a DATE file to a root file with the program "alimdc"
1514 const Int_t kDBSize = 2000000000;
1515 const Int_t kTagDBSize = 1000000000;
1516 const Bool_t kFilter = kFALSE;
1517 const Int_t kCompression = 1;
1519 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1521 AliError("the program alimdc was not found");
1522 if (fStopOnError) return kFALSE;
1527 AliInfo(Form("converting DATE file %s to root file %s",
1528 dateFileName, rootFileName));
1530 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1531 const char* tagDBFS = "/tmp/mdc1/tags";
1533 // User defined file system locations
1534 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1535 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1536 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1537 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1538 if (gSystem->Getenv("ALIMDC_TAGDB"))
1539 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1541 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1542 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1543 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1545 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1546 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1547 gSystem->Exec(Form("mkdir %s",tagDBFS));
1549 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1550 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1551 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1553 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1554 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1555 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1557 return (result == 0);
1561 //_____________________________________________________________________________
1562 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1564 // delete existing run loaders, open a new one and load gAlice
1566 delete AliRunLoader::Instance();
1567 AliRunLoader* runLoader =
1568 AliRunLoader::Open(fGAliceFileName.Data(),
1569 AliConfig::GetDefaultEventFolderName(), mode);
1571 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1574 runLoader->LoadgAlice();
1575 runLoader->LoadHeader();
1576 gAlice = runLoader->GetAliRun();
1578 AliError(Form("no gAlice object found in file %s",
1579 fGAliceFileName.Data()));
1585 //_____________________________________________________________________________
1586 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1588 // get or calculate the number of signal events per background event
1590 if (!fBkgrdFileNames) return 1;
1591 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1592 if (nBkgrdFiles == 0) return 1;
1594 // get the number of signal events
1596 AliRunLoader* runLoader =
1597 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1598 if (!runLoader) return 1;
1600 nEvents = runLoader->GetNumberOfEvents();
1605 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1606 // get the number of background events
1607 const char* fileName = ((TObjString*)
1608 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1609 AliRunLoader* runLoader =
1610 AliRunLoader::Open(fileName, "BKGRD");
1611 if (!runLoader) continue;
1612 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1615 // get or calculate the number of signal per background events
1616 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1617 if (nSignalPerBkgrd <= 0) {
1618 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1619 } else if (result && (result != nSignalPerBkgrd)) {
1620 AliInfo(Form("the number of signal events per background event "
1621 "will be changed from %d to %d for stream %d",
1622 nSignalPerBkgrd, result, iBkgrdFile+1));
1623 nSignalPerBkgrd = result;
1626 if (!result) result = nSignalPerBkgrd;
1627 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1628 AliWarning(Form("not enough background events (%d) for %d signal events "
1629 "using %d signal per background events for stream %d",
1630 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1637 //_____________________________________________________________________________
1638 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1640 // check whether detName is contained in detectors
1641 // if yes, it is removed from detectors
1643 // check if all detectors are selected
1644 if ((detectors.CompareTo("ALL") == 0) ||
1645 detectors.BeginsWith("ALL ") ||
1646 detectors.EndsWith(" ALL") ||
1647 detectors.Contains(" ALL ")) {
1652 // search for the given detector
1653 Bool_t result = kFALSE;
1654 if ((detectors.CompareTo(detName) == 0) ||
1655 detectors.BeginsWith(detName+" ") ||
1656 detectors.EndsWith(" "+detName) ||
1657 detectors.Contains(" "+detName+" ")) {
1658 detectors.ReplaceAll(detName, "");
1662 // clean up the detectors string
1663 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1664 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1665 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1670 //_____________________________________________________________________________
1671 Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName)
1674 // Steering routine to convert raw data in directory rawDirectory/ to fake SDigits.
1675 // These can be used for embedding of MC tracks into RAW data using the standard
1676 // merging procedure.
1678 // If an ESD file is given the reconstructed vertex is taken from it and stored in the event header.
1681 AliError("no gAlice object. Restart aliroot and try again.");
1684 if (gAlice->Modules()->GetEntries() > 0) {
1685 AliError("gAlice was already run. Restart aliroot and try again.");
1689 AliInfo(Form("initializing gAlice with config file %s",fConfigFileName.Data()));
1693 gROOT->LoadMacro(fConfigFileName.Data());
1694 gInterpreter->ProcessLine(gAlice->GetConfigFunction());
1696 if(AliCDBManager::Instance()->GetRun() >= 0) {
1697 SetRunNumber(AliCDBManager::Instance()->GetRun());
1699 AliWarning("Run number not initialized!!");
1702 AliRunLoader::Instance()->CdGAFile();
1704 AliPDG::AddParticlesToPdgDataBase();
1706 gAlice->GetMCApp()->Init();
1708 gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
1710 //Must be here because some MCs (G4) adds detectors here and not in Config.C
1711 gAlice->InitLoaders();
1712 AliRunLoader::Instance()->MakeTree("E");
1713 AliRunLoader::Instance()->LoadKinematics("RECREATE");
1714 AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
1715 AliRunLoader::Instance()->LoadHits("all","RECREATE");
1717 // Save stuff at the beginning of the file to avoid file corruption
1718 AliRunLoader::Instance()->CdGAFile();
1723 //AliCDBManager* man = AliCDBManager::Instance();
1724 //man->SetRun(0); // Should this come from rawdata header ?
1728 // Get the runloader
1729 AliRunLoader* runLoader = AliRunLoader::Instance();
1731 // Open esd file if available
1732 TFile* esdFile = TFile::Open(esdFileName);
1733 Bool_t esdOK = (esdFile != 0);
1734 AliESD* esd = new AliESD;
1737 treeESD = (TTree*) esdFile->Get("esdTree");
1739 AliWarning("No ESD tree found");
1742 treeESD->SetBranchAddress("ESD", &esd);
1746 // Create the RawReader
1747 TString fileName(rawDirectory);
1748 AliRawReader* rawReader = 0x0;
1749 if (fileName.EndsWith("/")) {
1750 rawReader = new AliRawReaderFile(fileName);
1751 } else if (fileName.EndsWith(".root")) {
1752 rawReader = new AliRawReaderRoot(fileName);
1753 } else if (!fileName.IsNull()) {
1754 rawReader = new AliRawReaderDate(fileName);
1756 // if (!fEquipIdMap.IsNull() && fRawReader)
1757 // fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1759 // Get list of detectors
1760 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1763 AliHeader* header = runLoader->GetHeader();
1765 TString detStr = fMakeSDigits;
1769 if (!(rawReader->NextEvent())) break;
1772 for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1773 AliModule* det = (AliModule*) detArray->At(iDet);
1774 if (!det || !det->IsActive()) continue;
1775 if (IsSelected(det->GetName(), detStr)) {
1776 AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
1777 det->Raw2SDigits(rawReader);
1784 // If ESD information available obtain reconstructed vertex and store in header.
1786 treeESD->GetEvent(nev);
1787 const AliESDVertex* esdVertex = esd->GetPrimaryVertex();
1788 Double_t position[3];
1789 esdVertex->GetXYZ(position);
1790 AliGenEventHeader* mcHeader = new AliGenEventHeader("ESD");
1793 for (Int_t i = 0; i < 3; i++) mcV[i] = position[i];
1794 mcHeader->SetPrimaryVertex(mcV);
1795 header->Reset(0,nev);
1796 header->SetGenEventHeader(mcHeader);
1797 printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
1802 runLoader->TreeE()->Fill();
1803 runLoader->SetNextEvent();
1809 runLoader->CdGAFile();
1810 runLoader->WriteHeader("OVERWRITE");
1811 runLoader->WriteRunLoader();
1816 //_____________________________________________________________________________
1817 void AliSimulation::FinishRun()
1820 // Called at the end of the run.
1825 AliDebug(1, "Finish Lego");
1826 AliRunLoader::Instance()->CdGAFile();
1830 // Clean detector information
1831 TIter next(gAlice->Modules());
1832 AliModule *detector;
1833 while((detector = dynamic_cast<AliModule*>(next()))) {
1834 AliDebug(2, Form("%s->FinishRun()", detector->GetName()));
1835 detector->FinishRun();
1838 AliDebug(1, "AliRunLoader::Instance()->WriteHeader(OVERWRITE)");
1839 AliRunLoader::Instance()->WriteHeader("OVERWRITE");
1841 // Write AliRun info and all detectors parameters
1842 AliRunLoader::Instance()->CdGAFile();
1843 gAlice->Write(0,TObject::kOverwrite);//write AliRun
1844 AliRunLoader::Instance()->Write(0,TObject::kOverwrite);//write RunLoader itself
1846 if(gAlice->GetMCApp()) gAlice->GetMCApp()->FinishRun();
1847 AliRunLoader::Instance()->Synchronize();
1850 //_____________________________________________________________________________
1851 Int_t AliSimulation::GetDetIndex(const char* detector)
1853 // return the detector index corresponding to detector
1855 for (index = 0; index < fgkNDetectors ; index++) {
1856 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
1862 //_____________________________________________________________________________
1863 Bool_t AliSimulation::RunHLT()
1865 // Run the HLT simulation
1866 // HLT simulation is implemented in HLT/sim/AliHLTSimulation
1867 // Disabled if fRunHLT is empty, default vaule is "default".
1868 // AliSimulation::SetRunHLT can be used to set the options for HLT simulation
1869 // The default simulation depends on the HLT component libraries and their
1870 // corresponding agents which define components and chains to run. See
1871 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/
1872 // http://web.ift.uib.no/~kjeks/doc/alice-hlt-current/classAliHLTModuleAgent.html
1874 // The libraries to be loaded can be specified as an option.
1876 // AliSimulation sim;
1877 // sim.SetRunHLT("libAliHLTSample.so");
1879 // will only load <tt>libAliHLTSample.so</tt>
1881 // Other available options:
1882 // \li loglevel=<i>level</i> <br>
1883 // logging level for this processing
1885 // disable redirection of log messages to AliLog class
1886 // \li config=<i>macro</i>
1887 // configuration macro
1888 // \li chains=<i>configuration</i>
1889 // comma separated list of configurations to be run during simulation
1890 // \li rawfile=<i>file</i>
1891 // source for the RawReader to be created, the default is <i>./</i> if
1892 // raw data is simulated
1895 AliRunLoader* pRunLoader = LoadRun("READ");
1896 if (!pRunLoader) return kFALSE;
1898 // initialize CDB storage, run number, set CDB lock
1900 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
1903 // load the library dynamically
1904 gSystem->Load(ALIHLTSIMULATION_LIBRARY);
1906 // check for the library version
1907 AliHLTSimulationGetLibraryVersion_t fctVersion=(AliHLTSimulationGetLibraryVersion_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_GET_LIBRARY_VERSION));
1909 AliError(Form("can not load library %s", ALIHLTSIMULATION_LIBRARY));
1912 if (fctVersion()!= ALIHLTSIMULATION_LIBRARY_VERSION) {
1913 AliError(Form("%s version does not match: compiled for version %d, loaded %d", ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_LIBRARY_VERSION, fctVersion()));
1917 // print compile info
1918 typedef void (*CompileInfo)( char*& date, char*& time);
1919 CompileInfo fctInfo=(CompileInfo)gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, "CompileInfo");
1923 (*fctInfo)(date,time);
1924 if (!date) {date=new Char_t[8]; strcpy(date,"unknown");}
1925 if (!time) {time=new Char_t[8]; strcpy(time,"unknown");}
1926 AliInfo(Form("%s build on %s (%s)", ALIHLTSIMULATION_LIBRARY, date, time));
1930 AliInfo(Form("no build info available for %s", ALIHLTSIMULATION_LIBRARY));
1933 // create instance of the HLT simulation
1934 AliHLTSimulationCreateInstance_t fctCreate=(AliHLTSimulationCreateInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_CREATE_INSTANCE));
1935 AliHLTSimulation* pHLT=NULL;
1936 if (fctCreate==NULL || (pHLT=(fctCreate()))==NULL) {
1937 AliError(Form("can not create instance of HLT simulation (creator %p)", fctCreate));
1941 // init the HLT simulation
1943 if (fRunHLT.CompareTo("default")!=0) options=fRunHLT;
1944 TString detStr = fWriteRawData;
1945 if (!IsSelected("HLT", detStr)) {
1946 options+=" writerawfiles=";
1948 options+=" writerawfiles=HLT";
1951 if (!detStr.IsNull() && !options.Contains("rawfile=")) {
1952 // as a matter of fact, HLT will run reconstruction and needs the RawReader
1953 // in order to get detector data. By default, RawReaderFile is used to read
1954 // the already simulated ddl files. Date and Root files from the raw data
1955 // are generated after the HLT simulation.
1956 options+=" rawfile=./";
1959 AliHLTSimulationInit_t fctInit=(AliHLTSimulationInit_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_INIT));
1960 if (fctInit==NULL || (iResult=(fctInit(pHLT, pRunLoader, options.Data())))<0) {
1961 AliError(Form("can not init HLT simulation: error %d (init %p)", iResult, fctInit));
1963 // run the HLT simulation
1964 AliHLTSimulationRun_t fctRun=(AliHLTSimulationRun_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_RUN));
1965 if (fctRun==NULL || (iResult=(fctRun(pHLT, pRunLoader)))<0) {
1966 AliError(Form("can not run HLT simulation: error %d (run %p)", iResult, fctRun));
1970 // delete the instance
1971 AliHLTSimulationDeleteInstance_t fctDelete=(AliHLTSimulationDeleteInstance_t)(gSystem->DynFindSymbol(ALIHLTSIMULATION_LIBRARY, ALIHLTSIMULATION_DELETE_INSTANCE));
1972 if (fctDelete==NULL || fctDelete(pHLT)<0) {
1973 AliError(Form("can not delete instance of HLT simulation (creator %p)", fctDelete));
1977 return iResult>=0?kTRUE:kFALSE;
1980 //_____________________________________________________________________________
1981 Bool_t AliSimulation::RunQA()
1983 // run the QA on summable hits, digits or digits
1985 if(!gAlice) return kFALSE;
1986 fQAManager->SetRunLoader(AliRunLoader::Instance()) ;
1988 TString detectorsw("") ;
1990 fQAManager->SetEventSpecie(fEventSpecie) ;
1991 detectorsw = fQAManager->Run(fQADetectors.Data()) ;
1992 if ( detectorsw.IsNull() )
1997 //_____________________________________________________________________________
1998 Bool_t AliSimulation::SetRunQA(TString detAndAction)
2000 // Allows to run QA for a selected set of detectors
2001 // and a selected set of tasks among HITS, SDIGITS and DIGITS
2002 // all selected detectors run the same selected tasks
2004 if (!detAndAction.Contains(":")) {
2005 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2009 Int_t colon = detAndAction.Index(":") ;
2010 fQADetectors = detAndAction(0, colon) ;
2011 if (fQADetectors.Contains("ALL") )
2012 fQADetectors = Form("%s %s", fMakeDigits.Data(), fMakeDigitsFromHits.Data()) ;
2013 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2014 if (fQATasks.Contains("ALL") ) {
2015 fQATasks = Form("%d %d %d", AliQA::kHITS, AliQA::kSDIGITS, AliQA::kDIGITS) ;
2017 fQATasks.ToUpper() ;
2019 if ( fQATasks.Contains("HIT") )
2020 tempo = Form("%d ", AliQA::kHITS) ;
2021 if ( fQATasks.Contains("SDIGIT") )
2022 tempo += Form("%d ", AliQA::kSDIGITS) ;
2023 if ( fQATasks.Contains("DIGIT") )
2024 tempo += Form("%d ", AliQA::kDIGITS) ;
2026 if (fQATasks.IsNull()) {
2027 AliInfo("No QA requested\n") ;
2032 TString tempo(fQATasks) ;
2033 tempo.ReplaceAll(Form("%d", AliQA::kHITS), AliQA::GetTaskName(AliQA::kHITS)) ;
2034 tempo.ReplaceAll(Form("%d", AliQA::kSDIGITS), AliQA::GetTaskName(AliQA::kSDIGITS)) ;
2035 tempo.ReplaceAll(Form("%d", AliQA::kDIGITS), AliQA::GetTaskName(AliQA::kDIGITS)) ;
2036 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2038 fQAManager->SetActiveDetectors(fQADetectors) ;
2039 fQAManager->SetTasks(fQATasks) ;
2040 for (Int_t det = 0 ; det < AliQA::kNDET ; det++)
2041 fQAManager->SetWriteExpert(AliQA::DETECTORINDEX_t(det)) ;
2046 //_____________________________________________________________________________
2047 void AliSimulation::ProcessEnvironmentVars()
2049 // Extract run number and random generator seed from env variables
2051 AliInfo("Processing environment variables");
2053 // Random Number seed
2055 // first check that seed is not already set
2057 if (gSystem->Getenv("CONFIG_SEED")) {
2058 fSeed = atoi(gSystem->Getenv("CONFIG_SEED"));
2061 if (gSystem->Getenv("CONFIG_SEED")) {
2062 AliInfo(Form("Seed for random number generation already set (%d)"
2063 ": CONFIG_SEED variable ignored!", fSeed));
2067 AliInfo(Form("Seed for random number generation = %d ", fSeed));
2071 // first check that run number is not already set
2073 if (gSystem->Getenv("DC_RUN")) {
2074 fRun = atoi(gSystem->Getenv("DC_RUN"));
2077 if (gSystem->Getenv("DC_RUN")) {
2078 AliInfo(Form("Run number already set (%d): DC_RUN variable ignored!", fRun));
2082 AliInfo(Form("Run number = %d", fRun));
2085 //---------------------------------------------------------------------
2086 void AliSimulation::WriteGRPEntry()
2088 // Get the necessary information from galice (generator, trigger etc) and
2089 // write a GRP entry corresponding to the settings in the Config.C used
2090 // note that Hall probes and Cavern and Surface Atmos pressures are not simulated.
2093 AliInfo("Writing global run parameters entry into the OCDB");
2095 AliGRPObject* grpObj = new AliGRPObject();
2097 grpObj->SetRunType("PHYSICS");
2098 grpObj->SetTimeStart(0);
2099 grpObj->SetTimeEnd(9999);
2101 const AliGenerator *gen = gAlice->GetMCApp()->Generator();
2103 grpObj->SetBeamEnergy(gen->GetEnergyCMS());
2106 gen->GetProjectile(projectile,a,z);
2108 gen->GetTarget(target,a,z);
2109 TString beamType = projectile + "-" + target;
2110 beamType.ReplaceAll(" ","");
2111 if (!beamType.CompareTo("-")) {
2112 grpObj->SetBeamType("UNKNOWN");
2115 grpObj->SetBeamType(beamType);
2116 // Heavy ion run, the event specie is set to kHighMult
2117 fEventSpecie = AliRecoParam::kHighMult;
2118 if ((strcmp(beamType,"p-p") == 0) ||
2119 (strcmp(beamType,"p-") == 0) ||
2120 (strcmp(beamType,"-p") == 0) ||
2121 (strcmp(beamType,"P-P") == 0) ||
2122 (strcmp(beamType,"P-") == 0) ||
2123 (strcmp(beamType,"-P") == 0)) {
2124 // Proton run, the event specie is set to kLowMult
2125 fEventSpecie = AliRecoParam::kLowMult;
2129 AliWarning("Unknown beam type and energy! Setting energy to 0");
2130 grpObj->SetBeamEnergy(0);
2131 grpObj->SetBeamType("UNKNOWN");
2134 UInt_t detectorPattern = 0;
2136 TObjArray *detArray = gAlice->Detectors();
2137 for (Int_t iDet = 0; iDet < AliDAQ::kNDetectors-1; iDet++) {
2138 if (detArray->FindObject(AliDAQ::OfflineModuleName(iDet))) {
2139 detectorPattern |= (1 << iDet);
2144 if (!fMakeTrigger.IsNull() || strcmp(gAlice->GetTriggerDescriptor(),""))
2145 detectorPattern |= (1 << AliDAQ::DetectorID("TRG"));
2148 if (!fRunHLT.IsNull())
2149 detectorPattern |= (1 << AliDAQ::kHLTId);
2151 grpObj->SetNumberOfDetectors((Char_t)nDets);
2152 grpObj->SetDetectorMask((Int_t)detectorPattern);
2153 grpObj->SetLHCPeriod("LHC08c");
2154 grpObj->SetLHCState("STABLE_BEAMS");
2155 grpObj->SetLHCLuminosity(0,(AliGRPObject::Stats)0);
2156 grpObj->SetBeamIntensity(0,(AliGRPObject::Stats)0);
2158 AliMagF *field = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
2159 Float_t solenoidField = field ? TMath::Abs(field->SolenoidField()) : 0;
2160 Float_t factor = field ? field->Factor() : 0;
2161 Float_t l3current = TMath::Abs(factor)*solenoidField*30000./5.;
2162 grpObj->SetL3Current(l3current,(AliGRPObject::Stats)0);
2165 grpObj->SetL3Polarity(0);
2166 grpObj->SetDipolePolarity(0);
2169 grpObj->SetL3Polarity(1);
2170 grpObj->SetDipolePolarity(1);
2173 if (TMath::Abs(factor) != 0)
2174 grpObj->SetDipoleCurrent(6000,(AliGRPObject::Stats)0);
2176 grpObj->SetDipoleCurrent(0,(AliGRPObject::Stats)0);
2178 grpObj->SetCavernTemperature(0,(AliGRPObject::Stats)0);
2180 //grpMap->Add(new TObjString("fCavernPressure"),new TObjString("0")); ---> not inserted in simulation with the new object, since it is now an AliDCSSensor
2182 // Now store the entry in OCDB
2183 AliCDBManager* man = AliCDBManager::Instance();
2185 AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
2186 AliCDBMetaData *metadata= new AliCDBMetaData();
2188 metadata->SetResponsible("alice-off@cern.ch");
2189 metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
2191 man->Put(grpObj,id,metadata);