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 fRunVertexFinder(kTRUE),
198 fRunVertexFinderTracks(kTRUE),
199 fRunHLTTracking(kFALSE),
200 fRunMuonTracking(kFALSE),
202 fRunCascadeFinder(kTRUE),
203 fStopOnError(kFALSE),
204 fWriteAlignmentData(kFALSE),
205 fWriteESDfriend(kFALSE),
206 fFillTriggerESD(kTRUE),
214 fRunLocalReconstruction("ALL"),
218 fUseTrackingErrorsForAlignment(""),
219 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),
257 fSameQACycle(kFALSE),
258 fInitQACalled(kFALSE),
259 fWriteQAExpertData(kTRUE),
260 fRunPlaneEff(kFALSE),
269 fIsNewRunLoader(kFALSE),
273 // create reconstruction object with default parameters
276 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
277 fReconstructor[iDet] = NULL;
278 fLoader[iDet] = NULL;
279 fTracker[iDet] = NULL;
281 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
282 fQACycles[iDet] = 999999 ;
283 fQAWriteExpert[iDet] = kFALSE ;
289 //_____________________________________________________________________________
290 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
292 fRunVertexFinder(rec.fRunVertexFinder),
293 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
294 fRunHLTTracking(rec.fRunHLTTracking),
295 fRunMuonTracking(rec.fRunMuonTracking),
296 fRunV0Finder(rec.fRunV0Finder),
297 fRunCascadeFinder(rec.fRunCascadeFinder),
298 fStopOnError(rec.fStopOnError),
299 fWriteAlignmentData(rec.fWriteAlignmentData),
300 fWriteESDfriend(rec.fWriteESDfriend),
301 fFillTriggerESD(rec.fFillTriggerESD),
303 fCleanESD(rec.fCleanESD),
304 fV0DCAmax(rec.fV0DCAmax),
305 fV0CsPmin(rec.fV0CsPmin),
309 fRunLocalReconstruction(rec.fRunLocalReconstruction),
310 fRunTracking(rec.fRunTracking),
311 fFillESD(rec.fFillESD),
312 fLoadCDB(rec.fLoadCDB),
313 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
314 fGAliceFileName(rec.fGAliceFileName),
315 fRawInput(rec.fRawInput),
316 fESDOutput(rec.fESDOutput),
317 fEquipIdMap(rec.fEquipIdMap),
318 fFirstEvent(rec.fFirstEvent),
319 fLastEvent(rec.fLastEvent),
320 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
322 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
323 fLoadAlignData(rec.fLoadAlignData),
324 fUseHLTData(rec.fUseHLTData),
330 fParentRawReader(NULL),
332 fRecoParam(rec.fRecoParam),
334 fSPDTrackleter(NULL),
336 fDiamondProfileSPD(rec.fDiamondProfileSPD),
337 fDiamondProfile(rec.fDiamondProfile),
338 fDiamondProfileTPC(rec.fDiamondProfileTPC),
342 fAlignObjArray(rec.fAlignObjArray),
343 fCDBUri(rec.fCDBUri),
344 fQARefUri(rec.fQARefUri),
346 fInitCDBCalled(rec.fInitCDBCalled),
347 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
348 fQADetectors(rec.fQADetectors),
349 fQATasks(rec.fQATasks),
351 fRunGlobalQA(rec.fRunGlobalQA),
352 fSameQACycle(rec.fSameQACycle),
353 fInitQACalled(rec.fInitQACalled),
354 fWriteQAExpertData(rec.fWriteQAExpertData),
355 fRunPlaneEff(rec.fRunPlaneEff),
364 fIsNewRunLoader(rec.fIsNewRunLoader),
370 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
371 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
373 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
374 fReconstructor[iDet] = NULL;
375 fLoader[iDet] = NULL;
376 fTracker[iDet] = NULL;
379 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
380 fQACycles[iDet] = rec.fQACycles[iDet];
381 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
384 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
385 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
389 //_____________________________________________________________________________
390 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
392 // assignment operator
393 // Used in PROOF mode
394 // Be very careful while modifing it!
395 // Simple rules to follow:
396 // for persistent data members - use their assignment operators
397 // for non-persistent ones - do nothing or take the default values from constructor
398 // TSelector members should not be touched
399 if(&rec == this) return *this;
401 fRunVertexFinder = rec.fRunVertexFinder;
402 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
403 fRunHLTTracking = rec.fRunHLTTracking;
404 fRunMuonTracking = rec.fRunMuonTracking;
405 fRunV0Finder = rec.fRunV0Finder;
406 fRunCascadeFinder = rec.fRunCascadeFinder;
407 fStopOnError = rec.fStopOnError;
408 fWriteAlignmentData = rec.fWriteAlignmentData;
409 fWriteESDfriend = rec.fWriteESDfriend;
410 fFillTriggerESD = rec.fFillTriggerESD;
412 fCleanESD = rec.fCleanESD;
413 fV0DCAmax = rec.fV0DCAmax;
414 fV0CsPmin = rec.fV0CsPmin;
418 fRunLocalReconstruction = rec.fRunLocalReconstruction;
419 fRunTracking = rec.fRunTracking;
420 fFillESD = rec.fFillESD;
421 fLoadCDB = rec.fLoadCDB;
422 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
423 fGAliceFileName = rec.fGAliceFileName;
424 fRawInput = rec.fRawInput;
425 fESDOutput = rec.fESDOutput;
426 fEquipIdMap = rec.fEquipIdMap;
427 fFirstEvent = rec.fFirstEvent;
428 fLastEvent = rec.fLastEvent;
429 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
431 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
432 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
435 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
436 fLoadAlignData = rec.fLoadAlignData;
437 fUseHLTData = rec.fUseHLTData;
439 delete fRunInfo; fRunInfo = NULL;
440 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
442 fEventInfo = rec.fEventInfo;
446 fParentRawReader = NULL;
448 fRecoParam = rec.fRecoParam;
450 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
451 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
452 delete fLoader[iDet]; fLoader[iDet] = NULL;
453 delete fTracker[iDet]; fTracker[iDet] = NULL;
456 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
457 fQACycles[iDet] = rec.fQACycles[iDet];
458 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
461 delete fSPDTrackleter; fSPDTrackleter = NULL;
463 delete fDiamondProfileSPD; fDiamondProfileSPD = NULL;
464 if (rec.fDiamondProfileSPD) fDiamondProfileSPD = new AliESDVertex(*rec.fDiamondProfileSPD);
465 delete fDiamondProfile; fDiamondProfile = NULL;
466 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
467 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
468 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
470 delete fGRPData; fGRPData = NULL;
471 // if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
472 if (rec.fGRPData) fGRPData = (AliGRPObject*)((rec.fGRPData)->Clone());
474 delete fAlignObjArray; fAlignObjArray = NULL;
477 fQARefUri = rec.fQARefUri;
478 fSpecCDBUri.Delete();
479 fInitCDBCalled = rec.fInitCDBCalled;
480 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
481 fQADetectors = rec.fQADetectors;
482 fQATasks = rec.fQATasks;
484 fRunGlobalQA = rec.fRunGlobalQA;
485 fSameQACycle = rec.fSameQACycle;
486 fInitQACalled = rec.fInitQACalled;
487 fWriteQAExpertData = rec.fWriteQAExpertData;
488 fRunPlaneEff = rec.fRunPlaneEff;
497 fIsNewRunLoader = rec.fIsNewRunLoader;
504 //_____________________________________________________________________________
505 AliReconstruction::~AliReconstruction()
512 if (fAlignObjArray) {
513 fAlignObjArray->Delete();
514 delete fAlignObjArray;
516 fSpecCDBUri.Delete();
518 AliCodeTimer::Instance()->Print();
521 //_____________________________________________________________________________
522 void AliReconstruction::InitQA()
524 //Initialize the QA and start of cycle
525 AliCodeTimerAuto("");
527 if (fInitQACalled) return;
528 fInitQACalled = kTRUE;
530 AliQAManager * qam = AliQAManager::QAManager("rec") ;
531 if (fWriteQAExpertData)
532 qam->SetWriteExpert() ;
534 if (qam->IsDefaultStorageSet()) {
535 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
536 AliWarning("Default QA reference storage has been already set !");
537 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fQARefUri.Data()));
538 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
539 fQARefUri = qam->GetDefaultStorage()->GetURI();
541 if (fQARefUri.Length() > 0) {
542 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
543 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
544 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
546 fQARefUri="local://$ALICE_ROOT/QAref";
547 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
548 AliWarning("Default QA refeference storage not yet set !!!!");
549 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
550 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
553 qam->SetDefaultStorage(fQARefUri);
557 qam->SetActiveDetectors(fQADetectors) ;
558 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++) {
559 qam->SetCycleLength(AliQAv1::DETECTORINDEX_t(det), fQACycles[det]) ;
560 qam->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
562 if (!fRawReader && !fInput && fQATasks.Contains(AliQAv1::kRAWS))
563 fQATasks.ReplaceAll(Form("%d",AliQAv1::kRAWS), "") ;
564 qam->SetTasks(fQATasks) ;
565 qam->InitQADataMaker(AliCDBManager::Instance()->GetRun()) ;
568 Bool_t sameCycle = kFALSE ;
569 AliQADataMaker *qadm = qam->GetQADataMaker(AliQAv1::kGLOBAL);
570 AliInfo(Form("Initializing the global QA data maker"));
571 if (fQATasks.Contains(Form("%d", AliQAv1::kRECPOINTS))) {
572 qadm->StartOfCycle(AliQAv1::kRECPOINTS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
573 TObjArray **arr=qadm->Init(AliQAv1::kRECPOINTS);
574 AliTracker::SetResidualsArray(arr);
577 if (fQATasks.Contains(Form("%d", AliQAv1::kESDS))) {
578 qadm->StartOfCycle(AliQAv1::kESDS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
579 qadm->Init(AliQAv1::kESDS);
582 AliSysInfo::AddStamp("InitQA") ;
585 //_____________________________________________________________________________
586 void AliReconstruction::MergeQA(const char *fileName)
588 //Initialize the QA and start of cycle
589 AliCodeTimerAuto("") ;
590 AliQAManager::QAManager()->Merge(AliCDBManager::Instance()->GetRun(),fileName) ;
591 AliSysInfo::AddStamp("MergeQA") ;
594 //_____________________________________________________________________________
595 void AliReconstruction::InitCDB()
597 // activate a default CDB storage
598 // First check if we have any CDB storage set, because it is used
599 // to retrieve the calibration and alignment constants
600 AliCodeTimerAuto("");
602 if (fInitCDBCalled) return;
603 fInitCDBCalled = kTRUE;
605 AliCDBManager* man = AliCDBManager::Instance();
606 if (man->IsDefaultStorageSet())
608 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
609 AliWarning("Default CDB storage has been already set !");
610 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
611 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
612 fCDBUri = man->GetDefaultStorage()->GetURI();
615 if (fCDBUri.Length() > 0)
617 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
618 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
619 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
621 fCDBUri="local://$ALICE_ROOT/OCDB";
622 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
623 AliWarning("Default CDB storage not yet set !!!!");
624 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
625 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
628 man->SetDefaultStorage(fCDBUri);
631 // Now activate the detector specific CDB storage locations
632 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
633 TObject* obj = fSpecCDBUri[i];
635 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
636 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
637 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
638 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
640 AliSysInfo::AddStamp("InitCDB");
643 //_____________________________________________________________________________
644 void AliReconstruction::SetDefaultStorage(const char* uri) {
645 // Store the desired default CDB storage location
646 // Activate it later within the Run() method
652 //_____________________________________________________________________________
653 void AliReconstruction::SetQARefDefaultStorage(const char* uri) {
654 // Store the desired default CDB storage location
655 // Activate it later within the Run() method
658 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
661 //_____________________________________________________________________________
662 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
663 // Store a detector-specific CDB storage location
664 // Activate it later within the Run() method
666 AliCDBPath aPath(calibType);
667 if(!aPath.IsValid()){
668 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
669 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
670 if(!strcmp(calibType, fgkDetectorName[iDet])) {
671 aPath.SetPath(Form("%s/*", calibType));
672 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
676 if(!aPath.IsValid()){
677 AliError(Form("Not a valid path or detector: %s", calibType));
682 // // check that calibType refers to a "valid" detector name
683 // Bool_t isDetector = kFALSE;
684 // for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
685 // TString detName = fgkDetectorName[iDet];
686 // if(aPath.GetLevel0() == detName) {
687 // isDetector = kTRUE;
693 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
697 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
698 if (obj) fSpecCDBUri.Remove(obj);
699 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
703 //_____________________________________________________________________________
704 Bool_t AliReconstruction::SetRunNumberFromData()
706 // The method is called in Run() in order
707 // to set a correct run number.
708 // In case of raw data reconstruction the
709 // run number is taken from the raw data header
711 if (fSetRunNumberFromDataCalled) return kTRUE;
712 fSetRunNumberFromDataCalled = kTRUE;
714 AliCDBManager* man = AliCDBManager::Instance();
717 if(fRawReader->NextEvent()) {
718 if(man->GetRun() > 0) {
719 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
721 man->SetRun(fRawReader->GetRunNumber());
722 fRawReader->RewindEvents();
725 if(man->GetRun() > 0) {
726 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
729 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
735 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
737 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
742 // read run number from gAlice
743 if(rl->GetHeader()) {
744 man->SetRun(rl->GetHeader()->GetRun());
749 AliError("Neither run-loader header nor RawReader objects are found !");
761 //_____________________________________________________________________________
762 void AliReconstruction::SetCDBLock() {
763 // Set CDB lock: from now on it is forbidden to reset the run number
764 // or the default storage or to activate any further storage!
766 AliCDBManager::Instance()->SetLock(1);
769 //_____________________________________________________________________________
770 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
772 // Read the alignment objects from CDB.
773 // Each detector is supposed to have the
774 // alignment objects in DET/Align/Data CDB path.
775 // All the detector objects are then collected,
776 // sorted by geometry level (starting from ALIC) and
777 // then applied to the TGeo geometry.
778 // Finally an overlaps check is performed.
780 // Load alignment data from CDB and fill fAlignObjArray
781 if(fLoadAlignFromCDB){
783 TString detStr = detectors;
784 TString loadAlObjsListOfDets = "";
786 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
787 if(!IsSelected(fgkDetectorName[iDet], detStr)) continue;
788 if(!strcmp(fgkDetectorName[iDet],"HLT")) continue;
790 if(AliGeomManager::GetNalignable(fgkDetectorName[iDet]) != 0)
792 loadAlObjsListOfDets += fgkDetectorName[iDet];
793 loadAlObjsListOfDets += " ";
795 } // end loop over detectors
797 if(AliGeomManager::GetNalignable("GRP") != 0)
798 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
799 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
800 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
802 // Check if the array with alignment objects was
803 // provided by the user. If yes, apply the objects
804 // to the present TGeo geometry
805 if (fAlignObjArray) {
806 if (gGeoManager && gGeoManager->IsClosed()) {
807 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
808 AliError("The misalignment of one or more volumes failed!"
809 "Compare the list of simulated detectors and the list of detector alignment data!");
814 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
820 if (fAlignObjArray) {
821 fAlignObjArray->Delete();
822 delete fAlignObjArray; fAlignObjArray=NULL;
828 //_____________________________________________________________________________
829 void AliReconstruction::SetGAliceFile(const char* fileName)
831 // set the name of the galice file
833 fGAliceFileName = fileName;
836 //_____________________________________________________________________________
837 void AliReconstruction::SetInput(const char* input)
839 // In case the input string starts with 'mem://', we run in an online mode
840 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
841 // file is assumed. One can give as an input:
842 // mem://: - events taken from DAQ monitoring libs online
844 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
845 if (input) fRawInput = input;
848 //_____________________________________________________________________________
849 void AliReconstruction::SetOutput(const char* output)
851 // Set the output ESD filename
852 // 'output' is a normalt ROOT url
853 // The method is used in case of raw-data reco with PROOF
854 if (output) fESDOutput.SetUrl(output);
857 //_____________________________________________________________________________
858 void AliReconstruction::SetOption(const char* detector, const char* option)
860 // set options for the reconstruction of a detector
862 TObject* obj = fOptions.FindObject(detector);
863 if (obj) fOptions.Remove(obj);
864 fOptions.Add(new TNamed(detector, option));
867 //_____________________________________________________________________________
868 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
870 // Set custom reconstruction parameters for a given detector
871 // Single set of parameters for all the events
873 // First check if the reco-params are global
874 if(!strcmp(detector, "GRP")) {
876 fRecoParam.AddDetRecoParam(kNDetectors,par);
880 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
881 if(!strcmp(detector, fgkDetectorName[iDet])) {
883 fRecoParam.AddDetRecoParam(iDet,par);
890 //_____________________________________________________________________________
891 Bool_t AliReconstruction::SetFieldMap(Float_t l3Cur, Float_t diCur, Float_t l3Pol,
892 Float_t diPol, Float_t beamenergy,
893 const Char_t *beamtype, const Char_t *path)
895 //------------------------------------------------
896 // The magnetic field map, defined externally...
897 // L3 current 30000 A -> 0.5 T
898 // L3 current 12000 A -> 0.2 T
899 // dipole current 6000 A
900 // The polarities must be the same
901 //------------------------------------------------
902 const Float_t l3NominalCurrent1=30000.; // (A)
903 const Float_t l3NominalCurrent2=12000.; // (A)
904 const Float_t diNominalCurrent =6000. ; // (A)
906 const Float_t tolerance=0.03; // relative current tolerance
907 const Float_t zero=77.; // "zero" current (A)
909 TString s=(l3Pol < 0) ? "L3: -" : "L3: +";
911 AliMagF::BMap_t map = AliMagF::k5kG;
915 l3Cur = TMath::Abs(l3Cur);
916 if (TMath::Abs(l3Cur-l3NominalCurrent1)/l3NominalCurrent1 < tolerance) {
917 fcL3 = l3Cur/l3NominalCurrent1;
920 } else if (TMath::Abs(l3Cur-l3NominalCurrent2)/l3NominalCurrent2 < tolerance) {
921 fcL3 = l3Cur/l3NominalCurrent2;
924 } else if (l3Cur <= zero) {
926 map = AliMagF::k5kGUniform;
929 AliError(Form("Wrong L3 current (%f A)!",l3Cur));
933 diCur = TMath::Abs(diCur);
934 if (TMath::Abs(diCur-diNominalCurrent)/diNominalCurrent < tolerance) {
935 // 3% current tolerance...
936 fcDip = diCur/diNominalCurrent;
938 } else if (diCur <= zero) { // some small current..
942 AliError(Form("Wrong dipole current (%f A)!",diCur));
946 if (l3Pol!=diPol && (map==AliMagF::k5kG || map==AliMagF::k2kG) && fcDip!=0) {
947 AliError("L3 and Dipole polarities must be the same");
951 if (l3Pol<0) fcL3 = -fcL3;
952 if (diPol<0) fcDip = -fcDip;
954 AliMagF::BeamType_t btype = AliMagF::kNoBeamField;
955 TString btypestr = beamtype;
957 TPRegexp protonBeam("(proton|p)\\s*-?\\s*\\1");
958 TPRegexp ionBeam("(lead|pb|ion|a)\\s*-?\\s*\\1");
959 if (btypestr.Contains(ionBeam)) btype = AliMagF::kBeamTypeAA;
960 else if (btypestr.Contains(protonBeam)) btype = AliMagF::kBeamTypepp;
962 AliInfo(Form("Cannot determine the beam type from %s, assume no LHC magnet field",beamtype));
965 AliMagF* fld = new AliMagF("MagneticFieldMap", s.Data(), 2, fcL3, fcDip, 10., map, path,
967 TGeoGlobalMagField::Instance()->SetField( fld );
968 TGeoGlobalMagField::Instance()->Lock();
974 Bool_t AliReconstruction::InitGRP() {
975 //------------------------------------
976 // Initialization of the GRP entry
977 //------------------------------------
978 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
982 TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
985 AliInfo("Found a TMap in GRP/GRP/Data, converting it into an AliGRPObject");
987 fGRPData = new AliGRPObject();
988 fGRPData->ReadValuesFromMap(m);
992 AliInfo("Found an AliGRPObject in GRP/GRP/Data, reading it");
993 fGRPData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
997 AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
1001 AliError("No GRP entry found in OCDB!");
1005 TString lhcState = fGRPData->GetLHCState();
1006 if (lhcState==AliGRPObject::GetInvalidString()) {
1007 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
1008 lhcState = "UNKNOWN";
1011 TString beamType = fGRPData->GetBeamType();
1012 if (beamType==AliGRPObject::GetInvalidString()) {
1013 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
1014 beamType = "UNKNOWN";
1017 Float_t beamEnergy = fGRPData->GetBeamEnergy();
1018 if (beamEnergy==AliGRPObject::GetInvalidFloat()) {
1019 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
1022 // energy is provided in MeV*120
1023 beamEnergy /= 120E3;
1025 TString runType = fGRPData->GetRunType();
1026 if (runType==AliGRPObject::GetInvalidString()) {
1027 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
1028 runType = "UNKNOWN";
1031 Int_t activeDetectors = fGRPData->GetDetectorMask();
1032 if (activeDetectors==AliGRPObject::GetInvalidUInt()) {
1033 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
1034 activeDetectors = 1074790399;
1037 fRunInfo = new AliRunInfo(lhcState, beamType, beamEnergy, runType, activeDetectors);
1041 // Process the list of active detectors
1042 if (activeDetectors) {
1043 UInt_t detMask = activeDetectors;
1044 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
1045 fRunTracking = MatchDetectorList(fRunTracking,detMask);
1046 fFillESD = MatchDetectorList(fFillESD,detMask);
1047 fQADetectors = MatchDetectorList(fQADetectors,detMask);
1048 fLoadCDB.Form("%s %s %s %s",
1049 fRunLocalReconstruction.Data(),
1050 fRunTracking.Data(),
1052 fQADetectors.Data());
1053 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
1054 if (!((detMask >> AliDAQ::DetectorID("ITSSPD")) & 0x1)) {
1055 // switch off the vertexer
1056 AliInfo("SPD is not in the list of active detectors. Vertexer switched off.");
1057 fRunVertexFinder = kFALSE;
1059 if (!((detMask >> AliDAQ::DetectorID("TRG")) & 0x1)) {
1060 // switch off the reading of CTP raw-data payload
1061 if (fFillTriggerESD) {
1062 AliInfo("CTP is not in the list of active detectors. CTP data reading switched off.");
1063 fFillTriggerESD = kFALSE;
1068 AliInfo("===================================================================================");
1069 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
1070 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
1071 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
1072 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
1073 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
1074 AliInfo("===================================================================================");
1076 //*** Dealing with the magnetic field map
1077 if ( TGeoGlobalMagField::Instance()->IsLocked() ) {AliInfo("Running with the externally locked B field !");}
1079 // Construct the field map out of the information retrieved from GRP.
1082 Float_t l3Current = fGRPData->GetL3Current((AliGRPObject::Stats)0);
1083 if (l3Current == AliGRPObject::GetInvalidFloat()) {
1084 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1088 Char_t l3Polarity = fGRPData->GetL3Polarity();
1089 if (l3Polarity == AliGRPObject::GetInvalidChar()) {
1090 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1095 Float_t diCurrent = fGRPData->GetDipoleCurrent((AliGRPObject::Stats)0);
1096 if (diCurrent == AliGRPObject::GetInvalidFloat()) {
1097 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1101 Char_t diPolarity = fGRPData->GetDipolePolarity();
1102 if (diPolarity == AliGRPObject::GetInvalidChar()) {
1103 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1108 TObjString *l3Current=
1109 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Current"));
1111 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1114 TObjString *l3Polarity=
1115 dynamic_cast<TObjString*>(fGRPData->GetValue("fL3Polarity"));
1117 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1122 TObjString *diCurrent=
1123 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipoleCurrent"));
1125 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1128 TObjString *diPolarity=
1129 dynamic_cast<TObjString*>(fGRPData->GetValue("fDipolePolarity"));
1131 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1137 if ( !SetFieldMap(l3Current, diCurrent, l3Polarity ? -1:1, diPolarity ? -1:1) )
1138 AliFatal("Failed to creat a B field map ! Exiting...");
1139 AliInfo("Running with the B field constructed out of GRP !");
1141 else AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
1145 //*** Get the diamond profiles from OCDB
1146 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1148 fDiamondProfileSPD = dynamic_cast<AliESDVertex*> (entry->GetObject());
1150 AliError("No SPD diamond profile found in OCDB!");
1153 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1155 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
1157 AliError("No diamond profile found in OCDB!");
1160 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
1162 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
1164 AliError("No TPC diamond profile found in OCDB!");
1170 //_____________________________________________________________________________
1171 Bool_t AliReconstruction::LoadCDB()
1173 AliCodeTimerAuto("");
1175 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1177 TString detStr = fLoadCDB;
1178 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1179 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1180 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1185 //_____________________________________________________________________________
1186 Bool_t AliReconstruction::Run(const char* input)
1189 AliCodeTimerAuto("");
1192 if (GetAbort() != TSelector::kContinue) return kFALSE;
1194 TChain *chain = NULL;
1195 if (fRawReader && (chain = fRawReader->GetChain())) {
1198 gProof->AddInput(this);
1199 if (!fESDOutput.IsValid()) {
1200 fESDOutput.SetProtocol("root",kTRUE);
1201 fESDOutput.SetHost(gSystem->HostName());
1202 fESDOutput.SetFile(Form("%s/AliESDs.root",gSystem->pwd()));
1204 AliInfo(Form("Output file with ESDs is %s",fESDOutput.GetUrl()));
1205 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",fESDOutput.GetUrl()));
1206 gProof->SetParameter("PROOF_MaxSlavesPerNode", 9999);
1208 chain->Process("AliReconstruction");
1211 chain->Process(this);
1216 if (GetAbort() != TSelector::kContinue) return kFALSE;
1218 if (GetAbort() != TSelector::kContinue) return kFALSE;
1219 //******* The loop over events
1220 AliInfo("Starting looping over events");
1222 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1223 (fRawReader && fRawReader->NextEvent())) {
1224 if (!ProcessEvent(iEvent)) {
1225 Abort("ProcessEvent",TSelector::kAbortFile);
1231 if (GetAbort() != TSelector::kContinue) return kFALSE;
1233 if (GetAbort() != TSelector::kContinue) return kFALSE;
1239 //_____________________________________________________________________________
1240 void AliReconstruction::InitRawReader(const char* input)
1242 AliCodeTimerAuto("");
1244 // Init raw-reader and
1245 // set the input in case of raw data
1246 if (input) fRawInput = input;
1247 fRawReader = AliRawReader::Create(fRawInput.Data());
1249 AliInfo("Reconstruction will run over digits");
1251 if (!fEquipIdMap.IsNull() && fRawReader)
1252 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1254 if (!fUseHLTData.IsNull()) {
1255 // create the RawReaderHLT which performs redirection of HLT input data for
1256 // the specified detectors
1257 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1259 fParentRawReader=fRawReader;
1260 fRawReader=pRawReader;
1262 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1265 AliSysInfo::AddStamp("CreateRawReader");
1268 //_____________________________________________________________________________
1269 void AliReconstruction::InitRun(const char* input)
1271 // Initialization of raw-reader,
1272 // run number, CDB etc.
1273 AliCodeTimerAuto("");
1274 AliSysInfo::AddStamp("Start");
1276 // Initialize raw-reader if any
1277 InitRawReader(input);
1279 // Initialize the CDB storage
1282 // Set run number in CDBManager (if it is not already set by the user)
1283 if (!SetRunNumberFromData()) {
1284 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1288 // Set CDB lock: from now on it is forbidden to reset the run number
1289 // or the default storage or to activate any further storage!
1294 //_____________________________________________________________________________
1295 void AliReconstruction::Begin(TTree *)
1297 // Initialize AlReconstruction before
1298 // going into the event loop
1299 // Should follow the TSelector convention
1300 // i.e. initialize only the object on the client side
1301 AliCodeTimerAuto("");
1303 AliReconstruction *reco = NULL;
1305 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1308 AliSysInfo::AddStamp("ReadInputInBegin");
1311 // Import ideal TGeo geometry and apply misalignment
1313 TString geom(gSystem->DirName(fGAliceFileName));
1314 geom += "/geometry.root";
1315 AliGeomManager::LoadGeometry(geom.Data());
1317 Abort("LoadGeometry", TSelector::kAbortProcess);
1320 AliSysInfo::AddStamp("LoadGeom");
1321 TString detsToCheck=fRunLocalReconstruction;
1322 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1323 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1326 AliSysInfo::AddStamp("CheckGeom");
1329 if (!MisalignGeometry(fLoadAlignData)) {
1330 Abort("MisalignGeometry", TSelector::kAbortProcess);
1333 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1334 AliSysInfo::AddStamp("MisalignGeom");
1337 Abort("InitGRP", TSelector::kAbortProcess);
1340 AliSysInfo::AddStamp("InitGRP");
1343 Abort("LoadCDB", TSelector::kAbortProcess);
1346 AliSysInfo::AddStamp("LoadCDB");
1348 // Read the reconstruction parameters from OCDB
1349 if (!InitRecoParams()) {
1350 AliWarning("Not all detectors have correct RecoParam objects initialized");
1352 AliSysInfo::AddStamp("InitRecoParams");
1354 if (fInput && gProof) {
1355 if (reco) *reco = *this;
1357 gProof->AddInputData(gGeoManager,kTRUE);
1359 gProof->AddInputData(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()),kTRUE);
1360 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1361 AliMagF *magFieldMap = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1362 magFieldMap->SetName("MagneticFieldMap");
1363 gProof->AddInputData(magFieldMap,kTRUE);
1368 //_____________________________________________________________________________
1369 void AliReconstruction::SlaveBegin(TTree*)
1371 // Initialization related to run-loader,
1372 // vertexer, trackers, recontructors
1373 // In proof mode it is executed on the slave
1374 AliCodeTimerAuto("");
1376 TProofOutputFile *outProofFile = NULL;
1378 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1381 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1383 AliGeomManager::SetGeometry(tgeo);
1385 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1386 Int_t runNumber = -1;
1387 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1388 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1389 man->SetCacheFlag(kTRUE);
1390 man->SetLock(kTRUE);
1394 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1395 TGeoGlobalMagField::Instance()->SetField(map);
1397 if (TNamed *outputFileName = (TNamed *) fInput->FindObject("PROOF_OUTPUTFILE")) {
1398 outProofFile = new TProofOutputFile(gSystem->BaseName(TUrl(outputFileName->GetTitle()).GetFile()));
1399 outProofFile->SetOutputFileName(outputFileName->GetTitle());
1400 fOutput->Add(outProofFile);
1402 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1405 // get the run loader
1406 if (!InitRunLoader()) {
1407 Abort("InitRunLoader", TSelector::kAbortProcess);
1410 AliSysInfo::AddStamp("LoadLoader");
1412 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1415 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1416 Abort("CreateTrackers", TSelector::kAbortProcess);
1419 AliSysInfo::AddStamp("CreateTrackers");
1421 // create the ESD output file and tree
1422 if (!outProofFile) {
1423 ffile = TFile::Open("AliESDs.root", "RECREATE");
1424 ffile->SetCompressionLevel(2);
1425 if (!ffile->IsOpen()) {
1426 Abort("OpenESDFile", TSelector::kAbortProcess);
1431 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1432 Abort(Form("Problems opening output PROOF file: %s/%s",
1433 outProofFile->GetDir(), outProofFile->GetFileName()),
1434 TSelector::kAbortProcess);
1439 ftree = new TTree("esdTree", "Tree with ESD objects");
1440 fesd = new AliESDEvent();
1441 fesd->CreateStdContent();
1443 fesd->WriteToTree(ftree);
1444 if (fWriteESDfriend) {
1446 // Since we add the branch manually we must
1447 // book and add it after WriteToTree
1448 // otherwise it is created twice,
1449 // once via writetotree and once here.
1450 // The case for AliESDfriend is now
1451 // caught also in AlIESDEvent::WriteToTree but
1452 // be careful when changing the name (AliESDfriend is not
1453 // a TNamed so we had to hardwire it)
1454 fesdf = new AliESDfriend();
1455 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1456 br->SetFile("AliESDfriends.root");
1457 fesd->AddObject(fesdf);
1459 ftree->GetUserInfo()->Add(fesd);
1461 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1462 fhltesd = new AliESDEvent();
1463 fhltesd->CreateStdContent();
1465 // read the ESD template from CDB
1466 // HLT is allowed to put non-std content to its ESD, the non-std
1467 // objects need to be created before invocation of WriteToTree in
1468 // order to create all branches. Initialization is done from an
1469 // ESD layout template in CDB
1470 AliCDBManager* man = AliCDBManager::Instance();
1471 AliCDBPath hltESDConfigPath("HLT/ConfigHLT/esdLayout");
1472 AliCDBEntry* hltESDConfig=NULL;
1473 if (man->GetId(hltESDConfigPath)!=NULL &&
1474 (hltESDConfig=man->Get(hltESDConfigPath))!=NULL) {
1475 AliESDEvent* pESDLayout=dynamic_cast<AliESDEvent*>(hltESDConfig->GetObject());
1477 // init all internal variables from the list of objects
1478 pESDLayout->GetStdContent();
1480 // copy content and create non-std objects
1481 *fhltesd=*pESDLayout;
1484 AliError(Form("error setting hltEsd layout from %s: invalid object type",
1485 hltESDConfigPath.GetPath().Data()));
1489 fhltesd->WriteToTree(fhlttree);
1490 fhlttree->GetUserInfo()->Add(fhltesd);
1492 ProcInfo_t procInfo;
1493 gSystem->GetProcInfo(&procInfo);
1494 AliInfo(Form("Current memory usage %d %d", procInfo.fMemResident, procInfo.fMemVirtual));
1497 //Initialize the QA and start of cycle
1498 if (fRunQA || fRunGlobalQA)
1501 //Initialize the Plane Efficiency framework
1502 if (fRunPlaneEff && !InitPlaneEff()) {
1503 Abort("InitPlaneEff", TSelector::kAbortProcess);
1507 if (strcmp(gProgName,"alieve") == 0)
1508 fRunAliEVE = InitAliEVE();
1513 //_____________________________________________________________________________
1514 Bool_t AliReconstruction::Process(Long64_t entry)
1516 // run the reconstruction over a single entry
1517 // from the chain with raw data
1518 AliCodeTimerAuto("");
1520 TTree *currTree = fChain->GetTree();
1521 AliRawVEvent *event = NULL;
1522 currTree->SetBranchAddress("rawevent",&event);
1523 currTree->GetEntry(entry);
1524 fRawReader = new AliRawReaderRoot(event);
1525 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1533 //_____________________________________________________________________________
1534 void AliReconstruction::Init(TTree *tree)
1537 AliError("The input tree is not found!");
1543 //_____________________________________________________________________________
1544 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1546 // run the reconstruction over a single event
1547 // The event loop is steered in Run method
1549 AliCodeTimerAuto("");
1551 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1552 fRunLoader->SetEventNumber(iEvent);
1553 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1555 fRunLoader->TreeE()->Fill();
1556 if (fRawReader && fRawReader->UseAutoSaveESD())
1557 fRunLoader->TreeE()->AutoSave("SaveSelf");
1560 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1564 AliInfo(Form("processing event %d", iEvent));
1566 fRunLoader->GetEvent(iEvent);
1568 // Fill Event-info object
1570 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo);
1571 AliInfo(Form("Current event specie: %s",fRecoParam.PrintEventSpecie()));
1573 // Set the reco-params
1575 TString detStr = fLoadCDB;
1576 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1577 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1578 AliReconstructor *reconstructor = GetReconstructor(iDet);
1579 if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
1580 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
1581 reconstructor->SetRecoParam(par);
1583 AliQAManager::QAManager()->SetRecoParam(iDet, par) ;
1591 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1592 AliQAManager::QAManager()->RunOneEvent(fRawReader) ;
1594 // local single event reconstruction
1595 if (!fRunLocalReconstruction.IsNull()) {
1596 TString detectors=fRunLocalReconstruction;
1597 // run HLT event reconstruction first
1598 // ;-( IsSelected changes the string
1599 if (IsSelected("HLT", detectors) &&
1600 !RunLocalEventReconstruction("HLT")) {
1601 if (fStopOnError) {CleanUp(); return kFALSE;}
1603 detectors=fRunLocalReconstruction;
1604 detectors.ReplaceAll("HLT", "");
1605 if (!RunLocalEventReconstruction(detectors)) {
1606 if (fStopOnError) {CleanUp(); return kFALSE;}
1610 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1611 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1612 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1613 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1615 // Set magnetic field from the tracker
1616 fesd->SetMagneticField(AliTracker::GetBz());
1617 fhltesd->SetMagneticField(AliTracker::GetBz());
1619 // Set most probable pt, for B=0 tracking
1620 // Get the global reco-params. They are atposition 16 inside the array of detectors in fRecoParam
1621 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(kNDetectors));
1622 if (grpRecoParam) AliExternalTrackParam::SetMostProbablePt(grpRecoParam->GetMostProbablePt());
1624 // Fill raw-data error log into the ESD
1625 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1628 if (fRunVertexFinder) {
1629 if (!RunVertexFinder(fesd)) {
1630 if (fStopOnError) {CleanUp(); return kFALSE;}
1634 // For Plane Efficiency: run the SPD trackleter
1635 if (fRunPlaneEff && fSPDTrackleter) {
1636 if (!RunSPDTrackleting(fesd)) {
1637 if (fStopOnError) {CleanUp(); return kFALSE;}
1642 if (!fRunTracking.IsNull()) {
1643 if (fRunMuonTracking) {
1644 if (!RunMuonTracking(fesd)) {
1645 if (fStopOnError) {CleanUp(); return kFALSE;}
1651 if (!fRunTracking.IsNull()) {
1652 if (!RunTracking(fesd)) {
1653 if (fStopOnError) {CleanUp(); return kFALSE;}
1658 if (!fFillESD.IsNull()) {
1659 TString detectors=fFillESD;
1660 // run HLT first and on hltesd
1661 // ;-( IsSelected changes the string
1662 if (IsSelected("HLT", detectors) &&
1663 !FillESD(fhltesd, "HLT")) {
1664 if (fStopOnError) {CleanUp(); return kFALSE;}
1667 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1668 if (detectors.Contains("ALL")) {
1670 for (Int_t idet=0; idet<kNDetectors; ++idet){
1671 detectors += fgkDetectorName[idet];
1675 detectors.ReplaceAll("HLT", "");
1676 if (!FillESD(fesd, detectors)) {
1677 if (fStopOnError) {CleanUp(); return kFALSE;}
1681 // fill Event header information from the RawEventHeader
1682 if (fRawReader){FillRawEventHeaderESD(fesd);}
1685 AliESDpid::MakePID(fesd);
1687 if (fFillTriggerESD) {
1688 if (!FillTriggerESD(fesd)) {
1689 if (fStopOnError) {CleanUp(); return kFALSE;}
1696 // Propagate track to the beam pipe (if not already done by ITS)
1698 const Int_t ntracks = fesd->GetNumberOfTracks();
1699 const Double_t kBz = fesd->GetMagneticField();
1700 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1703 UShort_t *selectedIdx=new UShort_t[ntracks];
1705 for (Int_t itrack=0; itrack<ntracks; itrack++){
1706 const Double_t kMaxStep = 1; //max step over the material
1709 AliESDtrack *track = fesd->GetTrack(itrack);
1710 if (!track) continue;
1712 AliExternalTrackParam *tpcTrack =
1713 (AliExternalTrackParam *)track->GetTPCInnerParam();
1717 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kFALSE);
1720 Int_t n=trkArray.GetEntriesFast();
1721 selectedIdx[n]=track->GetID();
1722 trkArray.AddLast(tpcTrack);
1725 //Tracks refitted by ITS should already be at the SPD vertex
1726 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1729 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kFALSE);
1730 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1735 // Improve the reconstructed primary vertex position using the tracks
1737 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
1738 if(fesd->GetPrimaryVertexSPD()) {
1739 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
1740 if(vtitle.Contains("cosmics")) {
1741 runVertexFinderTracks=kFALSE;
1745 if (runVertexFinderTracks) {
1746 // TPC + ITS primary vertex
1747 ftVertexer->SetITSMode();
1748 ftVertexer->SetConstraintOff();
1749 // get cuts for vertexer from AliGRPRecoParam
1751 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1752 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1753 grpRecoParam->GetVertexerTracksCutsITS(cutsVertexer);
1754 ftVertexer->SetCuts(cutsVertexer);
1755 delete [] cutsVertexer; cutsVertexer = NULL;
1756 if(fDiamondProfile && grpRecoParam->GetVertexerTracksConstraintITS())
1757 ftVertexer->SetVtxStart(fDiamondProfile);
1759 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1761 if (pvtx->GetStatus()) {
1762 fesd->SetPrimaryVertexTracks(pvtx);
1763 for (Int_t i=0; i<ntracks; i++) {
1764 AliESDtrack *t = fesd->GetTrack(i);
1765 t->RelateToVertex(pvtx, kBz, kVeryBig);
1770 // TPC-only primary vertex
1771 ftVertexer->SetTPCMode();
1772 ftVertexer->SetConstraintOff();
1773 // get cuts for vertexer from AliGRPRecoParam
1775 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1776 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1777 grpRecoParam->GetVertexerTracksCutsTPC(cutsVertexer);
1778 ftVertexer->SetCuts(cutsVertexer);
1779 delete [] cutsVertexer; cutsVertexer = NULL;
1780 if(fDiamondProfileTPC && grpRecoParam->GetVertexerTracksConstraintTPC())
1781 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1783 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1785 if (pvtx->GetStatus()) {
1786 fesd->SetPrimaryVertexTPC(pvtx);
1787 for (Int_t i=0; i<ntracks; i++) {
1788 AliESDtrack *t = fesd->GetTrack(i);
1789 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1795 delete[] selectedIdx;
1797 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1802 AliV0vertexer vtxer;
1803 vtxer.Tracks2V0vertices(fesd);
1805 if (fRunCascadeFinder) {
1807 AliCascadeVertexer cvtxer;
1808 cvtxer.V0sTracks2CascadeVertices(fesd);
1813 if (fCleanESD) CleanESD(fesd);
1816 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1817 AliQAManager::QAManager()->RunOneEvent(fesd) ;
1820 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
1821 qadm->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1822 if (qadm && fQATasks.Contains(Form("%d", AliQAv1::kESDS)))
1823 qadm->Exec(AliQAv1::kESDS, fesd);
1826 if (fWriteESDfriend) {
1827 // fesdf->~AliESDfriend();
1828 // new (fesdf) AliESDfriend(); // Reset...
1829 fesd->GetESDfriend(fesdf);
1833 // Auto-save the ESD tree in case of prompt reco @P2
1834 if (fRawReader && fRawReader->UseAutoSaveESD()) {
1835 ftree->AutoSave("SaveSelf");
1836 TFile *friendfile = (TFile *)(gROOT->GetListOfFiles()->FindObject("AliESDfriends.root"));
1837 if (friendfile) friendfile->Save();
1844 if (fRunAliEVE) RunAliEVE();
1848 if (fWriteESDfriend) {
1849 fesdf->~AliESDfriend();
1850 new (fesdf) AliESDfriend(); // Reset...
1853 ProcInfo_t procInfo;
1854 gSystem->GetProcInfo(&procInfo);
1855 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, procInfo.fMemResident, procInfo.fMemVirtual));
1858 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1859 if (fReconstructor[iDet])
1860 fReconstructor[iDet]->SetRecoParam(NULL);
1863 if (fRunQA || fRunGlobalQA)
1864 AliQAManager::QAManager()->Increment() ;
1869 //_____________________________________________________________________________
1870 void AliReconstruction::SlaveTerminate()
1872 // Finalize the run on the slave side
1873 // Called after the exit
1874 // from the event loop
1875 AliCodeTimerAuto("");
1877 if (fIsNewRunLoader) { // galice.root didn't exist
1878 fRunLoader->WriteHeader("OVERWRITE");
1879 fRunLoader->CdGAFile();
1880 fRunLoader->Write(0, TObject::kOverwrite);
1883 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1884 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1886 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1887 cdbMapCopy->SetOwner(1);
1888 cdbMapCopy->SetName("cdbMap");
1889 TIter iter(cdbMap->GetTable());
1892 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1893 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1894 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1895 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1898 TList *cdbListCopy = new TList();
1899 cdbListCopy->SetOwner(1);
1900 cdbListCopy->SetName("cdbList");
1902 TIter iter2(cdbList);
1905 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1906 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1909 ftree->GetUserInfo()->Add(cdbMapCopy);
1910 ftree->GetUserInfo()->Add(cdbListCopy);
1915 if (fWriteESDfriend)
1916 ftree->SetBranchStatus("ESDfriend*",0);
1917 // we want to have only one tree version number
1918 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1921 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1922 if (fRunPlaneEff && !FinishPlaneEff()) {
1923 AliWarning("Finish PlaneEff evaluation failed");
1926 // End of cycle for the in-loop
1928 AliQAManager::QAManager()->EndOfCycle() ;
1931 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
1933 if (fQATasks.Contains(Form("%d", AliQAv1::kRECPOINTS)))
1934 qadm->EndOfCycle(AliQAv1::kRECPOINTS);
1935 if (fQATasks.Contains(Form("%d", AliQAv1::kESDS)))
1936 qadm->EndOfCycle(AliQAv1::kESDS);
1941 if (fRunQA || fRunGlobalQA) {
1943 if (TNamed *outputFileName = (TNamed *) fInput->FindObject("PROOF_OUTPUTFILE")) {
1944 TString qaOutputFile = outputFileName->GetTitle();
1945 qaOutputFile.ReplaceAll(gSystem->BaseName(TUrl(outputFileName->GetTitle()).GetFile()),
1946 Form("Merged.%s.Data.root",AliQAv1::GetQADataFileName()));
1947 TProofOutputFile *qaProofFile = new TProofOutputFile(Form("Merged.%s.Data.root",AliQAv1::GetQADataFileName()));
1948 qaProofFile->SetOutputFileName(qaOutputFile.Data());
1949 fOutput->Add(qaProofFile);
1950 MergeQA(qaProofFile->GetFileName());
1962 //_____________________________________________________________________________
1963 void AliReconstruction::Terminate()
1965 // Create tags for the events in the ESD tree (the ESD tree is always present)
1966 // In case of empty events the tags will contain dummy values
1967 AliCodeTimerAuto("");
1969 // Do not call the ESD tag creator in case of PROOF-based reconstruction
1971 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1972 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData, AliQAv1::Instance()->GetQA(), AliQAv1::Instance()->GetEventSpecies(), AliQAv1::kNDET, AliRecoParam::kNSpecies);
1975 // Cleanup of CDB manager: cache and active storages!
1976 AliCDBManager::Instance()->ClearCache();
1979 //_____________________________________________________________________________
1980 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1982 // run the local reconstruction
1984 static Int_t eventNr=0;
1985 AliCodeTimerAuto("")
1987 TString detStr = detectors;
1988 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1989 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1990 AliReconstructor* reconstructor = GetReconstructor(iDet);
1991 if (!reconstructor) continue;
1992 AliLoader* loader = fLoader[iDet];
1993 // Matthias April 2008: temporary fix to run HLT reconstruction
1994 // although the HLT loader is missing
1995 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1997 reconstructor->Reconstruct(fRawReader, NULL);
2000 reconstructor->Reconstruct(dummy, NULL);
2005 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
2008 // conversion of digits
2009 if (fRawReader && reconstructor->HasDigitConversion()) {
2010 AliInfo(Form("converting raw data digits into root objects for %s",
2011 fgkDetectorName[iDet]));
2012 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
2013 // fgkDetectorName[iDet]));
2014 loader->LoadDigits("update");
2015 loader->CleanDigits();
2016 loader->MakeDigitsContainer();
2017 TTree* digitsTree = loader->TreeD();
2018 reconstructor->ConvertDigits(fRawReader, digitsTree);
2019 loader->WriteDigits("OVERWRITE");
2020 loader->UnloadDigits();
2022 // local reconstruction
2023 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2024 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2025 loader->LoadRecPoints("update");
2026 loader->CleanRecPoints();
2027 loader->MakeRecPointsContainer();
2028 TTree* clustersTree = loader->TreeR();
2029 if (fRawReader && !reconstructor->HasDigitConversion()) {
2030 reconstructor->Reconstruct(fRawReader, clustersTree);
2032 loader->LoadDigits("read");
2033 TTree* digitsTree = loader->TreeD();
2035 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
2036 if (fStopOnError) return kFALSE;
2038 reconstructor->Reconstruct(digitsTree, clustersTree);
2040 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2041 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, digitsTree) ;
2044 loader->UnloadDigits();
2047 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2048 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, clustersTree) ;
2050 loader->WriteRecPoints("OVERWRITE");
2051 loader->UnloadRecPoints();
2052 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
2054 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2055 AliError(Form("the following detectors were not found: %s",
2057 if (fStopOnError) return kFALSE;
2062 //_____________________________________________________________________________
2063 Bool_t AliReconstruction::RunSPDTrackleting(AliESDEvent*& esd)
2065 // run the SPD trackleting (for SPD efficiency purpouses)
2067 AliCodeTimerAuto("")
2069 Double_t vtxPos[3] = {0, 0, 0};
2070 Double_t vtxErr[3] = {0.0, 0.0, 0.0};
2072 TArrayF mcVertex(3);
2074 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
2075 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
2076 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
2079 const AliESDVertex *vertex = esd->GetVertex();
2081 AliWarning("Vertex not found");
2084 vertex->GetXYZ(vtxPos);
2085 vertex->GetSigmaXYZ(vtxErr);
2086 if (fSPDTrackleter) {
2087 AliInfo("running the SPD Trackleter for Plane Efficiency Evaluation");
2090 fLoader[0]->LoadRecPoints("read");
2091 TTree* tree = fLoader[0]->TreeR();
2093 AliError("Can't get the ITS cluster tree");
2096 fSPDTrackleter->LoadClusters(tree);
2097 fSPDTrackleter->SetVertex(vtxPos, vtxErr);
2099 if (fSPDTrackleter->Clusters2Tracks(esd) != 0) {
2100 AliError("AliITSTrackleterSPDEff Clusters2Tracks failed");
2101 // fLoader[0]->UnloadRecPoints();
2104 //fSPDTrackleter->UnloadRecPoints();
2106 AliWarning("SPDTrackleter not available");
2112 //_____________________________________________________________________________
2113 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
2115 // run the barrel tracking
2117 AliCodeTimerAuto("")
2119 AliVertexer *vertexer = CreateVertexer();
2120 if (!vertexer) return kFALSE;
2122 AliInfo("running the ITS vertex finder");
2123 AliESDVertex* vertex = NULL;
2125 fLoader[0]->LoadRecPoints();
2126 TTree* cltree = fLoader[0]->TreeR();
2128 if(fDiamondProfileSPD) vertexer->SetVtxStart(fDiamondProfileSPD);
2129 vertex = vertexer->FindVertexForCurrentEvent(cltree);
2132 AliError("Can't get the ITS cluster tree");
2134 fLoader[0]->UnloadRecPoints();
2137 AliError("Can't get the ITS loader");
2140 AliWarning("Vertex not found");
2141 vertex = new AliESDVertex();
2142 vertex->SetName("default");
2145 vertex->SetName("reconstructed");
2150 vertex->GetXYZ(vtxPos);
2151 vertex->GetSigmaXYZ(vtxErr);
2153 esd->SetPrimaryVertexSPD(vertex);
2154 AliESDVertex *vpileup = NULL;
2155 Int_t novertices = 0;
2156 vpileup = vertexer->GetAllVertices(novertices);
2158 for (Int_t kk=1; kk<novertices; kk++)esd->AddPileupVertexSPD(&vpileup[kk]);
2160 // if SPD multiplicity has been determined, it is stored in the ESD
2161 AliMultiplicity *mult = vertexer->GetMultiplicity();
2162 if(mult)esd->SetMultiplicity(mult);
2164 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2165 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
2174 //_____________________________________________________________________________
2175 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
2177 // run the HLT barrel tracking
2179 AliCodeTimerAuto("")
2182 AliError("Missing runLoader!");
2186 AliInfo("running HLT tracking");
2188 // Get a pointer to the HLT reconstructor
2189 AliReconstructor *reconstructor = GetReconstructor(kNDetectors-1);
2190 if (!reconstructor) return kFALSE;
2193 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2194 TString detName = fgkDetectorName[iDet];
2195 AliDebug(1, Form("%s HLT tracking", detName.Data()));
2196 reconstructor->SetOption(detName.Data());
2197 AliTracker *tracker = reconstructor->CreateTracker();
2199 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
2200 if (fStopOnError) return kFALSE;
2204 Double_t vtxErr[3]={0.005,0.005,0.010};
2205 const AliESDVertex *vertex = esd->GetVertex();
2206 vertex->GetXYZ(vtxPos);
2207 tracker->SetVertex(vtxPos,vtxErr);
2209 fLoader[iDet]->LoadRecPoints("read");
2210 TTree* tree = fLoader[iDet]->TreeR();
2212 AliError(Form("Can't get the %s cluster tree", detName.Data()));
2215 tracker->LoadClusters(tree);
2217 if (tracker->Clusters2Tracks(esd) != 0) {
2218 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
2222 tracker->UnloadClusters();
2230 //_____________________________________________________________________________
2231 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
2233 // run the muon spectrometer tracking
2235 AliCodeTimerAuto("")
2238 AliError("Missing runLoader!");
2241 Int_t iDet = 7; // for MUON
2243 AliInfo("is running...");
2245 // Get a pointer to the MUON reconstructor
2246 AliReconstructor *reconstructor = GetReconstructor(iDet);
2247 if (!reconstructor) return kFALSE;
2250 TString detName = fgkDetectorName[iDet];
2251 AliDebug(1, Form("%s tracking", detName.Data()));
2252 AliTracker *tracker = reconstructor->CreateTracker();
2254 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2259 fLoader[iDet]->LoadRecPoints("read");
2261 tracker->LoadClusters(fLoader[iDet]->TreeR());
2263 Int_t rv = tracker->Clusters2Tracks(esd);
2267 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2271 fLoader[iDet]->UnloadRecPoints();
2273 tracker->UnloadClusters();
2281 //_____________________________________________________________________________
2282 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2284 // run the barrel tracking
2285 static Int_t eventNr=0;
2286 AliCodeTimerAuto("")
2288 AliInfo("running tracking");
2291 //Fill the ESD with the T0 info (will be used by the TOF)
2292 if (fReconstructor[11] && fLoader[11]) {
2293 fLoader[11]->LoadRecPoints("READ");
2294 TTree *treeR = fLoader[11]->TreeR();
2296 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2300 // pass 1: TPC + ITS inwards
2301 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2302 if (!fTracker[iDet]) continue;
2303 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2306 fLoader[iDet]->LoadRecPoints("read");
2307 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2308 TTree* tree = fLoader[iDet]->TreeR();
2310 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2313 fTracker[iDet]->LoadClusters(tree);
2314 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2316 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2317 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2320 // preliminary PID in TPC needed by the ITS tracker
2322 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2323 AliESDpid::MakePID(esd);
2325 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2328 // pass 2: ALL backwards
2330 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2331 if (!fTracker[iDet]) continue;
2332 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2335 if (iDet > 1) { // all except ITS, TPC
2337 fLoader[iDet]->LoadRecPoints("read");
2338 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2339 tree = fLoader[iDet]->TreeR();
2341 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2344 fTracker[iDet]->LoadClusters(tree);
2345 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2349 if (iDet>1) // start filling residuals for the "outer" detectors
2351 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2352 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2354 TObjArray * elem = arr[fRecoParam.GetEventSpecie()];
2355 if ( elem && (! elem->At(0)) ) {
2356 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2357 if (qadm) qadm->InitRecPointsForTracker() ;
2361 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2362 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2367 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2368 fTracker[iDet]->UnloadClusters();
2369 fLoader[iDet]->UnloadRecPoints();
2371 // updated PID in TPC needed by the ITS tracker -MI
2373 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2374 AliESDpid::MakePID(esd);
2376 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2378 //stop filling residuals for the "outer" detectors
2379 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2381 // pass 3: TRD + TPC + ITS refit inwards
2383 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2384 if (!fTracker[iDet]) continue;
2385 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2388 if (iDet<2) // start filling residuals for TPC and ITS
2390 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2391 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2393 TObjArray * elem = arr[fRecoParam.GetEventSpecie()];
2394 if ( elem && (! elem->At(0)) ) {
2395 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2396 if (qadm) qadm->InitRecPointsForTracker() ;
2401 if (fTracker[iDet]->RefitInward(esd) != 0) {
2402 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2405 // run postprocessing
2406 if (fTracker[iDet]->PostProcess(esd) != 0) {
2407 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2410 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2413 // write space-points to the ESD in case alignment data output
2415 if (fWriteAlignmentData)
2416 WriteAlignmentData(esd);
2418 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2419 if (!fTracker[iDet]) continue;
2421 fTracker[iDet]->UnloadClusters();
2422 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2423 fLoader[iDet]->UnloadRecPoints();
2424 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2426 // stop filling residuals for TPC and ITS
2427 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2433 //_____________________________________________________________________________
2434 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2436 // Remove the data which are not needed for the physics analysis.
2439 Int_t nTracks=esd->GetNumberOfTracks();
2440 Int_t nV0s=esd->GetNumberOfV0s();
2442 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2444 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2445 Bool_t rc=esd->Clean(cleanPars);
2447 nTracks=esd->GetNumberOfTracks();
2448 nV0s=esd->GetNumberOfV0s();
2450 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2455 //_____________________________________________________________________________
2456 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2458 // fill the event summary data
2460 AliCodeTimerAuto("")
2461 static Int_t eventNr=0;
2462 TString detStr = detectors;
2464 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2465 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2466 AliReconstructor* reconstructor = GetReconstructor(iDet);
2467 if (!reconstructor) continue;
2468 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2469 TTree* clustersTree = NULL;
2470 if (fLoader[iDet]) {
2471 fLoader[iDet]->LoadRecPoints("read");
2472 clustersTree = fLoader[iDet]->TreeR();
2473 if (!clustersTree) {
2474 AliError(Form("Can't get the %s clusters tree",
2475 fgkDetectorName[iDet]));
2476 if (fStopOnError) return kFALSE;
2479 if (fRawReader && !reconstructor->HasDigitConversion()) {
2480 reconstructor->FillESD(fRawReader, clustersTree, esd);
2482 TTree* digitsTree = NULL;
2483 if (fLoader[iDet]) {
2484 fLoader[iDet]->LoadDigits("read");
2485 digitsTree = fLoader[iDet]->TreeD();
2487 AliError(Form("Can't get the %s digits tree",
2488 fgkDetectorName[iDet]));
2489 if (fStopOnError) return kFALSE;
2492 reconstructor->FillESD(digitsTree, clustersTree, esd);
2493 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2495 if (fLoader[iDet]) {
2496 fLoader[iDet]->UnloadRecPoints();
2500 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2501 AliError(Form("the following detectors were not found: %s",
2503 if (fStopOnError) return kFALSE;
2505 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2510 //_____________________________________________________________________________
2511 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2513 // Reads the trigger decision which is
2514 // stored in Trigger.root file and fills
2515 // the corresponding esd entries
2517 AliCodeTimerAuto("")
2519 AliInfo("Filling trigger information into the ESD");
2522 AliCTPRawStream input(fRawReader);
2523 if (!input.Next()) {
2524 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2527 if (esd->GetTriggerMask() != input.GetClassMask())
2528 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2529 input.GetClassMask(),esd->GetTriggerMask()));
2530 if (esd->GetOrbitNumber() != input.GetOrbitID())
2531 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2532 input.GetOrbitID(),esd->GetOrbitNumber()));
2533 if (esd->GetBunchCrossNumber() != input.GetBCID())
2534 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2535 input.GetBCID(),esd->GetBunchCrossNumber()));
2536 AliESDHeader* esdheader = esd->GetHeader();
2537 esdheader->SetL0TriggerInputs(input.GetL0Inputs());
2538 esdheader->SetL1TriggerInputs(input.GetL1Inputs());
2539 esdheader->SetL2TriggerInputs(input.GetL2Inputs());
2542 // Here one has to add the filling interaction records
2551 //_____________________________________________________________________________
2552 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2555 // Filling information from RawReader Header
2558 if (!fRawReader) return kFALSE;
2560 AliInfo("Filling information from RawReader Header");
2562 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2563 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2564 esd->SetPeriodNumber(fRawReader->GetPeriod());
2566 esd->SetTimeStamp(fRawReader->GetTimestamp());
2567 esd->SetEventType(fRawReader->GetType());
2573 //_____________________________________________________________________________
2574 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2576 // check whether detName is contained in detectors
2577 // if yes, it is removed from detectors
2579 // check if all detectors are selected
2580 if ((detectors.CompareTo("ALL") == 0) ||
2581 detectors.BeginsWith("ALL ") ||
2582 detectors.EndsWith(" ALL") ||
2583 detectors.Contains(" ALL ")) {
2588 // search for the given detector
2589 Bool_t result = kFALSE;
2590 if ((detectors.CompareTo(detName) == 0) ||
2591 detectors.BeginsWith(detName+" ") ||
2592 detectors.EndsWith(" "+detName) ||
2593 detectors.Contains(" "+detName+" ")) {
2594 detectors.ReplaceAll(detName, "");
2598 // clean up the detectors string
2599 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2600 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2601 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2606 //_____________________________________________________________________________
2607 Bool_t AliReconstruction::InitRunLoader()
2609 // get or create the run loader
2611 if (gAlice) delete gAlice;
2614 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2615 // load all base libraries to get the loader classes
2616 TString libs = gSystem->GetLibraries();
2617 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2618 TString detName = fgkDetectorName[iDet];
2619 if (detName == "HLT") continue;
2620 if (libs.Contains("lib" + detName + "base.so")) continue;
2621 gSystem->Load("lib" + detName + "base.so");
2623 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2625 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2630 fRunLoader->CdGAFile();
2631 fRunLoader->LoadgAlice();
2633 //PH This is a temporary fix to give access to the kinematics
2634 //PH that is needed for the labels of ITS clusters
2635 fRunLoader->LoadHeader();
2636 fRunLoader->LoadKinematics();
2638 } else { // galice.root does not exist
2640 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2642 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2643 AliConfig::GetDefaultEventFolderName(),
2646 AliError(Form("could not create run loader in file %s",
2647 fGAliceFileName.Data()));
2651 fIsNewRunLoader = kTRUE;
2652 fRunLoader->MakeTree("E");
2654 if (fNumberOfEventsPerFile > 0)
2655 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2657 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2663 //_____________________________________________________________________________
2664 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2666 // get the reconstructor object and the loader for a detector
2668 if (fReconstructor[iDet]) {
2669 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2670 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2671 fReconstructor[iDet]->SetRecoParam(par);
2673 return fReconstructor[iDet];
2676 // load the reconstructor object
2677 TPluginManager* pluginManager = gROOT->GetPluginManager();
2678 TString detName = fgkDetectorName[iDet];
2679 TString recName = "Ali" + detName + "Reconstructor";
2681 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2683 AliReconstructor* reconstructor = NULL;
2684 // first check if a plugin is defined for the reconstructor
2685 TPluginHandler* pluginHandler =
2686 pluginManager->FindHandler("AliReconstructor", detName);
2687 // if not, add a plugin for it
2688 if (!pluginHandler) {
2689 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2690 TString libs = gSystem->GetLibraries();
2691 if (libs.Contains("lib" + detName + "base.so") ||
2692 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2693 pluginManager->AddHandler("AliReconstructor", detName,
2694 recName, detName + "rec", recName + "()");
2696 pluginManager->AddHandler("AliReconstructor", detName,
2697 recName, detName, recName + "()");
2699 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2701 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2702 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2704 if (reconstructor) {
2705 TObject* obj = fOptions.FindObject(detName.Data());
2706 if (obj) reconstructor->SetOption(obj->GetTitle());
2707 reconstructor->Init();
2708 fReconstructor[iDet] = reconstructor;
2711 // get or create the loader
2712 if (detName != "HLT") {
2713 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2714 if (!fLoader[iDet]) {
2715 AliConfig::Instance()
2716 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2718 // first check if a plugin is defined for the loader
2720 pluginManager->FindHandler("AliLoader", detName);
2721 // if not, add a plugin for it
2722 if (!pluginHandler) {
2723 TString loaderName = "Ali" + detName + "Loader";
2724 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2725 pluginManager->AddHandler("AliLoader", detName,
2726 loaderName, detName + "base",
2727 loaderName + "(const char*, TFolder*)");
2728 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2730 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2732 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2733 fRunLoader->GetEventFolder());
2735 if (!fLoader[iDet]) { // use default loader
2736 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2738 if (!fLoader[iDet]) {
2739 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2740 if (fStopOnError) return NULL;
2742 fRunLoader->AddLoader(fLoader[iDet]);
2743 fRunLoader->CdGAFile();
2744 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2745 fRunLoader->Write(0, TObject::kOverwrite);
2750 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2751 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2752 reconstructor->SetRecoParam(par);
2754 return reconstructor;
2757 //_____________________________________________________________________________
2758 AliVertexer* AliReconstruction::CreateVertexer()
2760 // create the vertexer
2761 // Please note that the caller is the owner of the
2764 AliVertexer* vertexer = NULL;
2765 AliReconstructor* itsReconstructor = GetReconstructor(0);
2766 if (itsReconstructor && (fRunLocalReconstruction.Contains("ITS"))) {
2767 vertexer = itsReconstructor->CreateVertexer();
2770 AliWarning("couldn't create a vertexer for ITS");
2776 //_____________________________________________________________________________
2777 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2779 // create the trackers
2780 AliInfo("Creating trackers");
2782 TString detStr = detectors;
2783 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2784 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2785 AliReconstructor* reconstructor = GetReconstructor(iDet);
2786 if (!reconstructor) continue;
2787 TString detName = fgkDetectorName[iDet];
2788 if (detName == "HLT") {
2789 fRunHLTTracking = kTRUE;
2792 if (detName == "MUON") {
2793 fRunMuonTracking = kTRUE;
2798 fTracker[iDet] = reconstructor->CreateTracker();
2799 if (!fTracker[iDet] && (iDet < 7)) {
2800 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2801 if (fStopOnError) return kFALSE;
2803 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2809 //_____________________________________________________________________________
2810 void AliReconstruction::CleanUp()
2812 // delete trackers and the run loader and close and delete the file
2814 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2815 delete fReconstructor[iDet];
2816 fReconstructor[iDet] = NULL;
2817 fLoader[iDet] = NULL;
2818 delete fTracker[iDet];
2819 fTracker[iDet] = NULL;
2824 delete fSPDTrackleter;
2825 fSPDTrackleter = NULL;
2834 delete fParentRawReader;
2835 fParentRawReader=NULL;
2843 AliQAManager::Destroy() ;
2845 TGeoGlobalMagField::Instance()->SetField(NULL);
2848 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2850 // Write space-points which are then used in the alignment procedures
2851 // For the moment only ITS, TPC, TRD and TOF
2853 Int_t ntracks = esd->GetNumberOfTracks();
2854 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2856 AliESDtrack *track = esd->GetTrack(itrack);
2859 for (Int_t i=0; i<200; ++i) idx[i] = -1; //PH avoid uninitialized values
2860 for (Int_t iDet = 5; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2861 nsp += track->GetNcls(iDet);
2863 if (iDet==0) { // ITS "extra" clusters
2864 track->GetClusters(iDet,idx);
2865 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2870 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2871 track->SetTrackPointArray(sp);
2873 for (Int_t iDet = 5; iDet >= 0; iDet--) {
2874 AliTracker *tracker = fTracker[iDet];
2875 if (!tracker) continue;
2876 Int_t nspdet = track->GetClusters(iDet,idx);
2878 if (iDet==0) // ITS "extra" clusters
2879 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2881 if (nspdet <= 0) continue;
2885 while (isp2 < nspdet) {
2886 Bool_t isvalid=kTRUE;
2888 Int_t index=idx[isp++];
2889 if (index < 0) continue;
2891 TString dets = fgkDetectorName[iDet];
2892 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2893 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2894 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2895 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2896 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2898 isvalid = tracker->GetTrackPoint(index,p);
2901 if (!isvalid) continue;
2902 if (iDet==0 && (isp-1)>=6) p.SetExtra();
2903 sp->AddPoint(isptrack,&p); isptrack++;
2910 //_____________________________________________________________________________
2911 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2913 // The method reads the raw-data error log
2914 // accumulated within the rawReader.
2915 // It extracts the raw-data errors related to
2916 // the current event and stores them into
2917 // a TClonesArray inside the esd object.
2919 if (!fRawReader) return;
2921 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2923 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2925 if (iEvent != log->GetEventNumber()) continue;
2927 esd->AddRawDataErrorLog(log);
2932 //_____________________________________________________________________________
2933 void AliReconstruction::CheckQA()
2935 // check the QA of SIM for this run and remove the detectors
2936 // with status Fatal
2938 // TString newRunLocalReconstruction ;
2939 // TString newRunTracking ;
2940 // TString newFillESD ;
2942 // for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
2943 // TString detName(AliQAv1::GetDetName(iDet)) ;
2944 // AliQAv1 * qa = AliQAv1::Instance(AliQAv1::DETECTORINDEX_t(iDet)) ;
2945 // if ( qa->IsSet(AliQAv1::DETECTORINDEX_t(iDet), AliQAv1::kSIM, specie, AliQAv1::kFATAL)) {
2946 // AliInfo(Form("QA status for %s %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed",
2947 // detName.Data(), AliRecoParam::GetEventSpecieName(es))) ;
2949 // if ( fRunLocalReconstruction.Contains(AliQAv1::GetDetName(iDet)) ||
2950 // fRunLocalReconstruction.Contains("ALL") ) {
2951 // newRunLocalReconstruction += detName ;
2952 // newRunLocalReconstruction += " " ;
2954 // if ( fRunTracking.Contains(AliQAv1::GetDetName(iDet)) ||
2955 // fRunTracking.Contains("ALL") ) {
2956 // newRunTracking += detName ;
2957 // newRunTracking += " " ;
2959 // if ( fFillESD.Contains(AliQAv1::GetDetName(iDet)) ||
2960 // fFillESD.Contains("ALL") ) {
2961 // newFillESD += detName ;
2962 // newFillESD += " " ;
2966 // fRunLocalReconstruction = newRunLocalReconstruction ;
2967 // fRunTracking = newRunTracking ;
2968 // fFillESD = newFillESD ;
2971 //_____________________________________________________________________________
2972 Int_t AliReconstruction::GetDetIndex(const char* detector)
2974 // return the detector index corresponding to detector
2976 for (index = 0; index < kNDetectors ; index++) {
2977 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2982 //_____________________________________________________________________________
2983 Bool_t AliReconstruction::FinishPlaneEff() {
2985 // Here execute all the necessary operationis, at the end of the tracking phase,
2986 // in case that evaluation of PlaneEfficiencies was required for some detector.
2987 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2989 // This Preliminary version works only FOR ITS !!!!!
2990 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2993 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2996 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2997 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2998 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2999 if(fTracker[iDet] && fTracker[iDet]->GetPlaneEff()) {
3000 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
3001 TString name=planeeff->GetName();
3003 TFile* pefile = TFile::Open(name, "RECREATE");
3004 ret=(Bool_t)planeeff->Write();
3006 if(planeeff->GetCreateHistos()) {
3007 TString hname=planeeff->GetName();
3008 hname+="Histo.root";
3009 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
3012 if(fSPDTrackleter) {
3013 AliPlaneEff *planeeff=fSPDTrackleter->GetPlaneEff();
3014 TString name="AliITSPlaneEffSPDtracklet.root";
3015 TFile* pefile = TFile::Open(name, "RECREATE");
3016 ret=(Bool_t)planeeff->Write();
3018 AliESDEvent *dummy=NULL;
3019 ret=(Bool_t)fSPDTrackleter->PostProcess(dummy); // take care of writing other files
3024 //_____________________________________________________________________________
3025 Bool_t AliReconstruction::InitPlaneEff() {
3027 // Here execute all the necessary operations, before of the tracking phase,
3028 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
3029 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
3030 // which should be updated/recalculated.
3032 // This Preliminary version will work only FOR ITS !!!!!
3033 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3036 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3038 AliWarning(Form("Implementation of this method not yet completed !! Method return kTRUE"));
3040 fSPDTrackleter = NULL;
3041 AliReconstructor* itsReconstructor = GetReconstructor(0);
3042 if (itsReconstructor) {
3043 fSPDTrackleter = itsReconstructor->CreateTrackleter(); // this is NULL unless required in RecoParam
3045 if (fSPDTrackleter) {
3046 AliInfo("Trackleter for SPD has been created");
3052 //_____________________________________________________________________________
3053 Bool_t AliReconstruction::InitAliEVE()
3055 // This method should be called only in case
3056 // AliReconstruction is run
3057 // within the alieve environment.
3058 // It will initialize AliEVE in a way
3059 // so that it can visualize event processed
3060 // by AliReconstruction.
3061 // The return flag shows whenever the
3062 // AliEVE initialization was successful or not.
3065 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
3066 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
3067 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
3069 gROOT->ProcessLine("if (!AliEveEventManager::GetMaster()){new AliEveEventManager();AliEveEventManager::GetMaster()->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(AliEveEventManager::GetMaster());};");
3070 gROOT->ProcessLine("alieve_online_init()");
3075 //_____________________________________________________________________________
3076 void AliReconstruction::RunAliEVE()
3078 // Runs AliEVE visualisation of
3079 // the current event.
3080 // Should be executed only after
3081 // successful initialization of AliEVE.
3083 AliInfo("Running AliEVE...");
3084 gROOT->ProcessLine(Form("AliEveEventManager::GetMaster()->SetEvent((AliRunLoader*)0x%lx,(AliRawReader*)0x%lx,(AliESDEvent*)0x%lx,(AliESDfriend*)0x%lx);",fRunLoader,fRawReader,fesd,fesdf));
3088 //_____________________________________________________________________________
3089 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
3091 // Allows to run QA for a selected set of detectors
3092 // and a selected set of tasks among RAWS, DIGITSR, RECPOINTS and ESDS
3093 // all selected detectors run the same selected tasks
3095 if (!detAndAction.Contains(":")) {
3096 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
3100 Int_t colon = detAndAction.Index(":") ;
3101 fQADetectors = detAndAction(0, colon) ;
3102 if (fQADetectors.Contains("ALL") )
3103 fQADetectors = fFillESD ;
3104 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
3105 if (fQATasks.Contains("ALL") ) {
3106 fQATasks = Form("%d %d %d %d", AliQAv1::kRAWS, AliQAv1::kDIGITSR, AliQAv1::kRECPOINTS, AliQAv1::kESDS) ;
3108 fQATasks.ToUpper() ;
3110 if ( fQATasks.Contains("RAW") )
3111 tempo = Form("%d ", AliQAv1::kRAWS) ;
3112 if ( fQATasks.Contains("DIGIT") )
3113 tempo += Form("%d ", AliQAv1::kDIGITSR) ;
3114 if ( fQATasks.Contains("RECPOINT") )
3115 tempo += Form("%d ", AliQAv1::kRECPOINTS) ;
3116 if ( fQATasks.Contains("ESD") )
3117 tempo += Form("%d ", AliQAv1::kESDS) ;
3119 if (fQATasks.IsNull()) {
3120 AliInfo("No QA requested\n") ;
3125 TString tempo(fQATasks) ;
3126 tempo.ReplaceAll(Form("%d", AliQAv1::kRAWS), AliQAv1::GetTaskName(AliQAv1::kRAWS)) ;
3127 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITSR), AliQAv1::GetTaskName(AliQAv1::kDIGITSR)) ;
3128 tempo.ReplaceAll(Form("%d", AliQAv1::kRECPOINTS), AliQAv1::GetTaskName(AliQAv1::kRECPOINTS)) ;
3129 tempo.ReplaceAll(Form("%d", AliQAv1::kESDS), AliQAv1::GetTaskName(AliQAv1::kESDS)) ;
3130 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
3135 //_____________________________________________________________________________
3136 Bool_t AliReconstruction::InitRecoParams()
3138 // The method accesses OCDB and retrieves all
3139 // the available reco-param objects from there.
3141 Bool_t isOK = kTRUE;
3143 TString detStr = fLoadCDB;
3144 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3146 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3148 if (fRecoParam.GetDetRecoParamArray(iDet)) {
3149 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
3153 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
3155 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
3156 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3158 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
3162 TObject *recoParamObj = entry->GetObject();
3163 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3164 // The detector has a normal TobjArray of AliDetectorRecoParam objects
3165 // Registering them in AliRecoParam
3166 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
3168 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3169 // The detector has only onse set of reco parameters
3170 // Registering it in AliRecoParam
3171 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
3172 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3173 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3176 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
3180 AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
3184 if (AliDebugLevel() > 0) fRecoParam.Print();
3189 //_____________________________________________________________________________
3190 Bool_t AliReconstruction::GetEventInfo()
3192 // Fill the event info object
3194 AliCodeTimerAuto("")
3196 AliCentralTrigger *aCTP = NULL;
3198 fEventInfo.SetEventType(fRawReader->GetType());
3200 ULong64_t mask = fRawReader->GetClassMask();
3201 fEventInfo.SetTriggerMask(mask);
3202 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
3203 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
3205 aCTP = new AliCentralTrigger();
3206 TString configstr("");
3207 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
3208 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
3212 aCTP->SetClassMask(mask);
3213 aCTP->SetClusterMask(clmask);
3216 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
3218 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
3219 aCTP = fRunLoader->GetTrigger();
3220 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
3221 // get inputs from actp - just get
3222 AliESDHeader* esdheader = fesd->GetHeader();
3223 esdheader->SetL0TriggerInputs(aCTP->GetL0TriggerInputs());
3224 esdheader->SetL1TriggerInputs(aCTP->GetL1TriggerInputs());
3225 esdheader->SetL2TriggerInputs(aCTP->GetL2TriggerInputs());
3226 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
3229 AliWarning("No trigger can be loaded! The trigger information will not be used!");
3234 AliTriggerConfiguration *config = aCTP->GetConfiguration();
3236 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
3237 if (fRawReader) delete aCTP;
3241 UChar_t clustmask = 0;
3243 ULong64_t trmask = fEventInfo.GetTriggerMask();
3244 const TObjArray& classesArray = config->GetClasses();
3245 Int_t nclasses = classesArray.GetEntriesFast();
3246 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
3247 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
3249 Int_t trindex = TMath::Nint(TMath::Log2(trclass->GetMask()));
3250 fesd->SetTriggerClass(trclass->GetName(),trindex);
3251 if (fRawReader) fRawReader->LoadTriggerClass(trclass->GetName(),trindex);
3252 if (trmask & (1ull << trindex)) {
3254 trclasses += trclass->GetName();
3256 clustmask |= trclass->GetCluster()->GetClusterMask();
3260 fEventInfo.SetTriggerClasses(trclasses);
3262 // Set the information in ESD
3263 fesd->SetTriggerMask(trmask);
3264 fesd->SetTriggerCluster(clustmask);
3266 if (!aCTP->CheckTriggeredDetectors()) {
3267 if (fRawReader) delete aCTP;
3271 if (fRawReader) delete aCTP;
3273 // We have to fill also the HLT decision here!!
3279 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3281 // Match the detector list found in the rec.C or the default 'ALL'
3282 // to the list found in the GRP (stored there by the shuttle PP which
3283 // gets the information from ECS)
3284 static TString resultList;
3285 TString detList = detectorList;
3289 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3290 if ((detectorMask >> iDet) & 0x1) {
3291 TString det = AliDAQ::OfflineModuleName(iDet);
3292 if ((detList.CompareTo("ALL") == 0) ||
3293 ((detList.BeginsWith("ALL ") ||
3294 detList.EndsWith(" ALL") ||
3295 detList.Contains(" ALL ")) &&
3296 !(detList.BeginsWith("-"+det+" ") ||
3297 detList.EndsWith(" -"+det) ||
3298 detList.Contains(" -"+det+" "))) ||
3299 (detList.CompareTo(det) == 0) ||
3300 detList.BeginsWith(det+" ") ||
3301 detList.EndsWith(" "+det) ||
3302 detList.Contains( " "+det+" " )) {
3303 if (!resultList.EndsWith(det + " ")) {
3312 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3313 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3314 if ((detList.CompareTo("ALL") == 0) ||
3315 ((detList.BeginsWith("ALL ") ||
3316 detList.EndsWith(" ALL") ||
3317 detList.Contains(" ALL ")) &&
3318 !(detList.BeginsWith("-"+hltDet+" ") ||
3319 detList.EndsWith(" -"+hltDet) ||
3320 detList.Contains(" -"+hltDet+" "))) ||
3321 (detList.CompareTo(hltDet) == 0) ||
3322 detList.BeginsWith(hltDet+" ") ||
3323 detList.EndsWith(" "+hltDet) ||
3324 detList.Contains( " "+hltDet+" " )) {
3325 resultList += hltDet;
3329 return resultList.Data();
3333 //______________________________________________________________________________
3334 void AliReconstruction::Abort(const char *method, EAbort what)
3336 // Abort processing. If what = kAbortProcess, the Process() loop will be
3337 // aborted. If what = kAbortFile, the current file in a chain will be
3338 // aborted and the processing will continue with the next file, if there
3339 // is no next file then Process() will be aborted. Abort() can also be
3340 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3341 // the SlaveTerminate() and Terminate() are always called. The abort flag
3342 // can be checked in these methods using GetAbort().
3344 // The method is overwritten in AliReconstruction for better handling of
3345 // reco specific errors
3347 if (!fStopOnError) return;
3351 TString whyMess = method;
3352 whyMess += " failed! Aborting...";
3354 AliError(whyMess.Data());
3357 TString mess = "Abort";
3358 if (fAbort == kAbortProcess)
3359 mess = "AbortProcess";
3360 else if (fAbort == kAbortFile)
3363 Info(mess, whyMess.Data());