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 the reconstruction //
22 // Clusters and tracks are created for all detectors and all events by //
25 // AliReconstruction rec; //
28 // The Run method returns kTRUE in case of successful execution. //
30 // If the input to the reconstruction are not simulated digits but raw data, //
31 // this can be specified by an argument of the Run method or by the method //
33 // rec.SetInput("..."); //
35 // The input formats and the corresponding argument are: //
36 // - DDL raw data files: directory name, ends with "/" //
37 // - raw data root file: root file name, extension ".root" //
38 // - raw data DATE file: DATE file name, any other non-empty string //
39 // - MC root files : empty string, default //
41 // By default all events are reconstructed. The reconstruction can be //
42 // limited to a range of events by giving the index of the first and the //
43 // last event as an argument to the Run method or by calling //
45 // rec.SetEventRange(..., ...); //
47 // The index -1 (default) can be used for the last event to indicate no //
48 // upper limit of the event range. //
50 // In case of raw-data reconstruction the user can modify the default //
51 // number of events per digits/clusters/tracks file. In case the option //
52 // is not used the number is set 1. In case the user provides 0, than //
53 // the number of events is equal to the number of events inside the //
54 // raw-data file (i.e. one digits/clusters/tracks file): //
56 // rec.SetNumberOfEventsPerFile(...); //
59 // The name of the galice file can be changed from the default //
60 // "galice.root" by passing it as argument to the AliReconstruction //
61 // constructor or by //
63 // rec.SetGAliceFile("..."); //
65 // The local reconstruction can be switched on or off for individual //
68 // rec.SetRunLocalReconstruction("..."); //
70 // The argument is a (case sensitive) string with the names of the //
71 // detectors separated by a space. The special string "ALL" selects all //
72 // available detectors. This is the default. //
74 // The reconstruction of the primary vertex position can be switched off by //
76 // rec.SetRunVertexFinder(kFALSE); //
78 // The tracking and the creation of ESD tracks can be switched on for //
79 // selected detectors by //
81 // rec.SetRunTracking("..."); //
83 // Uniform/nonuniform field tracking switches (default: uniform field) //
85 // rec.SetUniformFieldTracking(); ( rec.SetUniformFieldTracking(kFALSE); ) //
87 // The filling of additional ESD information can be steered by //
89 // rec.SetFillESD("..."); //
91 // Again, for both methods the string specifies the list of detectors. //
92 // The default is "ALL". //
94 // The call of the shortcut method //
96 // rec.SetRunReconstruction("..."); //
98 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
99 // SetFillESD with the same detector selecting string as argument. //
101 // The reconstruction requires digits or raw data as input. For the creation //
102 // of digits and raw data have a look at the class AliSimulation. //
104 // The input data of a detector can be replaced by the corresponding HLT //
105 // data by calling (usual detector string) //
106 // SetUseHLTData("..."); //
109 ///////////////////////////////////////////////////////////////////////////////
116 #include <TPluginManager.h>
117 #include <TGeoManager.h>
118 #include <TLorentzVector.h>
121 #include <TObjArray.h>
124 #include "AliReconstruction.h"
125 #include "AliCodeTimer.h"
126 #include "AliReconstructor.h"
128 #include "AliRunLoader.h"
130 #include "AliRawReaderFile.h"
131 #include "AliRawReaderDate.h"
132 #include "AliRawReaderRoot.h"
133 #include "AliRawEventHeaderBase.h"
134 #include "AliESDEvent.h"
135 #include "AliESDMuonTrack.h"
136 #include "AliESDfriend.h"
137 #include "AliESDVertex.h"
138 #include "AliESDcascade.h"
139 #include "AliESDkink.h"
140 #include "AliESDtrack.h"
141 #include "AliESDCaloCluster.h"
142 #include "AliESDCaloCells.h"
143 #include "AliMultiplicity.h"
144 #include "AliTracker.h"
145 #include "AliVertexer.h"
146 #include "AliVertexerTracks.h"
147 #include "AliV0vertexer.h"
148 #include "AliCascadeVertexer.h"
149 #include "AliHeader.h"
150 #include "AliGenEventHeader.h"
152 #include "AliESDpid.h"
153 #include "AliESDtrack.h"
154 #include "AliESDPmdTrack.h"
156 #include "AliESDTagCreator.h"
157 #include "AliAODTagCreator.h"
159 #include "AliGeomManager.h"
160 #include "AliTrackPointArray.h"
161 #include "AliCDBManager.h"
162 #include "AliCDBStorage.h"
163 #include "AliCDBEntry.h"
164 #include "AliAlignObj.h"
166 #include "AliCentralTrigger.h"
167 #include "AliTriggerConfiguration.h"
168 #include "AliTriggerClass.h"
169 #include "AliTriggerCluster.h"
170 #include "AliCTPRawStream.h"
172 #include "AliQADataMakerRec.h"
173 #include "AliGlobalQADataMaker.h"
175 #include "AliQADataMakerSteer.h"
177 #include "AliPlaneEff.h"
179 #include "AliSysInfo.h" // memory snapshots
180 #include "AliRawHLTManager.h"
182 #include "AliMagWrapCheb.h"
184 #include "AliDetectorRecoParam.h"
185 #include "AliRunInfo.h"
186 #include "AliEventInfo.h"
190 ClassImp(AliReconstruction)
193 //_____________________________________________________________________________
194 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
196 //_____________________________________________________________________________
197 AliReconstruction::AliReconstruction(const char* gAliceFilename,
198 const char* name, const char* title) :
201 fUniformField(kFALSE),
202 fForcedFieldMap(0x0),
203 fRunVertexFinder(kTRUE),
204 fRunVertexFinderTracks(kTRUE),
205 fRunHLTTracking(kFALSE),
206 fRunMuonTracking(kFALSE),
208 fRunCascadeFinder(kTRUE),
209 fStopOnError(kFALSE),
210 fWriteAlignmentData(kFALSE),
211 fWriteESDfriend(kFALSE),
213 fFillTriggerESD(kTRUE),
221 fRunLocalReconstruction("ALL"),
224 fUseTrackingErrorsForAlignment(""),
225 fGAliceFileName(gAliceFilename),
230 fNumberOfEventsPerFile(1),
232 fLoadAlignFromCDB(kTRUE),
233 fLoadAlignData("ALL"),
241 fParentRawReader(NULL),
246 fDiamondProfile(NULL),
247 fDiamondProfileTPC(NULL),
248 fMeanVertexConstraint(kTRUE),
252 fAlignObjArray(NULL),
255 fInitCDBCalled(kFALSE),
256 fSetRunNumberFromDataCalled(kFALSE),
263 fSameQACycle(kFALSE),
265 fRunPlaneEff(kFALSE),
274 fIsNewRunLoader(kFALSE),
277 // create reconstruction object with default parameters
279 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
280 fReconstructor[iDet] = NULL;
281 fLoader[iDet] = NULL;
282 fTracker[iDet] = NULL;
283 fQACycles[iDet] = 999999;
285 fQASteer = new AliQADataMakerSteer("rec") ;
286 fQASteer->SetActiveDetectors(fQADetectors) ;
287 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
288 fQASteer->SetTasks(fQATasks) ;
292 //_____________________________________________________________________________
293 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
296 fUniformField(rec.fUniformField),
297 fForcedFieldMap(0x0),
298 fRunVertexFinder(rec.fRunVertexFinder),
299 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
300 fRunHLTTracking(rec.fRunHLTTracking),
301 fRunMuonTracking(rec.fRunMuonTracking),
302 fRunV0Finder(rec.fRunV0Finder),
303 fRunCascadeFinder(rec.fRunCascadeFinder),
304 fStopOnError(rec.fStopOnError),
305 fWriteAlignmentData(rec.fWriteAlignmentData),
306 fWriteESDfriend(rec.fWriteESDfriend),
307 fWriteAOD(rec.fWriteAOD),
308 fFillTriggerESD(rec.fFillTriggerESD),
310 fCleanESD(rec.fCleanESD),
311 fV0DCAmax(rec.fV0DCAmax),
312 fV0CsPmin(rec.fV0CsPmin),
316 fRunLocalReconstruction(rec.fRunLocalReconstruction),
317 fRunTracking(rec.fRunTracking),
318 fFillESD(rec.fFillESD),
319 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
320 fGAliceFileName(rec.fGAliceFileName),
322 fEquipIdMap(rec.fEquipIdMap),
323 fFirstEvent(rec.fFirstEvent),
324 fLastEvent(rec.fLastEvent),
325 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
327 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
328 fLoadAlignData(rec.fLoadAlignData),
329 fESDPar(rec.fESDPar),
330 fUseHLTData(rec.fUseHLTData),
336 fParentRawReader(NULL),
341 fDiamondProfile(NULL),
342 fDiamondProfileTPC(NULL),
343 fMeanVertexConstraint(rec.fMeanVertexConstraint),
347 fAlignObjArray(rec.fAlignObjArray),
348 fCDBUri(rec.fCDBUri),
350 fInitCDBCalled(rec.fInitCDBCalled),
351 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
352 fQADetectors(rec.fQADetectors),
353 fQASteer(rec.fQASteer),
354 fQATasks(rec.fQATasks),
356 fRunGlobalQA(rec.fRunGlobalQA),
357 fInLoopQA(rec.fInLoopQA),
358 fSameQACycle(rec.fSameQACycle),
359 fRunPlaneEff(rec.fRunPlaneEff),
368 fIsNewRunLoader(rec.fIsNewRunLoader),
373 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
374 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
376 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
377 fReconstructor[iDet] = NULL;
378 fLoader[iDet] = NULL;
379 fTracker[iDet] = NULL;
380 fQACycles[iDet] = rec.fQACycles[iDet];
382 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
383 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
386 fForcedFieldMap=new AliMagWrapCheb(*((AliMagWrapCheb*)rec.fForcedFieldMap));
389 //_____________________________________________________________________________
390 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
392 // assignment operator
394 this->~AliReconstruction();
395 new(this) AliReconstruction(rec);
399 //_____________________________________________________________________________
400 AliReconstruction::~AliReconstruction()
406 fSpecCDBUri.Delete();
407 delete fForcedFieldMap;
409 AliCodeTimer::Instance()->Print();
412 //_____________________________________________________________________________
413 void AliReconstruction::InitCDB()
415 // activate a default CDB storage
416 // First check if we have any CDB storage set, because it is used
417 // to retrieve the calibration and alignment constants
419 if (fInitCDBCalled) return;
420 fInitCDBCalled = kTRUE;
422 AliCDBManager* man = AliCDBManager::Instance();
423 if (man->IsDefaultStorageSet())
425 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
426 AliWarning("Default CDB storage has been already set !");
427 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
428 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
429 fCDBUri = man->GetDefaultStorage()->GetURI();
432 if (fCDBUri.Length() > 0)
434 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
435 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
436 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
438 fCDBUri="local://$ALICE_ROOT";
439 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
440 AliWarning("Default CDB storage not yet set !!!!");
441 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
442 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
445 man->SetDefaultStorage(fCDBUri);
448 // Now activate the detector specific CDB storage locations
449 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
450 TObject* obj = fSpecCDBUri[i];
452 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
453 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
454 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
455 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
460 //_____________________________________________________________________________
461 void AliReconstruction::SetDefaultStorage(const char* uri) {
462 // Store the desired default CDB storage location
463 // Activate it later within the Run() method
469 //_____________________________________________________________________________
470 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
471 // Store a detector-specific CDB storage location
472 // Activate it later within the Run() method
474 AliCDBPath aPath(calibType);
475 if(!aPath.IsValid()){
476 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
477 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
478 if(!strcmp(calibType, fgkDetectorName[iDet])) {
479 aPath.SetPath(Form("%s/*", calibType));
480 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
484 if(!aPath.IsValid()){
485 AliError(Form("Not a valid path or detector: %s", calibType));
490 // // check that calibType refers to a "valid" detector name
491 // Bool_t isDetector = kFALSE;
492 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
493 // TString detName = fgkDetectorName[iDet];
494 // if(aPath.GetLevel0() == detName) {
495 // isDetector = kTRUE;
501 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
505 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
506 if (obj) fSpecCDBUri.Remove(obj);
507 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
511 //_____________________________________________________________________________
512 Bool_t AliReconstruction::SetRunNumberFromData()
514 // The method is called in Run() in order
515 // to set a correct run number.
516 // In case of raw data reconstruction the
517 // run number is taken from the raw data header
519 if (fSetRunNumberFromDataCalled) return kTRUE;
520 fSetRunNumberFromDataCalled = kTRUE;
522 AliCDBManager* man = AliCDBManager::Instance();
524 if(man->GetRun() > 0) {
525 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
529 AliError("No run loader is found !");
532 // read run number from gAlice
533 if(fRunLoader->GetAliRun())
534 AliCDBManager::Instance()->SetRun(fRunLoader->GetHeader()->GetRun());
537 if(fRawReader->NextEvent()) {
538 AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
539 fRawReader->RewindEvents();
542 if(man->GetRun() > 0) {
543 AliWarning("No raw events is found ! Using settings in AliCDBManager !");
548 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
554 AliError("Neither gAlice nor RawReader objects are found !");
564 //_____________________________________________________________________________
565 void AliReconstruction::SetCDBLock() {
566 // Set CDB lock: from now on it is forbidden to reset the run number
567 // or the default storage or to activate any further storage!
569 AliCDBManager::Instance()->SetLock(1);
572 //_____________________________________________________________________________
573 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
575 // Read the alignment objects from CDB.
576 // Each detector is supposed to have the
577 // alignment objects in DET/Align/Data CDB path.
578 // All the detector objects are then collected,
579 // sorted by geometry level (starting from ALIC) and
580 // then applied to the TGeo geometry.
581 // Finally an overlaps check is performed.
583 // Load alignment data from CDB and fill fAlignObjArray
584 if(fLoadAlignFromCDB){
586 TString detStr = detectors;
587 TString loadAlObjsListOfDets = "";
589 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
590 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
591 loadAlObjsListOfDets += fgkDetectorName[iDet];
592 loadAlObjsListOfDets += " ";
593 } // end loop over detectors
594 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
595 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
597 // Check if the array with alignment objects was
598 // provided by the user. If yes, apply the objects
599 // to the present TGeo geometry
600 if (fAlignObjArray) {
601 if (gGeoManager && gGeoManager->IsClosed()) {
602 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
603 AliError("The misalignment of one or more volumes failed!"
604 "Compare the list of simulated detectors and the list of detector alignment data!");
609 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
615 delete fAlignObjArray; fAlignObjArray=0;
620 //_____________________________________________________________________________
621 void AliReconstruction::SetGAliceFile(const char* fileName)
623 // set the name of the galice file
625 fGAliceFileName = fileName;
628 //_____________________________________________________________________________
629 void AliReconstruction::SetInput(const char* input)
631 // In case the input string starts with 'mem://', we run in an online mode
632 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
633 // file is assumed. One can give as an input:
634 // mem://: - events taken from DAQ monitoring libs online
636 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
640 //_____________________________________________________________________________
641 void AliReconstruction::SetOption(const char* detector, const char* option)
643 // set options for the reconstruction of a detector
645 TObject* obj = fOptions.FindObject(detector);
646 if (obj) fOptions.Remove(obj);
647 fOptions.Add(new TNamed(detector, option));
650 //_____________________________________________________________________________
651 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
653 // Set custom reconstruction parameters for a given detector
654 // Single set of parameters for all the events
655 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
656 if(!strcmp(detector, fgkDetectorName[iDet])) {
658 fRecoParam.AddDetRecoParam(iDet,par);
665 //_____________________________________________________________________________
666 Bool_t AliReconstruction::SetFieldMap(Float_t l3Current, Float_t diCurrent, Float_t factor, const char *path) {
667 //------------------------------------------------
668 // The magnetic field map, defined externally...
669 // L3 current 30000 A -> 0.5 T
670 // L3 current 12000 A -> 0.2 T
671 // dipole current 6000 A
672 // The polarities must be the same
673 //------------------------------------------------
674 const Float_t l3NominalCurrent1=30000.; // (A)
675 const Float_t l3NominalCurrent2=12000.; // (A)
676 const Float_t diNominalCurrent =6000. ; // (A)
678 const Float_t tolerance=0.03; // relative current tolerance
679 const Float_t zero=77.; // "zero" current (A)
682 Bool_t dipoleON=kFALSE;
684 TString s=(factor < 0) ? "L3: -" : "L3: +";
686 l3Current = TMath::Abs(l3Current);
687 if (TMath::Abs(l3Current-l3NominalCurrent1)/l3NominalCurrent1 < tolerance) {
688 map=AliMagWrapCheb::k5kG;
691 if (TMath::Abs(l3Current-l3NominalCurrent2)/l3NominalCurrent2 < tolerance) {
692 map=AliMagWrapCheb::k2kG;
695 if (l3Current < zero) {
696 map=AliMagWrapCheb::k2kG;
698 factor=0.; // in fact, this is a global factor...
699 fUniformField=kTRUE; // track with the uniform (zero) B field
701 AliError(Form("Wrong L3 current (%f A)!",l3Current));
705 diCurrent = TMath::Abs(diCurrent);
706 if (TMath::Abs(diCurrent-diNominalCurrent)/diNominalCurrent < tolerance) {
707 // 3% current tolerance...
711 if (diCurrent < zero) { // some small current..
715 AliError(Form("Wrong dipole current (%f A)!",diCurrent));
719 delete fForcedFieldMap;
721 new AliMagWrapCheb("B field map ",s,2,factor,10.,map,dipoleON,path);
723 fForcedFieldMap->Print();
725 AliTracker::SetFieldMap(fForcedFieldMap,fUniformField);
731 Bool_t AliReconstruction::InitGRP() {
732 //------------------------------------
733 // Initialization of the GRP entry
734 //------------------------------------
735 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
737 if (entry) fGRPData = dynamic_cast<TMap*>(entry->GetObject());
740 AliError("No GRP entry found in OCDB!");
744 TObjString *lhcState=
745 dynamic_cast<TObjString*>(fGRPData->GetValue("fLHCState"));
747 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
750 TObjString *beamType=
751 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamType"));
753 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
756 TObjString *beamEnergyStr=
757 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamEnergy"));
758 if (!beamEnergyStr) {
759 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
763 dynamic_cast<TObjString*>(fGRPData->GetValue("fRunType"));
765 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
768 TObjString *activeDetectors=
769 dynamic_cast<TObjString*>(fGRPData->GetValue("fDetectorMask"));
770 if (!activeDetectors) {
771 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
774 fRunInfo = new AliRunInfo(lhcState ? lhcState->GetString().Data() : "UNKNOWN",
775 beamType ? beamType->GetString().Data() : "UNKNOWN",
776 beamEnergyStr ? beamEnergyStr->GetString().Atof() : 0,
777 runType ? runType->GetString().Data() : "UNKNOWN",
778 activeDetectors ? activeDetectors->GetString().Atoi() : 1074790399);
780 // Process the list of active detectors
781 if (activeDetectors && activeDetectors->GetString().IsDigit()) {
782 UInt_t detMask = activeDetectors->GetString().Atoi();
783 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
784 fRunTracking = MatchDetectorList(fRunTracking,detMask);
785 fFillESD = MatchDetectorList(fFillESD,detMask);
786 fQADetectors = MatchDetectorList(fQADetectors,detMask);
787 fQASteer->SetActiveDetectors(fQADetectors) ;
790 AliInfo("===================================================================================");
791 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
792 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
793 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
794 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
795 AliInfo("===================================================================================");
797 //*** Dealing with the magnetic field map
798 if (AliTracker::GetFieldMap()) {
799 AliInfo("Running with the externally set B field !");
801 // Construct the field map out of the information retrieved from GRP.
806 TObjString *l3Current=
807 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Current"));
809 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
812 TObjString *l3Polarity=
813 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Polarity"));
815 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
820 TObjString *diCurrent=
821 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipoleCurrent"));
823 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
826 TObjString *diPolarity=
827 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipolePolarity"));
829 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
834 Float_t l3Cur=TMath::Abs(atof(l3Current->GetName()));
835 Float_t diCur=TMath::Abs(atof(diCurrent->GetName()));
836 Float_t l3Pol=atof(l3Polarity->GetName());
838 if (l3Pol != 0.) factor=-1.;
841 if (!SetFieldMap(l3Cur, diCur, factor)) {
842 AliFatal("Failed to creat a B field map ! Exiting...");
844 AliInfo("Running with the B field constructed out of GRP !");
847 AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
853 //*** Get the diamond profile from OCDB
854 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
856 if (fMeanVertexConstraint)
857 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
859 AliError("No diamond profile found in OCDB!");
862 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
864 if (fMeanVertexConstraint)
865 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
867 AliError("No diamond profile found in OCDB!");
873 //_____________________________________________________________________________
874 Bool_t AliReconstruction::Run(const char* input)
877 AliCodeTimerAuto("");
879 if (!InitRun(input)) return kFALSE;
880 //******* The loop over events
882 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
883 (fRawReader && fRawReader->NextEvent())) {
884 if (!RunEvent(iEvent)) return kFALSE;
888 if (!FinishRun()) return kFALSE;
893 //_____________________________________________________________________________
894 Bool_t AliReconstruction::InitRun(const char* input)
896 // Initialize all the stuff before
897 // going into the event loop
898 // If the second argument is given, the first one is ignored and
899 // the reconstruction works in an online mode
900 AliCodeTimerAuto("");
902 // Overwrite the previous setting
903 if (input) fInput = input;
905 // set the input in case of raw data
906 fRawReader = AliRawReader::Create(fInput.Data());
908 AliInfo("Reconstruction will run over digits");
910 if (!fEquipIdMap.IsNull() && fRawReader)
911 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
913 if (!fUseHLTData.IsNull()) {
914 // create the RawReaderHLT which performs redirection of HLT input data for
915 // the specified detectors
916 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
918 fParentRawReader=fRawReader;
919 fRawReader=pRawReader;
921 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
925 AliSysInfo::AddStamp("Start");
926 // get the run loader
927 if (!InitRunLoader()) return kFALSE;
928 AliSysInfo::AddStamp("LoadLoader");
930 // Initialize the CDB storage
933 AliSysInfo::AddStamp("LoadCDB");
935 // Set run number in CDBManager (if it is not already set by the user)
936 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
938 // Set CDB lock: from now on it is forbidden to reset the run number
939 // or the default storage or to activate any further storage!
942 // Import ideal TGeo geometry and apply misalignment
944 TString geom(gSystem->DirName(fGAliceFileName));
945 geom += "/geometry.root";
946 AliGeomManager::LoadGeometry(geom.Data());
948 TString detsToCheck=fRunLocalReconstruction;
949 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data()))
950 AliFatalClass("Current loaded geometry differs in the definition of symbolic names!");
951 if (!gGeoManager) if (fStopOnError) return kFALSE;
954 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
955 AliSysInfo::AddStamp("LoadGeom");
958 if (!InitGRP()) return kFALSE;
960 // Read the reconstruction parameters from OCDB
961 if (!InitRecoParams()) {
967 AliSysInfo::AddStamp("ReadRecoParam");
969 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
970 if(fDiamondProfile && fMeanVertexConstraint) ftVertexer->SetVtxStart(fDiamondProfile);
973 if (fRunVertexFinder && !CreateVertexer()) {
979 AliSysInfo::AddStamp("Vertexer");
982 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
988 AliSysInfo::AddStamp("LoadTrackers");
990 // create the ESD output file and tree
991 ffile = TFile::Open("AliESDs.root", "RECREATE");
992 ffile->SetCompressionLevel(2);
993 if (!ffile->IsOpen()) {
994 AliError("opening AliESDs.root failed");
995 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
998 ftree = new TTree("esdTree", "Tree with ESD objects");
999 fesd = new AliESDEvent();
1000 fesd->CreateStdContent();
1001 fesd->WriteToTree(ftree);
1003 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1004 fhltesd = new AliESDEvent();
1005 fhltesd->CreateStdContent();
1006 fhltesd->WriteToTree(fhlttree);
1009 if (fWriteESDfriend) {
1010 fesdf = new AliESDfriend();
1011 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1012 br->SetFile("AliESDfriends.root");
1013 fesd->AddObject(fesdf);
1018 if (fRawReader) fRawReader->RewindEvents();
1020 ProcInfo_t ProcInfo;
1021 gSystem->GetProcInfo(&ProcInfo);
1022 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1026 if (fRunQA && fRawReader && fQATasks.Contains(Form("%d", AliQA::kRAWS))) {
1027 fQASteer->Run(fQADetectors, fRawReader) ;
1028 fSameQACycle = kTRUE ;
1032 //Initialize the QA and start of cycle for out-of-loop QA
1034 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam, fSameQACycle, !fInLoopQA) ;
1038 fSameQACycle = kFALSE;
1039 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1040 AliInfo(Form("Initializing the global QA data maker"));
1041 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1042 TObjArray *arr=qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun());
1043 AliTracker::SetResidualsArray(arr);
1045 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
1046 fSameQACycle = kTRUE;
1049 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1050 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
1052 qadm->StartOfCycle(AliQA::kESDS, fSameQACycle);
1053 fSameQACycle = kTRUE;
1058 //Initialize the Plane Efficiency framework
1059 if (fRunPlaneEff && !InitPlaneEff()) {
1060 if(fStopOnError) {CleanUp(ffile); return kFALSE;}
1063 if (strcmp(gProgName,"alieve") == 0)
1064 fRunAliEVE = InitAliEVE();
1069 //_____________________________________________________________________________
1070 Bool_t AliReconstruction::RunEvent(Int_t iEvent)
1072 // run the reconstruction over a single event
1073 // The event loop is steered in Run method
1075 AliCodeTimerAuto("");
1077 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1078 fRunLoader->SetEventNumber(iEvent);
1079 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1081 //?? fRunLoader->MakeTree("H");
1082 fRunLoader->TreeE()->Fill();
1085 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1089 AliInfo(Form("processing event %d", iEvent));
1091 // Fill Event-info object
1093 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo);
1095 //Start of cycle for the in-loop QA
1096 if (fInLoopQA && fRunQA) {
1097 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam, fSameQACycle, fInLoopQA) ;
1099 if (fInLoopQA && fRunGlobalQA) {
1100 fSameQACycle = kFALSE;
1101 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1102 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1103 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
1104 fSameQACycle = kTRUE;
1106 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1107 qadm->StartOfCycle(AliQA::kESDS, fSameQACycle);
1108 fSameQACycle = kTRUE;
1112 fRunLoader->GetEvent(iEvent);
1114 char aFileName[256];
1115 sprintf(aFileName, "ESD_%d.%d_final.root",
1116 fRunLoader->GetHeader()->GetRun(),
1117 fRunLoader->GetHeader()->GetEventNrInRun());
1118 if (!gSystem->AccessPathName(aFileName)) return kTRUE;
1121 if (fInLoopQA && fRunQA)
1122 fQASteer->RunOneEvent(fRawReader) ;
1124 // local single event reconstruction
1125 if (!fRunLocalReconstruction.IsNull()) {
1126 TString detectors=fRunLocalReconstruction;
1127 // run HLT event reconstruction first
1128 // ;-( IsSelected changes the string
1129 if (IsSelected("HLT", detectors) &&
1130 !RunLocalEventReconstruction("HLT")) {
1131 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1133 detectors=fRunLocalReconstruction;
1134 detectors.ReplaceAll("HLT", "");
1135 if (!RunLocalEventReconstruction(detectors)) {
1136 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1140 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1141 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1142 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1143 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1145 // Set magnetic field from the tracker
1146 fesd->SetMagneticField(AliTracker::GetBz());
1147 fhltesd->SetMagneticField(AliTracker::GetBz());
1151 // Fill raw-data error log into the ESD
1152 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1155 if (fRunVertexFinder) {
1156 if (!ReadESD(fesd, "vertex")) {
1157 if (!RunVertexFinder(fesd)) {
1158 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1164 if (!fRunTracking.IsNull()) {
1165 if (fRunMuonTracking) {
1166 if (!RunMuonTracking(fesd)) {
1167 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1173 if (!fRunTracking.IsNull()) {
1174 if (!ReadESD(fesd, "tracking")) {
1175 if (!RunTracking(fesd)) {
1176 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1182 if (!fFillESD.IsNull()) {
1183 TString detectors=fFillESD;
1184 // run HLT first and on hltesd
1185 // ;-( IsSelected changes the string
1186 if (IsSelected("HLT", detectors) &&
1187 !FillESD(fhltesd, "HLT")) {
1188 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1191 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1192 if (detectors.Contains("ALL")) {
1194 for (Int_t idet=0; idet<fgkNDetectors; ++idet){
1195 detectors += fgkDetectorName[idet];
1199 detectors.ReplaceAll("HLT", "");
1200 if (!FillESD(fesd, detectors)) {
1201 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1205 // fill Event header information from the RawEventHeader
1206 if (fRawReader){FillRawEventHeaderESD(fesd);}
1209 AliESDpid::MakePID(fesd);
1211 if (fFillTriggerESD) {
1212 if (!ReadESD(fesd, "trigger")) {
1213 if (!FillTriggerESD(fesd)) {
1214 if (fStopOnError) {CleanUp(ffile); return kFALSE;}
1222 // Propagate track to the beam pipe (if not already done by ITS)
1224 const Int_t ntracks = fesd->GetNumberOfTracks();
1225 const Double_t kBz = fesd->GetMagneticField();
1226 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1229 UShort_t *selectedIdx=new UShort_t[ntracks];
1231 for (Int_t itrack=0; itrack<ntracks; itrack++){
1232 const Double_t kMaxStep = 5; //max step over the material
1235 AliESDtrack *track = fesd->GetTrack(itrack);
1236 if (!track) continue;
1238 AliExternalTrackParam *tpcTrack =
1239 (AliExternalTrackParam *)track->GetTPCInnerParam();
1243 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kTRUE);
1246 Int_t n=trkArray.GetEntriesFast();
1247 selectedIdx[n]=track->GetID();
1248 trkArray.AddLast(tpcTrack);
1251 //Tracks refitted by ITS should already be at the SPD vertex
1252 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1255 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1256 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1261 // Improve the reconstructed primary vertex position using the tracks
1263 TObject *obj = fOptions.FindObject("ITS");
1265 TString optITS = obj->GetTitle();
1266 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
1267 fRunVertexFinderTracks=kFALSE;
1269 if (fRunVertexFinderTracks) {
1270 // TPC + ITS primary vertex
1271 ftVertexer->SetITSrefitRequired();
1272 if(fDiamondProfile && fMeanVertexConstraint) {
1273 ftVertexer->SetVtxStart(fDiamondProfile);
1275 ftVertexer->SetConstraintOff();
1277 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1279 if (pvtx->GetStatus()) {
1280 fesd->SetPrimaryVertex(pvtx);
1281 for (Int_t i=0; i<ntracks; i++) {
1282 AliESDtrack *t = fesd->GetTrack(i);
1283 t->RelateToVertex(pvtx, kBz, kVeryBig);
1288 // TPC-only primary vertex
1289 ftVertexer->SetITSrefitNotRequired();
1290 if(fDiamondProfileTPC && fMeanVertexConstraint) {
1291 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1293 ftVertexer->SetConstraintOff();
1295 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1297 if (pvtx->GetStatus()) {
1298 fesd->SetPrimaryVertexTPC(pvtx);
1299 for (Int_t i=0; i<ntracks; i++) {
1300 AliESDtrack *t = fesd->GetTrack(i);
1301 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1307 delete[] selectedIdx;
1309 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1314 AliV0vertexer vtxer;
1315 vtxer.Tracks2V0vertices(fesd);
1317 if (fRunCascadeFinder) {
1319 AliCascadeVertexer cvtxer;
1320 cvtxer.V0sTracks2CascadeVertices(fesd);
1325 if (fCleanESD) CleanESD(fesd);
1329 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1330 if (qadm && fQATasks.Contains(Form("%d", AliQA::kESDS)))
1331 qadm->Exec(AliQA::kESDS, fesd);
1335 if (fWriteESDfriend) {
1336 fesdf->~AliESDfriend();
1337 new (fesdf) AliESDfriend(); // Reset...
1338 fesd->GetESDfriend(fesdf);
1342 // Auto-save the ESD tree in case of prompt reco @P2
1343 if (fRawReader && fRawReader->UseAutoSaveESD())
1344 ftree->AutoSave("SaveSelf");
1350 if (fRunAliEVE) RunAliEVE();
1354 if (fWriteESDfriend) {
1355 fesdf->~AliESDfriend();
1356 new (fesdf) AliESDfriend(); // Reset...
1359 ProcInfo_t ProcInfo;
1360 gSystem->GetProcInfo(&ProcInfo);
1361 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1364 // End of cycle for the in-loop
1365 if (fInLoopQA && fRunQA) {
1366 fQASteer->RunOneEvent(fesd) ;
1367 fQASteer->EndOfCycle() ;
1370 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1372 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1373 qadm->EndOfCycle(AliQA::kRECPOINTS);
1374 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1375 qadm->EndOfCycle(AliQA::kESDS);
1381 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1382 if (fReconstructor[iDet])
1383 fReconstructor[iDet]->SetRecoParam(NULL);
1389 //_____________________________________________________________________________
1390 Bool_t AliReconstruction::FinishRun()
1393 // Called after the exit
1394 // from the event loop
1395 AliCodeTimerAuto("");
1397 if (fIsNewRunLoader) { // galice.root didn't exist
1398 fRunLoader->WriteHeader("OVERWRITE");
1399 fRunLoader->CdGAFile();
1400 fRunLoader->Write(0, TObject::kOverwrite);
1403 ftree->GetUserInfo()->Add(fesd);
1404 fhlttree->GetUserInfo()->Add(fhltesd);
1406 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1407 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1409 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1410 cdbMapCopy->SetOwner(1);
1411 cdbMapCopy->SetName("cdbMap");
1412 TIter iter(cdbMap->GetTable());
1415 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1416 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1417 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1418 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1421 TList *cdbListCopy = new TList();
1422 cdbListCopy->SetOwner(1);
1423 cdbListCopy->SetName("cdbList");
1425 TIter iter2(cdbList);
1428 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1429 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1432 ftree->GetUserInfo()->Add(cdbMapCopy);
1433 ftree->GetUserInfo()->Add(cdbListCopy);
1436 if(fESDPar.Contains("ESD.par")){
1437 AliInfo("Attaching ESD.par to Tree");
1438 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
1439 ftree->GetUserInfo()->Add(fn);
1445 if (fWriteESDfriend)
1446 ftree->SetBranchStatus("ESDfriend*",0);
1447 // we want to have only one tree version number
1448 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1451 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1452 if (fRunPlaneEff && !FinishPlaneEff()) {
1453 AliWarning("Finish PlaneEff evaluation failed");
1460 AliWarning("AOD creation not supported anymore during reconstruction. See ANALYSIS/AliAnalysisTaskESDfilter.cxx instead.");
1463 // Create tags for the events in the ESD tree (the ESD tree is always present)
1464 // In case of empty events the tags will contain dummy values
1465 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1466 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData);
1468 AliWarning("AOD tag creation not supported anymore during reconstruction.");
1471 //Finish QA and end of cycle for out-of-loop QA
1472 if (!fInLoopQA && fRunQA)
1473 fQASteer->Run(fRunLocalReconstruction.Data(), AliQA::kNULLTASKINDEX, fSameQACycle) ;
1474 if (!fInLoopQA && fRunGlobalQA) {
1475 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1477 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1478 qadm->EndOfCycle(AliQA::kRECPOINTS);
1479 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1480 qadm->EndOfCycle(AliQA::kESDS);
1485 // Cleanup of CDB manager: cache and active storages!
1486 AliCDBManager::Instance()->ClearCache();
1492 //_____________________________________________________________________________
1493 Bool_t AliReconstruction::RunLocalReconstruction(const TString& /*detectors*/)
1495 // run the local reconstruction
1496 static Int_t eventNr=0;
1497 AliCodeTimerAuto("")
1499 // AliCDBManager* man = AliCDBManager::Instance();
1500 // Bool_t origCache = man->GetCacheFlag();
1502 // TString detStr = detectors;
1503 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1504 // if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1505 // AliReconstructor* reconstructor = GetReconstructor(iDet);
1506 // if (!reconstructor) continue;
1507 // if (reconstructor->HasLocalReconstruction()) continue;
1509 // AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1510 // AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1512 // AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1513 // AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1515 // man->SetCacheFlag(kTRUE);
1516 // TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
1517 // man->GetAll(calibPath); // entries are cached!
1519 // AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1521 // if (fRawReader) {
1522 // fRawReader->RewindEvents();
1523 // reconstructor->Reconstruct(fRunLoader, fRawReader);
1525 // reconstructor->Reconstruct(fRunLoader);
1528 // AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1529 // AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
1531 // // unload calibration data
1532 // man->UnloadFromCache(calibPath);
1533 // //man->ClearCache();
1536 // man->SetCacheFlag(origCache);
1538 // if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1539 // AliError(Form("the following detectors were not found: %s",
1541 // if (fStopOnError) return kFALSE;
1548 //_____________________________________________________________________________
1549 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1551 // run the local reconstruction
1553 static Int_t eventNr=0;
1554 AliCodeTimerAuto("")
1556 TString detStr = detectors;
1557 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1558 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1559 AliReconstructor* reconstructor = GetReconstructor(iDet);
1560 if (!reconstructor) continue;
1561 AliLoader* loader = fLoader[iDet];
1562 // Matthias April 2008: temporary fix to run HLT reconstruction
1563 // although the HLT loader is missing
1564 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1566 reconstructor->Reconstruct(fRawReader, NULL);
1569 reconstructor->Reconstruct(dummy, NULL);
1574 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1577 // conversion of digits
1578 if (fRawReader && reconstructor->HasDigitConversion()) {
1579 AliInfo(Form("converting raw data digits into root objects for %s",
1580 fgkDetectorName[iDet]));
1581 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1582 // fgkDetectorName[iDet]));
1583 loader->LoadDigits("update");
1584 loader->CleanDigits();
1585 loader->MakeDigitsContainer();
1586 TTree* digitsTree = loader->TreeD();
1587 reconstructor->ConvertDigits(fRawReader, digitsTree);
1588 loader->WriteDigits("OVERWRITE");
1589 loader->UnloadDigits();
1591 // local reconstruction
1592 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1593 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1594 loader->LoadRecPoints("update");
1595 loader->CleanRecPoints();
1596 loader->MakeRecPointsContainer();
1597 TTree* clustersTree = loader->TreeR();
1598 if (fRawReader && !reconstructor->HasDigitConversion()) {
1599 reconstructor->Reconstruct(fRawReader, clustersTree);
1601 loader->LoadDigits("read");
1602 TTree* digitsTree = loader->TreeD();
1604 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1605 if (fStopOnError) return kFALSE;
1607 reconstructor->Reconstruct(digitsTree, clustersTree);
1609 loader->UnloadDigits();
1612 // In-loop QA for local reconstrucion
1613 TString detQAStr(fQADetectors) ;
1614 if (fRunQA && fInLoopQA)
1615 fQASteer->RunOneEventInOneDetector(iDet, clustersTree) ;
1617 loader->WriteRecPoints("OVERWRITE");
1618 loader->UnloadRecPoints();
1619 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1621 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1622 AliError(Form("the following detectors were not found: %s",
1624 if (fStopOnError) return kFALSE;
1630 //_____________________________________________________________________________
1631 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1633 // run the barrel tracking
1635 AliCodeTimerAuto("")
1637 AliESDVertex* vertex = NULL;
1638 Double_t vtxPos[3] = {0, 0, 0};
1639 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1640 TArrayF mcVertex(3);
1641 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1642 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1643 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1647 AliInfo("running the ITS vertex finder");
1649 fLoader[0]->LoadRecPoints();
1650 TTree* cltree = fLoader[0]->TreeR();
1652 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1653 vertex = fVertexer->FindVertexForCurrentEvent(cltree);
1656 AliError("Can't get the ITS cluster tree");
1658 fLoader[0]->UnloadRecPoints();
1661 AliError("Can't get the ITS loader");
1664 AliWarning("Vertex not found");
1665 vertex = new AliESDVertex();
1666 vertex->SetName("default");
1669 vertex->SetName("reconstructed");
1673 AliInfo("getting the primary vertex from MC");
1674 vertex = new AliESDVertex(vtxPos, vtxErr);
1678 vertex->GetXYZ(vtxPos);
1679 vertex->GetSigmaXYZ(vtxErr);
1681 AliWarning("no vertex reconstructed");
1682 vertex = new AliESDVertex(vtxPos, vtxErr);
1684 esd->SetPrimaryVertexSPD(vertex);
1685 // if SPD multiplicity has been determined, it is stored in the ESD
1686 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1687 if(mult)esd->SetMultiplicity(mult);
1689 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1690 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1697 //_____________________________________________________________________________
1698 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1700 // run the HLT barrel tracking
1702 AliCodeTimerAuto("")
1705 AliError("Missing runLoader!");
1709 AliInfo("running HLT tracking");
1711 // Get a pointer to the HLT reconstructor
1712 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1713 if (!reconstructor) return kFALSE;
1716 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1717 TString detName = fgkDetectorName[iDet];
1718 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1719 reconstructor->SetOption(detName.Data());
1720 AliTracker *tracker = reconstructor->CreateTracker();
1722 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1723 if (fStopOnError) return kFALSE;
1727 Double_t vtxErr[3]={0.005,0.005,0.010};
1728 const AliESDVertex *vertex = esd->GetVertex();
1729 vertex->GetXYZ(vtxPos);
1730 tracker->SetVertex(vtxPos,vtxErr);
1732 fLoader[iDet]->LoadRecPoints("read");
1733 TTree* tree = fLoader[iDet]->TreeR();
1735 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1738 tracker->LoadClusters(tree);
1740 if (tracker->Clusters2Tracks(esd) != 0) {
1741 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1745 tracker->UnloadClusters();
1753 //_____________________________________________________________________________
1754 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1756 // run the muon spectrometer tracking
1758 AliCodeTimerAuto("")
1761 AliError("Missing runLoader!");
1764 Int_t iDet = 7; // for MUON
1766 AliInfo("is running...");
1768 // Get a pointer to the MUON reconstructor
1769 AliReconstructor *reconstructor = GetReconstructor(iDet);
1770 if (!reconstructor) return kFALSE;
1773 TString detName = fgkDetectorName[iDet];
1774 AliDebug(1, Form("%s tracking", detName.Data()));
1775 AliTracker *tracker = reconstructor->CreateTracker();
1777 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1782 fLoader[iDet]->LoadRecPoints("read");
1784 tracker->LoadClusters(fLoader[iDet]->TreeR());
1786 Int_t rv = tracker->Clusters2Tracks(esd);
1790 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1794 fLoader[iDet]->UnloadRecPoints();
1796 tracker->UnloadClusters();
1804 //_____________________________________________________________________________
1805 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1807 // run the barrel tracking
1808 static Int_t eventNr=0;
1809 AliCodeTimerAuto("")
1811 AliInfo("running tracking");
1813 //Fill the ESD with the T0 info (will be used by the TOF)
1814 if (fReconstructor[11] && fLoader[11]) {
1815 fLoader[11]->LoadRecPoints("READ");
1816 TTree *treeR = fLoader[11]->TreeR();
1817 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
1820 // pass 1: TPC + ITS inwards
1821 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1822 if (!fTracker[iDet]) continue;
1823 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1826 fLoader[iDet]->LoadRecPoints("read");
1827 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
1828 TTree* tree = fLoader[iDet]->TreeR();
1830 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1833 fTracker[iDet]->LoadClusters(tree);
1834 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1836 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1837 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1840 // preliminary PID in TPC needed by the ITS tracker
1842 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1843 AliESDpid::MakePID(esd);
1845 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
1848 // pass 2: ALL backwards
1850 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1851 if (!fTracker[iDet]) continue;
1852 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1855 if (iDet > 1) { // all except ITS, TPC
1857 fLoader[iDet]->LoadRecPoints("read");
1858 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
1859 tree = fLoader[iDet]->TreeR();
1861 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1864 fTracker[iDet]->LoadClusters(tree);
1865 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1869 if (iDet>1) // start filling residuals for the "outer" detectors
1870 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
1872 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1873 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1878 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
1879 fTracker[iDet]->UnloadClusters();
1880 fLoader[iDet]->UnloadRecPoints();
1882 // updated PID in TPC needed by the ITS tracker -MI
1884 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1885 AliESDpid::MakePID(esd);
1887 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1889 //stop filling residuals for the "outer" detectors
1890 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
1892 // pass 3: TRD + TPC + ITS refit inwards
1894 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1895 if (!fTracker[iDet]) continue;
1896 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1899 if (iDet<2) // start filling residuals for TPC and ITS
1900 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
1902 if (fTracker[iDet]->RefitInward(esd) != 0) {
1903 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1906 // run postprocessing
1907 if (fTracker[iDet]->PostProcess(esd) != 0) {
1908 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
1911 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1914 // write space-points to the ESD in case alignment data output
1916 if (fWriteAlignmentData)
1917 WriteAlignmentData(esd);
1919 for (Int_t iDet = 3; iDet >= 0; iDet--) {
1920 if (!fTracker[iDet]) continue;
1922 fTracker[iDet]->UnloadClusters();
1923 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
1924 fLoader[iDet]->UnloadRecPoints();
1925 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
1927 // stop filling residuals for TPC and ITS
1928 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
1934 //_____________________________________________________________________________
1935 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1937 // Remove the data which are not needed for the physics analysis.
1940 Int_t nTracks=esd->GetNumberOfTracks();
1941 Int_t nV0s=esd->GetNumberOfV0s();
1943 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
1945 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
1946 Bool_t rc=esd->Clean(cleanPars);
1948 nTracks=esd->GetNumberOfTracks();
1949 nV0s=esd->GetNumberOfV0s();
1951 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
1956 //_____________________________________________________________________________
1957 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1959 // fill the event summary data
1961 AliCodeTimerAuto("")
1962 static Int_t eventNr=0;
1963 TString detStr = detectors;
1965 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1966 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1967 AliReconstructor* reconstructor = GetReconstructor(iDet);
1968 if (!reconstructor) continue;
1969 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1970 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1971 TTree* clustersTree = NULL;
1972 if (fLoader[iDet]) {
1973 fLoader[iDet]->LoadRecPoints("read");
1974 clustersTree = fLoader[iDet]->TreeR();
1975 if (!clustersTree) {
1976 AliError(Form("Can't get the %s clusters tree",
1977 fgkDetectorName[iDet]));
1978 if (fStopOnError) return kFALSE;
1981 if (fRawReader && !reconstructor->HasDigitConversion()) {
1982 reconstructor->FillESD(fRawReader, clustersTree, esd);
1984 TTree* digitsTree = NULL;
1985 if (fLoader[iDet]) {
1986 fLoader[iDet]->LoadDigits("read");
1987 digitsTree = fLoader[iDet]->TreeD();
1989 AliError(Form("Can't get the %s digits tree",
1990 fgkDetectorName[iDet]));
1991 if (fStopOnError) return kFALSE;
1994 reconstructor->FillESD(digitsTree, clustersTree, esd);
1995 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1997 if (fLoader[iDet]) {
1998 fLoader[iDet]->UnloadRecPoints();
2004 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2005 AliError(Form("the following detectors were not found: %s",
2007 if (fStopOnError) return kFALSE;
2009 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2014 //_____________________________________________________________________________
2015 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2017 // Reads the trigger decision which is
2018 // stored in Trigger.root file and fills
2019 // the corresponding esd entries
2021 AliCodeTimerAuto("")
2023 AliInfo("Filling trigger information into the ESD");
2026 AliCTPRawStream input(fRawReader);
2027 if (!input.Next()) {
2028 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2031 if (esd->GetTriggerMask() != input.GetClassMask())
2032 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2033 input.GetClassMask(),esd->GetTriggerMask()));
2034 if (esd->GetOrbitNumber() != input.GetOrbitID())
2035 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2036 input.GetOrbitID(),esd->GetOrbitNumber()));
2037 if (esd->GetBunchCrossNumber() != input.GetBCID())
2038 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2039 input.GetBCID(),esd->GetBunchCrossNumber()));
2042 // Here one has to add the filling of trigger inputs and
2043 // interaction records
2053 //_____________________________________________________________________________
2054 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2057 // Filling information from RawReader Header
2060 if (!fRawReader) return kFALSE;
2062 AliInfo("Filling information from RawReader Header");
2064 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2065 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2066 esd->SetPeriodNumber(fRawReader->GetPeriod());
2068 esd->SetTimeStamp(fRawReader->GetTimestamp());
2069 esd->SetEventType(fRawReader->GetType());
2075 //_____________________________________________________________________________
2076 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2078 // check whether detName is contained in detectors
2079 // if yes, it is removed from detectors
2081 // check if all detectors are selected
2082 if ((detectors.CompareTo("ALL") == 0) ||
2083 detectors.BeginsWith("ALL ") ||
2084 detectors.EndsWith(" ALL") ||
2085 detectors.Contains(" ALL ")) {
2090 // search for the given detector
2091 Bool_t result = kFALSE;
2092 if ((detectors.CompareTo(detName) == 0) ||
2093 detectors.BeginsWith(detName+" ") ||
2094 detectors.EndsWith(" "+detName) ||
2095 detectors.Contains(" "+detName+" ")) {
2096 detectors.ReplaceAll(detName, "");
2100 // clean up the detectors string
2101 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2102 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2103 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2108 //_____________________________________________________________________________
2109 Bool_t AliReconstruction::InitRunLoader()
2111 // get or create the run loader
2113 if (gAlice) delete gAlice;
2116 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2117 // load all base libraries to get the loader classes
2118 TString libs = gSystem->GetLibraries();
2119 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2120 TString detName = fgkDetectorName[iDet];
2121 if (detName == "HLT") continue;
2122 if (libs.Contains("lib" + detName + "base.so")) continue;
2123 gSystem->Load("lib" + detName + "base.so");
2125 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2127 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2132 fRunLoader->CdGAFile();
2133 fRunLoader->LoadgAlice();
2135 //PH This is a temporary fix to give access to the kinematics
2136 //PH that is needed for the labels of ITS clusters
2137 fRunLoader->LoadHeader();
2138 fRunLoader->LoadKinematics();
2140 } else { // galice.root does not exist
2142 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2146 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2147 AliConfig::GetDefaultEventFolderName(),
2150 AliError(Form("could not create run loader in file %s",
2151 fGAliceFileName.Data()));
2155 fIsNewRunLoader = kTRUE;
2156 fRunLoader->MakeTree("E");
2158 if (fNumberOfEventsPerFile > 0)
2159 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2161 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2167 //_____________________________________________________________________________
2168 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2170 // get the reconstructor object and the loader for a detector
2172 if (fReconstructor[iDet]) {
2173 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2174 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2175 fReconstructor[iDet]->SetRecoParam(par);
2177 return fReconstructor[iDet];
2180 // load the reconstructor object
2181 TPluginManager* pluginManager = gROOT->GetPluginManager();
2182 TString detName = fgkDetectorName[iDet];
2183 TString recName = "Ali" + detName + "Reconstructor";
2185 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2187 AliReconstructor* reconstructor = NULL;
2188 // first check if a plugin is defined for the reconstructor
2189 TPluginHandler* pluginHandler =
2190 pluginManager->FindHandler("AliReconstructor", detName);
2191 // if not, add a plugin for it
2192 if (!pluginHandler) {
2193 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2194 TString libs = gSystem->GetLibraries();
2195 if (libs.Contains("lib" + detName + "base.so") ||
2196 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2197 pluginManager->AddHandler("AliReconstructor", detName,
2198 recName, detName + "rec", recName + "()");
2200 pluginManager->AddHandler("AliReconstructor", detName,
2201 recName, detName, recName + "()");
2203 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2205 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2206 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2208 if (reconstructor) {
2209 TObject* obj = fOptions.FindObject(detName.Data());
2210 if (obj) reconstructor->SetOption(obj->GetTitle());
2211 reconstructor->Init();
2212 fReconstructor[iDet] = reconstructor;
2215 // get or create the loader
2216 if (detName != "HLT") {
2217 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2218 if (!fLoader[iDet]) {
2219 AliConfig::Instance()
2220 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2222 // first check if a plugin is defined for the loader
2224 pluginManager->FindHandler("AliLoader", detName);
2225 // if not, add a plugin for it
2226 if (!pluginHandler) {
2227 TString loaderName = "Ali" + detName + "Loader";
2228 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2229 pluginManager->AddHandler("AliLoader", detName,
2230 loaderName, detName + "base",
2231 loaderName + "(const char*, TFolder*)");
2232 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2234 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2236 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2237 fRunLoader->GetEventFolder());
2239 if (!fLoader[iDet]) { // use default loader
2240 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2242 if (!fLoader[iDet]) {
2243 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2244 if (fStopOnError) return NULL;
2246 fRunLoader->AddLoader(fLoader[iDet]);
2247 fRunLoader->CdGAFile();
2248 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2249 fRunLoader->Write(0, TObject::kOverwrite);
2254 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2255 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2256 reconstructor->SetRecoParam(par);
2258 return reconstructor;
2261 //_____________________________________________________________________________
2262 Bool_t AliReconstruction::CreateVertexer()
2264 // create the vertexer
2267 AliReconstructor* itsReconstructor = GetReconstructor(0);
2268 if (itsReconstructor) {
2269 fVertexer = itsReconstructor->CreateVertexer();
2272 AliWarning("couldn't create a vertexer for ITS");
2273 if (fStopOnError) return kFALSE;
2279 //_____________________________________________________________________________
2280 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2282 // create the trackers
2284 TString detStr = detectors;
2285 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2286 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2287 AliReconstructor* reconstructor = GetReconstructor(iDet);
2288 if (!reconstructor) continue;
2289 TString detName = fgkDetectorName[iDet];
2290 if (detName == "HLT") {
2291 fRunHLTTracking = kTRUE;
2294 if (detName == "MUON") {
2295 fRunMuonTracking = kTRUE;
2300 fTracker[iDet] = reconstructor->CreateTracker();
2301 if (!fTracker[iDet] && (iDet < 7)) {
2302 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2303 if (fStopOnError) return kFALSE;
2305 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2311 //_____________________________________________________________________________
2312 void AliReconstruction::CleanUp(TFile* file)
2314 // delete trackers and the run loader and close and delete the file
2316 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2317 delete fReconstructor[iDet];
2318 fReconstructor[iDet] = NULL;
2319 fLoader[iDet] = NULL;
2320 delete fTracker[iDet];
2321 fTracker[iDet] = NULL;
2323 if (fRunInfo) delete fRunInfo;
2329 if (ftVertexer) delete ftVertexer;
2332 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
2333 delete fDiamondProfile;
2334 fDiamondProfile = NULL;
2335 delete fDiamondProfileTPC;
2336 fDiamondProfileTPC = NULL;
2346 if (fParentRawReader) delete fParentRawReader;
2347 fParentRawReader=NULL;
2356 //_____________________________________________________________________________
2358 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
2360 // read the ESD event from a file
2362 if (!esd) return kFALSE;
2364 sprintf(fileName, "ESD_%d.%d_%s.root",
2365 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
2366 if (gSystem->AccessPathName(fileName)) return kFALSE;
2368 AliInfo(Form("reading ESD from file %s", fileName));
2369 AliDebug(1, Form("reading ESD from file %s", fileName));
2370 TFile* file = TFile::Open(fileName);
2371 if (!file || !file->IsOpen()) {
2372 AliError(Form("opening %s failed", fileName));
2379 esd = (AliESDEvent*) file->Get("ESD");
2388 //_____________________________________________________________________________
2389 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
2391 // write the ESD event to a file
2395 sprintf(fileName, "ESD_%d.%d_%s.root",
2396 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
2398 AliDebug(1, Form("writing ESD to file %s", fileName));
2399 TFile* file = TFile::Open(fileName, "recreate");
2400 if (!file || !file->IsOpen()) {
2401 AliError(Form("opening %s failed", fileName));
2410 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2412 // Write space-points which are then used in the alignment procedures
2413 // For the moment only ITS, TPC, TRD and TOF
2415 Int_t ntracks = esd->GetNumberOfTracks();
2416 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2418 AliESDtrack *track = esd->GetTrack(itrack);
2421 for (Int_t iDet = 3; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2422 nsp += track->GetNcls(iDet);
2424 if (iDet==0) { // ITS "extra" clusters
2425 track->GetClusters(iDet,idx);
2426 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2431 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2432 track->SetTrackPointArray(sp);
2434 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2435 AliTracker *tracker = fTracker[iDet];
2436 if (!tracker) continue;
2437 Int_t nspdet = track->GetClusters(iDet,idx);
2439 if (iDet==0) // ITS "extra" clusters
2440 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2442 if (nspdet <= 0) continue;
2446 while (isp2 < nspdet) {
2447 Bool_t isvalid=kTRUE;
2449 Int_t index=idx[isp++];
2450 if (index < 0) continue;
2452 TString dets = fgkDetectorName[iDet];
2453 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2454 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2455 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2456 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2457 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2459 isvalid = tracker->GetTrackPoint(index,p);
2462 if (!isvalid) continue;
2463 sp->AddPoint(isptrack,&p); isptrack++;
2470 //_____________________________________________________________________________
2471 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2473 // The method reads the raw-data error log
2474 // accumulated within the rawReader.
2475 // It extracts the raw-data errors related to
2476 // the current event and stores them into
2477 // a TClonesArray inside the esd object.
2479 if (!fRawReader) return;
2481 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2483 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2485 if (iEvent != log->GetEventNumber()) continue;
2487 esd->AddRawDataErrorLog(log);
2492 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString pName){
2493 // Dump a file content into a char in TNamed
2495 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2496 Int_t kBytes = (Int_t)in.tellg();
2497 printf("Size: %d \n",kBytes);
2500 char* memblock = new char [kBytes];
2501 in.seekg (0, ios::beg);
2502 in.read (memblock, kBytes);
2504 TString fData(memblock,kBytes);
2505 fn = new TNamed(pName,fData);
2506 printf("fData Size: %d \n",fData.Sizeof());
2507 printf("pName Size: %d \n",pName.Sizeof());
2508 printf("fn Size: %d \n",fn->Sizeof());
2512 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2518 void AliReconstruction::TNamedToFile(TTree* fTree, TString pName){
2519 // This is not really needed in AliReconstruction at the moment
2520 // but can serve as a template
2522 TList *fList = fTree->GetUserInfo();
2523 TNamed *fn = (TNamed*)fList->FindObject(pName.Data());
2524 printf("fn Size: %d \n",fn->Sizeof());
2526 TString fTmp(fn->GetName()); // to be 100% sure in principle pName also works
2527 const char* cdata = fn->GetTitle();
2528 printf("fTmp Size %d\n",fTmp.Sizeof());
2530 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2531 printf("calculated size %d\n",size);
2532 ofstream out(pName.Data(),ios::out | ios::binary);
2533 out.write(cdata,size);
2538 //_____________________________________________________________________________
2539 void AliReconstruction::CheckQA()
2541 // check the QA of SIM for this run and remove the detectors
2542 // with status Fatal
2544 TString newRunLocalReconstruction ;
2545 TString newRunTracking ;
2546 TString newFillESD ;
2548 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2549 TString detName(AliQA::GetDetName(iDet)) ;
2550 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX_t(iDet)) ;
2551 if ( qa->IsSet(AliQA::DETECTORINDEX_t(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2552 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2554 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2555 fRunLocalReconstruction.Contains("ALL") ) {
2556 newRunLocalReconstruction += detName ;
2557 newRunLocalReconstruction += " " ;
2559 if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) ||
2560 fRunTracking.Contains("ALL") ) {
2561 newRunTracking += detName ;
2562 newRunTracking += " " ;
2564 if ( fFillESD.Contains(AliQA::GetDetName(iDet)) ||
2565 fFillESD.Contains("ALL") ) {
2566 newFillESD += detName ;
2571 fRunLocalReconstruction = newRunLocalReconstruction ;
2572 fRunTracking = newRunTracking ;
2573 fFillESD = newFillESD ;
2576 //_____________________________________________________________________________
2577 Int_t AliReconstruction::GetDetIndex(const char* detector)
2579 // return the detector index corresponding to detector
2581 for (index = 0; index < fgkNDetectors ; index++) {
2582 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2587 //_____________________________________________________________________________
2588 Bool_t AliReconstruction::FinishPlaneEff() {
2590 // Here execute all the necessary operationis, at the end of the tracking phase,
2591 // in case that evaluation of PlaneEfficiencies was required for some detector.
2592 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2594 // This Preliminary version works only FOR ITS !!!!!
2595 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2598 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2601 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2602 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2603 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2604 if(fTracker[iDet]) {
2605 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
2606 TString name=planeeff->GetName();
2608 TFile* pefile = TFile::Open(name, "RECREATE");
2609 ret=(Bool_t)planeeff->Write();
2611 if(planeeff->GetCreateHistos()) {
2612 TString hname=planeeff->GetName();
2613 hname+="Histo.root";
2614 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
2620 //_____________________________________________________________________________
2621 Bool_t AliReconstruction::InitPlaneEff() {
2623 // Here execute all the necessary operations, before of the tracking phase,
2624 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
2625 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
2626 // which should be updated/recalculated.
2628 // This Preliminary version will work only FOR ITS !!!!!
2629 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2632 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2634 AliWarning(Form("Implementation of this method not yet done !! Method return kTRUE"));
2638 //_____________________________________________________________________________
2639 Bool_t AliReconstruction::InitAliEVE()
2641 // This method should be called only in case
2642 // AliReconstruction is run
2643 // within the alieve environment.
2644 // It will initialize AliEVE in a way
2645 // so that it can visualize event processed
2646 // by AliReconstruction.
2647 // The return flag shows whenever the
2648 // AliEVE initialization was successful or not.
2651 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
2652 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
2653 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
2655 gROOT->ProcessLine("if (!gAliEveEvent) {gAliEveEvent = new AliEveEventManager();gAliEveEvent->SetAutoLoad(kTRUE);gAliEveEvent->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(gAliEveEvent);};");
2656 gROOT->ProcessLine("alieve_online_init()");
2661 //_____________________________________________________________________________
2662 void AliReconstruction::RunAliEVE()
2664 // Runs AliEVE visualisation of
2665 // the current event.
2666 // Should be executed only after
2667 // successful initialization of AliEVE.
2669 AliInfo("Running AliEVE...");
2670 gROOT->ProcessLine(Form("gAliEveEvent->SetEvent((AliRunLoader*)%p,(AliRawReader*)%p,(AliESDEvent*)%p);",fRunLoader,fRawReader,fesd));
2671 gROOT->ProcessLine("gAliEveEvent->StartStopAutoLoadTimer();");
2675 //_____________________________________________________________________________
2676 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
2678 // Allows to run QA for a selected set of detectors
2679 // and a selected set of tasks among RAWS, RECPOINTS and ESDS
2680 // all selected detectors run the same selected tasks
2682 if (!detAndAction.Contains(":")) {
2683 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2687 Int_t colon = detAndAction.Index(":") ;
2688 fQADetectors = detAndAction(0, colon) ;
2689 if (fQADetectors.Contains("ALL") )
2690 fQADetectors = fFillESD ;
2691 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2692 if (fQATasks.Contains("ALL") ) {
2693 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
2695 fQATasks.ToUpper() ;
2697 if ( fQATasks.Contains("RAW") )
2698 tempo = Form("%d ", AliQA::kRAWS) ;
2699 if ( fQATasks.Contains("RECPOINT") )
2700 tempo += Form("%d ", AliQA::kRECPOINTS) ;
2701 if ( fQATasks.Contains("ESD") )
2702 tempo += Form("%d ", AliQA::kESDS) ;
2704 if (fQATasks.IsNull()) {
2705 AliInfo("No QA requested\n") ;
2710 TString tempo(fQATasks) ;
2711 tempo.ReplaceAll(Form("%d", AliQA::kRAWS), AliQA::GetTaskName(AliQA::kRAWS)) ;
2712 tempo.ReplaceAll(Form("%d", AliQA::kRECPOINTS), AliQA::GetTaskName(AliQA::kRECPOINTS)) ;
2713 tempo.ReplaceAll(Form("%d", AliQA::kESDS), AliQA::GetTaskName(AliQA::kESDS)) ;
2714 fQASteer->SetActiveDetectors(fQADetectors) ;
2715 fQASteer->SetTasks(fQATasks) ;
2716 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2721 //_____________________________________________________________________________
2722 Bool_t AliReconstruction::InitRecoParams()
2724 // The method accesses OCDB and retrieves all
2725 // the available reco-param objects from there.
2727 Bool_t isOK = kTRUE;
2729 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2731 if (fRecoParam.GetDetRecoParamArray(iDet)) {
2732 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
2736 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
2738 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
2739 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
2741 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
2745 TObject *recoParamObj = entry->GetObject();
2746 if (dynamic_cast<TObjArray*>(recoParamObj)) {
2747 // The detector has a normal TobjArray of AliDetectorRecoParam objects
2748 // Registering them in AliRecoParam
2749 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
2751 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
2752 // The detector has only onse set of reco parameters
2753 // Registering it in AliRecoParam
2754 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
2755 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
2756 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
2759 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
2771 //_____________________________________________________________________________
2772 Bool_t AliReconstruction::GetEventInfo()
2774 // Fill the event info object
2776 AliCodeTimerAuto("")
2778 AliCentralTrigger *aCTP = NULL;
2780 fEventInfo.SetEventType(fRawReader->GetType());
2782 ULong64_t mask = fRawReader->GetClassMask();
2783 fEventInfo.SetTriggerMask(mask);
2784 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
2785 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
2787 aCTP = new AliCentralTrigger();
2788 TString configstr("");
2789 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
2790 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
2794 aCTP->SetClassMask(mask);
2795 aCTP->SetClusterMask(clmask);
2798 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
2800 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
2801 aCTP = fRunLoader->GetTrigger();
2802 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
2803 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
2806 AliWarning("No trigger can be loaded! The trigger information will not be used!");
2811 AliTriggerConfiguration *config = aCTP->GetConfiguration();
2813 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
2814 if (fRawReader) delete aCTP;
2818 UChar_t clustmask = 0;
2820 ULong64_t trmask = fEventInfo.GetTriggerMask();
2821 const TObjArray& classesArray = config->GetClasses();
2822 Int_t nclasses = classesArray.GetEntriesFast();
2823 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
2824 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
2826 Int_t trindex = (Int_t)TMath::Log2(trclass->GetMask());
2827 fesd->SetTriggerClass(trclass->GetName(),trindex);
2828 if (trmask & (1 << trindex)) {
2830 trclasses += trclass->GetName();
2832 clustmask |= trclass->GetCluster()->GetClusterMask();
2836 fEventInfo.SetTriggerClasses(trclasses);
2838 // Set the information in ESD
2839 fesd->SetTriggerMask(trmask);
2840 fesd->SetTriggerCluster(clustmask);
2842 if (!aCTP->CheckTriggeredDetectors()) {
2843 if (fRawReader) delete aCTP;
2847 if (fRawReader) delete aCTP;
2849 // We have to fill also the HLT decision here!!
2855 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
2857 // Match the detector list found in the rec.C or the default 'ALL'
2858 // to the list found in the GRP (stored there by the shuttle PP which
2859 // gets the information from ECS)
2860 static TString resultList;
2861 TString detList = detectorList;
2865 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
2866 if ((detectorMask >> iDet) & 0x1) {
2867 TString det = AliDAQ::OfflineModuleName(iDet);
2868 if ((detList.CompareTo("ALL") == 0) ||
2869 detList.BeginsWith("ALL ") ||
2870 detList.EndsWith(" ALL") ||
2871 detList.Contains(" ALL ") ||
2872 (detList.CompareTo(det) == 0) ||
2873 detList.BeginsWith(det) ||
2874 detList.EndsWith(det) ||
2875 detList.Contains( " "+det+" " )) {
2876 if (!resultList.EndsWith(det + " ")) {
2885 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
2886 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
2887 if ((detList.CompareTo("ALL") == 0) ||
2888 detList.BeginsWith("ALL ") ||
2889 detList.EndsWith(" ALL") ||
2890 detList.Contains(" ALL ") ||
2891 (detList.CompareTo(hltDet) == 0) ||
2892 detList.BeginsWith(hltDet) ||
2893 detList.EndsWith(hltDet) ||
2894 detList.Contains( " "+hltDet+" " )) {
2895 resultList += hltDet;
2899 return resultList.Data();