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 <TProofOutputFile.h>
126 #include <TParameter.h>
128 #include "AliReconstruction.h"
129 #include "AliCodeTimer.h"
130 #include "AliReconstructor.h"
132 #include "AliRunLoader.h"
134 #include "AliRawReaderFile.h"
135 #include "AliRawReaderDate.h"
136 #include "AliRawReaderRoot.h"
137 #include "AliRawEventHeaderBase.h"
138 #include "AliRawEvent.h"
139 #include "AliESDEvent.h"
140 #include "AliESDMuonTrack.h"
141 #include "AliESDfriend.h"
142 #include "AliESDVertex.h"
143 #include "AliESDcascade.h"
144 #include "AliESDkink.h"
145 #include "AliESDtrack.h"
146 #include "AliESDCaloCluster.h"
147 #include "AliESDCaloCells.h"
148 #include "AliMultiplicity.h"
149 #include "AliTracker.h"
150 #include "AliVertexer.h"
151 #include "AliVertexerTracks.h"
152 #include "AliV0vertexer.h"
153 #include "AliCascadeVertexer.h"
154 #include "AliHeader.h"
155 #include "AliGenEventHeader.h"
157 #include "AliESDpid.h"
158 #include "AliESDtrack.h"
159 #include "AliESDPmdTrack.h"
161 #include "AliESDTagCreator.h"
163 #include "AliGeomManager.h"
164 #include "AliTrackPointArray.h"
165 #include "AliCDBManager.h"
166 #include "AliCDBStorage.h"
167 #include "AliCDBEntry.h"
168 #include "AliAlignObj.h"
170 #include "AliCentralTrigger.h"
171 #include "AliTriggerConfiguration.h"
172 #include "AliTriggerClass.h"
173 #include "AliTriggerCluster.h"
174 #include "AliCTPRawStream.h"
176 #include "AliQADataMakerRec.h"
177 #include "AliGlobalQADataMaker.h"
179 #include "AliQADataMakerSteer.h"
181 #include "AliPlaneEff.h"
183 #include "AliSysInfo.h" // memory snapshots
184 #include "AliRawHLTManager.h"
186 #include "AliMagWrapCheb.h"
188 #include "AliDetectorRecoParam.h"
189 #include "AliRunInfo.h"
190 #include "AliEventInfo.h"
194 ClassImp(AliReconstruction)
196 //_____________________________________________________________________________
197 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
199 //_____________________________________________________________________________
200 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
202 fUniformField(kFALSE),
203 fForcedFieldMap(NULL),
204 fRunVertexFinder(kTRUE),
205 fRunVertexFinderTracks(kTRUE),
206 fRunHLTTracking(kFALSE),
207 fRunMuonTracking(kFALSE),
209 fRunCascadeFinder(kTRUE),
210 fStopOnError(kFALSE),
211 fWriteAlignmentData(kFALSE),
212 fWriteESDfriend(kFALSE),
213 fFillTriggerESD(kTRUE),
221 fRunLocalReconstruction("ALL"),
225 fUseTrackingErrorsForAlignment(""),
226 fGAliceFileName(gAliceFilename),
231 fNumberOfEventsPerFile(1),
233 fLoadAlignFromCDB(kTRUE),
234 fLoadAlignData("ALL"),
242 fParentRawReader(NULL),
247 fDiamondProfile(NULL),
248 fDiamondProfileTPC(NULL),
249 fMeanVertexConstraint(kTRUE),
253 fAlignObjArray(NULL),
256 fInitCDBCalled(kFALSE),
257 fSetRunNumberFromDataCalled(kFALSE),
264 fSameQACycle(kFALSE),
266 fRunPlaneEff(kFALSE),
275 fIsNewRunLoader(kFALSE),
279 // create reconstruction object with default parameters
282 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
283 fReconstructor[iDet] = NULL;
284 fLoader[iDet] = NULL;
285 fTracker[iDet] = NULL;
286 fQACycles[iDet] = 999999;
288 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
292 //_____________________________________________________________________________
293 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
295 fUniformField(rec.fUniformField),
296 fForcedFieldMap(NULL),
297 fRunVertexFinder(rec.fRunVertexFinder),
298 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
299 fRunHLTTracking(rec.fRunHLTTracking),
300 fRunMuonTracking(rec.fRunMuonTracking),
301 fRunV0Finder(rec.fRunV0Finder),
302 fRunCascadeFinder(rec.fRunCascadeFinder),
303 fStopOnError(rec.fStopOnError),
304 fWriteAlignmentData(rec.fWriteAlignmentData),
305 fWriteESDfriend(rec.fWriteESDfriend),
306 fFillTriggerESD(rec.fFillTriggerESD),
308 fCleanESD(rec.fCleanESD),
309 fV0DCAmax(rec.fV0DCAmax),
310 fV0CsPmin(rec.fV0CsPmin),
314 fRunLocalReconstruction(rec.fRunLocalReconstruction),
315 fRunTracking(rec.fRunTracking),
316 fFillESD(rec.fFillESD),
317 fLoadCDB(rec.fLoadCDB),
318 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
319 fGAliceFileName(rec.fGAliceFileName),
320 fRawInput(rec.fRawInput),
321 fEquipIdMap(rec.fEquipIdMap),
322 fFirstEvent(rec.fFirstEvent),
323 fLastEvent(rec.fLastEvent),
324 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
326 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
327 fLoadAlignData(rec.fLoadAlignData),
328 fESDPar(rec.fESDPar),
329 fUseHLTData(rec.fUseHLTData),
335 fParentRawReader(NULL),
337 fRecoParam(rec.fRecoParam),
340 fDiamondProfile(rec.fDiamondProfile),
341 fDiamondProfileTPC(rec.fDiamondProfileTPC),
342 fMeanVertexConstraint(rec.fMeanVertexConstraint),
346 fAlignObjArray(rec.fAlignObjArray),
347 fCDBUri(rec.fCDBUri),
349 fInitCDBCalled(rec.fInitCDBCalled),
350 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
351 fQADetectors(rec.fQADetectors),
353 fQATasks(rec.fQATasks),
355 fRunGlobalQA(rec.fRunGlobalQA),
356 fInLoopQA(rec.fInLoopQA),
357 fSameQACycle(rec.fSameQACycle),
358 fRunPlaneEff(rec.fRunPlaneEff),
367 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());
388 //_____________________________________________________________________________
389 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
391 // assignment operator
392 // Used in PROOF mode
393 // Be very careful while modifing it!
394 // Simple rules to follow:
395 // for persistent data members - use their assignment operators
396 // for non-persistent ones - do nothing or take the default values from constructor
397 // TSelector members should not be touched
398 if(&rec == this) return *this;
400 fUniformField = rec.fUniformField;
401 fForcedFieldMap = NULL;
402 fRunVertexFinder = rec.fRunVertexFinder;
403 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
404 fRunHLTTracking = rec.fRunHLTTracking;
405 fRunMuonTracking = rec.fRunMuonTracking;
406 fRunV0Finder = rec.fRunV0Finder;
407 fRunCascadeFinder = rec.fRunCascadeFinder;
408 fStopOnError = rec.fStopOnError;
409 fWriteAlignmentData = rec.fWriteAlignmentData;
410 fWriteESDfriend = rec.fWriteESDfriend;
411 fFillTriggerESD = rec.fFillTriggerESD;
413 fCleanESD = rec.fCleanESD;
414 fV0DCAmax = rec.fV0DCAmax;
415 fV0CsPmin = rec.fV0CsPmin;
419 fRunLocalReconstruction = rec.fRunLocalReconstruction;
420 fRunTracking = rec.fRunTracking;
421 fFillESD = rec.fFillESD;
422 fLoadCDB = rec.fLoadCDB;
423 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
424 fGAliceFileName = rec.fGAliceFileName;
425 fRawInput = rec.fRawInput;
426 fEquipIdMap = rec.fEquipIdMap;
427 fFirstEvent = rec.fFirstEvent;
428 fLastEvent = rec.fLastEvent;
429 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
431 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
432 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
435 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
436 fLoadAlignData = rec.fLoadAlignData;
437 fESDPar = rec.fESDPar;
438 fUseHLTData = rec.fUseHLTData;
440 delete fRunInfo; fRunInfo = NULL;
441 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
443 fEventInfo = rec.fEventInfo;
447 fParentRawReader = NULL;
449 fRecoParam = rec.fRecoParam;
451 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
452 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
453 delete fLoader[iDet]; fLoader[iDet] = NULL;
454 delete fTracker[iDet]; fTracker[iDet] = NULL;
455 fQACycles[iDet] = rec.fQACycles[iDet];
459 delete fDiamondProfile; fDiamondProfile = NULL;
460 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
461 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
462 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
463 fMeanVertexConstraint = rec.fMeanVertexConstraint;
465 delete fGRPData; fGRPData = NULL;
466 if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
468 delete fAlignObjArray; fAlignObjArray = NULL;
471 fSpecCDBUri.Delete();
472 fInitCDBCalled = rec.fInitCDBCalled;
473 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
474 fQADetectors = rec.fQADetectors;
476 fQATasks = rec.fQATasks;
478 fRunGlobalQA = rec.fRunGlobalQA;
479 fInLoopQA = rec.fInLoopQA;
480 fSameQACycle = rec.fSameQACycle;
481 fRunPlaneEff = rec.fRunPlaneEff;
490 fIsNewRunLoader = rec.fIsNewRunLoader;
497 //_____________________________________________________________________________
498 AliReconstruction::~AliReconstruction()
503 delete fForcedFieldMap;
505 if (fAlignObjArray) {
506 fAlignObjArray->Delete();
507 delete fAlignObjArray;
509 fSpecCDBUri.Delete();
511 AliCodeTimer::Instance()->Print();
514 //_____________________________________________________________________________
515 void AliReconstruction::InitCDB()
517 // activate a default CDB storage
518 // First check if we have any CDB storage set, because it is used
519 // to retrieve the calibration and alignment constants
520 AliCodeTimerAuto("");
522 if (fInitCDBCalled) return;
523 fInitCDBCalled = kTRUE;
525 AliCDBManager* man = AliCDBManager::Instance();
526 if (man->IsDefaultStorageSet())
528 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
529 AliWarning("Default CDB storage has been already set !");
530 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
531 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
532 fCDBUri = man->GetDefaultStorage()->GetURI();
535 if (fCDBUri.Length() > 0)
537 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
538 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
539 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
541 fCDBUri="local://$ALICE_ROOT";
542 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
543 AliWarning("Default CDB storage not yet set !!!!");
544 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
545 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
548 man->SetDefaultStorage(fCDBUri);
551 // Now activate the detector specific CDB storage locations
552 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
553 TObject* obj = fSpecCDBUri[i];
555 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
556 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
557 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
558 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
560 AliSysInfo::AddStamp("InitCDB");
563 //_____________________________________________________________________________
564 void AliReconstruction::SetDefaultStorage(const char* uri) {
565 // Store the desired default CDB storage location
566 // Activate it later within the Run() method
572 //_____________________________________________________________________________
573 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
574 // Store a detector-specific CDB storage location
575 // Activate it later within the Run() method
577 AliCDBPath aPath(calibType);
578 if(!aPath.IsValid()){
579 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
580 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
581 if(!strcmp(calibType, fgkDetectorName[iDet])) {
582 aPath.SetPath(Form("%s/*", calibType));
583 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
587 if(!aPath.IsValid()){
588 AliError(Form("Not a valid path or detector: %s", calibType));
593 // // check that calibType refers to a "valid" detector name
594 // Bool_t isDetector = kFALSE;
595 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
596 // TString detName = fgkDetectorName[iDet];
597 // if(aPath.GetLevel0() == detName) {
598 // isDetector = kTRUE;
604 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
608 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
609 if (obj) fSpecCDBUri.Remove(obj);
610 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
614 //_____________________________________________________________________________
615 Bool_t AliReconstruction::SetRunNumberFromData()
617 // The method is called in Run() in order
618 // to set a correct run number.
619 // In case of raw data reconstruction the
620 // run number is taken from the raw data header
622 if (fSetRunNumberFromDataCalled) return kTRUE;
623 fSetRunNumberFromDataCalled = kTRUE;
625 AliCDBManager* man = AliCDBManager::Instance();
628 if(fRawReader->NextEvent()) {
629 if(man->GetRun() > 0) {
630 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
632 man->SetRun(fRawReader->GetRunNumber());
633 fRawReader->RewindEvents();
636 if(man->GetRun() > 0) {
637 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
640 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
646 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
648 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
653 // read run number from gAlice
654 if(rl->GetHeader()) {
655 man->SetRun(rl->GetHeader()->GetRun());
660 AliError("Neither run-loader header nor RawReader objects are found !");
672 //_____________________________________________________________________________
673 void AliReconstruction::SetCDBLock() {
674 // Set CDB lock: from now on it is forbidden to reset the run number
675 // or the default storage or to activate any further storage!
677 AliCDBManager::Instance()->SetLock(1);
680 //_____________________________________________________________________________
681 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
683 // Read the alignment objects from CDB.
684 // Each detector is supposed to have the
685 // alignment objects in DET/Align/Data CDB path.
686 // All the detector objects are then collected,
687 // sorted by geometry level (starting from ALIC) and
688 // then applied to the TGeo geometry.
689 // Finally an overlaps check is performed.
691 // Load alignment data from CDB and fill fAlignObjArray
692 if(fLoadAlignFromCDB){
694 TString detStr = detectors;
695 TString loadAlObjsListOfDets = "";
697 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
698 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
699 loadAlObjsListOfDets += fgkDetectorName[iDet];
700 loadAlObjsListOfDets += " ";
701 } // end loop over detectors
702 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
703 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
704 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
706 // Check if the array with alignment objects was
707 // provided by the user. If yes, apply the objects
708 // to the present TGeo geometry
709 if (fAlignObjArray) {
710 if (gGeoManager && gGeoManager->IsClosed()) {
711 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
712 AliError("The misalignment of one or more volumes failed!"
713 "Compare the list of simulated detectors and the list of detector alignment data!");
718 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
724 if (fAlignObjArray) {
725 fAlignObjArray->Delete();
726 delete fAlignObjArray; fAlignObjArray=NULL;
732 //_____________________________________________________________________________
733 void AliReconstruction::SetGAliceFile(const char* fileName)
735 // set the name of the galice file
737 fGAliceFileName = fileName;
740 //_____________________________________________________________________________
741 void AliReconstruction::SetInput(const char* input)
743 // In case the input string starts with 'mem://', we run in an online mode
744 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
745 // file is assumed. One can give as an input:
746 // mem://: - events taken from DAQ monitoring libs online
748 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
749 if (input) fRawInput = input;
752 //_____________________________________________________________________________
753 void AliReconstruction::SetOption(const char* detector, const char* option)
755 // set options for the reconstruction of a detector
757 TObject* obj = fOptions.FindObject(detector);
758 if (obj) fOptions.Remove(obj);
759 fOptions.Add(new TNamed(detector, option));
762 //_____________________________________________________________________________
763 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
765 // Set custom reconstruction parameters for a given detector
766 // Single set of parameters for all the events
767 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
768 if(!strcmp(detector, fgkDetectorName[iDet])) {
770 fRecoParam.AddDetRecoParam(iDet,par);
777 //_____________________________________________________________________________
778 Bool_t AliReconstruction::SetFieldMap(Float_t l3Current, Float_t diCurrent, Float_t factor, const char *path) {
779 //------------------------------------------------
780 // The magnetic field map, defined externally...
781 // L3 current 30000 A -> 0.5 T
782 // L3 current 12000 A -> 0.2 T
783 // dipole current 6000 A
784 // The polarities must be the same
785 //------------------------------------------------
786 const Float_t l3NominalCurrent1=30000.; // (A)
787 const Float_t l3NominalCurrent2=12000.; // (A)
788 const Float_t diNominalCurrent =6000. ; // (A)
790 const Float_t tolerance=0.03; // relative current tolerance
791 const Float_t zero=77.; // "zero" current (A)
794 Bool_t dipoleON=kFALSE;
796 TString s=(factor < 0) ? "L3: -" : "L3: +";
798 l3Current = TMath::Abs(l3Current);
799 if (TMath::Abs(l3Current-l3NominalCurrent1)/l3NominalCurrent1 < tolerance) {
800 map=AliMagWrapCheb::k5kG;
803 if (TMath::Abs(l3Current-l3NominalCurrent2)/l3NominalCurrent2 < tolerance) {
804 map=AliMagWrapCheb::k2kG;
807 if (l3Current < zero) {
808 map=AliMagWrapCheb::k2kG;
810 factor=0.; // in fact, this is a global factor...
811 fUniformField=kTRUE; // track with the uniform (zero) B field
813 AliError(Form("Wrong L3 current (%f A)!",l3Current));
817 diCurrent = TMath::Abs(diCurrent);
818 if (TMath::Abs(diCurrent-diNominalCurrent)/diNominalCurrent < tolerance) {
819 // 3% current tolerance...
823 if (diCurrent < zero) { // some small current..
827 AliError(Form("Wrong dipole current (%f A)!",diCurrent));
831 delete fForcedFieldMap;
833 new AliMagWrapCheb("B field map ",s,2,factor,10.,map,dipoleON,path);
835 fForcedFieldMap->Print();
837 AliTracker::SetFieldMap(fForcedFieldMap,fUniformField);
843 Bool_t AliReconstruction::InitGRP() {
844 //------------------------------------
845 // Initialization of the GRP entry
846 //------------------------------------
847 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
850 fGRPData = dynamic_cast<TMap*>(entry->GetObject());
852 AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
856 AliError("No GRP entry found in OCDB!");
860 TObjString *lhcState=
861 dynamic_cast<TObjString*>(fGRPData->GetValue("fLHCState"));
863 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
866 TObjString *beamType=
867 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamType"));
869 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
872 TObjString *beamEnergyStr=
873 dynamic_cast<TObjString*>(fGRPData->GetValue("fAliceBeamEnergy"));
874 if (!beamEnergyStr) {
875 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
879 dynamic_cast<TObjString*>(fGRPData->GetValue("fRunType"));
881 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
884 TObjString *activeDetectors=
885 dynamic_cast<TObjString*>(fGRPData->GetValue("fDetectorMask"));
886 if (!activeDetectors) {
887 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
890 fRunInfo = new AliRunInfo(lhcState ? lhcState->GetString().Data() : "UNKNOWN",
891 beamType ? beamType->GetString().Data() : "UNKNOWN",
892 beamEnergyStr ? beamEnergyStr->GetString().Atof() : 0,
893 runType ? runType->GetString().Data() : "UNKNOWN",
894 activeDetectors ? activeDetectors->GetString().Atoi() : 1074790399);
896 // Process the list of active detectors
897 if (activeDetectors && activeDetectors->GetString().IsDigit()) {
898 UInt_t detMask = activeDetectors->GetString().Atoi();
899 fLoadCDB.Form("%s %s %s %s",
900 fRunLocalReconstruction.Data(),
903 fQADetectors.Data());
904 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
905 fRunTracking = MatchDetectorList(fRunTracking,detMask);
906 fFillESD = MatchDetectorList(fFillESD,detMask);
907 fQADetectors = MatchDetectorList(fQADetectors,detMask);
908 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
911 AliInfo("===================================================================================");
912 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
913 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
914 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
915 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
916 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
917 AliInfo("===================================================================================");
919 //*** Dealing with the magnetic field map
920 if (AliTracker::GetFieldMap()) {
921 AliInfo("Running with the externally set B field !");
923 // Construct the field map out of the information retrieved from GRP.
928 TObjString *l3Current=
929 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Current"));
931 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
934 TObjString *l3Polarity=
935 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Polarity"));
937 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
942 TObjString *diCurrent=
943 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipoleCurrent"));
945 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
948 TObjString *diPolarity=
949 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipolePolarity"));
951 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
956 Float_t l3Cur=TMath::Abs(atof(l3Current->GetName()));
957 Float_t diCur=TMath::Abs(atof(diCurrent->GetName()));
958 Float_t l3Pol=atof(l3Polarity->GetName());
960 if (l3Pol != 0.) factor=-1.;
963 if (!SetFieldMap(l3Cur, diCur, factor)) {
964 AliFatal("Failed to creat a B field map ! Exiting...");
966 AliInfo("Running with the B field constructed out of GRP !");
969 AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
974 //*** Get the diamond profile from OCDB
975 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
977 if (fMeanVertexConstraint)
978 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
980 AliError("No diamond profile found in OCDB!");
983 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
985 if (fMeanVertexConstraint)
986 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
988 AliError("No diamond profile found in OCDB!");
994 //_____________________________________________________________________________
995 Bool_t AliReconstruction::LoadCDB()
997 AliCodeTimerAuto("");
999 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1001 TString detStr = fLoadCDB;
1002 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1003 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1004 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1009 //_____________________________________________________________________________
1010 Bool_t AliReconstruction::Run(const char* input)
1013 AliCodeTimerAuto("");
1016 if (GetAbort() != TSelector::kContinue) return kFALSE;
1018 TChain *chain = NULL;
1019 if (fRawReader && (chain = fRawReader->GetChain())) {
1022 gProof->AddInput(this);
1024 outputFile.SetProtocol("root",kTRUE);
1025 outputFile.SetHost(gSystem->HostName());
1026 outputFile.SetFile(Form("%s/AliESDs.root",gSystem->pwd()));
1027 AliInfo(Form("Output file with ESDs is %s",outputFile.GetUrl()));
1028 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",outputFile.GetUrl()));
1030 chain->Process("AliReconstruction");
1033 chain->Process(this);
1038 if (GetAbort() != TSelector::kContinue) return kFALSE;
1040 if (GetAbort() != TSelector::kContinue) return kFALSE;
1041 //******* The loop over events
1043 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1044 (fRawReader && fRawReader->NextEvent())) {
1045 if (!ProcessEvent(iEvent)) {
1046 Abort("ProcessEvent",TSelector::kAbortFile);
1052 if (GetAbort() != TSelector::kContinue) return kFALSE;
1054 if (GetAbort() != TSelector::kContinue) return kFALSE;
1060 //_____________________________________________________________________________
1061 void AliReconstruction::InitRawReader(const char* input)
1063 AliCodeTimerAuto("");
1065 // Init raw-reader and
1066 // set the input in case of raw data
1067 if (input) fRawInput = input;
1068 fRawReader = AliRawReader::Create(fRawInput.Data());
1070 AliInfo("Reconstruction will run over digits");
1072 if (!fEquipIdMap.IsNull() && fRawReader)
1073 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1075 if (!fUseHLTData.IsNull()) {
1076 // create the RawReaderHLT which performs redirection of HLT input data for
1077 // the specified detectors
1078 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1080 fParentRawReader=fRawReader;
1081 fRawReader=pRawReader;
1083 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1086 AliSysInfo::AddStamp("CreateRawReader");
1089 //_____________________________________________________________________________
1090 void AliReconstruction::InitRun(const char* input)
1092 // Initialization of raw-reader,
1093 // run number, CDB etc.
1094 AliCodeTimerAuto("");
1095 AliSysInfo::AddStamp("Start");
1097 // Initialize raw-reader if any
1098 InitRawReader(input);
1100 // Initialize the CDB storage
1103 // Set run number in CDBManager (if it is not already set by the user)
1104 if (!SetRunNumberFromData()) {
1105 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1109 // Set CDB lock: from now on it is forbidden to reset the run number
1110 // or the default storage or to activate any further storage!
1115 //_____________________________________________________________________________
1116 void AliReconstruction::Begin(TTree *)
1118 // Initialize AlReconstruction before
1119 // going into the event loop
1120 // Should follow the TSelector convention
1121 // i.e. initialize only the object on the client side
1122 AliCodeTimerAuto("");
1124 AliReconstruction *reco = NULL;
1126 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1129 AliSysInfo::AddStamp("ReadInputInBegin");
1132 // Import ideal TGeo geometry and apply misalignment
1134 TString geom(gSystem->DirName(fGAliceFileName));
1135 geom += "/geometry.root";
1136 AliGeomManager::LoadGeometry(geom.Data());
1138 Abort("LoadGeometry", TSelector::kAbortProcess);
1141 AliSysInfo::AddStamp("LoadGeom");
1142 TString detsToCheck=fRunLocalReconstruction;
1143 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1144 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1147 AliSysInfo::AddStamp("CheckGeom");
1150 if (!MisalignGeometry(fLoadAlignData)) {
1151 Abort("MisalignGeometry", TSelector::kAbortProcess);
1154 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1155 AliSysInfo::AddStamp("MisalignGeom");
1158 Abort("InitGRP", TSelector::kAbortProcess);
1161 AliSysInfo::AddStamp("InitGRP");
1164 Abort("LoadCDB", TSelector::kAbortProcess);
1167 AliSysInfo::AddStamp("LoadCDB");
1169 // Read the reconstruction parameters from OCDB
1170 if (!InitRecoParams()) {
1171 AliWarning("Not all detectors have correct RecoParam objects initialized");
1173 AliSysInfo::AddStamp("InitRecoParams");
1176 if (reco) *reco = *this;
1177 fInput->Add(gGeoManager);
1179 fInput->Add(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()));
1180 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1181 AliMagF *magFieldMap = (AliMagF*)AliTracker::GetFieldMap();
1182 magFieldMap->SetName("MagneticFieldMap");
1183 fInput->Add(magFieldMap);
1188 //_____________________________________________________________________________
1189 void AliReconstruction::SlaveBegin(TTree*)
1191 // Initialization related to run-loader,
1192 // vertexer, trackers, recontructors
1193 // In proof mode it is executed on the slave
1194 AliCodeTimerAuto("");
1196 TProofOutputFile *outProofFile = NULL;
1198 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1201 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1203 AliGeomManager::SetGeometry(tgeo);
1205 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1206 Int_t runNumber = -1;
1207 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1208 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1209 man->SetCacheFlag(kTRUE);
1210 man->SetLock(kTRUE);
1214 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1215 AliTracker::SetFieldMap(map,fUniformField);
1217 if (TNamed *outputFileName = (TNamed *) fInput->FindObject("PROOF_OUTPUTFILE")) {
1218 outProofFile = new TProofOutputFile(gSystem->BaseName(TUrl(outputFileName->GetTitle()).GetFile()));
1219 outProofFile->SetOutputFileName(outputFileName->GetTitle());
1220 fOutput->Add(outProofFile);
1222 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1225 // get the run loader
1226 if (!InitRunLoader()) {
1227 Abort("InitRunLoader", TSelector::kAbortProcess);
1230 AliSysInfo::AddStamp("LoadLoader");
1232 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1233 if(fDiamondProfile && fMeanVertexConstraint) ftVertexer->SetVtxStart(fDiamondProfile);
1236 if (fRunVertexFinder && !CreateVertexer()) {
1237 Abort("CreateVertexer", TSelector::kAbortProcess);
1240 AliSysInfo::AddStamp("CreateVertexer");
1243 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1244 Abort("CreateTrackers", TSelector::kAbortProcess);
1247 AliSysInfo::AddStamp("CreateTrackers");
1249 // create the ESD output file and tree
1250 if (!outProofFile) {
1251 ffile = TFile::Open("AliESDs.root", "RECREATE");
1252 ffile->SetCompressionLevel(2);
1253 if (!ffile->IsOpen()) {
1254 Abort("OpenESDFile", TSelector::kAbortProcess);
1259 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1260 Abort(Form("Problems opening output PROOF file: %s/%s",
1261 outProofFile->GetDir(), outProofFile->GetFileName()),
1262 TSelector::kAbortProcess);
1267 ftree = new TTree("esdTree", "Tree with ESD objects");
1268 fesd = new AliESDEvent();
1269 fesd->CreateStdContent();
1270 fesd->WriteToTree(ftree);
1272 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1273 fhltesd = new AliESDEvent();
1274 fhltesd->CreateStdContent();
1275 fhltesd->WriteToTree(fhlttree);
1278 if (fWriteESDfriend) {
1279 fesdf = new AliESDfriend();
1280 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1281 br->SetFile("AliESDfriends.root");
1282 fesd->AddObject(fesdf);
1285 ProcInfo_t ProcInfo;
1286 gSystem->GetProcInfo(&ProcInfo);
1287 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1290 fQASteer = new AliQADataMakerSteer("rec") ;
1291 fQASteer->SetActiveDetectors(fQADetectors) ;
1292 fQASteer->SetTasks(fQATasks) ;
1295 if (fRunQA && fRawReader && fQATasks.Contains(Form("%d", AliQA::kRAWS))) {
1296 fQASteer->Run(fQADetectors, fRawReader) ;
1297 fSameQACycle = kTRUE ;
1301 //Initialize the QA and start of cycle for out-of-loop QA
1303 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam, fSameQACycle, !fInLoopQA) ;
1307 fSameQACycle = kFALSE;
1308 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1309 AliInfo(Form("Initializing the global QA data maker"));
1310 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1311 TObjArray *arr=qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun());
1312 AliTracker::SetResidualsArray(arr);
1314 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
1315 fSameQACycle = kTRUE;
1318 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1319 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
1321 qadm->StartOfCycle(AliQA::kESDS, fSameQACycle);
1322 fSameQACycle = kTRUE;
1327 //Initialize the Plane Efficiency framework
1328 if (fRunPlaneEff && !InitPlaneEff()) {
1329 Abort("InitPlaneEff", TSelector::kAbortProcess);
1333 if (strcmp(gProgName,"alieve") == 0)
1334 fRunAliEVE = InitAliEVE();
1339 //_____________________________________________________________________________
1340 Bool_t AliReconstruction::Process(Long64_t entry)
1342 // run the reconstruction over a single entry
1343 // from the chain with raw data
1344 AliCodeTimerAuto("");
1346 TTree *currTree = fChain->GetTree();
1347 AliRawEvent *event = new AliRawEvent;
1348 currTree->SetBranchAddress("rawevent",&event);
1349 currTree->GetEntry(entry);
1350 fRawReader = new AliRawReaderRoot(event);
1351 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1359 //_____________________________________________________________________________
1360 void AliReconstruction::Init(TTree *tree)
1363 AliError("The input tree is not found!");
1369 //_____________________________________________________________________________
1370 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1372 // run the reconstruction over a single event
1373 // The event loop is steered in Run method
1375 AliCodeTimerAuto("");
1377 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1378 fRunLoader->SetEventNumber(iEvent);
1379 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1381 fRunLoader->TreeE()->Fill();
1382 if (fRawReader && fRawReader->UseAutoSaveESD())
1383 fRunLoader->TreeE()->AutoSave("SaveSelf");
1386 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1390 AliInfo(Form("processing event %d", iEvent));
1392 // Fill Event-info object
1394 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo);
1396 //Start of cycle for the in-loop QA
1397 if (fInLoopQA && fRunQA) {
1398 fQASteer->InitQADataMaker(AliCDBManager::Instance()->GetRun(), fRecoParam, fSameQACycle, fInLoopQA) ;
1400 if (fInLoopQA && fRunGlobalQA) {
1401 fSameQACycle = kFALSE;
1402 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1403 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS))) {
1404 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
1405 fSameQACycle = kTRUE;
1407 if (fQATasks.Contains(Form("%d", AliQA::kESDS))) {
1408 qadm->StartOfCycle(AliQA::kESDS, fSameQACycle);
1409 fSameQACycle = kTRUE;
1413 fRunLoader->GetEvent(iEvent);
1416 if (fInLoopQA && fRunQA)
1417 fQASteer->RunOneEvent(fRawReader) ;
1419 // local single event reconstruction
1420 if (!fRunLocalReconstruction.IsNull()) {
1421 TString detectors=fRunLocalReconstruction;
1422 // run HLT event reconstruction first
1423 // ;-( IsSelected changes the string
1424 if (IsSelected("HLT", detectors) &&
1425 !RunLocalEventReconstruction("HLT")) {
1426 if (fStopOnError) {CleanUp(); return kFALSE;}
1428 detectors=fRunLocalReconstruction;
1429 detectors.ReplaceAll("HLT", "");
1430 if (!RunLocalEventReconstruction(detectors)) {
1431 if (fStopOnError) {CleanUp(); return kFALSE;}
1435 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1436 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1437 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1438 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1440 // Set magnetic field from the tracker
1441 fesd->SetMagneticField(AliTracker::GetBz());
1442 fhltesd->SetMagneticField(AliTracker::GetBz());
1446 // Fill raw-data error log into the ESD
1447 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1450 if (fRunVertexFinder) {
1451 if (!RunVertexFinder(fesd)) {
1452 if (fStopOnError) {CleanUp(); return kFALSE;}
1457 if (!fRunTracking.IsNull()) {
1458 if (fRunMuonTracking) {
1459 if (!RunMuonTracking(fesd)) {
1460 if (fStopOnError) {CleanUp(); return kFALSE;}
1466 if (!fRunTracking.IsNull()) {
1467 if (!RunTracking(fesd)) {
1468 if (fStopOnError) {CleanUp(); return kFALSE;}
1473 if (!fFillESD.IsNull()) {
1474 TString detectors=fFillESD;
1475 // run HLT first and on hltesd
1476 // ;-( IsSelected changes the string
1477 if (IsSelected("HLT", detectors) &&
1478 !FillESD(fhltesd, "HLT")) {
1479 if (fStopOnError) {CleanUp(); return kFALSE;}
1482 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1483 if (detectors.Contains("ALL")) {
1485 for (Int_t idet=0; idet<fgkNDetectors; ++idet){
1486 detectors += fgkDetectorName[idet];
1490 detectors.ReplaceAll("HLT", "");
1491 if (!FillESD(fesd, detectors)) {
1492 if (fStopOnError) {CleanUp(); return kFALSE;}
1496 // fill Event header information from the RawEventHeader
1497 if (fRawReader){FillRawEventHeaderESD(fesd);}
1500 AliESDpid::MakePID(fesd);
1502 if (fFillTriggerESD) {
1503 if (!FillTriggerESD(fesd)) {
1504 if (fStopOnError) {CleanUp(); return kFALSE;}
1511 // Propagate track to the beam pipe (if not already done by ITS)
1513 const Int_t ntracks = fesd->GetNumberOfTracks();
1514 const Double_t kBz = fesd->GetMagneticField();
1515 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1518 UShort_t *selectedIdx=new UShort_t[ntracks];
1520 for (Int_t itrack=0; itrack<ntracks; itrack++){
1521 const Double_t kMaxStep = 5; //max step over the material
1524 AliESDtrack *track = fesd->GetTrack(itrack);
1525 if (!track) continue;
1527 AliExternalTrackParam *tpcTrack =
1528 (AliExternalTrackParam *)track->GetTPCInnerParam();
1532 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kTRUE);
1535 Int_t n=trkArray.GetEntriesFast();
1536 selectedIdx[n]=track->GetID();
1537 trkArray.AddLast(tpcTrack);
1540 //Tracks refitted by ITS should already be at the SPD vertex
1541 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1544 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1545 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1550 // Improve the reconstructed primary vertex position using the tracks
1552 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
1553 if(fesd->GetPrimaryVertexSPD()) {
1554 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
1555 if(vtitle.Contains("cosmics")) {
1556 runVertexFinderTracks=kFALSE;
1559 if (runVertexFinderTracks) {
1560 // TPC + ITS primary vertex
1561 ftVertexer->SetITSMode();
1562 if(fDiamondProfile && fMeanVertexConstraint) {
1563 ftVertexer->SetVtxStart(fDiamondProfile);
1565 ftVertexer->SetConstraintOff();
1567 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1569 if (pvtx->GetStatus()) {
1570 fesd->SetPrimaryVertex(pvtx);
1571 for (Int_t i=0; i<ntracks; i++) {
1572 AliESDtrack *t = fesd->GetTrack(i);
1573 t->RelateToVertex(pvtx, kBz, kVeryBig);
1578 // TPC-only primary vertex
1579 ftVertexer->SetTPCMode();
1580 if(fDiamondProfileTPC && fMeanVertexConstraint) {
1581 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1583 ftVertexer->SetConstraintOff();
1585 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1587 if (pvtx->GetStatus()) {
1588 fesd->SetPrimaryVertexTPC(pvtx);
1589 for (Int_t i=0; i<ntracks; i++) {
1590 AliESDtrack *t = fesd->GetTrack(i);
1591 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1597 delete[] selectedIdx;
1599 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1604 AliV0vertexer vtxer;
1605 vtxer.Tracks2V0vertices(fesd);
1607 if (fRunCascadeFinder) {
1609 AliCascadeVertexer cvtxer;
1610 cvtxer.V0sTracks2CascadeVertices(fesd);
1615 if (fCleanESD) CleanESD(fesd);
1618 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1619 if (qadm && fQATasks.Contains(Form("%d", AliQA::kESDS)))
1620 qadm->Exec(AliQA::kESDS, fesd);
1623 if (fWriteESDfriend) {
1624 fesdf->~AliESDfriend();
1625 new (fesdf) AliESDfriend(); // Reset...
1626 fesd->GetESDfriend(fesdf);
1630 // Auto-save the ESD tree in case of prompt reco @P2
1631 if (fRawReader && fRawReader->UseAutoSaveESD())
1632 ftree->AutoSave("SaveSelf");
1638 if (fRunAliEVE) RunAliEVE();
1642 if (fWriteESDfriend) {
1643 fesdf->~AliESDfriend();
1644 new (fesdf) AliESDfriend(); // Reset...
1647 ProcInfo_t ProcInfo;
1648 gSystem->GetProcInfo(&ProcInfo);
1649 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1652 // End of cycle for the in-loop
1653 if (fInLoopQA && fRunQA) {
1654 fQASteer->RunOneEvent(fesd) ;
1655 fQASteer->EndOfCycle() ;
1657 if (fInLoopQA && fRunGlobalQA) {
1658 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1660 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1661 qadm->EndOfCycle(AliQA::kRECPOINTS);
1662 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1663 qadm->EndOfCycle(AliQA::kESDS);
1669 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1670 if (fReconstructor[iDet])
1671 fReconstructor[iDet]->SetRecoParam(NULL);
1677 //_____________________________________________________________________________
1678 void AliReconstruction::SlaveTerminate()
1680 // Finalize the run on the slave side
1681 // Called after the exit
1682 // from the event loop
1683 AliCodeTimerAuto("");
1685 if (fIsNewRunLoader) { // galice.root didn't exist
1686 fRunLoader->WriteHeader("OVERWRITE");
1687 fRunLoader->CdGAFile();
1688 fRunLoader->Write(0, TObject::kOverwrite);
1691 ftree->GetUserInfo()->Add(fesd);
1692 fhlttree->GetUserInfo()->Add(fhltesd);
1694 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1695 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1697 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1698 cdbMapCopy->SetOwner(1);
1699 cdbMapCopy->SetName("cdbMap");
1700 TIter iter(cdbMap->GetTable());
1703 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1704 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1705 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1706 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1709 TList *cdbListCopy = new TList();
1710 cdbListCopy->SetOwner(1);
1711 cdbListCopy->SetName("cdbList");
1713 TIter iter2(cdbList);
1716 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1717 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1720 ftree->GetUserInfo()->Add(cdbMapCopy);
1721 ftree->GetUserInfo()->Add(cdbListCopy);
1724 if(fESDPar.Contains("ESD.par")){
1725 AliInfo("Attaching ESD.par to Tree");
1726 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
1727 ftree->GetUserInfo()->Add(fn);
1733 if (fWriteESDfriend)
1734 ftree->SetBranchStatus("ESDfriend*",0);
1735 // we want to have only one tree version number
1736 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1739 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1740 if (fRunPlaneEff && !FinishPlaneEff()) {
1741 AliWarning("Finish PlaneEff evaluation failed");
1744 //Finish QA and end of cycle for out-of-loop QA
1745 if (!fInLoopQA && fRunQA)
1746 fQASteer->Run(fRunLocalReconstruction.Data(), AliQA::kNULLTASKINDEX, fSameQACycle) ;
1747 if (!fInLoopQA && fRunGlobalQA) {
1748 AliQADataMaker *qadm = fQASteer->GetQADataMaker(AliQA::kGLOBAL);
1750 if (fQATasks.Contains(Form("%d", AliQA::kRECPOINTS)))
1751 qadm->EndOfCycle(AliQA::kRECPOINTS);
1752 if (fQATasks.Contains(Form("%d", AliQA::kESDS)))
1753 qadm->EndOfCycle(AliQA::kESDS);
1762 //_____________________________________________________________________________
1763 void AliReconstruction::Terminate()
1765 // Create tags for the events in the ESD tree (the ESD tree is always present)
1766 // In case of empty events the tags will contain dummy values
1767 AliCodeTimerAuto("");
1769 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1770 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData);
1772 // Cleanup of CDB manager: cache and active storages!
1773 AliCDBManager::Instance()->ClearCache();
1776 //_____________________________________________________________________________
1777 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1779 // run the local reconstruction
1781 static Int_t eventNr=0;
1782 AliCodeTimerAuto("")
1784 TString detStr = detectors;
1785 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1786 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1787 AliReconstructor* reconstructor = GetReconstructor(iDet);
1788 if (!reconstructor) continue;
1789 AliLoader* loader = fLoader[iDet];
1790 // Matthias April 2008: temporary fix to run HLT reconstruction
1791 // although the HLT loader is missing
1792 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1794 reconstructor->Reconstruct(fRawReader, NULL);
1797 reconstructor->Reconstruct(dummy, NULL);
1802 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1805 // conversion of digits
1806 if (fRawReader && reconstructor->HasDigitConversion()) {
1807 AliInfo(Form("converting raw data digits into root objects for %s",
1808 fgkDetectorName[iDet]));
1809 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1810 // fgkDetectorName[iDet]));
1811 loader->LoadDigits("update");
1812 loader->CleanDigits();
1813 loader->MakeDigitsContainer();
1814 TTree* digitsTree = loader->TreeD();
1815 reconstructor->ConvertDigits(fRawReader, digitsTree);
1816 loader->WriteDigits("OVERWRITE");
1817 loader->UnloadDigits();
1819 // local reconstruction
1820 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1821 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1822 loader->LoadRecPoints("update");
1823 loader->CleanRecPoints();
1824 loader->MakeRecPointsContainer();
1825 TTree* clustersTree = loader->TreeR();
1826 if (fRawReader && !reconstructor->HasDigitConversion()) {
1827 reconstructor->Reconstruct(fRawReader, clustersTree);
1829 loader->LoadDigits("read");
1830 TTree* digitsTree = loader->TreeD();
1832 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1833 if (fStopOnError) return kFALSE;
1835 reconstructor->Reconstruct(digitsTree, clustersTree);
1837 loader->UnloadDigits();
1840 // In-loop QA for local reconstrucion
1841 TString detQAStr(fQADetectors) ;
1842 if (fRunQA && fInLoopQA)
1843 fQASteer->RunOneEventInOneDetector(iDet, clustersTree) ;
1845 loader->WriteRecPoints("OVERWRITE");
1846 loader->UnloadRecPoints();
1847 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1849 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1850 AliError(Form("the following detectors were not found: %s",
1852 if (fStopOnError) return kFALSE;
1858 //_____________________________________________________________________________
1859 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1861 // run the barrel tracking
1863 AliCodeTimerAuto("")
1865 AliESDVertex* vertex = NULL;
1866 Double_t vtxPos[3] = {0, 0, 0};
1867 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1868 TArrayF mcVertex(3);
1869 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1870 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1871 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1875 AliInfo("running the ITS vertex finder");
1877 fLoader[0]->LoadRecPoints();
1878 TTree* cltree = fLoader[0]->TreeR();
1880 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1881 vertex = fVertexer->FindVertexForCurrentEvent(cltree);
1884 AliError("Can't get the ITS cluster tree");
1886 fLoader[0]->UnloadRecPoints();
1889 AliError("Can't get the ITS loader");
1892 AliWarning("Vertex not found");
1893 vertex = new AliESDVertex();
1894 vertex->SetName("default");
1897 vertex->SetName("reconstructed");
1901 AliInfo("getting the primary vertex from MC");
1902 vertex = new AliESDVertex(vtxPos, vtxErr);
1906 vertex->GetXYZ(vtxPos);
1907 vertex->GetSigmaXYZ(vtxErr);
1909 AliWarning("no vertex reconstructed");
1910 vertex = new AliESDVertex(vtxPos, vtxErr);
1912 esd->SetPrimaryVertexSPD(vertex);
1913 // if SPD multiplicity has been determined, it is stored in the ESD
1914 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1915 if(mult)esd->SetMultiplicity(mult);
1917 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1918 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1925 //_____________________________________________________________________________
1926 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1928 // run the HLT barrel tracking
1930 AliCodeTimerAuto("")
1933 AliError("Missing runLoader!");
1937 AliInfo("running HLT tracking");
1939 // Get a pointer to the HLT reconstructor
1940 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1941 if (!reconstructor) return kFALSE;
1944 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1945 TString detName = fgkDetectorName[iDet];
1946 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1947 reconstructor->SetOption(detName.Data());
1948 AliTracker *tracker = reconstructor->CreateTracker();
1950 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1951 if (fStopOnError) return kFALSE;
1955 Double_t vtxErr[3]={0.005,0.005,0.010};
1956 const AliESDVertex *vertex = esd->GetVertex();
1957 vertex->GetXYZ(vtxPos);
1958 tracker->SetVertex(vtxPos,vtxErr);
1960 fLoader[iDet]->LoadRecPoints("read");
1961 TTree* tree = fLoader[iDet]->TreeR();
1963 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1966 tracker->LoadClusters(tree);
1968 if (tracker->Clusters2Tracks(esd) != 0) {
1969 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1973 tracker->UnloadClusters();
1981 //_____________________________________________________________________________
1982 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1984 // run the muon spectrometer tracking
1986 AliCodeTimerAuto("")
1989 AliError("Missing runLoader!");
1992 Int_t iDet = 7; // for MUON
1994 AliInfo("is running...");
1996 // Get a pointer to the MUON reconstructor
1997 AliReconstructor *reconstructor = GetReconstructor(iDet);
1998 if (!reconstructor) return kFALSE;
2001 TString detName = fgkDetectorName[iDet];
2002 AliDebug(1, Form("%s tracking", detName.Data()));
2003 AliTracker *tracker = reconstructor->CreateTracker();
2005 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2010 fLoader[iDet]->LoadRecPoints("read");
2012 tracker->LoadClusters(fLoader[iDet]->TreeR());
2014 Int_t rv = tracker->Clusters2Tracks(esd);
2018 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2022 fLoader[iDet]->UnloadRecPoints();
2024 tracker->UnloadClusters();
2032 //_____________________________________________________________________________
2033 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2035 // run the barrel tracking
2036 static Int_t eventNr=0;
2037 AliCodeTimerAuto("")
2039 AliInfo("running tracking");
2041 //Fill the ESD with the T0 info (will be used by the TOF)
2042 if (fReconstructor[11] && fLoader[11]) {
2043 fLoader[11]->LoadRecPoints("READ");
2044 TTree *treeR = fLoader[11]->TreeR();
2045 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2048 // pass 1: TPC + ITS inwards
2049 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2050 if (!fTracker[iDet]) continue;
2051 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2054 fLoader[iDet]->LoadRecPoints("read");
2055 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2056 TTree* tree = fLoader[iDet]->TreeR();
2058 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2061 fTracker[iDet]->LoadClusters(tree);
2062 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2064 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2065 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2068 // preliminary PID in TPC needed by the ITS tracker
2070 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2071 AliESDpid::MakePID(esd);
2073 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2076 // pass 2: ALL backwards
2078 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2079 if (!fTracker[iDet]) continue;
2080 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2083 if (iDet > 1) { // all except ITS, TPC
2085 fLoader[iDet]->LoadRecPoints("read");
2086 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2087 tree = fLoader[iDet]->TreeR();
2089 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2092 fTracker[iDet]->LoadClusters(tree);
2093 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2097 if (iDet>1) // start filling residuals for the "outer" detectors
2098 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
2100 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2101 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2106 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2107 fTracker[iDet]->UnloadClusters();
2108 fLoader[iDet]->UnloadRecPoints();
2110 // updated PID in TPC needed by the ITS tracker -MI
2112 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2113 AliESDpid::MakePID(esd);
2115 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2117 //stop filling residuals for the "outer" detectors
2118 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
2120 // pass 3: TRD + TPC + ITS refit inwards
2122 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2123 if (!fTracker[iDet]) continue;
2124 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2127 if (iDet<2) // start filling residuals for TPC and ITS
2128 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
2130 if (fTracker[iDet]->RefitInward(esd) != 0) {
2131 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2134 // run postprocessing
2135 if (fTracker[iDet]->PostProcess(esd) != 0) {
2136 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2139 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2142 // write space-points to the ESD in case alignment data output
2144 if (fWriteAlignmentData)
2145 WriteAlignmentData(esd);
2147 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2148 if (!fTracker[iDet]) continue;
2150 fTracker[iDet]->UnloadClusters();
2151 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2152 fLoader[iDet]->UnloadRecPoints();
2153 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2155 // stop filling residuals for TPC and ITS
2156 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
2162 //_____________________________________________________________________________
2163 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2165 // Remove the data which are not needed for the physics analysis.
2168 Int_t nTracks=esd->GetNumberOfTracks();
2169 Int_t nV0s=esd->GetNumberOfV0s();
2171 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2173 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2174 Bool_t rc=esd->Clean(cleanPars);
2176 nTracks=esd->GetNumberOfTracks();
2177 nV0s=esd->GetNumberOfV0s();
2179 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2184 //_____________________________________________________________________________
2185 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2187 // fill the event summary data
2189 AliCodeTimerAuto("")
2190 static Int_t eventNr=0;
2191 TString detStr = detectors;
2193 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2194 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2195 AliReconstructor* reconstructor = GetReconstructor(iDet);
2196 if (!reconstructor) continue;
2197 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2198 TTree* clustersTree = NULL;
2199 if (fLoader[iDet]) {
2200 fLoader[iDet]->LoadRecPoints("read");
2201 clustersTree = fLoader[iDet]->TreeR();
2202 if (!clustersTree) {
2203 AliError(Form("Can't get the %s clusters tree",
2204 fgkDetectorName[iDet]));
2205 if (fStopOnError) return kFALSE;
2208 if (fRawReader && !reconstructor->HasDigitConversion()) {
2209 reconstructor->FillESD(fRawReader, clustersTree, esd);
2211 TTree* digitsTree = NULL;
2212 if (fLoader[iDet]) {
2213 fLoader[iDet]->LoadDigits("read");
2214 digitsTree = fLoader[iDet]->TreeD();
2216 AliError(Form("Can't get the %s digits tree",
2217 fgkDetectorName[iDet]));
2218 if (fStopOnError) return kFALSE;
2221 reconstructor->FillESD(digitsTree, clustersTree, esd);
2222 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2224 if (fLoader[iDet]) {
2225 fLoader[iDet]->UnloadRecPoints();
2229 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2230 AliError(Form("the following detectors were not found: %s",
2232 if (fStopOnError) return kFALSE;
2234 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2239 //_____________________________________________________________________________
2240 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2242 // Reads the trigger decision which is
2243 // stored in Trigger.root file and fills
2244 // the corresponding esd entries
2246 AliCodeTimerAuto("")
2248 AliInfo("Filling trigger information into the ESD");
2251 AliCTPRawStream input(fRawReader);
2252 if (!input.Next()) {
2253 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2256 if (esd->GetTriggerMask() != input.GetClassMask())
2257 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2258 input.GetClassMask(),esd->GetTriggerMask()));
2259 if (esd->GetOrbitNumber() != input.GetOrbitID())
2260 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2261 input.GetOrbitID(),esd->GetOrbitNumber()));
2262 if (esd->GetBunchCrossNumber() != input.GetBCID())
2263 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2264 input.GetBCID(),esd->GetBunchCrossNumber()));
2267 // Here one has to add the filling of trigger inputs and
2268 // interaction records
2278 //_____________________________________________________________________________
2279 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2282 // Filling information from RawReader Header
2285 if (!fRawReader) return kFALSE;
2287 AliInfo("Filling information from RawReader Header");
2289 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2290 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2291 esd->SetPeriodNumber(fRawReader->GetPeriod());
2293 esd->SetTimeStamp(fRawReader->GetTimestamp());
2294 esd->SetEventType(fRawReader->GetType());
2300 //_____________________________________________________________________________
2301 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2303 // check whether detName is contained in detectors
2304 // if yes, it is removed from detectors
2306 // check if all detectors are selected
2307 if ((detectors.CompareTo("ALL") == 0) ||
2308 detectors.BeginsWith("ALL ") ||
2309 detectors.EndsWith(" ALL") ||
2310 detectors.Contains(" ALL ")) {
2315 // search for the given detector
2316 Bool_t result = kFALSE;
2317 if ((detectors.CompareTo(detName) == 0) ||
2318 detectors.BeginsWith(detName+" ") ||
2319 detectors.EndsWith(" "+detName) ||
2320 detectors.Contains(" "+detName+" ")) {
2321 detectors.ReplaceAll(detName, "");
2325 // clean up the detectors string
2326 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2327 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2328 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2333 //_____________________________________________________________________________
2334 Bool_t AliReconstruction::InitRunLoader()
2336 // get or create the run loader
2338 if (gAlice) delete gAlice;
2341 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2342 // load all base libraries to get the loader classes
2343 TString libs = gSystem->GetLibraries();
2344 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2345 TString detName = fgkDetectorName[iDet];
2346 if (detName == "HLT") continue;
2347 if (libs.Contains("lib" + detName + "base.so")) continue;
2348 gSystem->Load("lib" + detName + "base.so");
2350 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2352 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2357 fRunLoader->CdGAFile();
2358 fRunLoader->LoadgAlice();
2360 //PH This is a temporary fix to give access to the kinematics
2361 //PH that is needed for the labels of ITS clusters
2362 fRunLoader->LoadHeader();
2363 fRunLoader->LoadKinematics();
2365 } else { // galice.root does not exist
2367 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2369 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2370 AliConfig::GetDefaultEventFolderName(),
2373 AliError(Form("could not create run loader in file %s",
2374 fGAliceFileName.Data()));
2378 fIsNewRunLoader = kTRUE;
2379 fRunLoader->MakeTree("E");
2381 if (fNumberOfEventsPerFile > 0)
2382 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2384 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2390 //_____________________________________________________________________________
2391 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2393 // get the reconstructor object and the loader for a detector
2395 if (fReconstructor[iDet]) {
2396 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2397 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2398 fReconstructor[iDet]->SetRecoParam(par);
2400 return fReconstructor[iDet];
2403 // load the reconstructor object
2404 TPluginManager* pluginManager = gROOT->GetPluginManager();
2405 TString detName = fgkDetectorName[iDet];
2406 TString recName = "Ali" + detName + "Reconstructor";
2408 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2410 AliReconstructor* reconstructor = NULL;
2411 // first check if a plugin is defined for the reconstructor
2412 TPluginHandler* pluginHandler =
2413 pluginManager->FindHandler("AliReconstructor", detName);
2414 // if not, add a plugin for it
2415 if (!pluginHandler) {
2416 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2417 TString libs = gSystem->GetLibraries();
2418 if (libs.Contains("lib" + detName + "base.so") ||
2419 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2420 pluginManager->AddHandler("AliReconstructor", detName,
2421 recName, detName + "rec", recName + "()");
2423 pluginManager->AddHandler("AliReconstructor", detName,
2424 recName, detName, recName + "()");
2426 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2428 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2429 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2431 if (reconstructor) {
2432 TObject* obj = fOptions.FindObject(detName.Data());
2433 if (obj) reconstructor->SetOption(obj->GetTitle());
2434 reconstructor->Init();
2435 fReconstructor[iDet] = reconstructor;
2438 // get or create the loader
2439 if (detName != "HLT") {
2440 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2441 if (!fLoader[iDet]) {
2442 AliConfig::Instance()
2443 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2445 // first check if a plugin is defined for the loader
2447 pluginManager->FindHandler("AliLoader", detName);
2448 // if not, add a plugin for it
2449 if (!pluginHandler) {
2450 TString loaderName = "Ali" + detName + "Loader";
2451 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2452 pluginManager->AddHandler("AliLoader", detName,
2453 loaderName, detName + "base",
2454 loaderName + "(const char*, TFolder*)");
2455 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2457 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2459 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2460 fRunLoader->GetEventFolder());
2462 if (!fLoader[iDet]) { // use default loader
2463 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2465 if (!fLoader[iDet]) {
2466 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2467 if (fStopOnError) return NULL;
2469 fRunLoader->AddLoader(fLoader[iDet]);
2470 fRunLoader->CdGAFile();
2471 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2472 fRunLoader->Write(0, TObject::kOverwrite);
2477 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2478 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2479 reconstructor->SetRecoParam(par);
2481 return reconstructor;
2484 //_____________________________________________________________________________
2485 Bool_t AliReconstruction::CreateVertexer()
2487 // create the vertexer
2490 AliReconstructor* itsReconstructor = GetReconstructor(0);
2491 if (itsReconstructor) {
2492 fVertexer = itsReconstructor->CreateVertexer();
2495 AliWarning("couldn't create a vertexer for ITS");
2496 if (fStopOnError) return kFALSE;
2502 //_____________________________________________________________________________
2503 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2505 // create the trackers
2507 TString detStr = detectors;
2508 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2509 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2510 AliReconstructor* reconstructor = GetReconstructor(iDet);
2511 if (!reconstructor) continue;
2512 TString detName = fgkDetectorName[iDet];
2513 if (detName == "HLT") {
2514 fRunHLTTracking = kTRUE;
2517 if (detName == "MUON") {
2518 fRunMuonTracking = kTRUE;
2523 fTracker[iDet] = reconstructor->CreateTracker();
2524 if (!fTracker[iDet] && (iDet < 7)) {
2525 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2526 if (fStopOnError) return kFALSE;
2528 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2534 //_____________________________________________________________________________
2535 void AliReconstruction::CleanUp()
2537 // delete trackers and the run loader and close and delete the file
2539 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2540 delete fReconstructor[iDet];
2541 fReconstructor[iDet] = NULL;
2542 fLoader[iDet] = NULL;
2543 delete fTracker[iDet];
2544 fTracker[iDet] = NULL;
2555 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
2556 delete fDiamondProfile;
2557 fDiamondProfile = NULL;
2558 delete fDiamondProfileTPC;
2559 fDiamondProfileTPC = NULL;
2568 delete fParentRawReader;
2569 fParentRawReader=NULL;
2578 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2580 // Write space-points which are then used in the alignment procedures
2581 // For the moment only ITS, TPC, TRD and TOF
2583 Int_t ntracks = esd->GetNumberOfTracks();
2584 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2586 AliESDtrack *track = esd->GetTrack(itrack);
2589 for (Int_t iDet = 3; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2590 nsp += track->GetNcls(iDet);
2592 if (iDet==0) { // ITS "extra" clusters
2593 track->GetClusters(iDet,idx);
2594 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2599 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2600 track->SetTrackPointArray(sp);
2602 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2603 AliTracker *tracker = fTracker[iDet];
2604 if (!tracker) continue;
2605 Int_t nspdet = track->GetClusters(iDet,idx);
2607 if (iDet==0) // ITS "extra" clusters
2608 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2610 if (nspdet <= 0) continue;
2614 while (isp2 < nspdet) {
2615 Bool_t isvalid=kTRUE;
2617 Int_t index=idx[isp++];
2618 if (index < 0) continue;
2620 TString dets = fgkDetectorName[iDet];
2621 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2622 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2623 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2624 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2625 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2627 isvalid = tracker->GetTrackPoint(index,p);
2630 if (!isvalid) continue;
2631 sp->AddPoint(isptrack,&p); isptrack++;
2638 //_____________________________________________________________________________
2639 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2641 // The method reads the raw-data error log
2642 // accumulated within the rawReader.
2643 // It extracts the raw-data errors related to
2644 // the current event and stores them into
2645 // a TClonesArray inside the esd object.
2647 if (!fRawReader) return;
2649 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2651 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2653 if (iEvent != log->GetEventNumber()) continue;
2655 esd->AddRawDataErrorLog(log);
2660 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString pName){
2661 // Dump a file content into a char in TNamed
2663 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2664 Int_t kBytes = (Int_t)in.tellg();
2665 printf("Size: %d \n",kBytes);
2668 char* memblock = new char [kBytes];
2669 in.seekg (0, ios::beg);
2670 in.read (memblock, kBytes);
2672 TString fData(memblock,kBytes);
2673 fn = new TNamed(pName,fData);
2674 printf("fData Size: %d \n",fData.Sizeof());
2675 printf("pName Size: %d \n",pName.Sizeof());
2676 printf("fn Size: %d \n",fn->Sizeof());
2680 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2686 void AliReconstruction::TNamedToFile(TTree* fTree, TString pName){
2687 // This is not really needed in AliReconstruction at the moment
2688 // but can serve as a template
2690 TList *fList = fTree->GetUserInfo();
2691 TNamed *fn = (TNamed*)fList->FindObject(pName.Data());
2692 printf("fn Size: %d \n",fn->Sizeof());
2694 TString fTmp(fn->GetName()); // to be 100% sure in principle pName also works
2695 const char* cdata = fn->GetTitle();
2696 printf("fTmp Size %d\n",fTmp.Sizeof());
2698 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2699 printf("calculated size %d\n",size);
2700 ofstream out(pName.Data(),ios::out | ios::binary);
2701 out.write(cdata,size);
2706 //_____________________________________________________________________________
2707 void AliReconstruction::CheckQA()
2709 // check the QA of SIM for this run and remove the detectors
2710 // with status Fatal
2712 TString newRunLocalReconstruction ;
2713 TString newRunTracking ;
2714 TString newFillESD ;
2716 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2717 TString detName(AliQA::GetDetName(iDet)) ;
2718 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX_t(iDet)) ;
2719 if ( qa->IsSet(AliQA::DETECTORINDEX_t(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2720 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2722 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2723 fRunLocalReconstruction.Contains("ALL") ) {
2724 newRunLocalReconstruction += detName ;
2725 newRunLocalReconstruction += " " ;
2727 if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) ||
2728 fRunTracking.Contains("ALL") ) {
2729 newRunTracking += detName ;
2730 newRunTracking += " " ;
2732 if ( fFillESD.Contains(AliQA::GetDetName(iDet)) ||
2733 fFillESD.Contains("ALL") ) {
2734 newFillESD += detName ;
2739 fRunLocalReconstruction = newRunLocalReconstruction ;
2740 fRunTracking = newRunTracking ;
2741 fFillESD = newFillESD ;
2744 //_____________________________________________________________________________
2745 Int_t AliReconstruction::GetDetIndex(const char* detector)
2747 // return the detector index corresponding to detector
2749 for (index = 0; index < fgkNDetectors ; index++) {
2750 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2755 //_____________________________________________________________________________
2756 Bool_t AliReconstruction::FinishPlaneEff() {
2758 // Here execute all the necessary operationis, at the end of the tracking phase,
2759 // in case that evaluation of PlaneEfficiencies was required for some detector.
2760 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2762 // This Preliminary version works only FOR ITS !!!!!
2763 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2766 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2769 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2770 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2771 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2772 if(fTracker[iDet]) {
2773 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
2774 TString name=planeeff->GetName();
2776 TFile* pefile = TFile::Open(name, "RECREATE");
2777 ret=(Bool_t)planeeff->Write();
2779 if(planeeff->GetCreateHistos()) {
2780 TString hname=planeeff->GetName();
2781 hname+="Histo.root";
2782 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
2788 //_____________________________________________________________________________
2789 Bool_t AliReconstruction::InitPlaneEff() {
2791 // Here execute all the necessary operations, before of the tracking phase,
2792 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
2793 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
2794 // which should be updated/recalculated.
2796 // This Preliminary version will work only FOR ITS !!!!!
2797 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2800 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2802 AliWarning(Form("Implementation of this method not yet done !! Method return kTRUE"));
2806 //_____________________________________________________________________________
2807 Bool_t AliReconstruction::InitAliEVE()
2809 // This method should be called only in case
2810 // AliReconstruction is run
2811 // within the alieve environment.
2812 // It will initialize AliEVE in a way
2813 // so that it can visualize event processed
2814 // by AliReconstruction.
2815 // The return flag shows whenever the
2816 // AliEVE initialization was successful or not.
2819 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
2820 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
2821 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
2823 gROOT->ProcessLine("if (!gAliEveEvent) {gAliEveEvent = new AliEveEventManager();gAliEveEvent->SetAutoLoad(kTRUE);gAliEveEvent->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(gAliEveEvent);};");
2824 gROOT->ProcessLine("alieve_online_init()");
2829 //_____________________________________________________________________________
2830 void AliReconstruction::RunAliEVE()
2832 // Runs AliEVE visualisation of
2833 // the current event.
2834 // Should be executed only after
2835 // successful initialization of AliEVE.
2837 AliInfo("Running AliEVE...");
2838 gROOT->ProcessLine(Form("gAliEveEvent->SetEvent((AliRunLoader*)%p,(AliRawReader*)%p,(AliESDEvent*)%p);",fRunLoader,fRawReader,fesd));
2839 gROOT->ProcessLine("gAliEveEvent->StartStopAutoLoadTimer();");
2843 //_____________________________________________________________________________
2844 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
2846 // Allows to run QA for a selected set of detectors
2847 // and a selected set of tasks among RAWS, RECPOINTS and ESDS
2848 // all selected detectors run the same selected tasks
2850 if (!detAndAction.Contains(":")) {
2851 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
2855 Int_t colon = detAndAction.Index(":") ;
2856 fQADetectors = detAndAction(0, colon) ;
2857 if (fQADetectors.Contains("ALL") )
2858 fQADetectors = fFillESD ;
2859 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
2860 if (fQATasks.Contains("ALL") ) {
2861 fQATasks = Form("%d %d %d", AliQA::kRAWS, AliQA::kRECPOINTS, AliQA::kESDS) ;
2863 fQATasks.ToUpper() ;
2865 if ( fQATasks.Contains("RAW") )
2866 tempo = Form("%d ", AliQA::kRAWS) ;
2867 if ( fQATasks.Contains("RECPOINT") )
2868 tempo += Form("%d ", AliQA::kRECPOINTS) ;
2869 if ( fQATasks.Contains("ESD") )
2870 tempo += Form("%d ", AliQA::kESDS) ;
2872 if (fQATasks.IsNull()) {
2873 AliInfo("No QA requested\n") ;
2878 TString tempo(fQATasks) ;
2879 tempo.ReplaceAll(Form("%d", AliQA::kRAWS), AliQA::GetTaskName(AliQA::kRAWS)) ;
2880 tempo.ReplaceAll(Form("%d", AliQA::kRECPOINTS), AliQA::GetTaskName(AliQA::kRECPOINTS)) ;
2881 tempo.ReplaceAll(Form("%d", AliQA::kESDS), AliQA::GetTaskName(AliQA::kESDS)) ;
2882 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
2887 //_____________________________________________________________________________
2888 Bool_t AliReconstruction::InitRecoParams()
2890 // The method accesses OCDB and retrieves all
2891 // the available reco-param objects from there.
2893 Bool_t isOK = kTRUE;
2895 TString detStr = fLoadCDB;
2896 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2898 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2900 if (fRecoParam.GetDetRecoParamArray(iDet)) {
2901 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
2905 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
2907 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
2908 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
2910 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
2914 TObject *recoParamObj = entry->GetObject();
2915 if (dynamic_cast<TObjArray*>(recoParamObj)) {
2916 // The detector has a normal TobjArray of AliDetectorRecoParam objects
2917 // Registering them in AliRecoParam
2918 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
2920 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
2921 // The detector has only onse set of reco parameters
2922 // Registering it in AliRecoParam
2923 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
2924 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
2925 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
2928 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
2932 AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
2936 if (AliDebugLevel() > 0) fRecoParam.Print();
2941 //_____________________________________________________________________________
2942 Bool_t AliReconstruction::GetEventInfo()
2944 // Fill the event info object
2946 AliCodeTimerAuto("")
2948 AliCentralTrigger *aCTP = NULL;
2950 fEventInfo.SetEventType(fRawReader->GetType());
2952 ULong64_t mask = fRawReader->GetClassMask();
2953 fEventInfo.SetTriggerMask(mask);
2954 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
2955 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
2957 aCTP = new AliCentralTrigger();
2958 TString configstr("");
2959 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
2960 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
2964 aCTP->SetClassMask(mask);
2965 aCTP->SetClusterMask(clmask);
2968 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
2970 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
2971 aCTP = fRunLoader->GetTrigger();
2972 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
2973 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
2976 AliWarning("No trigger can be loaded! The trigger information will not be used!");
2981 AliTriggerConfiguration *config = aCTP->GetConfiguration();
2983 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
2984 if (fRawReader) delete aCTP;
2988 UChar_t clustmask = 0;
2990 ULong64_t trmask = fEventInfo.GetTriggerMask();
2991 const TObjArray& classesArray = config->GetClasses();
2992 Int_t nclasses = classesArray.GetEntriesFast();
2993 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
2994 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
2996 Int_t trindex = (Int_t)TMath::Log2(trclass->GetMask());
2997 fesd->SetTriggerClass(trclass->GetName(),trindex);
2998 if (trmask & (1 << trindex)) {
3000 trclasses += trclass->GetName();
3002 clustmask |= trclass->GetCluster()->GetClusterMask();
3006 fEventInfo.SetTriggerClasses(trclasses);
3008 // Set the information in ESD
3009 fesd->SetTriggerMask(trmask);
3010 fesd->SetTriggerCluster(clustmask);
3012 if (!aCTP->CheckTriggeredDetectors()) {
3013 if (fRawReader) delete aCTP;
3017 if (fRawReader) delete aCTP;
3019 // We have to fill also the HLT decision here!!
3025 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3027 // Match the detector list found in the rec.C or the default 'ALL'
3028 // to the list found in the GRP (stored there by the shuttle PP which
3029 // gets the information from ECS)
3030 static TString resultList;
3031 TString detList = detectorList;
3035 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3036 if ((detectorMask >> iDet) & 0x1) {
3037 TString det = AliDAQ::OfflineModuleName(iDet);
3038 if ((detList.CompareTo("ALL") == 0) ||
3039 detList.BeginsWith("ALL ") ||
3040 detList.EndsWith(" ALL") ||
3041 detList.Contains(" ALL ") ||
3042 (detList.CompareTo(det) == 0) ||
3043 detList.BeginsWith(det) ||
3044 detList.EndsWith(det) ||
3045 detList.Contains( " "+det+" " )) {
3046 if (!resultList.EndsWith(det + " ")) {
3055 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3056 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3057 if ((detList.CompareTo("ALL") == 0) ||
3058 detList.BeginsWith("ALL ") ||
3059 detList.EndsWith(" ALL") ||
3060 detList.Contains(" ALL ") ||
3061 (detList.CompareTo(hltDet) == 0) ||
3062 detList.BeginsWith(hltDet) ||
3063 detList.EndsWith(hltDet) ||
3064 detList.Contains( " "+hltDet+" " )) {
3065 resultList += hltDet;
3069 return resultList.Data();
3073 //______________________________________________________________________________
3074 void AliReconstruction::Abort(const char *method, EAbort what)
3076 // Abort processing. If what = kAbortProcess, the Process() loop will be
3077 // aborted. If what = kAbortFile, the current file in a chain will be
3078 // aborted and the processing will continue with the next file, if there
3079 // is no next file then Process() will be aborted. Abort() can also be
3080 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3081 // the SlaveTerminate() and Terminate() are always called. The abort flag
3082 // can be checked in these methods using GetAbort().
3084 // The method is overwritten in AliReconstruction for better handling of
3085 // reco specific errors
3087 if (!fStopOnError) return;
3091 TString whyMess = method;
3092 whyMess += " failed! Aborting...";
3094 AliError(whyMess.Data());
3097 TString mess = "Abort";
3098 if (fAbort == kAbortProcess)
3099 mess = "AbortProcess";
3100 else if (fAbort == kAbortFile)
3103 Info(mess, whyMess.Data());