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 ///////////////////////////////////////////////////////////////////////////////
108 #include <TGeoManager.h>
109 #include <TObjString.h>
110 #include <TStopwatch.h>
114 #include "AliCDBStorage.h"
115 #include "AliCDBEntry.h"
116 #include "AliCDBManager.h"
117 #include "AliAlignObj.h"
118 #include "AliCentralTrigger.h"
119 #include "AliDAQConfig.h"
120 #include "AliDigitizer.h"
121 #include "AliGenerator.h"
123 #include "AliModule.h"
125 #include "AliRunDigitizer.h"
126 #include "AliRunLoader.h"
127 #include "AliSimulation.h"
128 #include "AliVertexGenFile.h"
129 #include "AliCentralTrigger.h"
130 #include "AliCTPRawData.h"
132 ClassImp(AliSimulation)
135 //_____________________________________________________________________________
136 AliSimulation::AliSimulation(const char* configFileName, const char* cdbUri,
137 const char* name, const char* title) :
140 fRunGeneration(kTRUE),
141 fRunSimulation(kTRUE),
142 fLoadAlignFromCDB(kTRUE),
143 fLoadAlignData("ALL"),
147 fMakeDigitsFromHits(""),
149 fRawDataFileName(""),
150 fDeleteIntermediateFiles(kFALSE),
151 fStopOnError(kFALSE),
154 fConfigFileName(configFileName),
155 fGAliceFileName("galice.root"),
157 fBkgrdFileNames(NULL),
158 fAlignObjArray(NULL),
159 fUseBkgrdVertex(kTRUE),
160 fRegionOfInterest(kFALSE),
163 // create simulation object with default parameters
165 SetGAliceFile("galice.root");
168 //_____________________________________________________________________________
169 AliSimulation::AliSimulation(const AliSimulation& sim) :
172 fRunGeneration(sim.fRunGeneration),
173 fRunSimulation(sim.fRunSimulation),
174 fLoadAlignFromCDB(sim.fLoadAlignFromCDB),
175 fLoadAlignData(sim.fLoadAlignData),
176 fMakeSDigits(sim.fMakeSDigits),
177 fMakeDigits(sim.fMakeDigits),
178 fMakeTrigger(sim.fMakeTrigger),
179 fMakeDigitsFromHits(sim.fMakeDigitsFromHits),
180 fWriteRawData(sim.fWriteRawData),
181 fRawDataFileName(""),
182 fDeleteIntermediateFiles(kFALSE),
183 fStopOnError(sim.fStopOnError),
185 fNEvents(sim.fNEvents),
186 fConfigFileName(sim.fConfigFileName),
187 fGAliceFileName(sim.fGAliceFileName),
189 fBkgrdFileNames(NULL),
190 fAlignObjArray(NULL),
191 fUseBkgrdVertex(sim.fUseBkgrdVertex),
192 fRegionOfInterest(sim.fRegionOfInterest),
197 for (Int_t i = 0; i < sim.fEventsPerFile.GetEntriesFast(); i++) {
198 if (!sim.fEventsPerFile[i]) continue;
199 fEventsPerFile.Add(sim.fEventsPerFile[i]->Clone());
202 fBkgrdFileNames = new TObjArray;
203 for (Int_t i = 0; i < sim.fBkgrdFileNames->GetEntriesFast(); i++) {
204 if (!sim.fBkgrdFileNames->At(i)) continue;
205 fBkgrdFileNames->Add(sim.fBkgrdFileNames->At(i)->Clone());
209 //_____________________________________________________________________________
210 AliSimulation& AliSimulation::operator = (const AliSimulation& sim)
212 // assignment operator
214 this->~AliSimulation();
215 new(this) AliSimulation(sim);
219 //_____________________________________________________________________________
220 AliSimulation::~AliSimulation()
224 fEventsPerFile.Delete();
225 // if(fAlignObjArray) fAlignObjArray->Delete(); // fAlignObjArray->RemoveAll() ???
226 // delete fAlignObjArray; fAlignObjArray=0;
228 if (fBkgrdFileNames) {
229 fBkgrdFileNames->Delete();
230 delete fBkgrdFileNames;
235 //_____________________________________________________________________________
236 void AliSimulation::SetNumberOfEvents(Int_t nEvents)
238 // set the number of events for one run
243 //_____________________________________________________________________________
244 void AliSimulation::InitCDBStorage(const char* uri)
246 // activate a default CDB storage
247 // First check if we have any CDB storage set, because it is used
248 // to retrieve the calibration and alignment constants
250 AliCDBManager* man = AliCDBManager::Instance();
251 if (!man->IsDefaultStorageSet())
253 AliWarningClass("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
254 AliWarningClass("Default CDB storage not yet set");
255 AliWarningClass(Form("Using default storage declared in AliSimulation: %s",uri));
256 AliWarningClass("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
257 SetDefaultStorage(uri);
262 //_____________________________________________________________________________
263 void AliSimulation::SetDefaultStorage(const char* uri) {
264 // activate a default CDB storage
266 AliCDBManager::Instance()->SetDefaultStorage(uri);
270 //_____________________________________________________________________________
271 void AliSimulation::SetSpecificStorage(const char* detName, const char* uri) {
272 // activate a detector-specific CDB storage
274 AliCDBManager::Instance()->SetSpecificStorage(detName, uri);
278 //_____________________________________________________________________________
279 void AliSimulation::SetConfigFile(const char* fileName)
281 // set the name of the config file
283 fConfigFileName = fileName;
286 //_____________________________________________________________________________
287 void AliSimulation::SetGAliceFile(const char* fileName)
289 // set the name of the galice file
290 // the path is converted to an absolute one if it is relative
292 fGAliceFileName = fileName;
293 if (!gSystem->IsAbsoluteFileName(fGAliceFileName)) {
294 char* absFileName = gSystem->ConcatFileName(gSystem->WorkingDirectory(),
296 fGAliceFileName = absFileName;
297 delete[] absFileName;
300 AliDebug(2, Form("galice file name set to %s", fileName));
303 //_____________________________________________________________________________
304 void AliSimulation::SetEventsPerFile(const char* detector, const char* type,
307 // set the number of events per file for the given detector and data type
308 // ("Hits", "Summable Digits", "Digits", "Reconstructed Points" or "Tracks")
310 TNamed* obj = new TNamed(detector, type);
311 obj->SetUniqueID(nEvents);
312 fEventsPerFile.Add(obj);
315 //_____________________________________________________________________________
316 Bool_t AliSimulation::ApplyAlignObjsToGeom(TObjArray* alObjArray)
318 // Read collection of alignment objects (AliAlignObj derived) saved
319 // in the TClonesArray ClArrayName and apply them to the geometry
320 // manager singleton.
323 Int_t nvols = alObjArray->GetEntriesFast();
327 for(Int_t j=0; j<nvols; j++)
329 AliAlignObj* alobj = (AliAlignObj*) alObjArray->UncheckedAt(j);
330 if (alobj->ApplyToGeometry() == kFALSE) flag = kFALSE;
333 if (AliDebugLevelClass() >= 1) {
334 gGeoManager->CheckOverlaps(20);
335 TObjArray* ovexlist = gGeoManager->GetListOfOverlaps();
336 if(ovexlist->GetEntriesFast()){
337 AliError("The application of alignment objects to the geometry caused huge overlaps/extrusions!");
345 //_____________________________________________________________________________
346 Bool_t AliSimulation::ApplyAlignObjsToGeom(const char* fileName, const char* clArrayName)
348 // read collection of alignment objects (AliAlignObj derived) saved
349 // in the TClonesArray ClArrayName in the file fileName and apply
350 // them to the TGeo geometry passed as argument
353 TFile* inFile = TFile::Open(fileName,"READ");
354 if (!inFile || !inFile->IsOpen()) {
355 AliErrorClass(Form("Could not open file %s !",fileName));
359 TClonesArray* alObjArray = ((TClonesArray*) inFile->Get(clArrayName));
362 AliErrorClass(Form("Could not get array (%s) from file (%s) !",clArrayName,fileName));
366 return ApplyAlignObjsToGeom(alObjArray);
370 //_____________________________________________________________________________
371 Bool_t AliSimulation::ApplyAlignObjsToGeom(AliCDBParam* param, AliCDBId& Id)
373 // read collection of alignment objects (AliAlignObj derived) saved
374 // in the TClonesArray ClArrayName in the AliCDBEntry identified by
375 // param (to get the AliCDBStorage) and Id; apply the alignment objects
376 // to the TGeo geometry passed as argument
379 AliCDBStorage* storage = AliCDBManager::Instance()->GetStorage(param);
380 AliCDBEntry* entry = storage->Get(Id);
381 TClonesArray* AlObjArray = ((TClonesArray*) entry->GetObject());
383 return ApplyAlignObjsToGeom(AlObjArray);
387 //_____________________________________________________________________________
388 Bool_t AliSimulation::ApplyAlignObjsToGeom(const char* uri, const char* path, Int_t runnum, Int_t version, Int_t sversion)
390 // read collection of alignment objects (AliAlignObj derived) saved
391 // in the TClonesArray ClArrayName in the AliCDBEntry identified by
392 // param (to get the AliCDBStorage) and Id; apply the alignment objects
393 // to the TGeo geometry passed as argument
396 AliCDBParam* param = AliCDBManager::Instance()->CreateParameter(uri);
397 AliCDBId id(path, runnum, runnum, version, sversion);
399 return ApplyAlignObjsToGeom(param, id);
403 //_____________________________________________________________________________
404 Bool_t AliSimulation::ApplyAlignObjsToGeom(const char* detName, Int_t runnum, Int_t version, Int_t sversion)
406 // read collection of alignment objects (AliAlignObj derived) saved
407 // in the TClonesArray ClArrayName in the AliCDBEntry identified by
408 // param (to get the AliCDBStorage) and Id; apply the alignment objects
409 // to the TGeo geometry passed as argument
412 InitCDBStorage("local://$ALICE_ROOT");
413 AliCDBPath path(detName,"Align","Data");
414 AliCDBEntry* entry = AliCDBManager::Instance()->Get(path.GetPath(),runnum,version,sversion);
416 if(!entry) return kFALSE;
417 TClonesArray* AlObjArray = ((TClonesArray*) entry->GetObject());
419 return ApplyAlignObjsToGeom(AlObjArray);
422 //_____________________________________________________________________________
423 Bool_t AliSimulation::SetAlignObjArraySingleDet(const char* detName)
425 // Fills array of single detector's alignable objects from CDB
427 AliDebug(2, Form("Loading alignment data for detector: %s",detName));
431 AliCDBPath path(detName,"Align","Data");
433 entry=AliCDBManager::Instance()->Get(path.GetPath());
435 AliDebug(2,Form("Couldn't load alignment data for detector %s",detName));
439 TClonesArray *alignArray = (TClonesArray*) entry->GetObject();
440 alignArray->SetOwner(0);
441 AliDebug(2,Form("Found %d alignment objects for %s",
442 alignArray->GetEntries(),detName));
444 AliAlignObj *alignObj=0;
445 TIter iter(alignArray);
447 // loop over align objects in detector
448 while( ( alignObj=(AliAlignObj *) iter.Next() ) ){
449 fAlignObjArray->Add(alignObj);
451 // delete entry --- Don't delete, it is cached!
453 AliDebug(2, Form("fAlignObjArray entries: %d",fAlignObjArray->GetEntries() ));
458 //_____________________________________________________________________________
459 Bool_t AliSimulation::MisalignGeometry(AliRunLoader *runLoader)
461 // Read the alignment objects from CDB.
462 // Each detector is supposed to have the
463 // alignment objects in DET/Align/Data CDB path.
464 // All the detector objects are then collected,
465 // sorted by geometry level (starting from ALIC) and
466 // then applied to the TGeo geometry.
467 // Finally an overlaps check is performed.
469 Bool_t delRunLoader = kFALSE;
471 runLoader = LoadRun("READ");
472 if (!runLoader) return kFALSE;
473 delRunLoader = kTRUE;
476 // Load alignment data from CDB and fill fAlignObjArray
477 if(fLoadAlignFromCDB){
478 if(!fAlignObjArray) fAlignObjArray = new TObjArray();
480 //fAlignObjArray->RemoveAll();
481 fAlignObjArray->Clear();
482 fAlignObjArray->SetOwner(0);
484 TString detStr = fLoadAlignData;
485 TString dataNotLoaded="";
486 TString dataLoaded="";
488 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
489 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
490 AliModule* det = (AliModule*) detArray->At(iDet);
491 if (!det || !det->IsActive()) continue;
492 if (IsSelected(det->GetName(), detStr)) {
493 if(!SetAlignObjArraySingleDet(det->GetName())){
494 dataNotLoaded += det->GetName();
495 dataNotLoaded += " ";
497 dataLoaded += det->GetName();
501 } // end loop over detectors
503 if ((detStr.CompareTo("ALL") == 0)) detStr = "";
504 dataNotLoaded += detStr;
505 AliInfo(Form("Alignment data loaded for: %s",
507 AliInfo(Form("Didn't/couldn't load alignment data for: %s",
508 dataNotLoaded.Data()));
509 } // fLoadAlignFromCDB flag
511 // Check if the array with alignment objects was
512 // provided by the user. If yes, apply the objects
513 // to the present TGeo geometry
514 if (fAlignObjArray) {
515 if (gGeoManager && gGeoManager->IsClosed()) {
516 if (ApplyAlignObjsToGeom(fAlignObjArray) == kFALSE) {
517 AliError("The misalignment of one or more volumes failed!"
518 "Compare the list of simulated detectors and the list of detector alignment data!");
519 if (delRunLoader) delete runLoader;
524 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
525 if (delRunLoader) delete runLoader;
530 if (delRunLoader) delete runLoader;
536 //_____________________________________________________________________________
537 Bool_t AliSimulation::SetRunNumber()
539 // Set the CDB manager run number
540 // The run number is retrieved from gAlice
542 if(AliCDBManager::Instance()->GetRun() < 0) {
543 AliRunLoader* runLoader = LoadRun("READ");
544 if (!runLoader) return kFALSE;
546 AliCDBManager::Instance()->SetRun(runLoader->GetAliRun()->GetRunNumber());
547 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
554 //_____________________________________________________________________________
555 void AliSimulation::MergeWith(const char* fileName, Int_t nSignalPerBkgrd)
557 // add a file with background events for merging
559 TObjString* fileNameStr = new TObjString(fileName);
560 fileNameStr->SetUniqueID(nSignalPerBkgrd);
561 if (!fBkgrdFileNames) fBkgrdFileNames = new TObjArray;
562 fBkgrdFileNames->Add(fileNameStr);
566 //_____________________________________________________________________________
567 Bool_t AliSimulation::Run(Int_t nEvents)
569 // run the generation, simulation and digitization
571 InitCDBStorage(fCDBUri);
573 if (nEvents > 0) fNEvents = nEvents;
575 // generation and simulation -> hits
576 if (fRunGeneration) {
577 if (!RunSimulation()) if (fStopOnError) return kFALSE;
580 // Set run number in CDBManager (if it is not already set in RunSimulation)
581 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
583 // Load and misalign the geometry
585 TGeoManager::Import("geometry.root");
586 if (!gGeoManager) if (fStopOnError) return kFALSE;
587 if (!MisalignGeometry()) if (fStopOnError) return kFALSE;
590 // hits -> summable digits
591 if (!fMakeSDigits.IsNull()) {
592 if (!RunSDigitization(fMakeSDigits)) if (fStopOnError) return kFALSE;
595 // summable digits -> digits
596 if (!fMakeDigits.IsNull()) {
597 if (!RunDigitization(fMakeDigits, fMakeDigitsFromHits)) {
598 if (fStopOnError) return kFALSE;
603 if (!fMakeDigitsFromHits.IsNull()) {
604 if (fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
605 AliWarning(Form("Merging and direct creation of digits from hits "
606 "was selected for some detectors. "
607 "No merging will be done for the following detectors: %s",
608 fMakeDigitsFromHits.Data()));
610 if (!RunHitsDigitization(fMakeDigitsFromHits)) {
611 if (fStopOnError) return kFALSE;
616 if (!fMakeTrigger.IsNull()) {
617 if (!RunTrigger(fMakeTrigger)) {
618 if (fStopOnError) return kFALSE;
622 // digits -> raw data
623 if (!fWriteRawData.IsNull()) {
624 if (!WriteRawData(fWriteRawData, fRawDataFileName,
625 fDeleteIntermediateFiles)) {
626 if (fStopOnError) return kFALSE;
633 //_____________________________________________________________________________
634 Bool_t AliSimulation::RunTrigger(const char* descriptors)
638 TStopwatch stopwatch;
641 AliRunLoader* runLoader = LoadRun("READ");
642 if (!runLoader) return kFALSE;
643 TString des = descriptors;
645 runLoader->MakeTree( "CT" );
646 AliCentralTrigger* aCTP = runLoader->GetTrigger();
648 aCTP->LoadDescriptor( des );
651 if( !aCTP->RunTrigger( runLoader ) ) {
658 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
659 stopwatch.RealTime(),stopwatch.CpuTime()));
666 //_____________________________________________________________________________
667 Bool_t AliSimulation::WriteTriggerRawData()
669 // Writes the CTP (trigger) DDL raw data
670 // Details of the format are given in the
671 // trigger TDR - pages 134 and 135.
672 AliCTPRawData writer;
678 //_____________________________________________________________________________
679 Bool_t AliSimulation::RunSimulation(Int_t nEvents)
681 // run the generation and simulation
683 TStopwatch stopwatch;
687 AliError("no gAlice object. Restart aliroot and try again.");
690 if (gAlice->Modules()->GetEntries() > 0) {
691 AliError("gAlice was already run. Restart aliroot and try again.");
695 AliInfo(Form("initializing gAlice with config file %s",
696 fConfigFileName.Data()));
697 StdoutToAliInfo(StderrToAliError(
698 gAlice->Init(fConfigFileName.Data());
701 // Set run number in CDBManager
702 AliCDBManager::Instance()->SetRun(gAlice->GetRunNumber());
703 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
705 AliRunLoader* runLoader = gAlice->GetRunLoader();
707 AliError(Form("gAlice has no run loader object. "
708 "Check your config file: %s", fConfigFileName.Data()));
711 SetGAliceFile(runLoader->GetFileName());
713 // Export ideal geometry
714 if (gGeoManager) gGeoManager->Export("geometry.root");
717 if (!MisalignGeometry(runLoader)) return kFALSE;
719 // Temporary fix by A.Gheata
720 // Could be removed with the next Root version (>5.11)
722 TIter next(gGeoManager->GetListOfVolumes());
724 while ((vol = (TGeoVolume *)next())) {
725 if (vol->GetVoxels()) {
726 if (vol->GetVoxels()->NeedRebuild()) {
727 vol->GetVoxels()->Voxelize();
734 // Export (mis)aligned geometry
735 if (gGeoManager) gGeoManager->Export("misaligned_geometry.root");
737 // AliRunLoader* runLoader = gAlice->GetRunLoader();
739 // AliError(Form("gAlice has no run loader object. "
740 // "Check your config file: %s", fConfigFileName.Data()));
743 // SetGAliceFile(runLoader->GetFileName());
745 if (!gAlice->Generator()) {
746 AliError(Form("gAlice has no generator object. "
747 "Check your config file: %s", fConfigFileName.Data()));
750 if (nEvents <= 0) nEvents = fNEvents;
752 // get vertex from background file in case of merging
753 if (fUseBkgrdVertex &&
754 fBkgrdFileNames && (fBkgrdFileNames->GetEntriesFast() > 0)) {
755 Int_t signalPerBkgrd = GetNSignalPerBkgrd(nEvents);
756 const char* fileName = ((TObjString*)
757 (fBkgrdFileNames->At(0)))->GetName();
758 AliInfo(Form("The vertex will be taken from the background "
759 "file %s with nSignalPerBackground = %d",
760 fileName, signalPerBkgrd));
761 AliVertexGenFile* vtxGen = new AliVertexGenFile(fileName, signalPerBkgrd);
762 gAlice->Generator()->SetVertexGenerator(vtxGen);
765 if (!fRunSimulation) {
766 gAlice->Generator()->SetTrackingFlag(0);
769 // set the number of events per file for given detectors and data types
770 for (Int_t i = 0; i < fEventsPerFile.GetEntriesFast(); i++) {
771 if (!fEventsPerFile[i]) continue;
772 const char* detName = fEventsPerFile[i]->GetName();
773 const char* typeName = fEventsPerFile[i]->GetTitle();
774 TString loaderName(detName);
775 loaderName += "Loader";
776 AliLoader* loader = runLoader->GetLoader(loaderName);
778 AliError(Form("RunSimulation", "no loader for %s found\n"
779 "Number of events per file not set for %s %s",
780 detName, typeName, detName));
783 AliDataLoader* dataLoader =
784 loader->GetDataLoader(typeName);
786 AliError(Form("no data loader for %s found\n"
787 "Number of events per file not set for %s %s",
788 typeName, detName, typeName));
791 dataLoader->SetNumberOfEventsPerFile(fEventsPerFile[i]->GetUniqueID());
792 AliDebug(1, Form("number of events per file set to %d for %s %s",
793 fEventsPerFile[i]->GetUniqueID(), detName, typeName));
796 AliInfo("running gAlice");
797 StdoutToAliInfo(StderrToAliError(
798 gAlice->Run(nEvents);
803 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
804 stopwatch.RealTime(),stopwatch.CpuTime()));
809 //_____________________________________________________________________________
810 Bool_t AliSimulation::RunSDigitization(const char* detectors)
812 // run the digitization and produce summable digits
814 TStopwatch stopwatch;
817 AliRunLoader* runLoader = LoadRun();
818 if (!runLoader) return kFALSE;
820 TString detStr = detectors;
821 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
822 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
823 AliModule* det = (AliModule*) detArray->At(iDet);
824 if (!det || !det->IsActive()) continue;
825 if (IsSelected(det->GetName(), detStr)) {
826 AliInfo(Form("creating summable digits for %s", det->GetName()));
827 TStopwatch stopwatchDet;
828 stopwatchDet.Start();
830 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
831 det->GetName(),stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
835 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
836 AliError(Form("the following detectors were not found: %s",
838 if (fStopOnError) return kFALSE;
843 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
844 stopwatch.RealTime(),stopwatch.CpuTime()));
850 //_____________________________________________________________________________
851 Bool_t AliSimulation::RunDigitization(const char* detectors,
852 const char* excludeDetectors)
854 // run the digitization and produce digits from sdigits
856 TStopwatch stopwatch;
859 while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
860 if (gAlice) delete gAlice;
864 if (fBkgrdFileNames) nStreams = fBkgrdFileNames->GetEntriesFast() + 1;
865 Int_t signalPerBkgrd = GetNSignalPerBkgrd();
866 AliRunDigitizer* manager = new AliRunDigitizer(nStreams, signalPerBkgrd);
867 manager->SetInputStream(0, fGAliceFileName.Data());
868 for (Int_t iStream = 1; iStream < nStreams; iStream++) {
869 const char* fileName = ((TObjString*)
870 (fBkgrdFileNames->At(iStream-1)))->GetName();
871 manager->SetInputStream(iStream, fileName);
874 TString detStr = detectors;
875 TString detExcl = excludeDetectors;
876 manager->GetInputStream(0)->ImportgAlice();
877 AliRunLoader* runLoader =
878 AliRunLoader::GetRunLoader(manager->GetInputStream(0)->GetFolderName());
879 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
880 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
881 AliModule* det = (AliModule*) detArray->At(iDet);
882 if (!det || !det->IsActive()) continue;
883 if (IsSelected(det->GetName(), detStr) &&
884 !IsSelected(det->GetName(), detExcl)) {
885 AliDigitizer* digitizer = det->CreateDigitizer(manager);
887 AliError(Form("no digitizer for %s", det->GetName()));
888 if (fStopOnError) return kFALSE;
890 digitizer->SetRegionOfInterest(fRegionOfInterest);
895 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
896 AliError(Form("the following detectors were not found: %s",
898 if (fStopOnError) return kFALSE;
901 if (!manager->GetListOfTasks()->IsEmpty()) {
902 AliInfo("executing digitization");
908 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
909 stopwatch.RealTime(),stopwatch.CpuTime()));
914 //_____________________________________________________________________________
915 Bool_t AliSimulation::RunHitsDigitization(const char* detectors)
917 // run the digitization and produce digits from hits
919 TStopwatch stopwatch;
922 AliRunLoader* runLoader = LoadRun("READ");
923 if (!runLoader) return kFALSE;
925 TString detStr = detectors;
926 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
927 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
928 AliModule* det = (AliModule*) detArray->At(iDet);
929 if (!det || !det->IsActive()) continue;
930 if (IsSelected(det->GetName(), detStr)) {
931 AliInfo(Form("creating digits from hits for %s", det->GetName()));
936 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
937 AliError(Form("the following detectors were not found: %s",
939 if (fStopOnError) return kFALSE;
943 //PH Temporary fix to avoid interference with the PHOS loder/getter
944 //PH The problem has to be solved in more general way 09/06/05
946 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
947 stopwatch.RealTime(),stopwatch.CpuTime()));
952 //_____________________________________________________________________________
953 Bool_t AliSimulation::WriteRawData(const char* detectors,
954 const char* fileName,
955 Bool_t deleteIntermediateFiles)
957 // convert the digits to raw data
958 // First DDL raw data files for the given detectors are created.
959 // If a file name is given, the DDL files are then converted to a DATE file.
960 // If deleteIntermediateFiles is true, the DDL raw files are deleted
962 // If the file name has the extension ".root", the DATE file is converted
964 // If deleteIntermediateFiles is true, the DATE file is deleted afterwards.
966 TStopwatch stopwatch;
969 if (!WriteRawFiles(detectors)) {
970 if (fStopOnError) return kFALSE;
973 TString dateFileName(fileName);
974 if (!dateFileName.IsNull()) {
975 Bool_t rootOutput = dateFileName.EndsWith(".root");
976 if (rootOutput) dateFileName += ".date";
977 if (!ConvertRawFilesToDate(dateFileName)) {
978 if (fStopOnError) return kFALSE;
980 if (deleteIntermediateFiles) {
981 AliRunLoader* runLoader = LoadRun("READ");
982 if (runLoader) for (Int_t iEvent = 0;
983 iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
985 sprintf(command, "rm -r raw%d", iEvent);
986 gSystem->Exec(command);
991 if (!ConvertDateToRoot(dateFileName, fileName)) {
992 if (fStopOnError) return kFALSE;
994 if (deleteIntermediateFiles) {
995 gSystem->Unlink(dateFileName);
1000 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
1001 stopwatch.RealTime(),stopwatch.CpuTime()));
1006 //_____________________________________________________________________________
1007 Bool_t AliSimulation::WriteRawFiles(const char* detectors)
1009 // convert the digits to raw data DDL files
1011 AliRunLoader* runLoader = LoadRun("READ");
1012 if (!runLoader) return kFALSE;
1014 // write raw data to DDL files
1015 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1016 AliInfo(Form("processing event %d", iEvent));
1017 runLoader->GetEvent(iEvent);
1018 TString baseDir = gSystem->WorkingDirectory();
1020 sprintf(dirName, "raw%d", iEvent);
1021 gSystem->MakeDirectory(dirName);
1022 if (!gSystem->ChangeDirectory(dirName)) {
1023 AliError(Form("couldn't change to directory %s", dirName));
1024 if (fStopOnError) return kFALSE; else continue;
1027 TString detStr = detectors;
1028 TObjArray* detArray = runLoader->GetAliRun()->Detectors();
1029 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
1030 AliModule* det = (AliModule*) detArray->At(iDet);
1031 if (!det || !det->IsActive()) continue;
1032 if (IsSelected(det->GetName(), detStr)) {
1033 AliInfo(Form("creating raw data from digits for %s", det->GetName()));
1038 if (!WriteTriggerRawData())
1039 if (fStopOnError) return kFALSE;
1041 gSystem->ChangeDirectory(baseDir);
1042 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1043 AliError(Form("the following detectors were not found: %s",
1045 if (fStopOnError) return kFALSE;
1053 //_____________________________________________________________________________
1054 Bool_t AliSimulation::ConvertRawFilesToDate(const char* dateFileName)
1056 // convert raw data DDL files to a DATE file with the program "dateStream"
1058 char* path = gSystem->Which(gSystem->Getenv("PATH"), "dateStream");
1060 AliError("the program dateStream was not found");
1061 if (fStopOnError) return kFALSE;
1066 AliRunLoader* runLoader = LoadRun("READ");
1067 if (!runLoader) return kFALSE;
1069 AliInfo(Form("converting raw data DDL files to DATE file %s", dateFileName));
1071 sprintf(command, "dateStream -o %s -# %d -C",
1072 dateFileName, runLoader->GetNumberOfEvents());
1073 FILE* pipe = gSystem->OpenPipe(command, "w");
1075 for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
1076 fprintf(pipe, "GDC\n");
1080 // loop over detectors and DDLs
1081 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1082 for (Int_t iDDL = 0; iDDL < kDetectorDDLs[iDet]; iDDL++) {
1084 Int_t ddlID = 0x100*iDet + iDDL;
1085 Int_t ldcID = Int_t(ldc + 0.0001);
1086 ldc += kDetectorLDCs[iDet] / kDetectorDDLs[iDet];
1088 char rawFileName[256];
1089 sprintf(rawFileName, "raw%d/%s_%d.ddl",
1090 iEvent, kDetectors[iDet], ddlID);
1092 // check existence and size of raw data file
1093 FILE* file = fopen(rawFileName, "rb");
1094 if (!file) continue;
1095 fseek(file, 0, SEEK_END);
1096 unsigned long size = ftell(file);
1098 if (!size) continue;
1100 if (ldcID != prevLDC) {
1101 fprintf(pipe, " LDC Id %d\n", ldcID);
1104 fprintf(pipe, " Equipment Id %d Payload %s\n", ddlID, rawFileName);
1109 Int_t result = gSystem->ClosePipe(pipe);
1112 return (result == 0);
1115 //_____________________________________________________________________________
1116 Bool_t AliSimulation::ConvertDateToRoot(const char* dateFileName,
1117 const char* rootFileName)
1119 // convert a DATE file to a root file with the program "alimdc"
1122 const Int_t kDBSize = 1000000000;
1123 const Int_t kTagDBSize = 1000000000;
1124 const Bool_t kFilter = kFALSE;
1125 const Int_t kCompression = 1;
1127 char* path = gSystem->Which(gSystem->Getenv("PATH"), "alimdc");
1129 AliError("the program alimdc was not found");
1130 if (fStopOnError) return kFALSE;
1135 AliInfo(Form("converting DATE file %s to root file %s",
1136 dateFileName, rootFileName));
1138 const char* rawDBFS[2] = { "/tmp/mdc1", "/tmp/mdc2" };
1139 const char* tagDBFS = "/tmp/mdc1/tags";
1140 const char* runDBFS = "/tmp/mdc1/meta";
1142 // User defined file system locations
1143 if (gSystem->Getenv("ALIMDC_RAWDB1"))
1144 rawDBFS[0] = gSystem->Getenv("ALIMDC_RAWDB1");
1145 if (gSystem->Getenv("ALIMDC_RAWDB2"))
1146 rawDBFS[1] = gSystem->Getenv("ALIMDC_RAWDB2");
1147 if (gSystem->Getenv("ALIMDC_TAGDB"))
1148 tagDBFS = gSystem->Getenv("ALIMDC_TAGDB");
1149 if (gSystem->Getenv("ALIMDC_RUNDB"))
1150 runDBFS = gSystem->Getenv("ALIMDC_RUNDB");
1152 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1153 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1154 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1155 gSystem->Exec(Form("rm -rf %s",runDBFS));
1157 gSystem->Exec(Form("mkdir %s",rawDBFS[0]));
1158 gSystem->Exec(Form("mkdir %s",rawDBFS[1]));
1159 gSystem->Exec(Form("mkdir %s",tagDBFS));
1160 gSystem->Exec(Form("mkdir %s",runDBFS));
1162 Int_t result = gSystem->Exec(Form("alimdc %d %d %d %d %s",
1163 kDBSize, kTagDBSize, kFilter, kCompression, dateFileName));
1164 gSystem->Exec(Form("mv %s/*.root %s", rawDBFS[0], rootFileName));
1166 gSystem->Exec(Form("rm -rf %s",rawDBFS[0]));
1167 gSystem->Exec(Form("rm -rf %s",rawDBFS[1]));
1168 gSystem->Exec(Form("rm -rf %s",tagDBFS));
1169 gSystem->Exec(Form("rm -rf %s",runDBFS));
1171 return (result == 0);
1175 //_____________________________________________________________________________
1176 AliRunLoader* AliSimulation::LoadRun(const char* mode) const
1178 // delete existing run loaders, open a new one and load gAlice
1180 while (AliRunLoader::GetRunLoader()) delete AliRunLoader::GetRunLoader();
1181 AliRunLoader* runLoader =
1182 AliRunLoader::Open(fGAliceFileName.Data(),
1183 AliConfig::GetDefaultEventFolderName(), mode);
1185 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1188 runLoader->LoadgAlice();
1189 gAlice = runLoader->GetAliRun();
1191 AliError(Form("no gAlice object found in file %s",
1192 fGAliceFileName.Data()));
1198 //_____________________________________________________________________________
1199 Int_t AliSimulation::GetNSignalPerBkgrd(Int_t nEvents) const
1201 // get or calculate the number of signal events per background event
1203 if (!fBkgrdFileNames) return 1;
1204 Int_t nBkgrdFiles = fBkgrdFileNames->GetEntriesFast();
1205 if (nBkgrdFiles == 0) return 1;
1207 // get the number of signal events
1209 AliRunLoader* runLoader =
1210 AliRunLoader::Open(fGAliceFileName.Data(), "SIGNAL");
1211 if (!runLoader) return 1;
1212 nEvents = runLoader->GetNumberOfEvents();
1217 for (Int_t iBkgrdFile = 0; iBkgrdFile < nBkgrdFiles; iBkgrdFile++) {
1218 // get the number of background events
1219 const char* fileName = ((TObjString*)
1220 (fBkgrdFileNames->At(iBkgrdFile)))->GetName();
1221 AliRunLoader* runLoader =
1222 AliRunLoader::Open(fileName, "BKGRD");
1223 if (!runLoader) continue;
1224 Int_t nBkgrdEvents = runLoader->GetNumberOfEvents();
1227 // get or calculate the number of signal per background events
1228 Int_t nSignalPerBkgrd = fBkgrdFileNames->At(iBkgrdFile)->GetUniqueID();
1229 if (nSignalPerBkgrd <= 0) {
1230 nSignalPerBkgrd = (nEvents-1) / nBkgrdEvents + 1;
1231 } else if (result && (result != nSignalPerBkgrd)) {
1232 AliInfo(Form("the number of signal events per background event "
1233 "will be changed from %d to %d for stream %d",
1234 nSignalPerBkgrd, result, iBkgrdFile+1));
1235 nSignalPerBkgrd = result;
1238 if (!result) result = nSignalPerBkgrd;
1239 if (nSignalPerBkgrd * nBkgrdEvents < nEvents) {
1240 AliWarning(Form("not enough background events (%d) for %d signal events "
1241 "using %d signal per background events for stream %d",
1242 nBkgrdEvents, nEvents, nSignalPerBkgrd, iBkgrdFile+1));
1249 //_____________________________________________________________________________
1250 Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
1252 // check whether detName is contained in detectors
1253 // if yes, it is removed from detectors
1255 // check if all detectors are selected
1256 if ((detectors.CompareTo("ALL") == 0) ||
1257 detectors.BeginsWith("ALL ") ||
1258 detectors.EndsWith(" ALL") ||
1259 detectors.Contains(" ALL ")) {
1264 // search for the given detector
1265 Bool_t result = kFALSE;
1266 if ((detectors.CompareTo(detName) == 0) ||
1267 detectors.BeginsWith(detName+" ") ||
1268 detectors.EndsWith(" "+detName) ||
1269 detectors.Contains(" "+detName+" ")) {
1270 detectors.ReplaceAll(detName, "");
1274 // clean up the detectors string
1275 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1276 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1277 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);