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>
125 #include <TParameter.h>
127 #include "AliReconstruction.h"
128 #include "AliCodeTimer.h"
129 #include "AliReconstructor.h"
131 #include "AliRunLoader.h"
133 #include "AliRawReaderFile.h"
134 #include "AliRawReaderDate.h"
135 #include "AliRawReaderRoot.h"
136 #include "AliRawEventHeaderBase.h"
137 #include "AliRawEvent.h"
138 #include "AliESDEvent.h"
139 #include "AliESDMuonTrack.h"
140 #include "AliESDfriend.h"
141 #include "AliESDVertex.h"
142 #include "AliESDcascade.h"
143 #include "AliESDkink.h"
144 #include "AliESDtrack.h"
145 #include "AliESDCaloCluster.h"
146 #include "AliESDCaloCells.h"
147 #include "AliMultiplicity.h"
148 #include "AliTracker.h"
149 #include "AliVertexer.h"
150 #include "AliVertexerTracks.h"
151 #include "AliV0vertexer.h"
152 #include "AliCascadeVertexer.h"
153 #include "AliHeader.h"
154 #include "AliGenEventHeader.h"
156 #include "AliESDpid.h"
157 #include "AliESDtrack.h"
158 #include "AliESDPmdTrack.h"
160 #include "AliESDTagCreator.h"
162 #include "AliGeomManager.h"
163 #include "AliTrackPointArray.h"
164 #include "AliCDBManager.h"
165 #include "AliCDBStorage.h"
166 #include "AliCDBEntry.h"
167 #include "AliAlignObj.h"
169 #include "AliCentralTrigger.h"
170 #include "AliTriggerConfiguration.h"
171 #include "AliTriggerClass.h"
172 #include "AliTriggerCluster.h"
173 #include "AliCTPRawStream.h"
175 #include "AliQADataMakerRec.h"
176 #include "AliGlobalQADataMaker.h"
178 #include "AliQADataMakerSteer.h"
180 #include "AliPlaneEff.h"
182 #include "AliSysInfo.h" // memory snapshots
183 #include "AliRawHLTManager.h"
185 #include "AliMagWrapCheb.h"
187 #include "AliDetectorRecoParam.h"
188 #include "AliRunInfo.h"
189 #include "AliEventInfo.h"
193 ClassImp(AliReconstruction)
195 //_____________________________________________________________________________
196 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
198 //_____________________________________________________________________________
199 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
201 fUniformField(kFALSE),
202 fForcedFieldMap(NULL),
203 fRunVertexFinder(kTRUE),
204 fRunVertexFinderTracks(kTRUE),
205 fRunHLTTracking(kFALSE),
206 fRunMuonTracking(kFALSE),
208 fRunCascadeFinder(kTRUE),
209 fStopOnError(kFALSE),
210 fWriteAlignmentData(kFALSE),
211 fWriteESDfriend(kFALSE),
212 fFillTriggerESD(kTRUE),
220 fRunLocalReconstruction("ALL"),
223 fUseTrackingErrorsForAlignment(""),
224 fGAliceFileName(gAliceFilename),
229 fNumberOfEventsPerFile(1),
231 fLoadAlignFromCDB(kTRUE),
232 fLoadAlignData("ALL"),
240 fParentRawReader(NULL),
245 fDiamondProfile(NULL),
246 fDiamondProfileTPC(NULL),
247 fMeanVertexConstraint(kTRUE),
251 fAlignObjArray(NULL),
254 fInitCDBCalled(kFALSE),
255 fSetRunNumberFromDataCalled(kFALSE),
262 fSameQACycle(kFALSE),
264 fRunPlaneEff(kFALSE),
273 fIsNewRunLoader(kFALSE),
277 // create reconstruction object with default parameters
280 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
281 fReconstructor[iDet] = NULL;
282 fLoader[iDet] = NULL;
283 fTracker[iDet] = NULL;
284 fQACycles[iDet] = 999999;
286 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
290 //_____________________________________________________________________________
291 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
293 fUniformField(rec.fUniformField),
294 fForcedFieldMap(NULL),
295 fRunVertexFinder(rec.fRunVertexFinder),
296 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
297 fRunHLTTracking(rec.fRunHLTTracking),
298 fRunMuonTracking(rec.fRunMuonTracking),
299 fRunV0Finder(rec.fRunV0Finder),
300 fRunCascadeFinder(rec.fRunCascadeFinder),
301 fStopOnError(rec.fStopOnError),
302 fWriteAlignmentData(rec.fWriteAlignmentData),
303 fWriteESDfriend(rec.fWriteESDfriend),
304 fFillTriggerESD(rec.fFillTriggerESD),
306 fCleanESD(rec.fCleanESD),
307 fV0DCAmax(rec.fV0DCAmax),
308 fV0CsPmin(rec.fV0CsPmin),
312 fRunLocalReconstruction(rec.fRunLocalReconstruction),
313 fRunTracking(rec.fRunTracking),
314 fFillESD(rec.fFillESD),
315 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
316 fGAliceFileName(rec.fGAliceFileName),
317 fRawInput(rec.fRawInput),
318 fEquipIdMap(rec.fEquipIdMap),
319 fFirstEvent(rec.fFirstEvent),
320 fLastEvent(rec.fLastEvent),
321 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
323 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
324 fLoadAlignData(rec.fLoadAlignData),
325 fESDPar(rec.fESDPar),
326 fUseHLTData(rec.fUseHLTData),
332 fParentRawReader(NULL),
334 fRecoParam(rec.fRecoParam),
337 fDiamondProfile(rec.fDiamondProfile),
338 fDiamondProfileTPC(rec.fDiamondProfileTPC),
339 fMeanVertexConstraint(rec.fMeanVertexConstraint),
343 fAlignObjArray(rec.fAlignObjArray),
344 fCDBUri(rec.fCDBUri),
346 fInitCDBCalled(rec.fInitCDBCalled),
347 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
348 fQADetectors(rec.fQADetectors),
350 fQATasks(rec.fQATasks),
352 fRunGlobalQA(rec.fRunGlobalQA),
353 fInLoopQA(rec.fInLoopQA),
354 fSameQACycle(rec.fSameQACycle),
355 fRunPlaneEff(rec.fRunPlaneEff),
364 fIsNewRunLoader(rec.fIsNewRunLoader),
370 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
371 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
373 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
374 fReconstructor[iDet] = NULL;
375 fLoader[iDet] = NULL;
376 fTracker[iDet] = NULL;
377 fQACycles[iDet] = rec.fQACycles[iDet];
379 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
380 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
385 //_____________________________________________________________________________
386 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
388 // assignment operator
389 // Used in PROOF mode
390 // Be very careful while modifing it!
391 // Simple rules to follow:
392 // for persistent data members - use their assignment operators
393 // for non-persistent ones - do nothing or take the default values from constructor
394 // TSelector members should not be touched
395 if(&rec == this) return *this;
397 fUniformField = rec.fUniformField;
398 fForcedFieldMap = NULL;
399 fRunVertexFinder = rec.fRunVertexFinder;
400 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
401 fRunHLTTracking = rec.fRunHLTTracking;
402 fRunMuonTracking = rec.fRunMuonTracking;
403 fRunV0Finder = rec.fRunV0Finder;
404 fRunCascadeFinder = rec.fRunCascadeFinder;
405 fStopOnError = rec.fStopOnError;
406 fWriteAlignmentData = rec.fWriteAlignmentData;
407 fWriteESDfriend = rec.fWriteESDfriend;
408 fFillTriggerESD = rec.fFillTriggerESD;
410 fCleanESD = rec.fCleanESD;
411 fV0DCAmax = rec.fV0DCAmax;
412 fV0CsPmin = rec.fV0CsPmin;
416 fRunLocalReconstruction = rec.fRunLocalReconstruction;
417 fRunTracking = rec.fRunTracking;
418 fFillESD = rec.fFillESD;
419 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
420 fGAliceFileName = rec.fGAliceFileName;
421 fRawInput = rec.fRawInput;
422 fEquipIdMap = rec.fEquipIdMap;
423 fFirstEvent = rec.fFirstEvent;
424 fLastEvent = rec.fLastEvent;
425 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
427 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
428 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
431 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
432 fLoadAlignData = rec.fLoadAlignData;
433 fESDPar = rec.fESDPar;
434 fUseHLTData = rec.fUseHLTData;
436 delete fRunInfo; fRunInfo = NULL;
437 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
439 fEventInfo = rec.fEventInfo;
443 fParentRawReader = NULL;
445 fRecoParam = rec.fRecoParam;
447 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
448 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
449 delete fLoader[iDet]; fLoader[iDet] = NULL;
450 delete fTracker[iDet]; fTracker[iDet] = NULL;
451 fQACycles[iDet] = rec.fQACycles[iDet];
455 delete fDiamondProfile; fDiamondProfile = NULL;
456 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
457 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
458 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
459 fMeanVertexConstraint = rec.fMeanVertexConstraint;
461 delete fGRPData; fGRPData = NULL;
462 if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
464 delete fAlignObjArray; fAlignObjArray = NULL;
467 fSpecCDBUri.Delete();
468 fInitCDBCalled = rec.fInitCDBCalled;
469 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
470 fQADetectors = rec.fQADetectors;
472 fQATasks = rec.fQATasks;
474 fRunGlobalQA = rec.fRunGlobalQA;
475 fInLoopQA = rec.fInLoopQA;
476 fSameQACycle = rec.fSameQACycle;
477 fRunPlaneEff = rec.fRunPlaneEff;
486 fIsNewRunLoader = rec.fIsNewRunLoader;
493 //_____________________________________________________________________________
494 AliReconstruction::~AliReconstruction()
499 delete fForcedFieldMap;
501 if (fAlignObjArray) {
502 fAlignObjArray->Delete();
503 delete fAlignObjArray;
505 fSpecCDBUri.Delete();
507 AliCodeTimer::Instance()->Print();
510 //_____________________________________________________________________________
511 void AliReconstruction::InitCDB()
513 // activate a default CDB storage
514 // First check if we have any CDB storage set, because it is used
515 // to retrieve the calibration and alignment constants
516 AliCodeTimerAuto("");
518 if (fInitCDBCalled) return;
519 fInitCDBCalled = kTRUE;
521 AliCDBManager* man = AliCDBManager::Instance();
522 if (man->IsDefaultStorageSet())
524 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
525 AliWarning("Default CDB storage has been already set !");
526 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
527 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
528 fCDBUri = man->GetDefaultStorage()->GetURI();
531 if (fCDBUri.Length() > 0)
533 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
534 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
535 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
537 fCDBUri="local://$ALICE_ROOT";
538 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
539 AliWarning("Default CDB storage not yet set !!!!");
540 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
541 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
544 man->SetDefaultStorage(fCDBUri);
547 // Now activate the detector specific CDB storage locations
548 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
549 TObject* obj = fSpecCDBUri[i];
551 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
552 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
553 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
554 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
556 AliSysInfo::AddStamp("InitCDB");
559 //_____________________________________________________________________________
560 void AliReconstruction::SetDefaultStorage(const char* uri) {
561 // Store the desired default CDB storage location
562 // Activate it later within the Run() method
568 //_____________________________________________________________________________
569 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
570 // Store a detector-specific CDB storage location
571 // Activate it later within the Run() method
573 AliCDBPath aPath(calibType);
574 if(!aPath.IsValid()){
575 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
576 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
577 if(!strcmp(calibType, fgkDetectorName[iDet])) {
578 aPath.SetPath(Form("%s/*", calibType));
579 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
583 if(!aPath.IsValid()){
584 AliError(Form("Not a valid path or detector: %s", calibType));
589 // // check that calibType refers to a "valid" detector name
590 // Bool_t isDetector = kFALSE;
591 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
592 // TString detName = fgkDetectorName[iDet];
593 // if(aPath.GetLevel0() == detName) {
594 // isDetector = kTRUE;
600 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
604 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
605 if (obj) fSpecCDBUri.Remove(obj);
606 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
610 //_____________________________________________________________________________
611 Bool_t AliReconstruction::SetRunNumberFromData()
613 // The method is called in Run() in order
614 // to set a correct run number.
615 // In case of raw data reconstruction the
616 // run number is taken from the raw data header
618 if (fSetRunNumberFromDataCalled) return kTRUE;
619 fSetRunNumberFromDataCalled = kTRUE;
621 AliCDBManager* man = AliCDBManager::Instance();
624 if(fRawReader->NextEvent()) {
625 if(man->GetRun() > 0) {
626 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
628 man->SetRun(fRawReader->GetRunNumber());
629 fRawReader->RewindEvents();
632 if(man->GetRun() > 0) {
633 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
636 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
642 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
644 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
649 // read run number from gAlice
650 if(rl->GetHeader()) {
651 man->SetRun(rl->GetHeader()->GetRun());
656 AliError("Neither run-loader header nor RawReader objects are found !");
668 //_____________________________________________________________________________
669 void AliReconstruction::SetCDBLock() {
670 // Set CDB lock: from now on it is forbidden to reset the run number
671 // or the default storage or to activate any further storage!
673 AliCDBManager::Instance()->SetLock(1);
676 //_____________________________________________________________________________
677 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
679 // Read the alignment objects from CDB.
680 // Each detector is supposed to have the
681 // alignment objects in DET/Align/Data CDB path.
682 // All the detector objects are then collected,
683 // sorted by geometry level (starting from ALIC) and
684 // then applied to the TGeo geometry.
685 // Finally an overlaps check is performed.
687 // Load alignment data from CDB and fill fAlignObjArray
688 if(fLoadAlignFromCDB){
690 TString detStr = detectors;
691 TString loadAlObjsListOfDets = "";
693 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
694 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
695 loadAlObjsListOfDets += fgkDetectorName[iDet];
696 loadAlObjsListOfDets += " ";
697 } // end loop over detectors
698 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
699 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
701 // Check if the array with alignment objects was
702 // provided by the user. If yes, apply the objects
703 // to the present TGeo geometry
704 if (fAlignObjArray) {
705 if (gGeoManager && gGeoManager->IsClosed()) {
706 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
707 AliError("The misalignment of one or more volumes failed!"
708 "Compare the list of simulated detectors and the list of detector alignment data!");
713 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
719 if (fAlignObjArray) {
720 fAlignObjArray->Delete();
721 delete fAlignObjArray; fAlignObjArray=NULL;
727 //_____________________________________________________________________________
728 void AliReconstruction::SetGAliceFile(const char* fileName)
730 // set the name of the galice file
732 fGAliceFileName = fileName;
735 //_____________________________________________________________________________
736 void AliReconstruction::SetInput(const char* input)
738 // In case the input string starts with 'mem://', we run in an online mode
739 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
740 // file is assumed. One can give as an input:
741 // mem://: - events taken from DAQ monitoring libs online
743 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
744 if (input) fRawInput = input;
747 //_____________________________________________________________________________
748 void AliReconstruction::SetOption(const char* detector, const char* option)
750 // set options for the reconstruction of a detector
752 TObject* obj = fOptions.FindObject(detector);
753 if (obj) fOptions.Remove(obj);
754 fOptions.Add(new TNamed(detector, option));
757 //_____________________________________________________________________________
758 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
760 // Set custom reconstruction parameters for a given detector
761 // Single set of parameters for all the events
762 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
763 if(!strcmp(detector, fgkDetectorName[iDet])) {
765 fRecoParam.AddDetRecoParam(iDet,par);
772 //_____________________________________________________________________________
773 Bool_t AliReconstruction::SetFieldMap(Float_t l3Current, Float_t diCurrent, Float_t factor, const char *path) {
774 //------------------------------------------------
775 // The magnetic field map, defined externally...
776 // L3 current 30000 A -> 0.5 T
777 // L3 current 12000 A -> 0.2 T
778 // dipole current 6000 A
779 // The polarities must be the same
780 //------------------------------------------------
781 const Float_t l3NominalCurrent1=30000.; // (A)
782 const Float_t l3NominalCurrent2=12000.; // (A)
783 const Float_t diNominalCurrent =6000. ; // (A)
785 const Float_t tolerance=0.03; // relative current tolerance
786 const Float_t zero=77.; // "zero" current (A)
789 Bool_t dipoleON=kFALSE;
791 TString s=(factor < 0) ? "L3: -" : "L3: +";
793 l3Current = TMath::Abs(l3Current);
794 if (TMath::Abs(l3Current-l3NominalCurrent1)/l3NominalCurrent1 < tolerance) {
795 map=AliMagWrapCheb::k5kG;
798 if (TMath::Abs(l3Current-l3NominalCurrent2)/l3NominalCurrent2 < tolerance) {
799 map=AliMagWrapCheb::k2kG;
802 if (l3Current < zero) {
803 map=AliMagWrapCheb::k2kG;
805 factor=0.; // in fact, this is a global factor...
806 fUniformField=kTRUE; // track with the uniform (zero) B field
808 AliError(Form("Wrong L3 current (%f A)!",l3Current));
812 diCurrent = TMath::Abs(diCurrent);
813 if (TMath::Abs(diCurrent-diNominalCurrent)/diNominalCurrent < tolerance) {
814 // 3% current tolerance...
818 if (diCurrent < zero) { // some small current..
822 AliError(Form("Wrong dipole current (%f A)!",diCurrent));
826 delete fForcedFieldMap;
828 new AliMagWrapCheb("B field map ",s,2,factor,10.,map,dipoleON,path);
830 fForcedFieldMap->Print();
832 AliTracker::SetFieldMap(fForcedFieldMap,fUniformField);
838 Bool_t AliReconstruction::InitGRP() {
839 //------------------------------------
840 // Initialization of the GRP entry
841 //------------------------------------
842 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
845 fGRPData = dynamic_cast<TMap*>(entry->GetObject());
847 AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
851 AliError("No GRP entry found in OCDB!");
855 TObjString *lhcState=
856 dynamic_cast<TObjString*>(fGRPData->GetValue("fLHCState"));
858 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
861 TObjString *beamType=
862 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamType"));
864 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
867 TObjString *beamEnergyStr=
868 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamEnergy"));
869 if (!beamEnergyStr) {
870 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
874 dynamic_cast<TObjString*>(fGRPData->GetValue("fRunType"));
876 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
879 TObjString *activeDetectors=
880 dynamic_cast<TObjString*>(fGRPData->GetValue("fDetectorMask"));
881 if (!activeDetectors) {
882 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
885 fRunInfo = new AliRunInfo(lhcState ? lhcState->GetString().Data() : "UNKNOWN",
886 beamType ? beamType->GetString().Data() : "UNKNOWN",
887 beamEnergyStr ? beamEnergyStr->GetString().Atof() : 0,
888 runType ? runType->GetString().Data() : "UNKNOWN",
889 activeDetectors ? activeDetectors->GetString().Atoi() : 1074790399);
891 // Process the list of active detectors
892 if (activeDetectors && activeDetectors->GetString().IsDigit()) {
893 UInt_t detMask = activeDetectors->GetString().Atoi();
894 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
895 fRunTracking = MatchDetectorList(fRunTracking,detMask);
896 fFillESD = MatchDetectorList(fFillESD,detMask);
897 fQADetectors = MatchDetectorList(fQADetectors,detMask);
900 AliInfo("===================================================================================");
901 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
902 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
903 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
904 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
905 AliInfo("===================================================================================");
907 //*** Dealing with the magnetic field map
908 if (AliTracker::GetFieldMap()) {
909 AliInfo("Running with the externally set B field !");
911 // Construct the field map out of the information retrieved from GRP.
916 TObjString *l3Current=
917 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Current"));
919 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
922 TObjString *l3Polarity=
923 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Polarity"));
925 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
930 TObjString *diCurrent=
931 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipoleCurrent"));
933 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
936 TObjString *diPolarity=
937 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipolePolarity"));
939 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
944 Float_t l3Cur=TMath::Abs(atof(l3Current->GetName()));
945 Float_t diCur=TMath::Abs(atof(diCurrent->GetName()));
946 Float_t l3Pol=atof(l3Polarity->GetName());
948 if (l3Pol != 0.) factor=-1.;
951 if (!SetFieldMap(l3Cur, diCur, factor)) {
952 AliFatal("Failed to creat a B field map ! Exiting...");
954 AliInfo("Running with the B field constructed out of GRP !");
957 AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
962 //*** Get the diamond profile from OCDB
963 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
965 if (fMeanVertexConstraint)
966 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
968 AliError("No diamond profile found in OCDB!");
971 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
973 if (fMeanVertexConstraint)
974 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
976 AliError("No diamond profile found in OCDB!");
982 //_____________________________________________________________________________
983 Bool_t AliReconstruction::LoadCDB()
985 AliCodeTimerAuto("");
987 AliCDBManager::Instance()->Get("GRP/CTP/Config");
989 TString detStr = fRunLocalReconstruction;
990 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
991 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
992 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
997 //_____________________________________________________________________________
998 Bool_t AliReconstruction::Run(const char* input)
1001 AliCodeTimerAuto("");
1004 if (GetAbort() != TSelector::kContinue) return kFALSE;
1006 TChain *chain = NULL;
1007 if (fRawReader && (chain = fRawReader->GetChain())) {
1010 gProof->AddInput(this);
1012 chain->Process("AliReconstruction");
1015 chain->Process(this);
1020 if (GetAbort() != TSelector::kContinue) return kFALSE;
1022 if (GetAbort() != TSelector::kContinue) return kFALSE;
1023 //******* The loop over events
1025 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1026 (fRawReader && fRawReader->NextEvent())) {
1027 if (!ProcessEvent(iEvent)) {
1028 Abort("ProcessEvent",TSelector::kAbortFile);
1034 if (GetAbort() != TSelector::kContinue) return kFALSE;
1036 if (GetAbort() != TSelector::kContinue) return kFALSE;
1042 //_____________________________________________________________________________
1043 void AliReconstruction::InitRawReader(const char* input)
1045 AliCodeTimerAuto("");
1047 // Init raw-reader and
1048 // set the input in case of raw data
1049 if (input) fRawInput = input;
1050 fRawReader = AliRawReader::Create(fRawInput.Data());
1052 AliInfo("Reconstruction will run over digits");
1054 if (!fEquipIdMap.IsNull() && fRawReader)
1055 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1057 if (!fUseHLTData.IsNull()) {
1058 // create the RawReaderHLT which performs redirection of HLT input data for
1059 // the specified detectors
1060 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1062 fParentRawReader=fRawReader;
1063 fRawReader=pRawReader;
1065 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1068 AliSysInfo::AddStamp("CreateRawReader");
1071 //_____________________________________________________________________________
1072 void AliReconstruction::InitRun(const char* input)
1074 // Initialization of raw-reader,
1075 // run number, CDB etc.
1076 AliCodeTimerAuto("");
1077 AliSysInfo::AddStamp("Start");
1079 // Initialize raw-reader if any
1080 InitRawReader(input);
1082 // Initialize the CDB storage
1085 // Set run number in CDBManager (if it is not already set by the user)
1086 if (!SetRunNumberFromData()) {
1087 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1091 // Set CDB lock: from now on it is forbidden to reset the run number
1092 // or the default storage or to activate any further storage!
1097 //_____________________________________________________________________________
1098 void AliReconstruction::Begin(TTree *)
1100 // Initialize AlReconstruction before
1101 // going into the event loop
1102 // Should follow the TSelector convention
1103 // i.e. initialize only the object on the client side
1105 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1107 fInput->Remove(reco);
1109 AliSysInfo::AddStamp("ReadInputInBegin");
1112 // Import ideal TGeo geometry and apply misalignment
1114 TString geom(gSystem->DirName(fGAliceFileName));
1115 geom += "/geometry.root";
1116 AliGeomManager::LoadGeometry(geom.Data());
1118 Abort("LoadGeometry", TSelector::kAbortProcess);
1121 AliSysInfo::AddStamp("LoadGeom");
1122 TString detsToCheck=fRunLocalReconstruction;
1123 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1124 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1127 AliSysInfo::AddStamp("CheckGeom");
1130 if (!MisalignGeometry(fLoadAlignData)) {
1131 Abort("MisalignGeometry", TSelector::kAbortProcess);
1134 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1135 AliSysInfo::AddStamp("MisalignGeom");
1138 Abort("InitGRP", TSelector::kAbortProcess);
1141 AliSysInfo::AddStamp("InitGRP");
1144 Abort("LoadCDB", TSelector::kAbortProcess);
1147 AliSysInfo::AddStamp("LoadCDB");
1149 // Read the reconstruction parameters from OCDB
1150 if (!InitRecoParams()) {
1151 AliWarning("Not all detectors have correct RecoParam objects initialized");
1153 AliSysInfo::AddStamp("InitRecoParams");
1156 fInput->Add(gGeoManager);
1158 fInput->Add(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()));
1159 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1160 fInput->Add((AliMagF*)AliTracker::GetFieldMap());
1166 //_____________________________________________________________________________
1167 void AliReconstruction::SlaveBegin(TTree*)
1169 // Initialization related to run-loader,
1170 // vertexer, trackers, recontructors
1171 // In proof mode it is executed on the slave
1172 AliCodeTimerAuto("");
1175 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1177 AliGeomManager::SetGeometry(tgeo);
1179 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1180 Int_t runNumber = -1;
1181 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1182 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1183 man->SetCacheFlag(kTRUE);
1184 man->SetLock(kTRUE);
1188 if (AliMagF *map = (AliMagF*)fInput->FindObject("Maps")) {
1189 AliTracker::SetFieldMap(map,fUniformField);
1191 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1194 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1197 // get the run loader
1198 if (!InitRunLoader()) {
1199 Abort("InitRunLoader", TSelector::kAbortProcess);
1202 AliSysInfo::AddStamp("LoadLoader");
1204 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1205 if(fDiamondProfile && fMeanVertexConstraint) ftVertexer->SetVtxStart(fDiamondProfile);
1208 if (fRunVertexFinder && !CreateVertexer()) {
1209 Abort("CreateVertexer", TSelector::kAbortProcess);
1212 AliSysInfo::AddStamp("CreateVertexer");
1215 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1216 Abort("CreateTrackers", TSelector::kAbortProcess);
1219 AliSysInfo::AddStamp("CreateTrackers");
1221 // create the ESD output file and tree
1222 ffile = TFile::Open("AliESDs.root", "RECREATE");
1223 ffile->SetCompressionLevel(2);
1224 if (!ffile->IsOpen()) {
1225 Abort("OpenESDFile", TSelector::kAbortProcess);
1229 ftree = new TTree("esdTree", "Tree with ESD objects");
1230 fesd = new AliESDEvent();
1231 fesd->CreateStdContent();
1232 fesd->WriteToTree(ftree);
1234 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1235 fhltesd = new AliESDEvent();
1236 fhltesd->CreateStdContent();
1237 fhltesd->WriteToTree(fhlttree);
1240 if (fWriteESDfriend) {
1241 fesdf = new AliESDfriend();
1242 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1243 br->SetFile("AliESDfriends.root");
1244 fesd->AddObject(fesdf);
1247 ProcInfo_t ProcInfo;
1248 gSystem->GetProcInfo(&ProcInfo);
1249 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1252 fQASteer = new AliQADataMakerSteer("rec") ;
1253 fQASteer->SetActiveDetectors(fQADetectors) ;
1254 fQASteer->SetTasks(fQATasks) ;
1257 if (fRunQA && fRawReader && fQATasks.Contains(Form("%d", AliQA::kRAWS))) {
1258 fQASteer->Run(fQADetectors, fRawReader) ;
1259 fSameQACycle = kTRUE ;
1263 //Initialize the QA and start of cycle for out-of-loop QA
1265 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam, fSameQACycle, !fInLoopQA) ;
1269 fSameQACycle = kFALSE;
1270 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1271 AliInfo(Form("Initializing the global QA data maker"));
1272 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1273 TObjArray *arr=qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun());
1274 AliTracker::SetResidualsArray(arr);
1276 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
1277 fSameQACycle = kTRUE;
1280 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1281 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
1283 qadm->StartOfCycle(AliQA::kESDS, fSameQACycle);
1284 fSameQACycle = kTRUE;
1289 //Initialize the Plane Efficiency framework
1290 if (fRunPlaneEff && !InitPlaneEff()) {
1291 Abort("InitPlaneEff", TSelector::kAbortProcess);
1295 if (strcmp(gProgName,"alieve") == 0)
1296 fRunAliEVE = InitAliEVE();
1301 //_____________________________________________________________________________
1302 Bool_t AliReconstruction::Process(Long64_t entry)
1304 // run the reconstruction over a single entry
1305 // from the chain with raw data
1306 AliCodeTimerAuto("");
1308 TTree *currTree = fChain->GetTree();
1309 AliRawEvent *event = new AliRawEvent;
1310 currTree->SetBranchAddress("rawevent",&event);
1311 currTree->GetEntry(entry);
1312 fRawReader = new AliRawReaderRoot(event);
1313 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1321 //_____________________________________________________________________________
1322 void AliReconstruction::Init(TTree *tree)
1325 AliError("The input tree is not found!");
1331 //_____________________________________________________________________________
1332 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1334 // run the reconstruction over a single event
1335 // The event loop is steered in Run method
1337 AliCodeTimerAuto("");
1339 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1340 fRunLoader->SetEventNumber(iEvent);
1341 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1343 fRunLoader->TreeE()->Fill();
1344 if (fRawReader && fRawReader->UseAutoSaveESD())
1345 fRunLoader->TreeE()->AutoSave("SaveSelf");
1348 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1352 AliInfo(Form("processing event %d", iEvent));
1354 // Fill Event-info object
1356 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo);
1358 //Start of cycle for the in-loop QA
1359 if (fInLoopQA && fRunQA) {
1360 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam, fSameQACycle, fInLoopQA) ;
1362 if (fInLoopQA && fRunGlobalQA) {
1363 fSameQACycle = kFALSE;
1364 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1365 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1366 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
1367 fSameQACycle = kTRUE;
1369 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1370 qadm->StartOfCycle(AliQA::kESDS, fSameQACycle);
1371 fSameQACycle = kTRUE;
1375 fRunLoader->GetEvent(iEvent);
1378 if (fInLoopQA && fRunQA)
1379 fQASteer->RunOneEvent(fRawReader) ;
1381 // local single event reconstruction
1382 if (!fRunLocalReconstruction.IsNull()) {
1383 TString detectors=fRunLocalReconstruction;
1384 // run HLT event reconstruction first
1385 // ;-( IsSelected changes the string
1386 if (IsSelected("HLT", detectors) &&
1387 !RunLocalEventReconstruction("HLT")) {
1388 if (fStopOnError) {CleanUp(); return kFALSE;}
1390 detectors=fRunLocalReconstruction;
1391 detectors.ReplaceAll("HLT", "");
1392 if (!RunLocalEventReconstruction(detectors)) {
1393 if (fStopOnError) {CleanUp(); return kFALSE;}
1397 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1398 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1399 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1400 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1402 // Set magnetic field from the tracker
1403 fesd->SetMagneticField(AliTracker::GetBz());
1404 fhltesd->SetMagneticField(AliTracker::GetBz());
1408 // Fill raw-data error log into the ESD
1409 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1412 if (fRunVertexFinder) {
1413 if (!RunVertexFinder(fesd)) {
1414 if (fStopOnError) {CleanUp(); return kFALSE;}
1419 if (!fRunTracking.IsNull()) {
1420 if (fRunMuonTracking) {
1421 if (!RunMuonTracking(fesd)) {
1422 if (fStopOnError) {CleanUp(); return kFALSE;}
1428 if (!fRunTracking.IsNull()) {
1429 if (!RunTracking(fesd)) {
1430 if (fStopOnError) {CleanUp(); return kFALSE;}
1435 if (!fFillESD.IsNull()) {
1436 TString detectors=fFillESD;
1437 // run HLT first and on hltesd
1438 // ;-( IsSelected changes the string
1439 if (IsSelected("HLT", detectors) &&
1440 !FillESD(fhltesd, "HLT")) {
1441 if (fStopOnError) {CleanUp(); return kFALSE;}
1444 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1445 if (detectors.Contains("ALL")) {
1447 for (Int_t idet=0; idet<fgkNDetectors; ++idet){
1448 detectors += fgkDetectorName[idet];
1452 detectors.ReplaceAll("HLT", "");
1453 if (!FillESD(fesd, detectors)) {
1454 if (fStopOnError) {CleanUp(); return kFALSE;}
1458 // fill Event header information from the RawEventHeader
1459 if (fRawReader){FillRawEventHeaderESD(fesd);}
1462 AliESDpid::MakePID(fesd);
1464 if (fFillTriggerESD) {
1465 if (!FillTriggerESD(fesd)) {
1466 if (fStopOnError) {CleanUp(); return kFALSE;}
1473 // Propagate track to the beam pipe (if not already done by ITS)
1475 const Int_t ntracks = fesd->GetNumberOfTracks();
1476 const Double_t kBz = fesd->GetMagneticField();
1477 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1480 UShort_t *selectedIdx=new UShort_t[ntracks];
1482 for (Int_t itrack=0; itrack<ntracks; itrack++){
1483 const Double_t kMaxStep = 5; //max step over the material
1486 AliESDtrack *track = fesd->GetTrack(itrack);
1487 if (!track) continue;
1489 AliExternalTrackParam *tpcTrack =
1490 (AliExternalTrackParam *)track->GetTPCInnerParam();
1494 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kTRUE);
1497 Int_t n=trkArray.GetEntriesFast();
1498 selectedIdx[n]=track->GetID();
1499 trkArray.AddLast(tpcTrack);
1502 //Tracks refitted by ITS should already be at the SPD vertex
1503 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1506 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1507 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1512 // Improve the reconstructed primary vertex position using the tracks
1514 TObject *obj = fOptions.FindObject("ITS");
1516 TString optITS = obj->GetTitle();
1517 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
1518 fRunVertexFinderTracks=kFALSE;
1520 if (fRunVertexFinderTracks) {
1521 // TPC + ITS primary vertex
1522 ftVertexer->SetITSrefitRequired();
1523 if(fDiamondProfile && fMeanVertexConstraint) {
1524 ftVertexer->SetVtxStart(fDiamondProfile);
1526 ftVertexer->SetConstraintOff();
1528 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1530 if (pvtx->GetStatus()) {
1531 fesd->SetPrimaryVertex(pvtx);
1532 for (Int_t i=0; i<ntracks; i++) {
1533 AliESDtrack *t = fesd->GetTrack(i);
1534 t->RelateToVertex(pvtx, kBz, kVeryBig);
1539 // TPC-only primary vertex
1540 ftVertexer->SetITSrefitNotRequired();
1541 if(fDiamondProfileTPC && fMeanVertexConstraint) {
1542 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1544 ftVertexer->SetConstraintOff();
1546 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1548 if (pvtx->GetStatus()) {
1549 fesd->SetPrimaryVertexTPC(pvtx);
1550 for (Int_t i=0; i<ntracks; i++) {
1551 AliESDtrack *t = fesd->GetTrack(i);
1552 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1558 delete[] selectedIdx;
1560 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1565 AliV0vertexer vtxer;
1566 vtxer.Tracks2V0vertices(fesd);
1568 if (fRunCascadeFinder) {
1570 AliCascadeVertexer cvtxer;
1571 cvtxer.V0sTracks2CascadeVertices(fesd);
1576 if (fCleanESD) CleanESD(fesd);
1579 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1580 if (qadm && fQATasks.Contains(Form("%d", AliQA::kESDS)))
1581 qadm->Exec(AliQA::kESDS, fesd);
1584 if (fWriteESDfriend) {
1585 fesdf->~AliESDfriend();
1586 new (fesdf) AliESDfriend(); // Reset...
1587 fesd->GetESDfriend(fesdf);
1591 // Auto-save the ESD tree in case of prompt reco @P2
1592 if (fRawReader && fRawReader->UseAutoSaveESD())
1593 ftree->AutoSave("SaveSelf");
1599 if (fRunAliEVE) RunAliEVE();
1603 if (fWriteESDfriend) {
1604 fesdf->~AliESDfriend();
1605 new (fesdf) AliESDfriend(); // Reset...
1608 ProcInfo_t ProcInfo;
1609 gSystem->GetProcInfo(&ProcInfo);
1610 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1613 // End of cycle for the in-loop
1614 if (fInLoopQA && fRunQA) {
1615 fQASteer->RunOneEvent(fesd) ;
1616 fQASteer->EndOfCycle() ;
1618 if (fInLoopQA && fRunGlobalQA) {
1619 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1621 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1622 qadm->EndOfCycle(AliQA::kRECPOINTS);
1623 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1624 qadm->EndOfCycle(AliQA::kESDS);
1630 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1631 if (fReconstructor[iDet])
1632 fReconstructor[iDet]->SetRecoParam(NULL);
1638 //_____________________________________________________________________________
1639 void AliReconstruction::SlaveTerminate()
1641 // Finalize the run on the slave side
1642 // Called after the exit
1643 // from the event loop
1644 AliCodeTimerAuto("");
1646 if (fIsNewRunLoader) { // galice.root didn't exist
1647 fRunLoader->WriteHeader("OVERWRITE");
1648 fRunLoader->CdGAFile();
1649 fRunLoader->Write(0, TObject::kOverwrite);
1652 ftree->GetUserInfo()->Add(fesd);
1653 fhlttree->GetUserInfo()->Add(fhltesd);
1655 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1656 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1658 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1659 cdbMapCopy->SetOwner(1);
1660 cdbMapCopy->SetName("cdbMap");
1661 TIter iter(cdbMap->GetTable());
1664 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1665 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1666 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1667 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1670 TList *cdbListCopy = new TList();
1671 cdbListCopy->SetOwner(1);
1672 cdbListCopy->SetName("cdbList");
1674 TIter iter2(cdbList);
1677 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1678 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1681 ftree->GetUserInfo()->Add(cdbMapCopy);
1682 ftree->GetUserInfo()->Add(cdbListCopy);
1685 if(fESDPar.Contains("ESD.par")){
1686 AliInfo("Attaching ESD.par to Tree");
1687 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
1688 ftree->GetUserInfo()->Add(fn);
1694 if (fWriteESDfriend)
1695 ftree->SetBranchStatus("ESDfriend*",0);
1696 // we want to have only one tree version number
1697 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1700 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1701 if (fRunPlaneEff && !FinishPlaneEff()) {
1702 AliWarning("Finish PlaneEff evaluation failed");
1705 //Finish QA and end of cycle for out-of-loop QA
1706 if (!fInLoopQA && fRunQA)
1707 fQASteer->Run(fRunLocalReconstruction.Data(), AliQA::kNULLTASKINDEX, fSameQACycle) ;
1708 if (!fInLoopQA && fRunGlobalQA) {
1709 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1711 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1712 qadm->EndOfCycle(AliQA::kRECPOINTS);
1713 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1714 qadm->EndOfCycle(AliQA::kESDS);
1723 //_____________________________________________________________________________
1724 void AliReconstruction::Terminate()
1726 // Create tags for the events in the ESD tree (the ESD tree is always present)
1727 // In case of empty events the tags will contain dummy values
1728 AliCodeTimerAuto("");
1730 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1731 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData);
1733 // Cleanup of CDB manager: cache and active storages!
1734 AliCDBManager::Instance()->ClearCache();
1737 //_____________________________________________________________________________
1738 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1740 // run the local reconstruction
1742 static Int_t eventNr=0;
1743 AliCodeTimerAuto("")
1745 TString detStr = detectors;
1746 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1747 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1748 AliReconstructor* reconstructor = GetReconstructor(iDet);
1749 if (!reconstructor) continue;
1750 AliLoader* loader = fLoader[iDet];
1751 // Matthias April 2008: temporary fix to run HLT reconstruction
1752 // although the HLT loader is missing
1753 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1755 reconstructor->Reconstruct(fRawReader, NULL);
1758 reconstructor->Reconstruct(dummy, NULL);
1763 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1766 // conversion of digits
1767 if (fRawReader && reconstructor->HasDigitConversion()) {
1768 AliInfo(Form("converting raw data digits into root objects for %s",
1769 fgkDetectorName[iDet]));
1770 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1771 // fgkDetectorName[iDet]));
1772 loader->LoadDigits("update");
1773 loader->CleanDigits();
1774 loader->MakeDigitsContainer();
1775 TTree* digitsTree = loader->TreeD();
1776 reconstructor->ConvertDigits(fRawReader, digitsTree);
1777 loader->WriteDigits("OVERWRITE");
1778 loader->UnloadDigits();
1780 // local reconstruction
1781 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1782 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1783 loader->LoadRecPoints("update");
1784 loader->CleanRecPoints();
1785 loader->MakeRecPointsContainer();
1786 TTree* clustersTree = loader->TreeR();
1787 if (fRawReader && !reconstructor->HasDigitConversion()) {
1788 reconstructor->Reconstruct(fRawReader, clustersTree);
1790 loader->LoadDigits("read");
1791 TTree* digitsTree = loader->TreeD();
1793 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1794 if (fStopOnError) return kFALSE;
1796 reconstructor->Reconstruct(digitsTree, clustersTree);
1798 loader->UnloadDigits();
1801 // In-loop QA for local reconstrucion
1802 TString detQAStr(fQADetectors) ;
1803 if (fRunQA && fInLoopQA)
1804 fQASteer->RunOneEventInOneDetector(iDet, clustersTree) ;
1806 loader->WriteRecPoints("OVERWRITE");
1807 loader->UnloadRecPoints();
1808 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1810 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1811 AliError(Form("the following detectors were not found: %s",
1813 if (fStopOnError) return kFALSE;
1819 //_____________________________________________________________________________
1820 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1822 // run the barrel tracking
1824 AliCodeTimerAuto("")
1826 AliESDVertex* vertex = NULL;
1827 Double_t vtxPos[3] = {0, 0, 0};
1828 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1829 TArrayF mcVertex(3);
1830 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1831 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1832 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1836 AliInfo("running the ITS vertex finder");
1838 fLoader[0]->LoadRecPoints();
1839 TTree* cltree = fLoader[0]->TreeR();
1841 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1842 vertex = fVertexer->FindVertexForCurrentEvent(cltree);
1845 AliError("Can't get the ITS cluster tree");
1847 fLoader[0]->UnloadRecPoints();
1850 AliError("Can't get the ITS loader");
1853 AliWarning("Vertex not found");
1854 vertex = new AliESDVertex();
1855 vertex->SetName("default");
1858 vertex->SetName("reconstructed");
1862 AliInfo("getting the primary vertex from MC");
1863 vertex = new AliESDVertex(vtxPos, vtxErr);
1867 vertex->GetXYZ(vtxPos);
1868 vertex->GetSigmaXYZ(vtxErr);
1870 AliWarning("no vertex reconstructed");
1871 vertex = new AliESDVertex(vtxPos, vtxErr);
1873 esd->SetPrimaryVertexSPD(vertex);
1874 // if SPD multiplicity has been determined, it is stored in the ESD
1875 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1876 if(mult)esd->SetMultiplicity(mult);
1878 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1879 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1886 //_____________________________________________________________________________
1887 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1889 // run the HLT barrel tracking
1891 AliCodeTimerAuto("")
1894 AliError("Missing runLoader!");
1898 AliInfo("running HLT tracking");
1900 // Get a pointer to the HLT reconstructor
1901 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1902 if (!reconstructor) return kFALSE;
1905 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1906 TString detName = fgkDetectorName[iDet];
1907 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1908 reconstructor->SetOption(detName.Data());
1909 AliTracker *tracker = reconstructor->CreateTracker();
1911 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1912 if (fStopOnError) return kFALSE;
1916 Double_t vtxErr[3]={0.005,0.005,0.010};
1917 const AliESDVertex *vertex = esd->GetVertex();
1918 vertex->GetXYZ(vtxPos);
1919 tracker->SetVertex(vtxPos,vtxErr);
1921 fLoader[iDet]->LoadRecPoints("read");
1922 TTree* tree = fLoader[iDet]->TreeR();
1924 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1927 tracker->LoadClusters(tree);
1929 if (tracker->Clusters2Tracks(esd) != 0) {
1930 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1934 tracker->UnloadClusters();
1942 //_____________________________________________________________________________
1943 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1945 // run the muon spectrometer tracking
1947 AliCodeTimerAuto("")
1950 AliError("Missing runLoader!");
1953 Int_t iDet = 7; // for MUON
1955 AliInfo("is running...");
1957 // Get a pointer to the MUON reconstructor
1958 AliReconstructor *reconstructor = GetReconstructor(iDet);
1959 if (!reconstructor) return kFALSE;
1962 TString detName = fgkDetectorName[iDet];
1963 AliDebug(1, Form("%s tracking", detName.Data()));
1964 AliTracker *tracker = reconstructor->CreateTracker();
1966 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1971 fLoader[iDet]->LoadRecPoints("read");
1973 tracker->LoadClusters(fLoader[iDet]->TreeR());
1975 Int_t rv = tracker->Clusters2Tracks(esd);
1979 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1983 fLoader[iDet]->UnloadRecPoints();
1985 tracker->UnloadClusters();
1993 //_____________________________________________________________________________
1994 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1996 // run the barrel tracking
1997 static Int_t eventNr=0;
1998 AliCodeTimerAuto("")
2000 AliInfo("running tracking");
2002 //Fill the ESD with the T0 info (will be used by the TOF)
2003 if (fReconstructor[11] && fLoader[11]) {
2004 fLoader[11]->LoadRecPoints("READ");
2005 TTree *treeR = fLoader[11]->TreeR();
2006 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2009 // pass 1: TPC + ITS inwards
2010 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2011 if (!fTracker[iDet]) continue;
2012 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2015 fLoader[iDet]->LoadRecPoints("read");
2016 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2017 TTree* tree = fLoader[iDet]->TreeR();
2019 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2022 fTracker[iDet]->LoadClusters(tree);
2023 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2025 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2026 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2029 // preliminary PID in TPC needed by the ITS tracker
2031 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2032 AliESDpid::MakePID(esd);
2034 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2037 // pass 2: ALL backwards
2039 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2040 if (!fTracker[iDet]) continue;
2041 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2044 if (iDet > 1) { // all except ITS, TPC
2046 fLoader[iDet]->LoadRecPoints("read");
2047 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2048 tree = fLoader[iDet]->TreeR();
2050 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2053 fTracker[iDet]->LoadClusters(tree);
2054 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2058 if (iDet>1) // start filling residuals for the "outer" detectors
2059 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
2061 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2062 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2067 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2068 fTracker[iDet]->UnloadClusters();
2069 fLoader[iDet]->UnloadRecPoints();
2071 // updated PID in TPC needed by the ITS tracker -MI
2073 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2074 AliESDpid::MakePID(esd);
2076 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2078 //stop filling residuals for the "outer" detectors
2079 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
2081 // pass 3: TRD + TPC + ITS refit inwards
2083 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2084 if (!fTracker[iDet]) continue;
2085 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2088 if (iDet<2) // start filling residuals for TPC and ITS
2089 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
2091 if (fTracker[iDet]->RefitInward(esd) != 0) {
2092 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2095 // run postprocessing
2096 if (fTracker[iDet]->PostProcess(esd) != 0) {
2097 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2100 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2103 // write space-points to the ESD in case alignment data output
2105 if (fWriteAlignmentData)
2106 WriteAlignmentData(esd);
2108 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2109 if (!fTracker[iDet]) continue;
2111 fTracker[iDet]->UnloadClusters();
2112 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2113 fLoader[iDet]->UnloadRecPoints();
2114 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2116 // stop filling residuals for TPC and ITS
2117 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
2123 //_____________________________________________________________________________
2124 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2126 // Remove the data which are not needed for the physics analysis.
2129 Int_t nTracks=esd->GetNumberOfTracks();
2130 Int_t nV0s=esd->GetNumberOfV0s();
2132 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2134 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2135 Bool_t rc=esd->Clean(cleanPars);
2137 nTracks=esd->GetNumberOfTracks();
2138 nV0s=esd->GetNumberOfV0s();
2140 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2145 //_____________________________________________________________________________
2146 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2148 // fill the event summary data
2150 AliCodeTimerAuto("")
2151 static Int_t eventNr=0;
2152 TString detStr = detectors;
2154 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2155 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2156 AliReconstructor* reconstructor = GetReconstructor(iDet);
2157 if (!reconstructor) continue;
2158 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2159 TTree* clustersTree = NULL;
2160 if (fLoader[iDet]) {
2161 fLoader[iDet]->LoadRecPoints("read");
2162 clustersTree = fLoader[iDet]->TreeR();
2163 if (!clustersTree) {
2164 AliError(Form("Can't get the %s clusters tree",
2165 fgkDetectorName[iDet]));
2166 if (fStopOnError) return kFALSE;
2169 if (fRawReader && !reconstructor->HasDigitConversion()) {
2170 reconstructor->FillESD(fRawReader, clustersTree, esd);
2172 TTree* digitsTree = NULL;
2173 if (fLoader[iDet]) {
2174 fLoader[iDet]->LoadDigits("read");
2175 digitsTree = fLoader[iDet]->TreeD();
2177 AliError(Form("Can't get the %s digits tree",
2178 fgkDetectorName[iDet]));
2179 if (fStopOnError) return kFALSE;
2182 reconstructor->FillESD(digitsTree, clustersTree, esd);
2183 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2185 if (fLoader[iDet]) {
2186 fLoader[iDet]->UnloadRecPoints();
2190 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2191 AliError(Form("the following detectors were not found: %s",
2193 if (fStopOnError) return kFALSE;
2195 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2200 //_____________________________________________________________________________
2201 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2203 // Reads the trigger decision which is
2204 // stored in Trigger.root file and fills
2205 // the corresponding esd entries
2207 AliCodeTimerAuto("")
2209 AliInfo("Filling trigger information into the ESD");
2212 AliCTPRawStream input(fRawReader);
2213 if (!input.Next()) {
2214 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2217 if (esd->GetTriggerMask() != input.GetClassMask())
2218 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2219 input.GetClassMask(),esd->GetTriggerMask()));
2220 if (esd->GetOrbitNumber() != input.GetOrbitID())
2221 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2222 input.GetOrbitID(),esd->GetOrbitNumber()));
2223 if (esd->GetBunchCrossNumber() != input.GetBCID())
2224 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2225 input.GetBCID(),esd->GetBunchCrossNumber()));
2228 // Here one has to add the filling of trigger inputs and
2229 // interaction records
2239 //_____________________________________________________________________________
2240 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2243 // Filling information from RawReader Header
2246 if (!fRawReader) return kFALSE;
2248 AliInfo("Filling information from RawReader Header");
2250 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2251 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2252 esd->SetPeriodNumber(fRawReader->GetPeriod());
2254 esd->SetTimeStamp(fRawReader->GetTimestamp());
2255 esd->SetEventType(fRawReader->GetType());
2261 //_____________________________________________________________________________
2262 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2264 // check whether detName is contained in detectors
2265 // if yes, it is removed from detectors
2267 // check if all detectors are selected
2268 if ((detectors.CompareTo("ALL") == 0) ||
2269 detectors.BeginsWith("ALL ") ||
2270 detectors.EndsWith(" ALL") ||
2271 detectors.Contains(" ALL ")) {
2276 // search for the given detector
2277 Bool_t result = kFALSE;
2278 if ((detectors.CompareTo(detName) == 0) ||
2279 detectors.BeginsWith(detName+" ") ||
2280 detectors.EndsWith(" "+detName) ||
2281 detectors.Contains(" "+detName+" ")) {
2282 detectors.ReplaceAll(detName, "");
2286 // clean up the detectors string
2287 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2288 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2289 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2294 //_____________________________________________________________________________
2295 Bool_t AliReconstruction::InitRunLoader()
2297 // get or create the run loader
2299 if (gAlice) delete gAlice;
2302 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2303 // load all base libraries to get the loader classes
2304 TString libs = gSystem->GetLibraries();
2305 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2306 TString detName = fgkDetectorName[iDet];
2307 if (detName == "HLT") continue;
2308 if (libs.Contains("lib" + detName + "base.so")) continue;
2309 gSystem->Load("lib" + detName + "base.so");
2311 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2313 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2318 fRunLoader->CdGAFile();
2319 fRunLoader->LoadgAlice();
2321 //PH This is a temporary fix to give access to the kinematics
2322 //PH that is needed for the labels of ITS clusters
2323 fRunLoader->LoadHeader();
2324 fRunLoader->LoadKinematics();
2326 } else { // galice.root does not exist
2328 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2330 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2331 AliConfig::GetDefaultEventFolderName(),
2334 AliError(Form("could not create run loader in file %s",
2335 fGAliceFileName.Data()));
2339 fIsNewRunLoader = kTRUE;
2340 fRunLoader->MakeTree("E");
2342 if (fNumberOfEventsPerFile > 0)
2343 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2345 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2351 //_____________________________________________________________________________
2352 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2354 // get the reconstructor object and the loader for a detector
2356 if (fReconstructor[iDet]) {
2357 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2358 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2359 fReconstructor[iDet]->SetRecoParam(par);
2361 return fReconstructor[iDet];
2364 // load the reconstructor object
2365 TPluginManager* pluginManager = gROOT->GetPluginManager();
2366 TString detName = fgkDetectorName[iDet];
2367 TString recName = "Ali" + detName + "Reconstructor";
2369 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2371 AliReconstructor* reconstructor = NULL;
2372 // first check if a plugin is defined for the reconstructor
2373 TPluginHandler* pluginHandler =
2374 pluginManager->FindHandler("AliReconstructor", detName);
2375 // if not, add a plugin for it
2376 if (!pluginHandler) {
2377 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2378 TString libs = gSystem->GetLibraries();
2379 if (libs.Contains("lib" + detName + "base.so") ||
2380 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2381 pluginManager->AddHandler("AliReconstructor", detName,
2382 recName, detName + "rec", recName + "()");
2384 pluginManager->AddHandler("AliReconstructor", detName,
2385 recName, detName, recName + "()");
2387 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2389 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2390 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2392 if (reconstructor) {
2393 TObject* obj = fOptions.FindObject(detName.Data());
2394 if (obj) reconstructor->SetOption(obj->GetTitle());
2395 reconstructor->Init();
2396 fReconstructor[iDet] = reconstructor;
2399 // get or create the loader
2400 if (detName != "HLT") {
2401 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2402 if (!fLoader[iDet]) {
2403 AliConfig::Instance()
2404 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2406 // first check if a plugin is defined for the loader
2408 pluginManager->FindHandler("AliLoader", detName);
2409 // if not, add a plugin for it
2410 if (!pluginHandler) {
2411 TString loaderName = "Ali" + detName + "Loader";
2412 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2413 pluginManager->AddHandler("AliLoader", detName,
2414 loaderName, detName + "base",
2415 loaderName + "(const char*, TFolder*)");
2416 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2418 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2420 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2421 fRunLoader->GetEventFolder());
2423 if (!fLoader[iDet]) { // use default loader
2424 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2426 if (!fLoader[iDet]) {
2427 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2428 if (fStopOnError) return NULL;
2430 fRunLoader->AddLoader(fLoader[iDet]);
2431 fRunLoader->CdGAFile();
2432 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2433 fRunLoader->Write(0, TObject::kOverwrite);
2438 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2439 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2440 reconstructor->SetRecoParam(par);
2442 return reconstructor;
2445 //_____________________________________________________________________________
2446 Bool_t AliReconstruction::CreateVertexer()
2448 // create the vertexer
2451 AliReconstructor* itsReconstructor = GetReconstructor(0);
2452 if (itsReconstructor) {
2453 fVertexer = itsReconstructor->CreateVertexer();
2456 AliWarning("couldn't create a vertexer for ITS");
2457 if (fStopOnError) return kFALSE;
2463 //_____________________________________________________________________________
2464 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2466 // create the trackers
2468 TString detStr = detectors;
2469 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2470 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2471 AliReconstructor* reconstructor = GetReconstructor(iDet);
2472 if (!reconstructor) continue;
2473 TString detName = fgkDetectorName[iDet];
2474 if (detName == "HLT") {
2475 fRunHLTTracking = kTRUE;
2478 if (detName == "MUON") {
2479 fRunMuonTracking = kTRUE;
2484 fTracker[iDet] = reconstructor->CreateTracker();
2485 if (!fTracker[iDet] && (iDet < 7)) {
2486 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2487 if (fStopOnError) return kFALSE;
2489 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2495 //_____________________________________________________________________________
2496 void AliReconstruction::CleanUp()
2498 // delete trackers and the run loader and close and delete the file
2500 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2501 delete fReconstructor[iDet];
2502 fReconstructor[iDet] = NULL;
2503 fLoader[iDet] = NULL;
2504 delete fTracker[iDet];
2505 fTracker[iDet] = NULL;
2516 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
2517 delete fDiamondProfile;
2518 fDiamondProfile = NULL;
2519 delete fDiamondProfileTPC;
2520 fDiamondProfileTPC = NULL;
2529 delete fParentRawReader;
2530 fParentRawReader=NULL;
2539 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2541 // Write space-points which are then used in the alignment procedures
2542 // For the moment only ITS, TPC, TRD and TOF
2544 Int_t ntracks = esd->GetNumberOfTracks();
2545 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2547 AliESDtrack *track = esd->GetTrack(itrack);
2550 for (Int_t iDet = 3; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2551 nsp += track->GetNcls(iDet);
2553 if (iDet==0) { // ITS "extra" clusters
2554 track->GetClusters(iDet,idx);
2555 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2560 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2561 track->SetTrackPointArray(sp);
2563 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2564 AliTracker *tracker = fTracker[iDet];
2565 if (!tracker) continue;
2566 Int_t nspdet = track->GetClusters(iDet,idx);
2568 if (iDet==0) // ITS "extra" clusters
2569 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2571 if (nspdet <= 0) continue;
2575 while (isp2 < nspdet) {
2576 Bool_t isvalid=kTRUE;
2578 Int_t index=idx[isp++];
2579 if (index < 0) continue;
2581 TString dets = fgkDetectorName[iDet];
2582 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2583 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2584 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2585 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2586 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2588 isvalid = tracker->GetTrackPoint(index,p);
2591 if (!isvalid) continue;
2592 sp->AddPoint(isptrack,&p); isptrack++;
2599 //_____________________________________________________________________________
2600 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2602 // The method reads the raw-data error log
2603 // accumulated within the rawReader.
2604 // It extracts the raw-data errors related to
2605 // the current event and stores them into
2606 // a TClonesArray inside the esd object.
2608 if (!fRawReader) return;
2610 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2612 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2614 if (iEvent != log->GetEventNumber()) continue;
2616 esd->AddRawDataErrorLog(log);
2621 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString pName){
2622 // Dump a file content into a char in TNamed
2624 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2625 Int_t kBytes = (Int_t)in.tellg();
2626 printf("Size: %d \n",kBytes);
2629 char* memblock = new char [kBytes];
2630 in.seekg (0, ios::beg);
2631 in.read (memblock, kBytes);
2633 TString fData(memblock,kBytes);
2634 fn = new TNamed(pName,fData);
2635 printf("fData Size: %d \n",fData.Sizeof());
2636 printf("pName Size: %d \n",pName.Sizeof());
2637 printf("fn Size: %d \n",fn->Sizeof());
2641 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2647 void AliReconstruction::TNamedToFile(TTree* fTree, TString pName){
2648 // This is not really needed in AliReconstruction at the moment
2649 // but can serve as a template
2651 TList *fList = fTree->GetUserInfo();
2652 TNamed *fn = (TNamed*)fList->FindObject(pName.Data());
2653 printf("fn Size: %d \n",fn->Sizeof());
2655 TString fTmp(fn->GetName()); // to be 100% sure in principle pName also works
2656 const char* cdata = fn->GetTitle();
2657 printf("fTmp Size %d\n",fTmp.Sizeof());
2659 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2660 printf("calculated size %d\n",size);
2661 ofstream out(pName.Data(),ios::out | ios::binary);
2662 out.write(cdata,size);
2667 //_____________________________________________________________________________
2668 void AliReconstruction::CheckQA()
2670 // check the QA of SIM for this run and remove the detectors
2671 // with status Fatal
2673 TString newRunLocalReconstruction ;
2674 TString newRunTracking ;
2675 TString newFillESD ;
2677 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2678 TString detName(AliQA::GetDetName(iDet)) ;
2679 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX_t(iDet)) ;
2680 if ( qa->IsSet(AliQA::DETECTORINDEX_t(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2681 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2683 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2684 fRunLocalReconstruction.Contains("ALL") ) {
2685 newRunLocalReconstruction += detName ;
2686 newRunLocalReconstruction += " " ;
2688 if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) ||
2689 fRunTracking.Contains("ALL") ) {
2690 newRunTracking += detName ;
2691 newRunTracking += " " ;
2693 if ( fFillESD.Contains(AliQA::GetDetName(iDet)) ||
2694 fFillESD.Contains("ALL") ) {
2695 newFillESD += detName ;
2700 fRunLocalReconstruction = newRunLocalReconstruction ;
2701 fRunTracking = newRunTracking ;
2702 fFillESD = newFillESD ;
2705 //_____________________________________________________________________________
2706 Int_t AliReconstruction::GetDetIndex(const char* detector)
2708 // return the detector index corresponding to detector
2710 for (index = 0; index < fgkNDetectors ; index++) {
2711 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2716 //_____________________________________________________________________________
2717 Bool_t AliReconstruction::FinishPlaneEff() {
2719 // Here execute all the necessary operationis, at the end of the tracking phase,
2720 // in case that evaluation of PlaneEfficiencies was required for some detector.
2721 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2723 // This Preliminary version works only FOR ITS !!!!!
2724 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2727 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2730 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2731 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2732 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2733 if(fTracker[iDet]) {
2734 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
2735 TString name=planeeff->GetName();
2737 TFile* pefile = TFile::Open(name, "RECREATE");
2738 ret=(Bool_t)planeeff->Write();
2740 if(planeeff->GetCreateHistos()) {
2741 TString hname=planeeff->GetName();
2742 hname+="Histo.root";
2743 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
2749 //_____________________________________________________________________________
2750 Bool_t AliReconstruction::InitPlaneEff() {
2752 // Here execute all the necessary operations, before of the tracking phase,
2753 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
2754 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
2755 // which should be updated/recalculated.
2757 // This Preliminary version will work only FOR ITS !!!!!
2758 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2761 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2763 AliWarning(Form("Implementation of this method not yet done !! Method return kTRUE"));
2767 //_____________________________________________________________________________
2768 Bool_t AliReconstruction::InitAliEVE()
2770 // This method should be called only in case
2771 // AliReconstruction is run
2772 // within the alieve environment.
2773 // It will initialize AliEVE in a way
2774 // so that it can visualize event processed
2775 // by AliReconstruction.
2776 // The return flag shows whenever the
2777 // AliEVE initialization was successful or not.
2780 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
2781 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
2782 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
2784 gROOT->ProcessLine("if (!gAliEveEvent) {gAliEveEvent = new AliEveEventManager();gAliEveEvent->SetAutoLoad(kTRUE);gAliEveEvent->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(gAliEveEvent);};");
2785 gROOT->ProcessLine("alieve_online_init()");
2790 //_____________________________________________________________________________
2791 void AliReconstruction::RunAliEVE()
2793 // Runs AliEVE visualisation of
2794 // the current event.
2795 // Should be executed only after
2796 // successful initialization of AliEVE.
2798 AliInfo("Running AliEVE...");
2799 gROOT->ProcessLine(Form("gAliEveEvent->SetEvent((AliRunLoader*)%p,(AliRawReader*)%p,(AliESDEvent*)%p);",fRunLoader,fRawReader,fesd));
2800 gROOT->ProcessLine("gAliEveEvent->StartStopAutoLoadTimer();");
2804 //_____________________________________________________________________________
2805 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
2807 // Allows to run QA for a selected set of detectors
2808 // and a selected set of tasks among RAWS, RECPOINTS and ESDS
2809 // all selected detectors run the same selected tasks
2811 if (!detAndAction.Contains(":")) {
2812 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2816 Int_t colon = detAndAction.Index(":") ;
2817 fQADetectors = detAndAction(0, colon) ;
2818 if (fQADetectors.Contains("ALL") )
2819 fQADetectors = fFillESD ;
2820 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2821 if (fQATasks.Contains("ALL") ) {
2822 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
2824 fQATasks.ToUpper() ;
2826 if ( fQATasks.Contains("RAW") )
2827 tempo = Form("%d ", AliQA::kRAWS) ;
2828 if ( fQATasks.Contains("RECPOINT") )
2829 tempo += Form("%d ", AliQA::kRECPOINTS) ;
2830 if ( fQATasks.Contains("ESD") )
2831 tempo += Form("%d ", AliQA::kESDS) ;
2833 if (fQATasks.IsNull()) {
2834 AliInfo("No QA requested\n") ;
2839 TString tempo(fQATasks) ;
2840 tempo.ReplaceAll(Form("%d", AliQA::kRAWS), AliQA::GetTaskName(AliQA::kRAWS)) ;
2841 tempo.ReplaceAll(Form("%d", AliQA::kRECPOINTS), AliQA::GetTaskName(AliQA::kRECPOINTS)) ;
2842 tempo.ReplaceAll(Form("%d", AliQA::kESDS), AliQA::GetTaskName(AliQA::kESDS)) ;
2843 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2848 //_____________________________________________________________________________
2849 Bool_t AliReconstruction::InitRecoParams()
2851 // The method accesses OCDB and retrieves all
2852 // the available reco-param objects from there.
2854 Bool_t isOK = kTRUE;
2856 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2858 if (fRecoParam.GetDetRecoParamArray(iDet)) {
2859 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
2863 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
2865 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
2866 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
2868 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
2872 TObject *recoParamObj = entry->GetObject();
2873 if (dynamic_cast<TObjArray*>(recoParamObj)) {
2874 // The detector has a normal TobjArray of AliDetectorRecoParam objects
2875 // Registering them in AliRecoParam
2876 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
2878 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
2879 // The detector has only onse set of reco parameters
2880 // Registering it in AliRecoParam
2881 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
2882 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
2883 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
2886 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
2890 AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
2899 //_____________________________________________________________________________
2900 Bool_t AliReconstruction::GetEventInfo()
2902 // Fill the event info object
2904 AliCodeTimerAuto("")
2906 AliCentralTrigger *aCTP = NULL;
2908 fEventInfo.SetEventType(fRawReader->GetType());
2910 ULong64_t mask = fRawReader->GetClassMask();
2911 fEventInfo.SetTriggerMask(mask);
2912 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
2913 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
2915 aCTP = new AliCentralTrigger();
2916 TString configstr("");
2917 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
2918 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
2922 aCTP->SetClassMask(mask);
2923 aCTP->SetClusterMask(clmask);
2926 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
2928 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
2929 aCTP = fRunLoader->GetTrigger();
2930 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
2931 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
2934 AliWarning("No trigger can be loaded! The trigger information will not be used!");
2939 AliTriggerConfiguration *config = aCTP->GetConfiguration();
2941 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
2942 if (fRawReader) delete aCTP;
2946 UChar_t clustmask = 0;
2948 ULong64_t trmask = fEventInfo.GetTriggerMask();
2949 const TObjArray& classesArray = config->GetClasses();
2950 Int_t nclasses = classesArray.GetEntriesFast();
2951 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
2952 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
2954 Int_t trindex = (Int_t)TMath::Log2(trclass->GetMask());
2955 fesd->SetTriggerClass(trclass->GetName(),trindex);
2956 if (trmask & (1 << trindex)) {
2958 trclasses += trclass->GetName();
2960 clustmask |= trclass->GetCluster()->GetClusterMask();
2964 fEventInfo.SetTriggerClasses(trclasses);
2966 // Set the information in ESD
2967 fesd->SetTriggerMask(trmask);
2968 fesd->SetTriggerCluster(clustmask);
2970 if (!aCTP->CheckTriggeredDetectors()) {
2971 if (fRawReader) delete aCTP;
2975 if (fRawReader) delete aCTP;
2977 // We have to fill also the HLT decision here!!
2983 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
2985 // Match the detector list found in the rec.C or the default 'ALL'
2986 // to the list found in the GRP (stored there by the shuttle PP which
2987 // gets the information from ECS)
2988 static TString resultList;
2989 TString detList = detectorList;
2993 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
2994 if ((detectorMask >> iDet) & 0x1) {
2995 TString det = AliDAQ::OfflineModuleName(iDet);
2996 if ((detList.CompareTo("ALL") == 0) ||
2997 detList.BeginsWith("ALL ") ||
2998 detList.EndsWith(" ALL") ||
2999 detList.Contains(" ALL ") ||
3000 (detList.CompareTo(det) == 0) ||
3001 detList.BeginsWith(det) ||
3002 detList.EndsWith(det) ||
3003 detList.Contains( " "+det+" " )) {
3004 if (!resultList.EndsWith(det + " ")) {
3013 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3014 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3015 if ((detList.CompareTo("ALL") == 0) ||
3016 detList.BeginsWith("ALL ") ||
3017 detList.EndsWith(" ALL") ||
3018 detList.Contains(" ALL ") ||
3019 (detList.CompareTo(hltDet) == 0) ||
3020 detList.BeginsWith(hltDet) ||
3021 detList.EndsWith(hltDet) ||
3022 detList.Contains( " "+hltDet+" " )) {
3023 resultList += hltDet;
3027 return resultList.Data();
3031 //______________________________________________________________________________
3032 void AliReconstruction::Abort(const char *method, EAbort what)
3034 // Abort processing. If what = kAbortProcess, the Process() loop will be
3035 // aborted. If what = kAbortFile, the current file in a chain will be
3036 // aborted and the processing will continue with the next file, if there
3037 // is no next file then Process() will be aborted. Abort() can also be
3038 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3039 // the SlaveTerminate() and Terminate() are always called. The abort flag
3040 // can be checked in these methods using GetAbort().
3042 // The method is overwritten in AliReconstruction for better handling of
3043 // reco specific errors
3045 if (!fStopOnError) return;
3049 TString whyMess = method;
3050 whyMess += " failed! Aborting...";
3052 AliError(whyMess.Data());
3055 TString mess = "Abort";
3056 if (fAbort == kAbortProcess)
3057 mess = "AbortProcess";
3058 else if (fAbort == kAbortFile)
3061 Info(mess, whyMess.Data());