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 <TGeoGlobalMagField.h>
117 #include <TGeoManager.h>
119 #include <TLorentzVector.h>
121 #include <TObjArray.h>
122 #include <TPRegexp.h>
123 #include <TParameter.h>
124 #include <TPluginManager.h>
126 #include <TProofOutputFile.h>
130 #include "AliAlignObj.h"
131 #include "AliCDBEntry.h"
132 #include "AliCDBManager.h"
133 #include "AliCDBStorage.h"
134 #include "AliCTPRawStream.h"
135 #include "AliCascadeVertexer.h"
136 #include "AliCentralTrigger.h"
137 #include "AliCodeTimer.h"
139 #include "AliDetectorRecoParam.h"
140 #include "AliESDCaloCells.h"
141 #include "AliESDCaloCluster.h"
142 #include "AliESDEvent.h"
143 #include "AliESDMuonTrack.h"
144 #include "AliESDPmdTrack.h"
145 #include "AliESDTagCreator.h"
146 #include "AliESDVertex.h"
147 #include "AliESDcascade.h"
148 #include "AliESDfriend.h"
149 #include "AliESDkink.h"
150 #include "AliESDpid.h"
151 #include "AliESDtrack.h"
152 #include "AliESDtrack.h"
153 #include "AliEventInfo.h"
154 #include "AliGRPObject.h"
155 #include "AliGRPRecoParam.h"
156 #include "AliGenEventHeader.h"
157 #include "AliGeomManager.h"
158 #include "AliGlobalQADataMaker.h"
159 #include "AliHeader.h"
162 #include "AliMultiplicity.h"
164 #include "AliPlaneEff.h"
166 #include "AliQADataMakerRec.h"
167 #include "AliQAManager.h"
168 #include "AliRawVEvent.h"
169 #include "AliRawEventHeaderBase.h"
170 #include "AliRawHLTManager.h"
171 #include "AliRawReaderDate.h"
172 #include "AliRawReaderFile.h"
173 #include "AliRawReaderRoot.h"
174 #include "AliReconstruction.h"
175 #include "AliReconstructor.h"
177 #include "AliRunInfo.h"
178 #include "AliRunLoader.h"
179 #include "AliSysInfo.h" // memory snapshots
180 #include "AliTrackPointArray.h"
181 #include "AliTracker.h"
182 #include "AliTriggerClass.h"
183 #include "AliTriggerCluster.h"
184 #include "AliTriggerConfiguration.h"
185 #include "AliV0vertexer.h"
186 #include "AliVertexer.h"
187 #include "AliVertexerTracks.h"
189 ClassImp(AliReconstruction)
191 //_____________________________________________________________________________
192 const char* AliReconstruction::fgkDetectorName[AliReconstruction::kNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
194 //_____________________________________________________________________________
195 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
197 fUniformField(kFALSE),
198 fRunVertexFinder(kTRUE),
199 fRunVertexFinderTracks(kTRUE),
200 fRunHLTTracking(kFALSE),
201 fRunMuonTracking(kFALSE),
203 fRunCascadeFinder(kTRUE),
204 fStopOnError(kFALSE),
205 fWriteAlignmentData(kFALSE),
206 fWriteESDfriend(kFALSE),
207 fFillTriggerESD(kTRUE),
215 fRunLocalReconstruction("ALL"),
219 fUseTrackingErrorsForAlignment(""),
220 fGAliceFileName(gAliceFilename),
225 fNumberOfEventsPerFile((UInt_t)-1),
227 fLoadAlignFromCDB(kTRUE),
228 fLoadAlignData("ALL"),
235 fParentRawReader(NULL),
239 fSPDTrackleter(NULL),
241 fDiamondProfileSPD(NULL),
242 fDiamondProfile(NULL),
243 fDiamondProfileTPC(NULL),
247 fAlignObjArray(NULL),
251 fInitCDBCalled(kFALSE),
252 fSetRunNumberFromDataCalled(kFALSE),
258 fSameQACycle(kFALSE),
259 fInitQACalled(kFALSE),
260 fWriteQAExpertData(kTRUE),
261 fRunPlaneEff(kFALSE),
270 fIsNewRunLoader(kFALSE),
274 // create reconstruction object with default parameters
277 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
278 fReconstructor[iDet] = NULL;
279 fLoader[iDet] = NULL;
280 fTracker[iDet] = NULL;
282 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
283 fQACycles[iDet] = 999999 ;
284 fQAWriteExpert[iDet] = kFALSE ;
290 //_____________________________________________________________________________
291 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
293 fUniformField(rec.fUniformField),
294 fRunVertexFinder(rec.fRunVertexFinder),
295 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
296 fRunHLTTracking(rec.fRunHLTTracking),
297 fRunMuonTracking(rec.fRunMuonTracking),
298 fRunV0Finder(rec.fRunV0Finder),
299 fRunCascadeFinder(rec.fRunCascadeFinder),
300 fStopOnError(rec.fStopOnError),
301 fWriteAlignmentData(rec.fWriteAlignmentData),
302 fWriteESDfriend(rec.fWriteESDfriend),
303 fFillTriggerESD(rec.fFillTriggerESD),
305 fCleanESD(rec.fCleanESD),
306 fV0DCAmax(rec.fV0DCAmax),
307 fV0CsPmin(rec.fV0CsPmin),
311 fRunLocalReconstruction(rec.fRunLocalReconstruction),
312 fRunTracking(rec.fRunTracking),
313 fFillESD(rec.fFillESD),
314 fLoadCDB(rec.fLoadCDB),
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 fUseHLTData(rec.fUseHLTData),
331 fParentRawReader(NULL),
333 fRecoParam(rec.fRecoParam),
335 fSPDTrackleter(NULL),
337 fDiamondProfileSPD(rec.fDiamondProfileSPD),
338 fDiamondProfile(rec.fDiamondProfile),
339 fDiamondProfileTPC(rec.fDiamondProfileTPC),
343 fAlignObjArray(rec.fAlignObjArray),
344 fCDBUri(rec.fCDBUri),
345 fQARefUri(rec.fQARefUri),
347 fInitCDBCalled(rec.fInitCDBCalled),
348 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
349 fQADetectors(rec.fQADetectors),
351 fQATasks(rec.fQATasks),
353 fRunGlobalQA(rec.fRunGlobalQA),
354 fSameQACycle(rec.fSameQACycle),
355 fInitQACalled(rec.fInitQACalled),
356 fWriteQAExpertData(rec.fWriteQAExpertData),
357 fRunPlaneEff(rec.fRunPlaneEff),
366 fIsNewRunLoader(rec.fIsNewRunLoader),
372 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
373 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
375 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
376 fReconstructor[iDet] = NULL;
377 fLoader[iDet] = NULL;
378 fTracker[iDet] = NULL;
381 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
382 fQACycles[iDet] = rec.fQACycles[iDet];
383 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
386 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
387 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
391 //_____________________________________________________________________________
392 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
394 // assignment operator
395 // Used in PROOF mode
396 // Be very careful while modifing it!
397 // Simple rules to follow:
398 // for persistent data members - use their assignment operators
399 // for non-persistent ones - do nothing or take the default values from constructor
400 // TSelector members should not be touched
401 if(&rec == this) return *this;
403 fUniformField = rec.fUniformField;
404 fRunVertexFinder = rec.fRunVertexFinder;
405 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
406 fRunHLTTracking = rec.fRunHLTTracking;
407 fRunMuonTracking = rec.fRunMuonTracking;
408 fRunV0Finder = rec.fRunV0Finder;
409 fRunCascadeFinder = rec.fRunCascadeFinder;
410 fStopOnError = rec.fStopOnError;
411 fWriteAlignmentData = rec.fWriteAlignmentData;
412 fWriteESDfriend = rec.fWriteESDfriend;
413 fFillTriggerESD = rec.fFillTriggerESD;
415 fCleanESD = rec.fCleanESD;
416 fV0DCAmax = rec.fV0DCAmax;
417 fV0CsPmin = rec.fV0CsPmin;
421 fRunLocalReconstruction = rec.fRunLocalReconstruction;
422 fRunTracking = rec.fRunTracking;
423 fFillESD = rec.fFillESD;
424 fLoadCDB = rec.fLoadCDB;
425 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
426 fGAliceFileName = rec.fGAliceFileName;
427 fRawInput = rec.fRawInput;
428 fEquipIdMap = rec.fEquipIdMap;
429 fFirstEvent = rec.fFirstEvent;
430 fLastEvent = rec.fLastEvent;
431 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
433 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
434 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
437 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
438 fLoadAlignData = rec.fLoadAlignData;
439 fUseHLTData = rec.fUseHLTData;
441 delete fRunInfo; fRunInfo = NULL;
442 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
444 fEventInfo = rec.fEventInfo;
448 fParentRawReader = NULL;
450 fRecoParam = rec.fRecoParam;
452 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
453 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
454 delete fLoader[iDet]; fLoader[iDet] = NULL;
455 delete fTracker[iDet]; fTracker[iDet] = NULL;
458 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
459 fQACycles[iDet] = rec.fQACycles[iDet];
460 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
463 delete fSPDTrackleter; fSPDTrackleter = NULL;
465 delete fDiamondProfileSPD; fDiamondProfileSPD = NULL;
466 if (rec.fDiamondProfileSPD) fDiamondProfileSPD = new AliESDVertex(*rec.fDiamondProfileSPD);
467 delete fDiamondProfile; fDiamondProfile = NULL;
468 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
469 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
470 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
472 delete fGRPData; fGRPData = NULL;
473 // if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
474 if (rec.fGRPData) fGRPData = (AliGRPObject*)((rec.fGRPData)->Clone());
476 delete fAlignObjArray; fAlignObjArray = NULL;
479 fQARefUri = rec.fQARefUri;
480 fSpecCDBUri.Delete();
481 fInitCDBCalled = rec.fInitCDBCalled;
482 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
483 fQADetectors = rec.fQADetectors;
485 fQATasks = rec.fQATasks;
487 fRunGlobalQA = rec.fRunGlobalQA;
488 fSameQACycle = rec.fSameQACycle;
489 fInitQACalled = rec.fInitQACalled;
490 fWriteQAExpertData = rec.fWriteQAExpertData;
491 fRunPlaneEff = rec.fRunPlaneEff;
500 fIsNewRunLoader = rec.fIsNewRunLoader;
507 //_____________________________________________________________________________
508 AliReconstruction::~AliReconstruction()
515 if (fAlignObjArray) {
516 fAlignObjArray->Delete();
517 delete fAlignObjArray;
519 fSpecCDBUri.Delete();
521 AliCodeTimer::Instance()->Print();
524 //_____________________________________________________________________________
525 void AliReconstruction::InitQA()
527 //Initialize the QA and start of cycle
528 AliCodeTimerAuto("");
530 if (fInitQACalled) return;
531 fInitQACalled = kTRUE;
533 fQAManager = AliQAManager::QAManager("rec") ;
534 if (fWriteQAExpertData)
535 fQAManager->SetWriteExpert() ;
537 if (fQAManager->IsDefaultStorageSet()) {
538 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
539 AliWarning("Default QA reference storage has been already set !");
540 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fQARefUri.Data()));
541 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
542 fQARefUri = fQAManager->GetDefaultStorage()->GetURI();
544 if (fQARefUri.Length() > 0) {
545 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
546 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
547 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
549 fQARefUri="local://$ALICE_ROOT/QAref";
550 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
551 AliWarning("Default QA refeference storage not yet set !!!!");
552 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
553 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
556 fQAManager->SetDefaultStorage(fQARefUri);
560 fQAManager->SetActiveDetectors(fQADetectors) ;
561 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++) {
562 fQAManager->SetCycleLength(AliQAv1::DETECTORINDEX_t(det), fQACycles[det]) ;
563 fQAManager->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
565 if (!fRawReader && fQATasks.Contains(AliQAv1::kRAWS))
566 fQATasks.ReplaceAll(Form("%d",AliQAv1::kRAWS), "") ;
567 fQAManager->SetTasks(fQATasks) ;
568 fQAManager->InitQADataMaker(AliCDBManager::Instance()->GetRun()) ;
571 Bool_t sameCycle = kFALSE ;
572 AliQADataMaker *qadm = fQAManager->GetQADataMaker(AliQAv1::kGLOBAL);
573 AliInfo(Form("Initializing the global QA data maker"));
574 if (fQATasks.Contains(Form("%d", AliQAv1::kRECPOINTS))) {
575 qadm->StartOfCycle(AliQAv1::kRECPOINTS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
576 TObjArray **arr=qadm->Init(AliQAv1::kRECPOINTS);
577 AliTracker::SetResidualsArray(arr);
580 if (fQATasks.Contains(Form("%d", AliQAv1::kESDS))) {
581 qadm->StartOfCycle(AliQAv1::kESDS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
582 qadm->Init(AliQAv1::kESDS);
585 AliSysInfo::AddStamp("InitQA") ;
588 //_____________________________________________________________________________
589 void AliReconstruction::MergeQA()
591 //Initialize the QA and start of cycle
592 AliCodeTimerAuto("") ;
593 if ( ! fQAManager ) {
594 AliFatal("Hum... this should not happen") ;
596 fQAManager->Merge(AliCDBManager::Instance()->GetRun()) ;
598 AliSysInfo::AddStamp("MergeQA") ;
601 //_____________________________________________________________________________
602 void AliReconstruction::InitCDB()
604 // activate a default CDB storage
605 // First check if we have any CDB storage set, because it is used
606 // to retrieve the calibration and alignment constants
607 AliCodeTimerAuto("");
609 if (fInitCDBCalled) return;
610 fInitCDBCalled = kTRUE;
612 AliCDBManager* man = AliCDBManager::Instance();
613 if (man->IsDefaultStorageSet())
615 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
616 AliWarning("Default CDB storage has been already set !");
617 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
618 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
619 fCDBUri = man->GetDefaultStorage()->GetURI();
622 if (fCDBUri.Length() > 0)
624 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
625 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
626 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
628 fCDBUri="local://$ALICE_ROOT/OCDB";
629 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
630 AliWarning("Default CDB storage not yet set !!!!");
631 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
632 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
635 man->SetDefaultStorage(fCDBUri);
638 // Now activate the detector specific CDB storage locations
639 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
640 TObject* obj = fSpecCDBUri[i];
642 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
643 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
644 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
645 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
647 AliSysInfo::AddStamp("InitCDB");
650 //_____________________________________________________________________________
651 void AliReconstruction::SetDefaultStorage(const char* uri) {
652 // Store the desired default CDB storage location
653 // Activate it later within the Run() method
659 //_____________________________________________________________________________
660 void AliReconstruction::SetQARefDefaultStorage(const char* uri) {
661 // Store the desired default CDB storage location
662 // Activate it later within the Run() method
665 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
668 //_____________________________________________________________________________
669 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
670 // Store a detector-specific CDB storage location
671 // Activate it later within the Run() method
673 AliCDBPath aPath(calibType);
674 if(!aPath.IsValid()){
675 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
676 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
677 if(!strcmp(calibType, fgkDetectorName[iDet])) {
678 aPath.SetPath(Form("%s/*", calibType));
679 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
683 if(!aPath.IsValid()){
684 AliError(Form("Not a valid path or detector: %s", calibType));
689 // // check that calibType refers to a "valid" detector name
690 // Bool_t isDetector = kFALSE;
691 // for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
692 // TString detName = fgkDetectorName[iDet];
693 // if(aPath.GetLevel0() == detName) {
694 // isDetector = kTRUE;
700 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
704 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
705 if (obj) fSpecCDBUri.Remove(obj);
706 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
710 //_____________________________________________________________________________
711 Bool_t AliReconstruction::SetRunNumberFromData()
713 // The method is called in Run() in order
714 // to set a correct run number.
715 // In case of raw data reconstruction the
716 // run number is taken from the raw data header
718 if (fSetRunNumberFromDataCalled) return kTRUE;
719 fSetRunNumberFromDataCalled = kTRUE;
721 AliCDBManager* man = AliCDBManager::Instance();
724 if(fRawReader->NextEvent()) {
725 if(man->GetRun() > 0) {
726 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
728 man->SetRun(fRawReader->GetRunNumber());
729 fRawReader->RewindEvents();
732 if(man->GetRun() > 0) {
733 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
736 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
742 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
744 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
749 // read run number from gAlice
750 if(rl->GetHeader()) {
751 man->SetRun(rl->GetHeader()->GetRun());
756 AliError("Neither run-loader header nor RawReader objects are found !");
768 //_____________________________________________________________________________
769 void AliReconstruction::SetCDBLock() {
770 // Set CDB lock: from now on it is forbidden to reset the run number
771 // or the default storage or to activate any further storage!
773 AliCDBManager::Instance()->SetLock(1);
776 //_____________________________________________________________________________
777 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
779 // Read the alignment objects from CDB.
780 // Each detector is supposed to have the
781 // alignment objects in DET/Align/Data CDB path.
782 // All the detector objects are then collected,
783 // sorted by geometry level (starting from ALIC) and
784 // then applied to the TGeo geometry.
785 // Finally an overlaps check is performed.
787 // Load alignment data from CDB and fill fAlignObjArray
788 if(fLoadAlignFromCDB){
790 TString detStr = detectors;
791 TString loadAlObjsListOfDets = "";
793 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
794 if(!IsSelected(fgkDetectorName[iDet], detStr)) continue;
795 if(!strcmp(fgkDetectorName[iDet],"HLT")) continue;
797 if(AliGeomManager::GetNalignable(fgkDetectorName[iDet]) != 0)
799 loadAlObjsListOfDets += fgkDetectorName[iDet];
800 loadAlObjsListOfDets += " ";
802 } // end loop over detectors
804 if(AliGeomManager::GetNalignable("GRP") != 0)
805 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
806 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
807 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
809 // Check if the array with alignment objects was
810 // provided by the user. If yes, apply the objects
811 // to the present TGeo geometry
812 if (fAlignObjArray) {
813 if (gGeoManager && gGeoManager->IsClosed()) {
814 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
815 AliError("The misalignment of one or more volumes failed!"
816 "Compare the list of simulated detectors and the list of detector alignment data!");
821 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
827 if (fAlignObjArray) {
828 fAlignObjArray->Delete();
829 delete fAlignObjArray; fAlignObjArray=NULL;
835 //_____________________________________________________________________________
836 void AliReconstruction::SetGAliceFile(const char* fileName)
838 // set the name of the galice file
840 fGAliceFileName = fileName;
843 //_____________________________________________________________________________
844 void AliReconstruction::SetInput(const char* input)
846 // In case the input string starts with 'mem://', we run in an online mode
847 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
848 // file is assumed. One can give as an input:
849 // mem://: - events taken from DAQ monitoring libs online
851 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
852 if (input) fRawInput = input;
855 //_____________________________________________________________________________
856 void AliReconstruction::SetOption(const char* detector, const char* option)
858 // set options for the reconstruction of a detector
860 TObject* obj = fOptions.FindObject(detector);
861 if (obj) fOptions.Remove(obj);
862 fOptions.Add(new TNamed(detector, option));
865 //_____________________________________________________________________________
866 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
868 // Set custom reconstruction parameters for a given detector
869 // Single set of parameters for all the events
871 // First check if the reco-params are global
872 if(!strcmp(detector, "GRP")) {
874 fRecoParam.AddDetRecoParam(kNDetectors,par);
878 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
879 if(!strcmp(detector, fgkDetectorName[iDet])) {
881 fRecoParam.AddDetRecoParam(iDet,par);
888 //_____________________________________________________________________________
889 Bool_t AliReconstruction::SetFieldMap(Float_t l3Cur, Float_t diCur, Float_t l3Pol,
890 Float_t diPol, Float_t beamenergy,
891 const Char_t *beamtype, const Char_t *path)
893 //------------------------------------------------
894 // The magnetic field map, defined externally...
895 // L3 current 30000 A -> 0.5 T
896 // L3 current 12000 A -> 0.2 T
897 // dipole current 6000 A
898 // The polarities must be the same
899 //------------------------------------------------
900 const Float_t l3NominalCurrent1=30000.; // (A)
901 const Float_t l3NominalCurrent2=12000.; // (A)
902 const Float_t diNominalCurrent =6000. ; // (A)
904 const Float_t tolerance=0.03; // relative current tolerance
905 const Float_t zero=77.; // "zero" current (A)
907 TString s=(l3Pol < 0) ? "L3: -" : "L3: +";
909 AliMagF::BMap_t map = AliMagF::k5kG;
913 l3Cur = TMath::Abs(l3Cur);
914 if (TMath::Abs(l3Cur-l3NominalCurrent1)/l3NominalCurrent1 < tolerance) {
915 fcL3 = l3Cur/l3NominalCurrent1;
918 } else if (TMath::Abs(l3Cur-l3NominalCurrent2)/l3NominalCurrent2 < tolerance) {
919 fcL3 = l3Cur/l3NominalCurrent2;
922 } else if (l3Cur <= zero) {
924 map = AliMagF::k5kGUniform;
926 fUniformField=kTRUE; // track with the uniform (zero) B field
928 AliError(Form("Wrong L3 current (%f A)!",l3Cur));
932 diCur = TMath::Abs(diCur);
933 if (TMath::Abs(diCur-diNominalCurrent)/diNominalCurrent < tolerance) {
934 // 3% current tolerance...
935 fcDip = diCur/diNominalCurrent;
937 } else if (diCur <= zero) { // some small current..
941 AliError(Form("Wrong dipole current (%f A)!",diCur));
945 if (l3Pol!=diPol && (map==AliMagF::k5kG || map==AliMagF::k2kG) && fcDip!=0) {
946 AliError("L3 and Dipole polarities must be the same");
950 if (l3Pol<0) fcL3 = -fcL3;
951 if (diPol<0) fcDip = -fcDip;
953 AliMagF::BeamType_t btype = AliMagF::kNoBeamField;
954 TString btypestr = beamtype;
956 TPRegexp protonBeam("(proton|p)\\s*-?\\s*\\1");
957 TPRegexp ionBeam("(lead|pb|ion|a)\\s*-?\\s*\\1");
958 if (btypestr.Contains(ionBeam)) btype = AliMagF::kBeamTypeAA;
959 else if (btypestr.Contains(protonBeam)) btype = AliMagF::kBeamTypepp;
961 AliInfo(Form("Cannot determine the beam type from %s, assume no LHC magnet field",beamtype));
964 AliMagF* fld = new AliMagF("MagneticFieldMap", s.Data(), 2, fcL3, fcDip, 10., map, path,
966 TGeoGlobalMagField::Instance()->SetField( fld );
967 TGeoGlobalMagField::Instance()->Lock();
973 Bool_t AliReconstruction::InitGRP() {
974 //------------------------------------
975 // Initialization of the GRP entry
976 //------------------------------------
977 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
981 TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
984 AliInfo("Found a TMap in GRP/GRP/Data, converting it into an AliGRPObject");
986 fGRPData = new AliGRPObject();
987 fGRPData->ReadValuesFromMap(m);
991 AliInfo("Found an AliGRPObject in GRP/GRP/Data, reading it");
992 fGRPData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
996 AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
1000 AliError("No GRP entry found in OCDB!");
1004 TString lhcState = fGRPData->GetLHCState();
1005 if (lhcState==AliGRPObject::GetInvalidString()) {
1006 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
1007 lhcState = "UNKNOWN";
1010 TString beamType = fGRPData->GetBeamType();
1011 if (beamType==AliGRPObject::GetInvalidString()) {
1012 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
1013 beamType = "UNKNOWN";
1016 Float_t beamEnergy = fGRPData->GetBeamEnergy();
1017 if (beamEnergy==AliGRPObject::GetInvalidFloat()) {
1018 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
1021 // energy is provided in MeV*120
1022 beamEnergy /= 120E3;
1024 TString runType = fGRPData->GetRunType();
1025 if (runType==AliGRPObject::GetInvalidString()) {
1026 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
1027 runType = "UNKNOWN";
1030 Int_t activeDetectors = fGRPData->GetDetectorMask();
1031 if (activeDetectors==AliGRPObject::GetInvalidUInt()) {
1032 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
1033 activeDetectors = 1074790399;
1036 fRunInfo = new AliRunInfo(lhcState, beamType, beamEnergy, runType, activeDetectors);
1040 // Process the list of active detectors
1041 if (activeDetectors) {
1042 UInt_t detMask = activeDetectors;
1043 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
1044 fRunTracking = MatchDetectorList(fRunTracking,detMask);
1045 fFillESD = MatchDetectorList(fFillESD,detMask);
1046 fQADetectors = MatchDetectorList(fQADetectors,detMask);
1047 fLoadCDB.Form("%s %s %s %s",
1048 fRunLocalReconstruction.Data(),
1049 fRunTracking.Data(),
1051 fQADetectors.Data());
1052 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
1053 if (!((detMask >> AliDAQ::DetectorID("ITSSPD")) & 0x1)) {
1054 // switch off the vertexer
1055 AliInfo("SPD is not in the list of active detectors. Vertexer switched off.");
1056 fRunVertexFinder = kFALSE;
1058 if (!((detMask >> AliDAQ::DetectorID("TRG")) & 0x1)) {
1059 // switch off the reading of CTP raw-data payload
1060 if (fFillTriggerESD) {
1061 AliInfo("CTP is not in the list of active detectors. CTP data reading switched off.");
1062 fFillTriggerESD = kFALSE;
1067 AliInfo("===================================================================================");
1068 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
1069 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
1070 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
1071 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
1072 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
1073 AliInfo("===================================================================================");
1075 //*** Dealing with the magnetic field map
1076 if ( TGeoGlobalMagField::Instance()->IsLocked() ) {AliInfo("Running with the externally locked B field !");}
1078 // Construct the field map out of the information retrieved from GRP.
1081 Float_t l3Current = fGRPData->GetL3Current((AliGRPObject::Stats)0);
1082 if (l3Current == AliGRPObject::GetInvalidFloat()) {
1083 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1087 Char_t l3Polarity = fGRPData->GetL3Polarity();
1088 if (l3Polarity == AliGRPObject::GetInvalidChar()) {
1089 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1094 Float_t diCurrent = fGRPData->GetDipoleCurrent((AliGRPObject::Stats)0);
1095 if (diCurrent == AliGRPObject::GetInvalidFloat()) {
1096 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1100 Char_t diPolarity = fGRPData->GetDipolePolarity();
1101 if (diPolarity == AliGRPObject::GetInvalidChar()) {
1102 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1107 TObjString *l3Current=
1108 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Current"));
1110 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1113 TObjString *l3Polarity=
1114 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Polarity"));
1116 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1121 TObjString *diCurrent=
1122 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipoleCurrent"));
1124 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1127 TObjString *diPolarity=
1128 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipolePolarity"));
1130 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1136 if ( !SetFieldMap(l3Current, diCurrent, l3Polarity ? -1:1, diPolarity ? -1:1) )
1137 AliFatal("Failed to creat a B field map ! Exiting...");
1138 AliInfo("Running with the B field constructed out of GRP !");
1140 else AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
1144 //*** Get the diamond profiles from OCDB
1145 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1147 fDiamondProfileSPD = dynamic_cast<AliESDVertex*> (entry->GetObject());
1149 AliError("No SPD diamond profile found in OCDB!");
1152 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1154 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
1156 AliError("No diamond profile found in OCDB!");
1159 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
1161 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
1163 AliError("No TPC diamond profile found in OCDB!");
1169 //_____________________________________________________________________________
1170 Bool_t AliReconstruction::LoadCDB()
1172 AliCodeTimerAuto("");
1174 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1176 TString detStr = fLoadCDB;
1177 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1178 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1179 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1184 //_____________________________________________________________________________
1185 Bool_t AliReconstruction::Run(const char* input)
1188 AliCodeTimerAuto("");
1191 if (GetAbort() != TSelector::kContinue) return kFALSE;
1193 TChain *chain = NULL;
1194 if (fRawReader && (chain = fRawReader->GetChain())) {
1197 gProof->AddInput(this);
1199 outputFile.SetProtocol("root",kTRUE);
1200 outputFile.SetHost(gSystem->HostName());
1201 outputFile.SetFile(Form("%s/AliESDs.root",gSystem->pwd()));
1202 AliInfo(Form("Output file with ESDs is %s",outputFile.GetUrl()));
1203 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",outputFile.GetUrl()));
1205 chain->Process("AliReconstruction");
1208 chain->Process(this);
1213 if (GetAbort() != TSelector::kContinue) return kFALSE;
1215 if (GetAbort() != TSelector::kContinue) return kFALSE;
1216 //******* The loop over events
1217 AliInfo("Starting looping over events");
1219 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1220 (fRawReader && fRawReader->NextEvent())) {
1221 if (!ProcessEvent(iEvent)) {
1222 Abort("ProcessEvent",TSelector::kAbortFile);
1228 if (GetAbort() != TSelector::kContinue) return kFALSE;
1230 if (GetAbort() != TSelector::kContinue) return kFALSE;
1236 //_____________________________________________________________________________
1237 void AliReconstruction::InitRawReader(const char* input)
1239 AliCodeTimerAuto("");
1241 // Init raw-reader and
1242 // set the input in case of raw data
1243 if (input) fRawInput = input;
1244 fRawReader = AliRawReader::Create(fRawInput.Data());
1246 AliInfo("Reconstruction will run over digits");
1248 if (!fEquipIdMap.IsNull() && fRawReader)
1249 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1251 if (!fUseHLTData.IsNull()) {
1252 // create the RawReaderHLT which performs redirection of HLT input data for
1253 // the specified detectors
1254 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1256 fParentRawReader=fRawReader;
1257 fRawReader=pRawReader;
1259 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1262 AliSysInfo::AddStamp("CreateRawReader");
1265 //_____________________________________________________________________________
1266 void AliReconstruction::InitRun(const char* input)
1268 // Initialization of raw-reader,
1269 // run number, CDB etc.
1270 AliCodeTimerAuto("");
1271 AliSysInfo::AddStamp("Start");
1273 // Initialize raw-reader if any
1274 InitRawReader(input);
1276 // Initialize the CDB storage
1279 // Set run number in CDBManager (if it is not already set by the user)
1280 if (!SetRunNumberFromData()) {
1281 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1285 // Set CDB lock: from now on it is forbidden to reset the run number
1286 // or the default storage or to activate any further storage!
1291 //_____________________________________________________________________________
1292 void AliReconstruction::Begin(TTree *)
1294 // Initialize AlReconstruction before
1295 // going into the event loop
1296 // Should follow the TSelector convention
1297 // i.e. initialize only the object on the client side
1298 AliCodeTimerAuto("");
1300 AliReconstruction *reco = NULL;
1302 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1305 AliSysInfo::AddStamp("ReadInputInBegin");
1308 // Import ideal TGeo geometry and apply misalignment
1310 TString geom(gSystem->DirName(fGAliceFileName));
1311 geom += "/geometry.root";
1312 AliGeomManager::LoadGeometry(geom.Data());
1314 Abort("LoadGeometry", TSelector::kAbortProcess);
1317 AliSysInfo::AddStamp("LoadGeom");
1318 TString detsToCheck=fRunLocalReconstruction;
1319 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1320 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1323 AliSysInfo::AddStamp("CheckGeom");
1326 if (!MisalignGeometry(fLoadAlignData)) {
1327 Abort("MisalignGeometry", TSelector::kAbortProcess);
1330 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1331 AliSysInfo::AddStamp("MisalignGeom");
1334 Abort("InitGRP", TSelector::kAbortProcess);
1337 AliSysInfo::AddStamp("InitGRP");
1340 Abort("LoadCDB", TSelector::kAbortProcess);
1343 AliSysInfo::AddStamp("LoadCDB");
1345 // Read the reconstruction parameters from OCDB
1346 if (!InitRecoParams()) {
1347 AliWarning("Not all detectors have correct RecoParam objects initialized");
1349 AliSysInfo::AddStamp("InitRecoParams");
1351 if (fInput && gProof) {
1352 if (reco) *reco = *this;
1354 gProof->AddInputData(gGeoManager,kTRUE);
1356 gProof->AddInputData(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()),kTRUE);
1357 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1358 AliMagF *magFieldMap = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1359 magFieldMap->SetName("MagneticFieldMap");
1360 gProof->AddInputData(magFieldMap,kTRUE);
1365 //_____________________________________________________________________________
1366 void AliReconstruction::SlaveBegin(TTree*)
1368 // Initialization related to run-loader,
1369 // vertexer, trackers, recontructors
1370 // In proof mode it is executed on the slave
1371 AliCodeTimerAuto("");
1373 TProofOutputFile *outProofFile = NULL;
1375 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1378 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1380 AliGeomManager::SetGeometry(tgeo);
1382 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1383 Int_t runNumber = -1;
1384 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1385 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1386 man->SetCacheFlag(kTRUE);
1387 man->SetLock(kTRUE);
1391 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1392 TGeoGlobalMagField::Instance()->SetField(map);
1394 if (TNamed *outputFileName = (TNamed *) fInput->FindObject("PROOF_OUTPUTFILE")) {
1395 outProofFile = new TProofOutputFile(gSystem->BaseName(TUrl(outputFileName->GetTitle()).GetFile()));
1396 outProofFile->SetOutputFileName(outputFileName->GetTitle());
1397 fOutput->Add(outProofFile);
1399 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1402 // get the run loader
1403 if (!InitRunLoader()) {
1404 Abort("InitRunLoader", TSelector::kAbortProcess);
1407 AliSysInfo::AddStamp("LoadLoader");
1409 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1412 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1413 Abort("CreateTrackers", TSelector::kAbortProcess);
1416 AliSysInfo::AddStamp("CreateTrackers");
1418 // create the ESD output file and tree
1419 if (!outProofFile) {
1420 ffile = TFile::Open("AliESDs.root", "RECREATE");
1421 ffile->SetCompressionLevel(2);
1422 if (!ffile->IsOpen()) {
1423 Abort("OpenESDFile", TSelector::kAbortProcess);
1428 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1429 Abort(Form("Problems opening output PROOF file: %s/%s",
1430 outProofFile->GetDir(), outProofFile->GetFileName()),
1431 TSelector::kAbortProcess);
1436 ftree = new TTree("esdTree", "Tree with ESD objects");
1437 fesd = new AliESDEvent();
1438 fesd->CreateStdContent();
1440 fesd->WriteToTree(ftree);
1441 if (fWriteESDfriend) {
1443 // Since we add the branch manually we must
1444 // book and add it after WriteToTree
1445 // otherwise it is created twice,
1446 // once via writetotree and once here.
1447 // The case for AliESDfriend is now
1448 // caught also in AlIESDEvent::WriteToTree but
1449 // be careful when changing the name (AliESDfriend is not
1450 // a TNamed so we had to hardwire it)
1451 fesdf = new AliESDfriend();
1452 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1453 br->SetFile("AliESDfriends.root");
1454 fesd->AddObject(fesdf);
1456 ftree->GetUserInfo()->Add(fesd);
1458 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1459 fhltesd = new AliESDEvent();
1460 fhltesd->CreateStdContent();
1462 // read the ESD template from CDB
1463 // HLT is allowed to put non-std content to its ESD, the non-std
1464 // objects need to be created before invocation of WriteToTree in
1465 // order to create all branches. Initialization is done from an
1466 // ESD layout template in CDB
1467 AliCDBManager* man = AliCDBManager::Instance();
1468 AliCDBPath hltESDConfigPath("HLT/ConfigHLT/esdLayout");
1469 AliCDBEntry* hltESDConfig=NULL;
1470 if (man->GetId(hltESDConfigPath)!=NULL &&
1471 (hltESDConfig=man->Get(hltESDConfigPath))!=NULL) {
1472 AliESDEvent* pESDLayout=dynamic_cast<AliESDEvent*>(hltESDConfig->GetObject());
1474 // init all internal variables from the list of objects
1475 pESDLayout->GetStdContent();
1477 // copy content and create non-std objects
1478 *fhltesd=*pESDLayout;
1481 AliError(Form("error setting hltEsd layout from %s: invalid object type",
1482 hltESDConfigPath.GetPath().Data()));
1486 fhltesd->WriteToTree(fhlttree);
1487 fhlttree->GetUserInfo()->Add(fhltesd);
1489 ProcInfo_t procInfo;
1490 gSystem->GetProcInfo(&procInfo);
1491 AliInfo(Form("Current memory usage %d %d", procInfo.fMemResident, procInfo.fMemVirtual));
1494 //Initialize the QA and start of cycle
1495 if (fRunQA || fRunGlobalQA)
1498 //Initialize the Plane Efficiency framework
1499 if (fRunPlaneEff && !InitPlaneEff()) {
1500 Abort("InitPlaneEff", TSelector::kAbortProcess);
1504 if (strcmp(gProgName,"alieve") == 0)
1505 fRunAliEVE = InitAliEVE();
1510 //_____________________________________________________________________________
1511 Bool_t AliReconstruction::Process(Long64_t entry)
1513 // run the reconstruction over a single entry
1514 // from the chain with raw data
1515 AliCodeTimerAuto("");
1517 TTree *currTree = fChain->GetTree();
1518 AliRawVEvent *event = NULL;
1519 currTree->SetBranchAddress("rawevent",&event);
1520 currTree->GetEntry(entry);
1521 fRawReader = new AliRawReaderRoot(event);
1522 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1530 //_____________________________________________________________________________
1531 void AliReconstruction::Init(TTree *tree)
1534 AliError("The input tree is not found!");
1540 //_____________________________________________________________________________
1541 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1543 // run the reconstruction over a single event
1544 // The event loop is steered in Run method
1546 AliCodeTimerAuto("");
1548 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1549 fRunLoader->SetEventNumber(iEvent);
1550 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1552 fRunLoader->TreeE()->Fill();
1553 if (fRawReader && fRawReader->UseAutoSaveESD())
1554 fRunLoader->TreeE()->AutoSave("SaveSelf");
1557 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1561 AliInfo(Form("processing event %d", iEvent));
1563 fRunLoader->GetEvent(iEvent);
1565 // Fill Event-info object
1567 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo);
1568 AliInfo(Form("Current event specie: %s",fRecoParam.PrintEventSpecie()));
1570 // Set the reco-params
1572 TString detStr = fLoadCDB;
1573 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1574 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1575 AliReconstructor *reconstructor = GetReconstructor(iDet);
1576 if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
1577 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
1578 reconstructor->SetRecoParam(par);
1580 fQAManager->SetRecoParam(iDet, par) ;
1588 fQAManager->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1589 fQAManager->RunOneEvent(fRawReader) ;
1591 // local single event reconstruction
1592 if (!fRunLocalReconstruction.IsNull()) {
1593 TString detectors=fRunLocalReconstruction;
1594 // run HLT event reconstruction first
1595 // ;-( IsSelected changes the string
1596 if (IsSelected("HLT", detectors) &&
1597 !RunLocalEventReconstruction("HLT")) {
1598 if (fStopOnError) {CleanUp(); return kFALSE;}
1600 detectors=fRunLocalReconstruction;
1601 detectors.ReplaceAll("HLT", "");
1602 if (!RunLocalEventReconstruction(detectors)) {
1603 if (fStopOnError) {CleanUp(); return kFALSE;}
1607 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1608 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1609 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1610 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1612 // Set magnetic field from the tracker
1613 fesd->SetMagneticField(AliTracker::GetBz());
1614 fhltesd->SetMagneticField(AliTracker::GetBz());
1616 // Set most probable pt, for B=0 tracking
1617 // Get the global reco-params. They are atposition 16 inside the array of detectors in fRecoParam
1618 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(kNDetectors));
1619 if (grpRecoParam) AliExternalTrackParam::SetMostProbablePt(grpRecoParam->GetMostProbablePt());
1621 // Fill raw-data error log into the ESD
1622 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1625 if (fRunVertexFinder) {
1626 if (!RunVertexFinder(fesd)) {
1627 if (fStopOnError) {CleanUp(); return kFALSE;}
1631 // For Plane Efficiency: run the SPD trackleter
1632 if (fRunPlaneEff && fSPDTrackleter) {
1633 if (!RunSPDTrackleting(fesd)) {
1634 if (fStopOnError) {CleanUp(); return kFALSE;}
1639 if (!fRunTracking.IsNull()) {
1640 if (fRunMuonTracking) {
1641 if (!RunMuonTracking(fesd)) {
1642 if (fStopOnError) {CleanUp(); return kFALSE;}
1648 if (!fRunTracking.IsNull()) {
1649 if (!RunTracking(fesd)) {
1650 if (fStopOnError) {CleanUp(); return kFALSE;}
1655 if (!fFillESD.IsNull()) {
1656 TString detectors=fFillESD;
1657 // run HLT first and on hltesd
1658 // ;-( IsSelected changes the string
1659 if (IsSelected("HLT", detectors) &&
1660 !FillESD(fhltesd, "HLT")) {
1661 if (fStopOnError) {CleanUp(); return kFALSE;}
1664 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1665 if (detectors.Contains("ALL")) {
1667 for (Int_t idet=0; idet<kNDetectors; ++idet){
1668 detectors += fgkDetectorName[idet];
1672 detectors.ReplaceAll("HLT", "");
1673 if (!FillESD(fesd, detectors)) {
1674 if (fStopOnError) {CleanUp(); return kFALSE;}
1678 // fill Event header information from the RawEventHeader
1679 if (fRawReader){FillRawEventHeaderESD(fesd);}
1682 AliESDpid::MakePID(fesd);
1684 if (fFillTriggerESD) {
1685 if (!FillTriggerESD(fesd)) {
1686 if (fStopOnError) {CleanUp(); return kFALSE;}
1693 // Propagate track to the beam pipe (if not already done by ITS)
1695 const Int_t ntracks = fesd->GetNumberOfTracks();
1696 const Double_t kBz = fesd->GetMagneticField();
1697 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1700 UShort_t *selectedIdx=new UShort_t[ntracks];
1702 for (Int_t itrack=0; itrack<ntracks; itrack++){
1703 const Double_t kMaxStep = 1; //max step over the material
1706 AliESDtrack *track = fesd->GetTrack(itrack);
1707 if (!track) continue;
1709 AliExternalTrackParam *tpcTrack =
1710 (AliExternalTrackParam *)track->GetTPCInnerParam();
1714 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kFALSE);
1717 Int_t n=trkArray.GetEntriesFast();
1718 selectedIdx[n]=track->GetID();
1719 trkArray.AddLast(tpcTrack);
1722 //Tracks refitted by ITS should already be at the SPD vertex
1723 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1726 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kFALSE);
1727 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1732 // Improve the reconstructed primary vertex position using the tracks
1734 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
1735 if(fesd->GetPrimaryVertexSPD()) {
1736 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
1737 if(vtitle.Contains("cosmics")) {
1738 runVertexFinderTracks=kFALSE;
1742 if (runVertexFinderTracks) {
1743 // TPC + ITS primary vertex
1744 ftVertexer->SetITSMode();
1745 ftVertexer->SetConstraintOff();
1746 // get cuts for vertexer from AliGRPRecoParam
1748 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1749 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1750 grpRecoParam->GetVertexerTracksCutsITS(cutsVertexer);
1751 ftVertexer->SetCuts(cutsVertexer);
1752 delete [] cutsVertexer; cutsVertexer = NULL;
1753 if(fDiamondProfile && grpRecoParam->GetVertexerTracksConstraintITS())
1754 ftVertexer->SetVtxStart(fDiamondProfile);
1756 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1758 if (pvtx->GetStatus()) {
1759 fesd->SetPrimaryVertexTracks(pvtx);
1760 for (Int_t i=0; i<ntracks; i++) {
1761 AliESDtrack *t = fesd->GetTrack(i);
1762 t->RelateToVertex(pvtx, kBz, kVeryBig);
1767 // TPC-only primary vertex
1768 ftVertexer->SetTPCMode();
1769 ftVertexer->SetConstraintOff();
1770 // get cuts for vertexer from AliGRPRecoParam
1772 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1773 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1774 grpRecoParam->GetVertexerTracksCutsTPC(cutsVertexer);
1775 ftVertexer->SetCuts(cutsVertexer);
1776 delete [] cutsVertexer; cutsVertexer = NULL;
1777 if(fDiamondProfileTPC && grpRecoParam->GetVertexerTracksConstraintTPC())
1778 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1780 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1782 if (pvtx->GetStatus()) {
1783 fesd->SetPrimaryVertexTPC(pvtx);
1784 for (Int_t i=0; i<ntracks; i++) {
1785 AliESDtrack *t = fesd->GetTrack(i);
1786 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1792 delete[] selectedIdx;
1794 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1799 AliV0vertexer vtxer;
1800 vtxer.Tracks2V0vertices(fesd);
1802 if (fRunCascadeFinder) {
1804 AliCascadeVertexer cvtxer;
1805 cvtxer.V0sTracks2CascadeVertices(fesd);
1810 if (fCleanESD) CleanESD(fesd);
1813 fQAManager->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1814 fQAManager->RunOneEvent(fesd) ;
1817 AliQADataMaker *qadm = fQAManager->GetQADataMaker(AliQAv1::kGLOBAL);
1818 qadm->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1819 if (qadm && fQATasks.Contains(Form("%d", AliQAv1::kESDS)))
1820 qadm->Exec(AliQAv1::kESDS, fesd);
1823 if (fWriteESDfriend) {
1824 // fesdf->~AliESDfriend();
1825 // new (fesdf) AliESDfriend(); // Reset...
1826 fesd->GetESDfriend(fesdf);
1830 // Auto-save the ESD tree in case of prompt reco @P2
1831 if (fRawReader && fRawReader->UseAutoSaveESD()) {
1832 ftree->AutoSave("SaveSelf");
1833 TFile *friendfile = (TFile *)(gROOT->GetListOfFiles()->FindObject("AliESDfriends.root"));
1834 if (friendfile) friendfile->Save();
1841 if (fRunAliEVE) RunAliEVE();
1845 if (fWriteESDfriend) {
1846 fesdf->~AliESDfriend();
1847 new (fesdf) AliESDfriend(); // Reset...
1850 ProcInfo_t procInfo;
1851 gSystem->GetProcInfo(&procInfo);
1852 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, procInfo.fMemResident, procInfo.fMemVirtual));
1855 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1856 if (fReconstructor[iDet])
1857 fReconstructor[iDet]->SetRecoParam(NULL);
1860 if (fRunQA || fRunGlobalQA)
1861 fQAManager->Increment() ;
1866 //_____________________________________________________________________________
1867 void AliReconstruction::SlaveTerminate()
1869 // Finalize the run on the slave side
1870 // Called after the exit
1871 // from the event loop
1872 AliCodeTimerAuto("");
1874 if (fIsNewRunLoader) { // galice.root didn't exist
1875 fRunLoader->WriteHeader("OVERWRITE");
1876 fRunLoader->CdGAFile();
1877 fRunLoader->Write(0, TObject::kOverwrite);
1880 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1881 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1883 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1884 cdbMapCopy->SetOwner(1);
1885 cdbMapCopy->SetName("cdbMap");
1886 TIter iter(cdbMap->GetTable());
1889 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1890 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1891 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1892 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1895 TList *cdbListCopy = new TList();
1896 cdbListCopy->SetOwner(1);
1897 cdbListCopy->SetName("cdbList");
1899 TIter iter2(cdbList);
1902 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1903 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1906 ftree->GetUserInfo()->Add(cdbMapCopy);
1907 ftree->GetUserInfo()->Add(cdbListCopy);
1912 if (fWriteESDfriend)
1913 ftree->SetBranchStatus("ESDfriend*",0);
1914 // we want to have only one tree version number
1915 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1918 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1919 if (fRunPlaneEff && !FinishPlaneEff()) {
1920 AliWarning("Finish PlaneEff evaluation failed");
1923 // End of cycle for the in-loop
1925 fQAManager->EndOfCycle() ;
1928 AliQADataMaker *qadm = fQAManager->GetQADataMaker(AliQAv1::kGLOBAL);
1930 if (fQATasks.Contains(Form("%d", AliQAv1::kRECPOINTS)))
1931 qadm->EndOfCycle(AliQAv1::kRECPOINTS);
1932 if (fQATasks.Contains(Form("%d", AliQAv1::kESDS)))
1933 qadm->EndOfCycle(AliQAv1::kESDS);
1941 //_____________________________________________________________________________
1942 void AliReconstruction::Terminate()
1944 // Create tags for the events in the ESD tree (the ESD tree is always present)
1945 // In case of empty events the tags will contain dummy values
1946 AliCodeTimerAuto("");
1948 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1949 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData, AliQAv1::Instance()->GetQA(), AliQAv1::Instance()->GetEventSpecies(), AliQAv1::kNDET, AliRecoParam::kNSpecies);
1951 // Cleanup of CDB manager: cache and active storages!
1952 AliCDBManager::Instance()->ClearCache();
1955 //_____________________________________________________________________________
1956 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1958 // run the local reconstruction
1960 static Int_t eventNr=0;
1961 AliCodeTimerAuto("")
1963 TString detStr = detectors;
1964 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1965 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1966 AliReconstructor* reconstructor = GetReconstructor(iDet);
1967 if (!reconstructor) continue;
1968 AliLoader* loader = fLoader[iDet];
1969 // Matthias April 2008: temporary fix to run HLT reconstruction
1970 // although the HLT loader is missing
1971 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1973 reconstructor->Reconstruct(fRawReader, NULL);
1976 reconstructor->Reconstruct(dummy, NULL);
1981 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1984 // conversion of digits
1985 if (fRawReader && reconstructor->HasDigitConversion()) {
1986 AliInfo(Form("converting raw data digits into root objects for %s",
1987 fgkDetectorName[iDet]));
1988 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1989 // fgkDetectorName[iDet]));
1990 loader->LoadDigits("update");
1991 loader->CleanDigits();
1992 loader->MakeDigitsContainer();
1993 TTree* digitsTree = loader->TreeD();
1994 reconstructor->ConvertDigits(fRawReader, digitsTree);
1995 loader->WriteDigits("OVERWRITE");
1996 loader->UnloadDigits();
1998 // local reconstruction
1999 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2000 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2001 loader->LoadRecPoints("update");
2002 loader->CleanRecPoints();
2003 loader->MakeRecPointsContainer();
2004 TTree* clustersTree = loader->TreeR();
2005 if (fRawReader && !reconstructor->HasDigitConversion()) {
2006 reconstructor->Reconstruct(fRawReader, clustersTree);
2008 loader->LoadDigits("read");
2009 TTree* digitsTree = loader->TreeD();
2011 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
2012 if (fStopOnError) return kFALSE;
2014 reconstructor->Reconstruct(digitsTree, clustersTree);
2016 loader->UnloadDigits();
2019 TString detQAStr(fQADetectors) ;
2021 fQAManager->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2022 fQAManager->RunOneEventInOneDetector(iDet, clustersTree) ;
2024 loader->WriteRecPoints("OVERWRITE");
2025 loader->UnloadRecPoints();
2026 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
2028 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2029 AliError(Form("the following detectors were not found: %s",
2031 if (fStopOnError) return kFALSE;
2036 //_____________________________________________________________________________
2037 Bool_t AliReconstruction::RunSPDTrackleting(AliESDEvent*& esd)
2039 // run the SPD trackleting (for SPD efficiency purpouses)
2041 AliCodeTimerAuto("")
2043 Double_t vtxPos[3] = {0, 0, 0};
2044 Double_t vtxErr[3] = {0.0, 0.0, 0.0};
2046 TArrayF mcVertex(3);
2048 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
2049 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
2050 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
2053 const AliESDVertex *vertex = esd->GetVertex();
2055 AliWarning("Vertex not found");
2058 vertex->GetXYZ(vtxPos);
2059 vertex->GetSigmaXYZ(vtxErr);
2060 if (fSPDTrackleter) {
2061 AliInfo("running the SPD Trackleter for Plane Efficiency Evaluation");
2064 fLoader[0]->LoadRecPoints("read");
2065 TTree* tree = fLoader[0]->TreeR();
2067 AliError("Can't get the ITS cluster tree");
2070 fSPDTrackleter->LoadClusters(tree);
2071 fSPDTrackleter->SetVertex(vtxPos, vtxErr);
2073 if (fSPDTrackleter->Clusters2Tracks(esd) != 0) {
2074 AliError("AliITSTrackleterSPDEff Clusters2Tracks failed");
2075 // fLoader[0]->UnloadRecPoints();
2078 //fSPDTrackleter->UnloadRecPoints();
2080 AliWarning("SPDTrackleter not available");
2086 //_____________________________________________________________________________
2087 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
2089 // run the barrel tracking
2091 AliCodeTimerAuto("")
2093 AliVertexer *vertexer = CreateVertexer();
2094 if (!vertexer) return kFALSE;
2096 AliInfo("running the ITS vertex finder");
2097 AliESDVertex* vertex = NULL;
2099 fLoader[0]->LoadRecPoints();
2100 TTree* cltree = fLoader[0]->TreeR();
2102 if(fDiamondProfileSPD) vertexer->SetVtxStart(fDiamondProfileSPD);
2103 vertex = vertexer->FindVertexForCurrentEvent(cltree);
2106 AliError("Can't get the ITS cluster tree");
2108 fLoader[0]->UnloadRecPoints();
2111 AliError("Can't get the ITS loader");
2114 AliWarning("Vertex not found");
2115 vertex = new AliESDVertex();
2116 vertex->SetName("default");
2119 vertex->SetName("reconstructed");
2124 vertex->GetXYZ(vtxPos);
2125 vertex->GetSigmaXYZ(vtxErr);
2127 esd->SetPrimaryVertexSPD(vertex);
2128 AliESDVertex *vpileup = NULL;
2129 Int_t novertices = 0;
2130 vpileup = vertexer->GetAllVertices(novertices);
2132 for (Int_t kk=1; kk<novertices; kk++)esd->AddPileupVertexSPD(&vpileup[kk]);
2134 // if SPD multiplicity has been determined, it is stored in the ESD
2135 AliMultiplicity *mult = vertexer->GetMultiplicity();
2136 if(mult)esd->SetMultiplicity(mult);
2138 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2139 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
2148 //_____________________________________________________________________________
2149 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
2151 // run the HLT barrel tracking
2153 AliCodeTimerAuto("")
2156 AliError("Missing runLoader!");
2160 AliInfo("running HLT tracking");
2162 // Get a pointer to the HLT reconstructor
2163 AliReconstructor *reconstructor = GetReconstructor(kNDetectors-1);
2164 if (!reconstructor) return kFALSE;
2167 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2168 TString detName = fgkDetectorName[iDet];
2169 AliDebug(1, Form("%s HLT tracking", detName.Data()));
2170 reconstructor->SetOption(detName.Data());
2171 AliTracker *tracker = reconstructor->CreateTracker();
2173 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
2174 if (fStopOnError) return kFALSE;
2178 Double_t vtxErr[3]={0.005,0.005,0.010};
2179 const AliESDVertex *vertex = esd->GetVertex();
2180 vertex->GetXYZ(vtxPos);
2181 tracker->SetVertex(vtxPos,vtxErr);
2183 fLoader[iDet]->LoadRecPoints("read");
2184 TTree* tree = fLoader[iDet]->TreeR();
2186 AliError(Form("Can't get the %s cluster tree", detName.Data()));
2189 tracker->LoadClusters(tree);
2191 if (tracker->Clusters2Tracks(esd) != 0) {
2192 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
2196 tracker->UnloadClusters();
2204 //_____________________________________________________________________________
2205 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
2207 // run the muon spectrometer tracking
2209 AliCodeTimerAuto("")
2212 AliError("Missing runLoader!");
2215 Int_t iDet = 7; // for MUON
2217 AliInfo("is running...");
2219 // Get a pointer to the MUON reconstructor
2220 AliReconstructor *reconstructor = GetReconstructor(iDet);
2221 if (!reconstructor) return kFALSE;
2224 TString detName = fgkDetectorName[iDet];
2225 AliDebug(1, Form("%s tracking", detName.Data()));
2226 AliTracker *tracker = reconstructor->CreateTracker();
2228 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2233 fLoader[iDet]->LoadRecPoints("read");
2235 tracker->LoadClusters(fLoader[iDet]->TreeR());
2237 Int_t rv = tracker->Clusters2Tracks(esd);
2241 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2245 fLoader[iDet]->UnloadRecPoints();
2247 tracker->UnloadClusters();
2255 //_____________________________________________________________________________
2256 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2258 // run the barrel tracking
2259 static Int_t eventNr=0;
2260 AliCodeTimerAuto("")
2262 AliInfo("running tracking");
2264 //Fill the ESD with the T0 info (will be used by the TOF)
2265 if (fReconstructor[11] && fLoader[11]) {
2266 fLoader[11]->LoadRecPoints("READ");
2267 TTree *treeR = fLoader[11]->TreeR();
2269 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2273 // pass 1: TPC + ITS inwards
2274 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2275 if (!fTracker[iDet]) continue;
2276 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2279 fLoader[iDet]->LoadRecPoints("read");
2280 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2281 TTree* tree = fLoader[iDet]->TreeR();
2283 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2286 fTracker[iDet]->LoadClusters(tree);
2287 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2289 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2290 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2293 // preliminary PID in TPC needed by the ITS tracker
2295 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2296 AliESDpid::MakePID(esd);
2298 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2301 // pass 2: ALL backwards
2303 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2304 if (!fTracker[iDet]) continue;
2305 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2308 if (iDet > 1) { // all except ITS, TPC
2310 fLoader[iDet]->LoadRecPoints("read");
2311 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2312 tree = fLoader[iDet]->TreeR();
2314 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2317 fTracker[iDet]->LoadClusters(tree);
2318 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2322 if (iDet>1) // start filling residuals for the "outer" detectors
2323 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2325 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2326 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2331 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2332 fTracker[iDet]->UnloadClusters();
2333 fLoader[iDet]->UnloadRecPoints();
2335 // updated PID in TPC needed by the ITS tracker -MI
2337 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2338 AliESDpid::MakePID(esd);
2340 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2342 //stop filling residuals for the "outer" detectors
2343 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2345 // pass 3: TRD + TPC + ITS refit inwards
2347 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2348 if (!fTracker[iDet]) continue;
2349 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2352 if (iDet<2) // start filling residuals for TPC and ITS
2353 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2355 if (fTracker[iDet]->RefitInward(esd) != 0) {
2356 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2359 // run postprocessing
2360 if (fTracker[iDet]->PostProcess(esd) != 0) {
2361 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2364 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2367 // write space-points to the ESD in case alignment data output
2369 if (fWriteAlignmentData)
2370 WriteAlignmentData(esd);
2372 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2373 if (!fTracker[iDet]) continue;
2375 fTracker[iDet]->UnloadClusters();
2376 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2377 fLoader[iDet]->UnloadRecPoints();
2378 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2380 // stop filling residuals for TPC and ITS
2381 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2387 //_____________________________________________________________________________
2388 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2390 // Remove the data which are not needed for the physics analysis.
2393 Int_t nTracks=esd->GetNumberOfTracks();
2394 Int_t nV0s=esd->GetNumberOfV0s();
2396 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2398 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2399 Bool_t rc=esd->Clean(cleanPars);
2401 nTracks=esd->GetNumberOfTracks();
2402 nV0s=esd->GetNumberOfV0s();
2404 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2409 //_____________________________________________________________________________
2410 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2412 // fill the event summary data
2414 AliCodeTimerAuto("")
2415 static Int_t eventNr=0;
2416 TString detStr = detectors;
2418 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2419 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2420 AliReconstructor* reconstructor = GetReconstructor(iDet);
2421 if (!reconstructor) continue;
2422 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2423 TTree* clustersTree = NULL;
2424 if (fLoader[iDet]) {
2425 fLoader[iDet]->LoadRecPoints("read");
2426 clustersTree = fLoader[iDet]->TreeR();
2427 if (!clustersTree) {
2428 AliError(Form("Can't get the %s clusters tree",
2429 fgkDetectorName[iDet]));
2430 if (fStopOnError) return kFALSE;
2433 if (fRawReader && !reconstructor->HasDigitConversion()) {
2434 reconstructor->FillESD(fRawReader, clustersTree, esd);
2436 TTree* digitsTree = NULL;
2437 if (fLoader[iDet]) {
2438 fLoader[iDet]->LoadDigits("read");
2439 digitsTree = fLoader[iDet]->TreeD();
2441 AliError(Form("Can't get the %s digits tree",
2442 fgkDetectorName[iDet]));
2443 if (fStopOnError) return kFALSE;
2446 reconstructor->FillESD(digitsTree, clustersTree, esd);
2447 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2449 if (fLoader[iDet]) {
2450 fLoader[iDet]->UnloadRecPoints();
2454 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2455 AliError(Form("the following detectors were not found: %s",
2457 if (fStopOnError) return kFALSE;
2459 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2464 //_____________________________________________________________________________
2465 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2467 // Reads the trigger decision which is
2468 // stored in Trigger.root file and fills
2469 // the corresponding esd entries
2471 AliCodeTimerAuto("")
2473 AliInfo("Filling trigger information into the ESD");
2476 AliCTPRawStream input(fRawReader);
2477 if (!input.Next()) {
2478 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2481 if (esd->GetTriggerMask() != input.GetClassMask())
2482 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2483 input.GetClassMask(),esd->GetTriggerMask()));
2484 if (esd->GetOrbitNumber() != input.GetOrbitID())
2485 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2486 input.GetOrbitID(),esd->GetOrbitNumber()));
2487 if (esd->GetBunchCrossNumber() != input.GetBCID())
2488 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2489 input.GetBCID(),esd->GetBunchCrossNumber()));
2492 // Here one has to add the filling of trigger inputs and
2493 // interaction records
2503 //_____________________________________________________________________________
2504 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2507 // Filling information from RawReader Header
2510 if (!fRawReader) return kFALSE;
2512 AliInfo("Filling information from RawReader Header");
2514 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2515 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2516 esd->SetPeriodNumber(fRawReader->GetPeriod());
2518 esd->SetTimeStamp(fRawReader->GetTimestamp());
2519 esd->SetEventType(fRawReader->GetType());
2525 //_____________________________________________________________________________
2526 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2528 // check whether detName is contained in detectors
2529 // if yes, it is removed from detectors
2531 // check if all detectors are selected
2532 if ((detectors.CompareTo("ALL") == 0) ||
2533 detectors.BeginsWith("ALL ") ||
2534 detectors.EndsWith(" ALL") ||
2535 detectors.Contains(" ALL ")) {
2540 // search for the given detector
2541 Bool_t result = kFALSE;
2542 if ((detectors.CompareTo(detName) == 0) ||
2543 detectors.BeginsWith(detName+" ") ||
2544 detectors.EndsWith(" "+detName) ||
2545 detectors.Contains(" "+detName+" ")) {
2546 detectors.ReplaceAll(detName, "");
2550 // clean up the detectors string
2551 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2552 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2553 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2558 //_____________________________________________________________________________
2559 Bool_t AliReconstruction::InitRunLoader()
2561 // get or create the run loader
2563 if (gAlice) delete gAlice;
2566 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2567 // load all base libraries to get the loader classes
2568 TString libs = gSystem->GetLibraries();
2569 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2570 TString detName = fgkDetectorName[iDet];
2571 if (detName == "HLT") continue;
2572 if (libs.Contains("lib" + detName + "base.so")) continue;
2573 gSystem->Load("lib" + detName + "base.so");
2575 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2577 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2582 fRunLoader->CdGAFile();
2583 fRunLoader->LoadgAlice();
2585 //PH This is a temporary fix to give access to the kinematics
2586 //PH that is needed for the labels of ITS clusters
2587 fRunLoader->LoadHeader();
2588 fRunLoader->LoadKinematics();
2590 } else { // galice.root does not exist
2592 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2594 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2595 AliConfig::GetDefaultEventFolderName(),
2598 AliError(Form("could not create run loader in file %s",
2599 fGAliceFileName.Data()));
2603 fIsNewRunLoader = kTRUE;
2604 fRunLoader->MakeTree("E");
2606 if (fNumberOfEventsPerFile > 0)
2607 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2609 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2615 //_____________________________________________________________________________
2616 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2618 // get the reconstructor object and the loader for a detector
2620 if (fReconstructor[iDet]) {
2621 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2622 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2623 fReconstructor[iDet]->SetRecoParam(par);
2625 return fReconstructor[iDet];
2628 // load the reconstructor object
2629 TPluginManager* pluginManager = gROOT->GetPluginManager();
2630 TString detName = fgkDetectorName[iDet];
2631 TString recName = "Ali" + detName + "Reconstructor";
2633 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2635 AliReconstructor* reconstructor = NULL;
2636 // first check if a plugin is defined for the reconstructor
2637 TPluginHandler* pluginHandler =
2638 pluginManager->FindHandler("AliReconstructor", detName);
2639 // if not, add a plugin for it
2640 if (!pluginHandler) {
2641 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2642 TString libs = gSystem->GetLibraries();
2643 if (libs.Contains("lib" + detName + "base.so") ||
2644 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2645 pluginManager->AddHandler("AliReconstructor", detName,
2646 recName, detName + "rec", recName + "()");
2648 pluginManager->AddHandler("AliReconstructor", detName,
2649 recName, detName, recName + "()");
2651 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2653 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2654 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2656 if (reconstructor) {
2657 TObject* obj = fOptions.FindObject(detName.Data());
2658 if (obj) reconstructor->SetOption(obj->GetTitle());
2659 reconstructor->Init();
2660 fReconstructor[iDet] = reconstructor;
2663 // get or create the loader
2664 if (detName != "HLT") {
2665 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2666 if (!fLoader[iDet]) {
2667 AliConfig::Instance()
2668 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2670 // first check if a plugin is defined for the loader
2672 pluginManager->FindHandler("AliLoader", detName);
2673 // if not, add a plugin for it
2674 if (!pluginHandler) {
2675 TString loaderName = "Ali" + detName + "Loader";
2676 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2677 pluginManager->AddHandler("AliLoader", detName,
2678 loaderName, detName + "base",
2679 loaderName + "(const char*, TFolder*)");
2680 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2682 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2684 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2685 fRunLoader->GetEventFolder());
2687 if (!fLoader[iDet]) { // use default loader
2688 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2690 if (!fLoader[iDet]) {
2691 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2692 if (fStopOnError) return NULL;
2694 fRunLoader->AddLoader(fLoader[iDet]);
2695 fRunLoader->CdGAFile();
2696 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2697 fRunLoader->Write(0, TObject::kOverwrite);
2702 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2703 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2704 reconstructor->SetRecoParam(par);
2706 return reconstructor;
2709 //_____________________________________________________________________________
2710 AliVertexer* AliReconstruction::CreateVertexer()
2712 // create the vertexer
2713 // Please note that the caller is the owner of the
2716 AliVertexer* vertexer = NULL;
2717 AliReconstructor* itsReconstructor = GetReconstructor(0);
2718 if (itsReconstructor) {
2719 vertexer = itsReconstructor->CreateVertexer();
2722 AliWarning("couldn't create a vertexer for ITS");
2728 //_____________________________________________________________________________
2729 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2731 // create the trackers
2732 AliInfo("Creating trackers");
2734 TString detStr = detectors;
2735 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2736 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2737 AliReconstructor* reconstructor = GetReconstructor(iDet);
2738 if (!reconstructor) continue;
2739 TString detName = fgkDetectorName[iDet];
2740 if (detName == "HLT") {
2741 fRunHLTTracking = kTRUE;
2744 if (detName == "MUON") {
2745 fRunMuonTracking = kTRUE;
2750 fTracker[iDet] = reconstructor->CreateTracker();
2751 if (!fTracker[iDet] && (iDet < 7)) {
2752 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2753 if (fStopOnError) return kFALSE;
2755 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2761 //_____________________________________________________________________________
2762 void AliReconstruction::CleanUp()
2764 // delete trackers and the run loader and close and delete the file
2766 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2767 delete fReconstructor[iDet];
2768 fReconstructor[iDet] = NULL;
2769 fLoader[iDet] = NULL;
2770 delete fTracker[iDet];
2771 fTracker[iDet] = NULL;
2776 delete fSPDTrackleter;
2777 fSPDTrackleter = NULL;
2782 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
2783 delete fDiamondProfileSPD;
2784 fDiamondProfileSPD = NULL;
2785 delete fDiamondProfile;
2786 fDiamondProfile = NULL;
2787 delete fDiamondProfileTPC;
2788 fDiamondProfileTPC = NULL;
2794 delete fParentRawReader;
2795 fParentRawReader=NULL;
2804 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2806 // Write space-points which are then used in the alignment procedures
2807 // For the moment only ITS, TPC, TRD and TOF
2809 Int_t ntracks = esd->GetNumberOfTracks();
2810 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2812 AliESDtrack *track = esd->GetTrack(itrack);
2815 for (Int_t i=0; i<200; ++i) idx[i] = -1; //PH avoid uninitialized values
2816 for (Int_t iDet = 5; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2817 nsp += track->GetNcls(iDet);
2819 if (iDet==0) { // ITS "extra" clusters
2820 track->GetClusters(iDet,idx);
2821 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2826 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2827 track->SetTrackPointArray(sp);
2829 for (Int_t iDet = 5; iDet >= 0; iDet--) {
2830 AliTracker *tracker = fTracker[iDet];
2831 if (!tracker) continue;
2832 Int_t nspdet = track->GetClusters(iDet,idx);
2834 if (iDet==0) // ITS "extra" clusters
2835 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2837 if (nspdet <= 0) continue;
2841 while (isp2 < nspdet) {
2842 Bool_t isvalid=kTRUE;
2844 Int_t index=idx[isp++];
2845 if (index < 0) continue;
2847 TString dets = fgkDetectorName[iDet];
2848 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2849 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2850 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2851 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2852 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2854 isvalid = tracker->GetTrackPoint(index,p);
2857 if (!isvalid) continue;
2858 if (iDet==0 && (isp-1)>=6) p.SetExtra();
2859 sp->AddPoint(isptrack,&p); isptrack++;
2866 //_____________________________________________________________________________
2867 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2869 // The method reads the raw-data error log
2870 // accumulated within the rawReader.
2871 // It extracts the raw-data errors related to
2872 // the current event and stores them into
2873 // a TClonesArray inside the esd object.
2875 if (!fRawReader) return;
2877 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2879 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2881 if (iEvent != log->GetEventNumber()) continue;
2883 esd->AddRawDataErrorLog(log);
2888 //_____________________________________________________________________________
2889 void AliReconstruction::CheckQA()
2891 // check the QA of SIM for this run and remove the detectors
2892 // with status Fatal
2894 // TString newRunLocalReconstruction ;
2895 // TString newRunTracking ;
2896 // TString newFillESD ;
2898 // for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
2899 // TString detName(AliQAv1::GetDetName(iDet)) ;
2900 // AliQAv1 * qa = AliQAv1::Instance(AliQAv1::DETECTORINDEX_t(iDet)) ;
2901 // if ( qa->IsSet(AliQAv1::DETECTORINDEX_t(iDet), AliQAv1::kSIM, specie, AliQAv1::kFATAL)) {
2902 // AliInfo(Form("QA status for %s %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed",
2903 // detName.Data(), AliRecoParam::GetEventSpecieName(es))) ;
2905 // if ( fRunLocalReconstruction.Contains(AliQAv1::GetDetName(iDet)) ||
2906 // fRunLocalReconstruction.Contains("ALL") ) {
2907 // newRunLocalReconstruction += detName ;
2908 // newRunLocalReconstruction += " " ;
2910 // if ( fRunTracking.Contains(AliQAv1::GetDetName(iDet)) ||
2911 // fRunTracking.Contains("ALL") ) {
2912 // newRunTracking += detName ;
2913 // newRunTracking += " " ;
2915 // if ( fFillESD.Contains(AliQAv1::GetDetName(iDet)) ||
2916 // fFillESD.Contains("ALL") ) {
2917 // newFillESD += detName ;
2918 // newFillESD += " " ;
2922 // fRunLocalReconstruction = newRunLocalReconstruction ;
2923 // fRunTracking = newRunTracking ;
2924 // fFillESD = newFillESD ;
2927 //_____________________________________________________________________________
2928 Int_t AliReconstruction::GetDetIndex(const char* detector)
2930 // return the detector index corresponding to detector
2932 for (index = 0; index < kNDetectors ; index++) {
2933 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2938 //_____________________________________________________________________________
2939 Bool_t AliReconstruction::FinishPlaneEff() {
2941 // Here execute all the necessary operationis, at the end of the tracking phase,
2942 // in case that evaluation of PlaneEfficiencies was required for some detector.
2943 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2945 // This Preliminary version works only FOR ITS !!!!!
2946 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2949 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2952 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2953 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2954 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2955 if(fTracker[iDet] && fTracker[iDet]->GetPlaneEff()) {
2956 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
2957 TString name=planeeff->GetName();
2959 TFile* pefile = TFile::Open(name, "RECREATE");
2960 ret=(Bool_t)planeeff->Write();
2962 if(planeeff->GetCreateHistos()) {
2963 TString hname=planeeff->GetName();
2964 hname+="Histo.root";
2965 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
2968 if(fSPDTrackleter) {
2969 AliPlaneEff *planeeff=fSPDTrackleter->GetPlaneEff();
2970 TString name="AliITSPlaneEffSPDtracklet.root";
2971 TFile* pefile = TFile::Open(name, "RECREATE");
2972 ret=(Bool_t)planeeff->Write();
2974 AliESDEvent *dummy=NULL;
2975 ret=(Bool_t)fSPDTrackleter->PostProcess(dummy); // take care of writing other files
2980 //_____________________________________________________________________________
2981 Bool_t AliReconstruction::InitPlaneEff() {
2983 // Here execute all the necessary operations, before of the tracking phase,
2984 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
2985 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
2986 // which should be updated/recalculated.
2988 // This Preliminary version will work only FOR ITS !!!!!
2989 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2992 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2994 AliWarning(Form("Implementation of this method not yet completed !! Method return kTRUE"));
2996 fSPDTrackleter = NULL;
2997 AliReconstructor* itsReconstructor = GetReconstructor(0);
2998 if (itsReconstructor) {
2999 fSPDTrackleter = itsReconstructor->CreateTrackleter(); // this is NULL unless required in RecoParam
3001 if (fSPDTrackleter) {
3002 AliInfo("Trackleter for SPD has been created");
3008 //_____________________________________________________________________________
3009 Bool_t AliReconstruction::InitAliEVE()
3011 // This method should be called only in case
3012 // AliReconstruction is run
3013 // within the alieve environment.
3014 // It will initialize AliEVE in a way
3015 // so that it can visualize event processed
3016 // by AliReconstruction.
3017 // The return flag shows whenever the
3018 // AliEVE initialization was successful or not.
3021 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
3022 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
3023 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
3025 gROOT->ProcessLine("if (!AliEveEventManager::GetMaster()){new AliEveEventManager();AliEveEventManager::GetMaster()->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(AliEveEventManager::GetMaster());};");
3026 gROOT->ProcessLine("alieve_online_init()");
3031 //_____________________________________________________________________________
3032 void AliReconstruction::RunAliEVE()
3034 // Runs AliEVE visualisation of
3035 // the current event.
3036 // Should be executed only after
3037 // successful initialization of AliEVE.
3039 AliInfo("Running AliEVE...");
3040 gROOT->ProcessLine(Form("AliEveEventManager::GetMaster()->SetEvent((AliRunLoader*)0x%lx,(AliRawReader*)0x%lx,(AliESDEvent*)0x%lx,(AliESDfriend*)0x%lx);",fRunLoader,fRawReader,fesd,fesdf));
3044 //_____________________________________________________________________________
3045 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
3047 // Allows to run QA for a selected set of detectors
3048 // and a selected set of tasks among RAWS, RECPOINTS and ESDS
3049 // all selected detectors run the same selected tasks
3051 if (!detAndAction.Contains(":")) {
3052 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
3056 Int_t colon = detAndAction.Index(":") ;
3057 fQADetectors = detAndAction(0, colon) ;
3058 if (fQADetectors.Contains("ALL") )
3059 fQADetectors = fFillESD ;
3060 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
3061 if (fQATasks.Contains("ALL") ) {
3062 fQATasks = Form("%d %d %d", AliQAv1::kRAWS, AliQAv1::kRECPOINTS, AliQAv1::kESDS) ;
3064 fQATasks.ToUpper() ;
3066 if ( fQATasks.Contains("RAW") )
3067 tempo = Form("%d ", AliQAv1::kRAWS) ;
3068 if ( fQATasks.Contains("RECPOINT") )
3069 tempo += Form("%d ", AliQAv1::kRECPOINTS) ;
3070 if ( fQATasks.Contains("ESD") )
3071 tempo += Form("%d ", AliQAv1::kESDS) ;
3073 if (fQATasks.IsNull()) {
3074 AliInfo("No QA requested\n") ;
3079 TString tempo(fQATasks) ;
3080 tempo.ReplaceAll(Form("%d", AliQAv1::kRAWS), AliQAv1::GetTaskName(AliQAv1::kRAWS)) ;
3081 tempo.ReplaceAll(Form("%d", AliQAv1::kRECPOINTS), AliQAv1::GetTaskName(AliQAv1::kRECPOINTS)) ;
3082 tempo.ReplaceAll(Form("%d", AliQAv1::kESDS), AliQAv1::GetTaskName(AliQAv1::kESDS)) ;
3083 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
3088 //_____________________________________________________________________________
3089 Bool_t AliReconstruction::InitRecoParams()
3091 // The method accesses OCDB and retrieves all
3092 // the available reco-param objects from there.
3094 Bool_t isOK = kTRUE;
3096 TString detStr = fLoadCDB;
3097 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3099 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3101 if (fRecoParam.GetDetRecoParamArray(iDet)) {
3102 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
3106 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
3108 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
3109 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3111 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
3115 TObject *recoParamObj = entry->GetObject();
3116 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3117 // The detector has a normal TobjArray of AliDetectorRecoParam objects
3118 // Registering them in AliRecoParam
3119 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
3121 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3122 // The detector has only onse set of reco parameters
3123 // Registering it in AliRecoParam
3124 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
3125 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3126 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3129 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
3133 AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
3137 if (AliDebugLevel() > 0) fRecoParam.Print();
3142 //_____________________________________________________________________________
3143 Bool_t AliReconstruction::GetEventInfo()
3145 // Fill the event info object
3147 AliCodeTimerAuto("")
3149 AliCentralTrigger *aCTP = NULL;
3151 fEventInfo.SetEventType(fRawReader->GetType());
3153 ULong64_t mask = fRawReader->GetClassMask();
3154 fEventInfo.SetTriggerMask(mask);
3155 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
3156 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
3158 aCTP = new AliCentralTrigger();
3159 TString configstr("");
3160 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
3161 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
3165 aCTP->SetClassMask(mask);
3166 aCTP->SetClusterMask(clmask);
3169 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
3171 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
3172 aCTP = fRunLoader->GetTrigger();
3173 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
3174 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
3177 AliWarning("No trigger can be loaded! The trigger information will not be used!");
3182 AliTriggerConfiguration *config = aCTP->GetConfiguration();
3184 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
3185 if (fRawReader) delete aCTP;
3189 UChar_t clustmask = 0;
3191 ULong64_t trmask = fEventInfo.GetTriggerMask();
3192 const TObjArray& classesArray = config->GetClasses();
3193 Int_t nclasses = classesArray.GetEntriesFast();
3194 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
3195 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
3197 Int_t trindex = TMath::Nint(TMath::Log2(trclass->GetMask()));
3198 fesd->SetTriggerClass(trclass->GetName(),trindex);
3199 if (fRawReader) fRawReader->LoadTriggerClass(trclass->GetName(),trindex);
3200 if (trmask & (1 << trindex)) {
3202 trclasses += trclass->GetName();
3204 clustmask |= trclass->GetCluster()->GetClusterMask();
3208 fEventInfo.SetTriggerClasses(trclasses);
3210 // Set the information in ESD
3211 fesd->SetTriggerMask(trmask);
3212 fesd->SetTriggerCluster(clustmask);
3214 if (!aCTP->CheckTriggeredDetectors()) {
3215 if (fRawReader) delete aCTP;
3219 if (fRawReader) delete aCTP;
3221 // We have to fill also the HLT decision here!!
3227 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3229 // Match the detector list found in the rec.C or the default 'ALL'
3230 // to the list found in the GRP (stored there by the shuttle PP which
3231 // gets the information from ECS)
3232 static TString resultList;
3233 TString detList = detectorList;
3237 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3238 if ((detectorMask >> iDet) & 0x1) {
3239 TString det = AliDAQ::OfflineModuleName(iDet);
3240 if ((detList.CompareTo("ALL") == 0) ||
3241 ((detList.BeginsWith("ALL ") ||
3242 detList.EndsWith(" ALL") ||
3243 detList.Contains(" ALL ")) &&
3244 !(detList.BeginsWith("-"+det+" ") ||
3245 detList.EndsWith(" -"+det) ||
3246 detList.Contains(" -"+det+" "))) ||
3247 (detList.CompareTo(det) == 0) ||
3248 detList.BeginsWith(det+" ") ||
3249 detList.EndsWith(" "+det) ||
3250 detList.Contains( " "+det+" " )) {
3251 if (!resultList.EndsWith(det + " ")) {
3260 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3261 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3262 if ((detList.CompareTo("ALL") == 0) ||
3263 ((detList.BeginsWith("ALL ") ||
3264 detList.EndsWith(" ALL") ||
3265 detList.Contains(" ALL ")) &&
3266 !(detList.BeginsWith("-"+hltDet+" ") ||
3267 detList.EndsWith(" -"+hltDet) ||
3268 detList.Contains(" -"+hltDet+" "))) ||
3269 (detList.CompareTo(hltDet) == 0) ||
3270 detList.BeginsWith(hltDet+" ") ||
3271 detList.EndsWith(" "+hltDet) ||
3272 detList.Contains( " "+hltDet+" " )) {
3273 resultList += hltDet;
3277 return resultList.Data();
3281 //______________________________________________________________________________
3282 void AliReconstruction::Abort(const char *method, EAbort what)
3284 // Abort processing. If what = kAbortProcess, the Process() loop will be
3285 // aborted. If what = kAbortFile, the current file in a chain will be
3286 // aborted and the processing will continue with the next file, if there
3287 // is no next file then Process() will be aborted. Abort() can also be
3288 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3289 // the SlaveTerminate() and Terminate() are always called. The abort flag
3290 // can be checked in these methods using GetAbort().
3292 // The method is overwritten in AliReconstruction for better handling of
3293 // reco specific errors
3295 if (!fStopOnError) return;
3299 TString whyMess = method;
3300 whyMess += " failed! Aborting...";
3302 AliError(whyMess.Data());
3305 TString mess = "Abort";
3306 if (fAbort == kAbortProcess)
3307 mess = "AbortProcess";
3308 else if (fAbort == kAbortFile)
3311 Info(mess, whyMess.Data());