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>
129 #include <THashTable.h>
131 #include <TMessage.h>
133 #include "AliAlignObj.h"
134 #include "AliCDBEntry.h"
135 #include "AliCDBManager.h"
136 #include "AliCDBStorage.h"
137 #include "AliCTPRawStream.h"
138 #include "AliCascadeVertexer.h"
139 #include "AliCentralTrigger.h"
140 #include "AliCodeTimer.h"
142 #include "AliDetectorRecoParam.h"
143 #include "AliESDCaloCells.h"
144 #include "AliESDCaloCluster.h"
145 #include "AliESDEvent.h"
146 #include "AliESDMuonTrack.h"
147 #include "AliESDPmdTrack.h"
148 #include "AliESDTagCreator.h"
149 #include "AliESDVertex.h"
150 #include "AliESDcascade.h"
151 #include "AliESDfriend.h"
152 #include "AliESDkink.h"
153 #include "AliESDpid.h"
154 #include "AliESDtrack.h"
155 #include "AliESDtrack.h"
156 #include "AliEventInfo.h"
157 #include "AliGRPObject.h"
158 #include "AliGRPRecoParam.h"
159 #include "AliGenEventHeader.h"
160 #include "AliGeomManager.h"
161 #include "AliGlobalQADataMaker.h"
162 #include "AliHeader.h"
165 #include "AliMultiplicity.h"
167 #include "AliPlaneEff.h"
169 #include "AliQADataMakerRec.h"
170 #include "AliQAManager.h"
171 #include "AliRawVEvent.h"
172 #include "AliRawEventHeaderBase.h"
173 #include "AliRawHLTManager.h"
174 #include "AliRawReaderDate.h"
175 #include "AliRawReaderFile.h"
176 #include "AliRawReaderRoot.h"
177 #include "AliReconstruction.h"
178 #include "AliReconstructor.h"
180 #include "AliRunInfo.h"
181 #include "AliRunLoader.h"
182 #include "AliSysInfo.h" // memory snapshots
183 #include "AliTrackPointArray.h"
184 #include "AliTracker.h"
185 #include "AliTriggerClass.h"
186 #include "AliTriggerCluster.h"
187 #include "AliTriggerIR.h"
188 #include "AliTriggerConfiguration.h"
189 #include "AliV0vertexer.h"
190 #include "AliVertexer.h"
191 #include "AliVertexerTracks.h"
192 #include "AliTriggerRunScalers.h"
193 #include "AliCTPTimeParams.h"
195 ClassImp(AliReconstruction)
197 //_____________________________________________________________________________
198 const char* AliReconstruction::fgkDetectorName[AliReconstruction::kNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
200 //_____________________________________________________________________________
201 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
203 fRunVertexFinder(kTRUE),
204 fRunVertexFinderTracks(kTRUE),
205 fRunHLTTracking(kFALSE),
206 fRunMuonTracking(kFALSE),
208 fRunCascadeFinder(kTRUE),
209 fStopOnError(kFALSE),
210 fWriteAlignmentData(kFALSE),
211 fWriteESDfriend(kFALSE),
212 fFillTriggerESD(kTRUE),
220 fRunLocalReconstruction("ALL"),
224 fUseTrackingErrorsForAlignment(""),
225 fGAliceFileName(gAliceFilename),
228 fProofOutputFileName(""),
229 fProofOutputLocation(""),
230 fProofOutputDataset(kFALSE),
231 fProofOutputArchive(""),
235 fNumberOfEventsPerFile((UInt_t)-1),
237 fLoadAlignFromCDB(kTRUE),
238 fLoadAlignData("ALL"),
243 fCTPTimeParams(NULL),
247 fParentRawReader(NULL),
251 fSPDTrackleter(NULL),
253 fDiamondProfileSPD(NULL),
254 fDiamondProfile(NULL),
255 fDiamondProfileTPC(NULL),
256 fListOfCosmicTriggers(NULL),
260 fAlignObjArray(NULL),
264 fInitCDBCalled(kFALSE),
265 fSetRunNumberFromDataCalled(kFALSE),
270 fSameQACycle(kFALSE),
271 fInitQACalled(kFALSE),
272 fWriteQAExpertData(kTRUE),
273 fRunPlaneEff(kFALSE),
282 fIsNewRunLoader(kFALSE),
286 // create reconstruction object with default parameters
289 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
290 fReconstructor[iDet] = NULL;
291 fLoader[iDet] = NULL;
292 fTracker[iDet] = NULL;
294 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
295 fQACycles[iDet] = 999999 ;
296 fQAWriteExpert[iDet] = kFALSE ;
302 //_____________________________________________________________________________
303 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
305 fRunVertexFinder(rec.fRunVertexFinder),
306 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
307 fRunHLTTracking(rec.fRunHLTTracking),
308 fRunMuonTracking(rec.fRunMuonTracking),
309 fRunV0Finder(rec.fRunV0Finder),
310 fRunCascadeFinder(rec.fRunCascadeFinder),
311 fStopOnError(rec.fStopOnError),
312 fWriteAlignmentData(rec.fWriteAlignmentData),
313 fWriteESDfriend(rec.fWriteESDfriend),
314 fFillTriggerESD(rec.fFillTriggerESD),
316 fCleanESD(rec.fCleanESD),
317 fV0DCAmax(rec.fV0DCAmax),
318 fV0CsPmin(rec.fV0CsPmin),
322 fRunLocalReconstruction(rec.fRunLocalReconstruction),
323 fRunTracking(rec.fRunTracking),
324 fFillESD(rec.fFillESD),
325 fLoadCDB(rec.fLoadCDB),
326 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
327 fGAliceFileName(rec.fGAliceFileName),
328 fRawInput(rec.fRawInput),
329 fESDOutput(rec.fESDOutput),
330 fProofOutputFileName(rec.fProofOutputFileName),
331 fProofOutputLocation(rec.fProofOutputLocation),
332 fProofOutputDataset(rec.fProofOutputDataset),
333 fProofOutputArchive(rec.fProofOutputArchive),
334 fEquipIdMap(rec.fEquipIdMap),
335 fFirstEvent(rec.fFirstEvent),
336 fLastEvent(rec.fLastEvent),
337 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
339 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
340 fLoadAlignData(rec.fLoadAlignData),
341 fUseHLTData(rec.fUseHLTData),
345 fCTPTimeParams(NULL),
349 fParentRawReader(NULL),
351 fRecoParam(rec.fRecoParam),
353 fSPDTrackleter(NULL),
355 fDiamondProfileSPD(rec.fDiamondProfileSPD),
356 fDiamondProfile(rec.fDiamondProfile),
357 fDiamondProfileTPC(rec.fDiamondProfileTPC),
358 fListOfCosmicTriggers(NULL),
362 fAlignObjArray(rec.fAlignObjArray),
363 fCDBUri(rec.fCDBUri),
364 fQARefUri(rec.fQARefUri),
366 fInitCDBCalled(rec.fInitCDBCalled),
367 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
368 fQADetectors(rec.fQADetectors),
369 fQATasks(rec.fQATasks),
371 fRunGlobalQA(rec.fRunGlobalQA),
372 fSameQACycle(rec.fSameQACycle),
373 fInitQACalled(rec.fInitQACalled),
374 fWriteQAExpertData(rec.fWriteQAExpertData),
375 fRunPlaneEff(rec.fRunPlaneEff),
384 fIsNewRunLoader(rec.fIsNewRunLoader),
390 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
391 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
393 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
394 fReconstructor[iDet] = NULL;
395 fLoader[iDet] = NULL;
396 fTracker[iDet] = NULL;
399 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
400 fQACycles[iDet] = rec.fQACycles[iDet];
401 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
404 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
405 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
409 //_____________________________________________________________________________
410 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
412 // assignment operator
413 // Used in PROOF mode
414 // Be very careful while modifing it!
415 // Simple rules to follow:
416 // for persistent data members - use their assignment operators
417 // for non-persistent ones - do nothing or take the default values from constructor
418 // TSelector members should not be touched
419 if(&rec == this) return *this;
421 fRunVertexFinder = rec.fRunVertexFinder;
422 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
423 fRunHLTTracking = rec.fRunHLTTracking;
424 fRunMuonTracking = rec.fRunMuonTracking;
425 fRunV0Finder = rec.fRunV0Finder;
426 fRunCascadeFinder = rec.fRunCascadeFinder;
427 fStopOnError = rec.fStopOnError;
428 fWriteAlignmentData = rec.fWriteAlignmentData;
429 fWriteESDfriend = rec.fWriteESDfriend;
430 fFillTriggerESD = rec.fFillTriggerESD;
432 fCleanESD = rec.fCleanESD;
433 fV0DCAmax = rec.fV0DCAmax;
434 fV0CsPmin = rec.fV0CsPmin;
438 fRunLocalReconstruction = rec.fRunLocalReconstruction;
439 fRunTracking = rec.fRunTracking;
440 fFillESD = rec.fFillESD;
441 fLoadCDB = rec.fLoadCDB;
442 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
443 fGAliceFileName = rec.fGAliceFileName;
444 fRawInput = rec.fRawInput;
445 fESDOutput = rec.fESDOutput;
446 fProofOutputFileName = rec.fProofOutputFileName;
447 fProofOutputLocation = rec.fProofOutputLocation;
448 fProofOutputDataset = rec.fProofOutputDataset;
449 fProofOutputArchive = rec.fProofOutputArchive;
450 fEquipIdMap = rec.fEquipIdMap;
451 fFirstEvent = rec.fFirstEvent;
452 fLastEvent = rec.fLastEvent;
453 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
455 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
456 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
459 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
460 fLoadAlignData = rec.fLoadAlignData;
461 fUseHLTData = rec.fUseHLTData;
463 delete fRunInfo; fRunInfo = NULL;
464 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
466 fEventInfo = rec.fEventInfo;
468 delete fRunScalers; fRunScalers = NULL;
469 if (rec.fRunScalers) fRunScalers = new AliTriggerRunScalers(*rec.fRunScalers);
471 delete fCTPTimeParams; fCTPTimeParams = NULL;
472 if (rec.fCTPTimeParams) fCTPTimeParams = new AliCTPTimeParams(*rec.fCTPTimeParams);
476 fParentRawReader = NULL;
478 fRecoParam = rec.fRecoParam;
480 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
481 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
482 delete fLoader[iDet]; fLoader[iDet] = NULL;
483 delete fTracker[iDet]; fTracker[iDet] = NULL;
486 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
487 fQACycles[iDet] = rec.fQACycles[iDet];
488 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
491 delete fSPDTrackleter; fSPDTrackleter = NULL;
493 delete fDiamondProfileSPD; fDiamondProfileSPD = NULL;
494 if (rec.fDiamondProfileSPD) fDiamondProfileSPD = new AliESDVertex(*rec.fDiamondProfileSPD);
495 delete fDiamondProfile; fDiamondProfile = NULL;
496 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
497 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
498 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
500 delete fListOfCosmicTriggers; fListOfCosmicTriggers = NULL;
501 if (rec.fListOfCosmicTriggers) fListOfCosmicTriggers = (THashTable*)((rec.fListOfCosmicTriggers)->Clone());
503 delete fGRPData; fGRPData = NULL;
504 // if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
505 if (rec.fGRPData) fGRPData = (AliGRPObject*)((rec.fGRPData)->Clone());
507 delete fAlignObjArray; fAlignObjArray = NULL;
510 fQARefUri = rec.fQARefUri;
511 fSpecCDBUri.Delete();
512 fInitCDBCalled = rec.fInitCDBCalled;
513 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
514 fQADetectors = rec.fQADetectors;
515 fQATasks = rec.fQATasks;
517 fRunGlobalQA = rec.fRunGlobalQA;
518 fSameQACycle = rec.fSameQACycle;
519 fInitQACalled = rec.fInitQACalled;
520 fWriteQAExpertData = rec.fWriteQAExpertData;
521 fRunPlaneEff = rec.fRunPlaneEff;
530 fIsNewRunLoader = rec.fIsNewRunLoader;
537 //_____________________________________________________________________________
538 AliReconstruction::~AliReconstruction()
543 if (fListOfCosmicTriggers) {
544 fListOfCosmicTriggers->Delete();
545 delete fListOfCosmicTriggers;
549 delete fCTPTimeParams;
551 if (fAlignObjArray) {
552 fAlignObjArray->Delete();
553 delete fAlignObjArray;
555 fSpecCDBUri.Delete();
557 AliCodeTimer::Instance()->Print();
560 //_____________________________________________________________________________
561 void AliReconstruction::InitQA()
563 //Initialize the QA and start of cycle
564 AliCodeTimerAuto("",0);
566 if (fInitQACalled) return;
567 fInitQACalled = kTRUE;
569 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kRECMODE) ;
570 if (fWriteQAExpertData)
571 qam->SetWriteExpert() ;
573 if (qam->IsDefaultStorageSet()) {
574 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
575 AliWarning("Default QA reference storage has been already set !");
576 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fQARefUri.Data()));
577 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
578 fQARefUri = qam->GetDefaultStorage()->GetURI();
580 if (fQARefUri.Length() > 0) {
581 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
582 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
583 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
585 fQARefUri="local://$ALICE_ROOT/QAref";
586 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
587 AliWarning("Default QA refeference storage not yet set !!!!");
588 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
589 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
592 qam->SetDefaultStorage(fQARefUri);
596 qam->SetActiveDetectors(fQADetectors) ;
597 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++) {
598 qam->SetCycleLength(AliQAv1::DETECTORINDEX_t(det), fQACycles[det]) ;
599 qam->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
601 if (!fRawReader && !fInput && IsInTasks(AliQAv1::kRAWS))
602 fQATasks.ReplaceAll(Form("%d",AliQAv1::kRAWS), "") ;
603 qam->SetTasks(fQATasks) ;
604 qam->InitQADataMaker(AliCDBManager::Instance()->GetRun()) ;
607 Bool_t sameCycle = kFALSE ;
608 AliQADataMaker *qadm = qam->GetQADataMaker(AliQAv1::kGLOBAL);
609 AliInfo(Form("Initializing the global QA data maker"));
610 if (IsInTasks(AliQAv1::kRECPOINTS)) {
611 qadm->StartOfCycle(AliQAv1::kRECPOINTS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
612 TObjArray **arr=qadm->Init(AliQAv1::kRECPOINTS);
613 AliTracker::SetResidualsArray(arr);
616 if (IsInTasks(AliQAv1::kESDS)) {
617 qadm->StartOfCycle(AliQAv1::kESDS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
618 qadm->Init(AliQAv1::kESDS);
621 AliSysInfo::AddStamp("InitQA") ;
624 //_____________________________________________________________________________
625 void AliReconstruction::MergeQA(const char *fileName)
627 //Initialize the QA and start of cycle
628 AliCodeTimerAuto("",0) ;
629 AliQAManager::QAManager()->Merge(AliCDBManager::Instance()->GetRun(),fileName) ;
630 AliSysInfo::AddStamp("MergeQA") ;
633 //_____________________________________________________________________________
634 void AliReconstruction::InitCDB()
636 // activate a default CDB storage
637 // First check if we have any CDB storage set, because it is used
638 // to retrieve the calibration and alignment constants
639 AliCodeTimerAuto("",0);
641 if (fInitCDBCalled) return;
642 fInitCDBCalled = kTRUE;
644 AliCDBManager* man = AliCDBManager::Instance();
645 if (man->IsDefaultStorageSet())
647 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
648 AliWarning("Default CDB storage has been already set !");
649 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
650 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
651 fCDBUri = man->GetDefaultStorage()->GetURI();
654 if (fCDBUri.Length() > 0)
656 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
657 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
658 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
659 man->SetDefaultStorage(fCDBUri);
661 else if (!man->GetRaw()){
662 fCDBUri="local://$ALICE_ROOT/OCDB";
663 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
664 AliWarning("Default CDB storage not yet set !!!!");
665 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
666 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
667 man->SetDefaultStorage(fCDBUri);
670 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
671 AliWarning("Default storage will be set after setting the Run Number!!!");
672 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
676 // Now activate the detector specific CDB storage locations
677 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
678 TObject* obj = fSpecCDBUri[i];
680 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
681 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
682 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
683 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
685 AliSysInfo::AddStamp("InitCDB");
688 //_____________________________________________________________________________
689 void AliReconstruction::SetDefaultStorage(const char* uri) {
690 // Store the desired default CDB storage location
691 // Activate it later within the Run() method
697 //_____________________________________________________________________________
698 void AliReconstruction::SetQARefDefaultStorage(const char* uri) {
699 // Store the desired default CDB storage location
700 // Activate it later within the Run() method
703 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
706 //_____________________________________________________________________________
707 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
708 // Store a detector-specific CDB storage location
709 // Activate it later within the Run() method
711 AliCDBPath aPath(calibType);
712 if(!aPath.IsValid()){
713 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
714 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
715 if(!strcmp(calibType, fgkDetectorName[iDet])) {
716 aPath.SetPath(Form("%s/*", calibType));
717 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
721 if(!aPath.IsValid()){
722 AliError(Form("Not a valid path or detector: %s", calibType));
727 // // check that calibType refers to a "valid" detector name
728 // Bool_t isDetector = kFALSE;
729 // for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
730 // TString detName = fgkDetectorName[iDet];
731 // if(aPath.GetLevel0() == detName) {
732 // isDetector = kTRUE;
738 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
742 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
743 if (obj) fSpecCDBUri.Remove(obj);
744 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
748 //_____________________________________________________________________________
749 Bool_t AliReconstruction::SetRunNumberFromData()
751 // The method is called in Run() in order
752 // to set a correct run number.
753 // In case of raw data reconstruction the
754 // run number is taken from the raw data header
756 if (fSetRunNumberFromDataCalled) return kTRUE;
757 fSetRunNumberFromDataCalled = kTRUE;
759 AliCDBManager* man = AliCDBManager::Instance();
762 if(fRawReader->NextEvent()) {
763 if(man->GetRun() > 0) {
764 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
766 man->SetRun(fRawReader->GetRunNumber());
767 fRawReader->RewindEvents();
770 if(man->GetRun() > 0) {
771 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
774 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
780 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
782 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
787 // read run number from gAlice
788 if(rl->GetHeader()) {
789 man->SetRun(rl->GetHeader()->GetRun());
794 AliError("Neither run-loader header nor RawReader objects are found !");
806 //_____________________________________________________________________________
807 void AliReconstruction::SetCDBLock() {
808 // Set CDB lock: from now on it is forbidden to reset the run number
809 // or the default storage or to activate any further storage!
811 AliCDBManager::Instance()->SetLock(1);
814 //_____________________________________________________________________________
815 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
817 // Read the alignment objects from CDB.
818 // Each detector is supposed to have the
819 // alignment objects in DET/Align/Data CDB path.
820 // All the detector objects are then collected,
821 // sorted by geometry level (starting from ALIC) and
822 // then applied to the TGeo geometry.
823 // Finally an overlaps check is performed.
825 // Load alignment data from CDB and fill fAlignObjArray
826 if(fLoadAlignFromCDB){
828 TString detStr = detectors;
829 TString loadAlObjsListOfDets = "";
831 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
832 if(!IsSelected(fgkDetectorName[iDet], detStr)) continue;
833 if(!strcmp(fgkDetectorName[iDet],"HLT")) continue;
835 if(AliGeomManager::GetNalignable(fgkDetectorName[iDet]) != 0)
837 loadAlObjsListOfDets += fgkDetectorName[iDet];
838 loadAlObjsListOfDets += " ";
840 } // end loop over detectors
842 if(AliGeomManager::GetNalignable("GRP") != 0)
843 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
844 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
845 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
847 // Check if the array with alignment objects was
848 // provided by the user. If yes, apply the objects
849 // to the present TGeo geometry
850 if (fAlignObjArray) {
851 if (gGeoManager && gGeoManager->IsClosed()) {
852 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
853 AliError("The misalignment of one or more volumes failed!"
854 "Compare the list of simulated detectors and the list of detector alignment data!");
859 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
865 if (fAlignObjArray) {
866 fAlignObjArray->Delete();
867 delete fAlignObjArray; fAlignObjArray=NULL;
873 //_____________________________________________________________________________
874 void AliReconstruction::SetGAliceFile(const char* fileName)
876 // set the name of the galice file
878 fGAliceFileName = fileName;
881 //_____________________________________________________________________________
882 void AliReconstruction::SetInput(const char* input)
884 // In case the input string starts with 'mem://', we run in an online mode
885 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
886 // file is assumed. One can give as an input:
887 // mem://: - events taken from DAQ monitoring libs online
889 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
890 if (input) fRawInput = input;
893 //_____________________________________________________________________________
894 void AliReconstruction::SetOutput(const char* output)
896 // Set the output ESD filename
897 // 'output' is a normalt ROOT url
898 // The method is used in case of raw-data reco with PROOF
899 if (output) fESDOutput = output;
902 //_____________________________________________________________________________
903 void AliReconstruction::SetOption(const char* detector, const char* option)
905 // set options for the reconstruction of a detector
907 TObject* obj = fOptions.FindObject(detector);
908 if (obj) fOptions.Remove(obj);
909 fOptions.Add(new TNamed(detector, option));
912 //_____________________________________________________________________________
913 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
915 // Set custom reconstruction parameters for a given detector
916 // Single set of parameters for all the events
918 // First check if the reco-params are global
919 if(!strcmp(detector, "GRP")) {
921 fRecoParam.AddDetRecoParam(kNDetectors,par);
925 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
926 if(!strcmp(detector, fgkDetectorName[iDet])) {
928 fRecoParam.AddDetRecoParam(iDet,par);
935 //_____________________________________________________________________________
936 Bool_t AliReconstruction::InitGRP() {
937 //------------------------------------
938 // Initialization of the GRP entry
939 //------------------------------------
940 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
944 TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
947 AliInfo("Found a TMap in GRP/GRP/Data, converting it into an AliGRPObject");
949 fGRPData = new AliGRPObject();
950 fGRPData->ReadValuesFromMap(m);
954 AliInfo("Found an AliGRPObject in GRP/GRP/Data, reading it");
955 fGRPData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
959 // FIX ME: The unloading of GRP entry is temporarily disabled
960 // because ZDC and VZERO are using it in order to initialize
961 // their reconstructor objects. In the future one has to think
962 // of propagating AliRunInfo to the reconstructors.
963 // AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
967 AliError("No GRP entry found in OCDB!");
971 TString lhcState = fGRPData->GetLHCState();
972 if (lhcState==AliGRPObject::GetInvalidString()) {
973 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
974 lhcState = "UNKNOWN";
977 TString beamType = fGRPData->GetBeamType();
978 if (beamType==AliGRPObject::GetInvalidString()) {
979 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
980 beamType = "UNKNOWN";
983 Float_t beamEnergy = fGRPData->GetBeamEnergy();
984 if (beamEnergy==AliGRPObject::GetInvalidFloat()) {
985 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
988 // LHC: "multiply by 120 to get the energy in MeV"
991 TString runType = fGRPData->GetRunType();
992 if (runType==AliGRPObject::GetInvalidString()) {
993 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
997 Int_t activeDetectors = fGRPData->GetDetectorMask();
998 if (activeDetectors==AliGRPObject::GetInvalidUInt()) {
999 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
1000 activeDetectors = 1074790399;
1003 fRunInfo = new AliRunInfo(lhcState, beamType, beamEnergy, runType, activeDetectors);
1007 // Process the list of active detectors
1008 if (activeDetectors) {
1009 UInt_t detMask = activeDetectors;
1010 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
1011 fRunTracking = MatchDetectorList(fRunTracking,detMask);
1012 fFillESD = MatchDetectorList(fFillESD,detMask);
1013 fQADetectors = MatchDetectorList(fQADetectors,detMask);
1014 fLoadCDB.Form("%s %s %s %s",
1015 fRunLocalReconstruction.Data(),
1016 fRunTracking.Data(),
1018 fQADetectors.Data());
1019 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
1020 if (!((detMask >> AliDAQ::DetectorID("ITSSPD")) & 0x1) &&
1021 !((detMask >> AliDAQ::DetectorID("ITSSDD")) & 0x1) &&
1022 !((detMask >> AliDAQ::DetectorID("ITSSSD")) & 0x1) ) {
1023 // switch off the vertexer
1024 AliInfo("SPD,SDD,SSD is not in the list of active detectors. Vertexer switched off.");
1025 fRunVertexFinder = kFALSE;
1027 if (!((detMask >> AliDAQ::DetectorID("TRG")) & 0x1)) {
1028 // switch off the reading of CTP raw-data payload
1029 if (fFillTriggerESD) {
1030 AliInfo("CTP is not in the list of active detectors. CTP data reading switched off.");
1031 fFillTriggerESD = kFALSE;
1036 AliInfo("===================================================================================");
1037 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
1038 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
1039 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
1040 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
1041 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
1042 AliInfo("===================================================================================");
1044 //*** Dealing with the magnetic field map
1045 if ( TGeoGlobalMagField::Instance()->IsLocked() ) {
1046 if (TGeoGlobalMagField::Instance()->GetField()->TestBit(AliMagF::kOverrideGRP)) {
1047 AliInfo("ExpertMode!!! GRP information will be ignored !");
1048 AliInfo("ExpertMode!!! Running with the externally locked B field !");
1051 AliInfo("Destroying existing B field instance!");
1052 delete TGeoGlobalMagField::Instance();
1055 if ( !TGeoGlobalMagField::Instance()->IsLocked() ) {
1056 // Construct the field map out of the information retrieved from GRP.
1059 Float_t l3Current = fGRPData->GetL3Current((AliGRPObject::Stats)0);
1060 if (l3Current == AliGRPObject::GetInvalidFloat()) {
1061 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1065 Char_t l3Polarity = fGRPData->GetL3Polarity();
1066 if (l3Polarity == AliGRPObject::GetInvalidChar()) {
1067 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1072 Float_t diCurrent = fGRPData->GetDipoleCurrent((AliGRPObject::Stats)0);
1073 if (diCurrent == AliGRPObject::GetInvalidFloat()) {
1074 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1078 Char_t diPolarity = fGRPData->GetDipolePolarity();
1079 if (diPolarity == AliGRPObject::GetInvalidChar()) {
1080 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1084 // read special bits for the polarity convention and map type
1085 Int_t polConvention = fGRPData->IsPolarityConventionLHC() ? AliMagF::kConvLHC : AliMagF::kConvDCS2008;
1086 Bool_t uniformB = fGRPData->IsUniformBMap();
1089 AliMagF* fld = AliMagF::CreateFieldMap(TMath::Abs(l3Current) * (l3Polarity ? -1:1),
1090 TMath::Abs(diCurrent) * (diPolarity ? -1:1),
1091 polConvention,uniformB,beamEnergy, beamType.Data());
1093 TGeoGlobalMagField::Instance()->SetField( fld );
1094 TGeoGlobalMagField::Instance()->Lock();
1095 AliInfo("Running with the B field constructed out of GRP !");
1097 else AliFatal("Failed to create a B field map !");
1099 else AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
1102 //*** Get the diamond profiles from OCDB
1103 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1105 fDiamondProfileSPD = dynamic_cast<AliESDVertex*> (entry->GetObject());
1107 AliError("No SPD diamond profile found in OCDB!");
1110 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1112 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
1114 AliError("No diamond profile found in OCDB!");
1117 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
1119 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
1121 AliError("No TPC diamond profile found in OCDB!");
1124 entry = AliCDBManager::Instance()->Get("GRP/Calib/CosmicTriggers");
1126 fListOfCosmicTriggers = dynamic_cast<THashTable*>(entry->GetObject());
1128 AliCDBManager::Instance()->UnloadFromCache("GRP/Calib/CosmicTriggers");
1131 if (!fListOfCosmicTriggers) {
1132 AliWarning("Can not get list of cosmic triggers from OCDB! Cosmic event specie will be effectively disabled!");
1138 //_____________________________________________________________________________
1139 Bool_t AliReconstruction::LoadCDB()
1141 AliCodeTimerAuto("",0);
1143 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1145 TString detStr = fLoadCDB;
1146 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1147 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1148 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1151 // Temporary fix - one has to define the correct policy in order
1152 // to load the trigger OCDB entries only for the detectors that
1153 // in the trigger or that are needed in order to put correct
1154 // information in ESD
1155 AliCDBManager::Instance()->GetAll("TRIGGER/*/*");
1159 //_____________________________________________________________________________
1160 Bool_t AliReconstruction::LoadTriggerScalersCDB()
1162 AliCodeTimerAuto("",0);
1164 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/CTP/Scalers");
1168 AliInfo("Found an AliTriggerRunScalers in GRP/CTP/Scalers, reading it");
1169 fRunScalers = dynamic_cast<AliTriggerRunScalers*> (entry->GetObject());
1171 if (fRunScalers->CorrectScalersOverflow() == 0) AliInfo("32bit Trigger counters corrected for overflow");
1176 //_____________________________________________________________________________
1177 Bool_t AliReconstruction::LoadCTPTimeParamsCDB()
1179 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/CTP/CTPtiming");
1183 AliInfo("Found an AliCTPTimeParams in GRP/CTP/CTPtiming, reading it");
1184 fCTPTimeParams = dynamic_cast<AliCTPTimeParams*> (entry->GetObject());
1191 //_____________________________________________________________________________
1192 Bool_t AliReconstruction::Run(const char* input)
1195 AliCodeTimerAuto("",0);
1198 if (GetAbort() != TSelector::kContinue) return kFALSE;
1200 TChain *chain = NULL;
1201 if (fRawReader && (chain = fRawReader->GetChain())) {
1202 Long64_t nEntries = (fLastEvent < 0) ? (TChain::kBigNumber) : (fLastEvent - fFirstEvent + 1);
1205 // Temporary fix for long raw-data runs (until socket timeout handling in PROOF is revised)
1206 gProof->Exec("gEnv->SetValue(\"Proof.SocketActivityTimeout\",-1)", kTRUE);
1209 gProof->Exec("TGrid::Connect(\"alien://\")",kTRUE);
1211 TMessage::EnableSchemaEvolutionForAll(kTRUE);
1212 gProof->Exec("TMessage::EnableSchemaEvolutionForAll(kTRUE)",kTRUE);
1214 gProof->AddInput(this);
1216 if (!ParseOutput()) return kFALSE;
1218 gProof->SetParameter("PROOF_MaxSlavesPerNode", 9999);
1220 chain->Process("AliReconstruction","",nEntries,fFirstEvent);
1223 chain->Process(this,"",nEntries,fFirstEvent);
1228 if (GetAbort() != TSelector::kContinue) return kFALSE;
1230 if (GetAbort() != TSelector::kContinue) return kFALSE;
1231 //******* The loop over events
1232 AliInfo("Starting looping over events");
1234 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1235 (fRawReader && fRawReader->NextEvent())) {
1236 if (!ProcessEvent(iEvent)) {
1237 Abort("ProcessEvent",TSelector::kAbortFile);
1243 if (GetAbort() != TSelector::kContinue) return kFALSE;
1245 if (GetAbort() != TSelector::kContinue) return kFALSE;
1251 //_____________________________________________________________________________
1252 void AliReconstruction::InitRawReader(const char* input)
1254 AliCodeTimerAuto("",0);
1256 // Init raw-reader and
1257 // set the input in case of raw data
1258 if (input) fRawInput = input;
1259 fRawReader = AliRawReader::Create(fRawInput.Data());
1261 if (fRawInput.IsNull()) {
1262 AliInfo("Reconstruction will run over digits");
1265 AliFatal("Can not create raw-data reader ! Exiting...");
1269 if (!fEquipIdMap.IsNull() && fRawReader)
1270 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1272 if (!fUseHLTData.IsNull()) {
1273 // create the RawReaderHLT which performs redirection of HLT input data for
1274 // the specified detectors
1275 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1277 fParentRawReader=fRawReader;
1278 fRawReader=pRawReader;
1280 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1283 AliSysInfo::AddStamp("CreateRawReader");
1286 //_____________________________________________________________________________
1287 void AliReconstruction::InitRun(const char* input)
1289 // Initialization of raw-reader,
1290 // run number, CDB etc.
1291 AliCodeTimerAuto("",0);
1292 AliSysInfo::AddStamp("Start");
1294 // Initialize raw-reader if any
1295 InitRawReader(input);
1297 // Initialize the CDB storage
1300 // Set run number in CDBManager (if it is not already set by the user)
1301 if (!SetRunNumberFromData()) {
1302 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1306 // Set CDB lock: from now on it is forbidden to reset the run number
1307 // or the default storage or to activate any further storage!
1312 //_____________________________________________________________________________
1313 void AliReconstruction::Begin(TTree *)
1315 // Initialize AlReconstruction before
1316 // going into the event loop
1317 // Should follow the TSelector convention
1318 // i.e. initialize only the object on the client side
1319 AliCodeTimerAuto("",0);
1321 AliReconstruction *reco = NULL;
1323 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1326 AliSysInfo::AddStamp("ReadInputInBegin");
1329 // Import ideal TGeo geometry and apply misalignment
1331 TString geom(gSystem->DirName(fGAliceFileName));
1332 geom += "/geometry.root";
1333 AliGeomManager::LoadGeometry(geom.Data());
1335 Abort("LoadGeometry", TSelector::kAbortProcess);
1338 AliSysInfo::AddStamp("LoadGeom");
1339 TString detsToCheck=fRunLocalReconstruction;
1340 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1341 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1344 AliSysInfo::AddStamp("CheckGeom");
1347 if (!MisalignGeometry(fLoadAlignData)) {
1348 Abort("MisalignGeometry", TSelector::kAbortProcess);
1351 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1352 AliSysInfo::AddStamp("MisalignGeom");
1355 Abort("InitGRP", TSelector::kAbortProcess);
1358 AliSysInfo::AddStamp("InitGRP");
1361 Abort("LoadCDB", TSelector::kAbortProcess);
1364 AliSysInfo::AddStamp("LoadCDB");
1366 if (!LoadTriggerScalersCDB()) {
1367 Abort("LoadTriggerScalersCDB", TSelector::kAbortProcess);
1370 AliSysInfo::AddStamp("LoadTriggerScalersCDB");
1372 if (!LoadCTPTimeParamsCDB()) {
1373 Abort("LoadCTPTimeParamsCDB", TSelector::kAbortProcess);
1376 AliSysInfo::AddStamp("LoadCTPTimeParamsCDB");
1378 // Read the reconstruction parameters from OCDB
1379 if (!InitRecoParams()) {
1380 AliWarning("Not all detectors have correct RecoParam objects initialized");
1382 AliSysInfo::AddStamp("InitRecoParams");
1384 if (fInput && gProof) {
1385 if (reco) *reco = *this;
1387 gGeoManager->SetName("Geometry");
1388 gProof->AddInputData(gGeoManager,kTRUE);
1390 gProof->AddInputData(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()),kTRUE);
1391 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1392 AliMagF *magFieldMap = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1393 magFieldMap->SetName("MagneticFieldMap");
1394 gProof->AddInputData(magFieldMap,kTRUE);
1399 //_____________________________________________________________________________
1400 void AliReconstruction::SlaveBegin(TTree*)
1402 // Initialization related to run-loader,
1403 // vertexer, trackers, recontructors
1404 // In proof mode it is executed on the slave
1405 AliCodeTimerAuto("",0);
1407 TProofOutputFile *outProofFile = NULL;
1409 if (AliDebugLevel() > 0) fInput->Print();
1410 if (AliDebugLevel() > 10) fInput->Dump();
1411 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1414 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1416 AliGeomManager::SetGeometry(tgeo);
1418 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1419 Int_t runNumber = -1;
1420 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1421 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1422 man->SetCacheFlag(kTRUE);
1423 man->SetLock(kTRUE);
1427 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1428 AliMagF *newMap = new AliMagF(*map);
1429 if (!newMap->LoadParameterization()) {
1430 Abort("AliMagF::LoadParameterization", TSelector::kAbortProcess);
1433 TGeoGlobalMagField::Instance()->SetField(newMap);
1434 TGeoGlobalMagField::Instance()->Lock();
1436 if (TNamed *outputFileName = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE"))
1437 fProofOutputFileName = outputFileName->GetTitle();
1438 if (TNamed *outputLocation = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE_LOCATION"))
1439 fProofOutputLocation = outputLocation->GetTitle();
1440 if (fInput->FindObject("PROOF_OUTPUTFILE_DATASET"))
1441 fProofOutputDataset = kTRUE;
1442 if (TNamed *archiveList = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE_ARCHIVE"))
1443 fProofOutputArchive = archiveList->GetTitle();
1444 if (!fProofOutputFileName.IsNull() &&
1445 !fProofOutputLocation.IsNull() &&
1446 fProofOutputArchive.IsNull()) {
1447 if (!fProofOutputDataset) {
1448 outProofFile = new TProofOutputFile(fProofOutputFileName.Data(),"M");
1449 outProofFile->SetOutputFileName(Form("%s%s",fProofOutputLocation.Data(),fProofOutputFileName.Data()));
1452 outProofFile = new TProofOutputFile(fProofOutputFileName.Data(),"DROV",fProofOutputLocation.Data());
1454 if (AliDebugLevel() > 0) outProofFile->Dump();
1455 fOutput->Add(outProofFile);
1457 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1460 // get the run loader
1461 if (!InitRunLoader()) {
1462 Abort("InitRunLoader", TSelector::kAbortProcess);
1465 AliSysInfo::AddStamp("LoadLoader");
1467 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1470 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1471 Abort("CreateTrackers", TSelector::kAbortProcess);
1474 AliSysInfo::AddStamp("CreateTrackers");
1476 // create the ESD output file and tree
1477 if (!outProofFile) {
1478 ffile = TFile::Open("AliESDs.root", "RECREATE");
1479 ffile->SetCompressionLevel(2);
1480 if (!ffile->IsOpen()) {
1481 Abort("OpenESDFile", TSelector::kAbortProcess);
1486 AliInfo(Form("Opening output PROOF file: %s/%s",
1487 outProofFile->GetDir(), outProofFile->GetFileName()));
1488 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1489 Abort(Form("Problems opening output PROOF file: %s/%s",
1490 outProofFile->GetDir(), outProofFile->GetFileName()),
1491 TSelector::kAbortProcess);
1496 ftree = new TTree("esdTree", "Tree with ESD objects");
1497 fesd = new AliESDEvent();
1498 fesd->CreateStdContent();
1500 fesd->WriteToTree(ftree);
1501 if (fWriteESDfriend) {
1503 // Since we add the branch manually we must
1504 // book and add it after WriteToTree
1505 // otherwise it is created twice,
1506 // once via writetotree and once here.
1507 // The case for AliESDfriend is now
1508 // caught also in AlIESDEvent::WriteToTree but
1509 // be careful when changing the name (AliESDfriend is not
1510 // a TNamed so we had to hardwire it)
1511 fesdf = new AliESDfriend();
1512 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1513 br->SetFile("AliESDfriends.root");
1514 fesd->AddObject(fesdf);
1516 ftree->GetUserInfo()->Add(fesd);
1518 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1519 fhltesd = new AliESDEvent();
1520 fhltesd->CreateStdContent();
1522 // read the ESD template from CDB
1523 // HLT is allowed to put non-std content to its ESD, the non-std
1524 // objects need to be created before invocation of WriteToTree in
1525 // order to create all branches. Initialization is done from an
1526 // ESD layout template in CDB
1527 AliCDBManager* man = AliCDBManager::Instance();
1528 AliCDBPath hltESDConfigPath("HLT/ConfigHLT/esdLayout");
1529 AliCDBEntry* hltESDConfig=NULL;
1530 if (man->GetId(hltESDConfigPath)!=NULL &&
1531 (hltESDConfig=man->Get(hltESDConfigPath))!=NULL) {
1532 AliESDEvent* pESDLayout=dynamic_cast<AliESDEvent*>(hltESDConfig->GetObject());
1534 // init all internal variables from the list of objects
1535 pESDLayout->GetStdContent();
1537 // copy content and create non-std objects
1538 *fhltesd=*pESDLayout;
1541 AliError(Form("error setting hltEsd layout from %s: invalid object type",
1542 hltESDConfigPath.GetPath().Data()));
1546 fhltesd->WriteToTree(fhlttree);
1547 fhlttree->GetUserInfo()->Add(fhltesd);
1549 ProcInfo_t procInfo;
1550 gSystem->GetProcInfo(&procInfo);
1551 AliInfo(Form("Current memory usage %d %d", procInfo.fMemResident, procInfo.fMemVirtual));
1554 //Initialize the QA and start of cycle
1555 if (fRunQA || fRunGlobalQA)
1558 //Initialize the Plane Efficiency framework
1559 if (fRunPlaneEff && !InitPlaneEff()) {
1560 Abort("InitPlaneEff", TSelector::kAbortProcess);
1564 if (strcmp(gProgName,"alieve") == 0)
1565 fRunAliEVE = InitAliEVE();
1570 //_____________________________________________________________________________
1571 Bool_t AliReconstruction::Process(Long64_t entry)
1573 // run the reconstruction over a single entry
1574 // from the chain with raw data
1575 AliCodeTimerAuto("",0);
1577 TTree *currTree = fChain->GetTree();
1578 AliRawVEvent *event = NULL;
1579 currTree->SetBranchAddress("rawevent",&event);
1580 currTree->GetEntry(entry);
1581 fRawReader = new AliRawReaderRoot(event);
1582 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1590 //_____________________________________________________________________________
1591 void AliReconstruction::Init(TTree *tree)
1594 AliError("The input tree is not found!");
1600 //_____________________________________________________________________________
1601 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1603 // run the reconstruction over a single event
1604 // The event loop is steered in Run method
1606 AliCodeTimerAuto("",0);
1608 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1609 fRunLoader->SetEventNumber(iEvent);
1610 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1612 fRunLoader->TreeE()->Fill();
1613 if (fRawReader && fRawReader->UseAutoSaveESD())
1614 fRunLoader->TreeE()->AutoSave("SaveSelf");
1617 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1621 AliInfo(Form("processing event %d", iEvent));
1623 fRunLoader->GetEvent(iEvent);
1625 // Fill Event-info object
1627 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo,fListOfCosmicTriggers);
1628 AliInfo(Form("Current event specie: %s",fRecoParam.PrintEventSpecie()));
1630 // Set the reco-params
1632 TString detStr = fLoadCDB;
1633 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1634 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1635 AliReconstructor *reconstructor = GetReconstructor(iDet);
1636 if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
1637 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
1638 reconstructor->SetRecoParam(par);
1639 reconstructor->SetEventInfo(&fEventInfo);
1641 AliQAManager::QAManager()->SetRecoParam(iDet, par) ;
1642 AliQAManager::QAManager()->SetEventSpecie(AliRecoParam::Convert(par->GetEventSpecie())) ;
1649 if (fRunQA && IsInTasks(AliQAv1::kRAWS)) {
1650 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1651 AliQAManager::QAManager()->RunOneEvent(fRawReader) ;
1653 // local single event reconstruction
1654 if (!fRunLocalReconstruction.IsNull()) {
1655 TString detectors=fRunLocalReconstruction;
1656 // run HLT event reconstruction first
1657 // ;-( IsSelected changes the string
1658 if (IsSelected("HLT", detectors) &&
1659 !RunLocalEventReconstruction("HLT")) {
1660 if (fStopOnError) {CleanUp(); return kFALSE;}
1662 detectors=fRunLocalReconstruction;
1663 detectors.ReplaceAll("HLT", "");
1664 if (!RunLocalEventReconstruction(detectors)) {
1665 if (fStopOnError) {CleanUp(); return kFALSE;}
1669 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1670 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1671 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1672 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1674 fesd->SetEventSpecie(fRecoParam.GetEventSpecie());
1675 fhltesd->SetEventSpecie(fRecoParam.GetEventSpecie());
1677 // Set magnetic field from the tracker
1678 fesd->SetMagneticField(AliTracker::GetBz());
1679 fhltesd->SetMagneticField(AliTracker::GetBz());
1681 AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1682 if (fld) { // set info needed for field initialization
1683 fesd->SetCurrentL3(fld->GetCurrentSol());
1684 fesd->SetCurrentDip(fld->GetCurrentDip());
1685 fesd->SetBeamEnergy(fld->GetBeamEnergy());
1686 fesd->SetBeamType(fld->GetBeamTypeText());
1687 fesd->SetUniformBMap(fld->IsUniform());
1688 fesd->SetBInfoStored();
1690 fhltesd->SetCurrentL3(fld->GetCurrentSol());
1691 fhltesd->SetCurrentDip(fld->GetCurrentDip());
1692 fhltesd->SetBeamEnergy(fld->GetBeamEnergy());
1693 fhltesd->SetBeamType(fld->GetBeamTypeText());
1694 fhltesd->SetUniformBMap(fld->IsUniform());
1695 fhltesd->SetBInfoStored();
1698 // Set most probable pt, for B=0 tracking
1699 // Get the global reco-params. They are atposition 16 inside the array of detectors in fRecoParam
1700 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(kNDetectors));
1701 if (grpRecoParam) AliExternalTrackParam::SetMostProbablePt(grpRecoParam->GetMostProbablePt());
1703 // Fill raw-data error log into the ESD
1704 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1707 if (fRunVertexFinder) {
1708 if (!RunVertexFinder(fesd)) {
1709 if (fStopOnError) {CleanUp(); return kFALSE;}
1713 // For Plane Efficiency: run the SPD trackleter
1714 if (fRunPlaneEff && fSPDTrackleter) {
1715 if (!RunSPDTrackleting(fesd)) {
1716 if (fStopOnError) {CleanUp(); return kFALSE;}
1721 if (!fRunTracking.IsNull()) {
1722 if (fRunMuonTracking) {
1723 if (!RunMuonTracking(fesd)) {
1724 if (fStopOnError) {CleanUp(); return kFALSE;}
1730 if (!fRunTracking.IsNull()) {
1731 if (!RunTracking(fesd)) {
1732 if (fStopOnError) {CleanUp(); return kFALSE;}
1737 if (!fFillESD.IsNull()) {
1738 TString detectors=fFillESD;
1739 // run HLT first and on hltesd
1740 // ;-( IsSelected changes the string
1741 if (IsSelected("HLT", detectors) &&
1742 !FillESD(fhltesd, "HLT")) {
1743 if (fStopOnError) {CleanUp(); return kFALSE;}
1746 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1747 if (detectors.Contains("ALL")) {
1749 for (Int_t idet=0; idet<kNDetectors; ++idet){
1750 detectors += fgkDetectorName[idet];
1754 detectors.ReplaceAll("HLT", "");
1755 if (!FillESD(fesd, detectors)) {
1756 if (fStopOnError) {CleanUp(); return kFALSE;}
1760 // fill Event header information from the RawEventHeader
1761 if (fRawReader){FillRawEventHeaderESD(fesd);}
1762 if (fRawReader){FillRawEventHeaderESD(fhltesd);}
1765 AliESDpid::MakePID(fesd);
1767 if (fFillTriggerESD) {
1768 if (!FillTriggerESD(fesd)) {
1769 if (fStopOnError) {CleanUp(); return kFALSE;}
1772 // Always fill scalers
1773 if (!FillTriggerScalers(fesd)) {
1774 if (fStopOnError) {CleanUp(); return kFALSE;}
1781 // Propagate track to the beam pipe (if not already done by ITS)
1783 const Int_t ntracks = fesd->GetNumberOfTracks();
1784 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1787 UShort_t *selectedIdx=new UShort_t[ntracks];
1789 for (Int_t itrack=0; itrack<ntracks; itrack++){
1790 const Double_t kMaxStep = 1; //max step over the material
1793 AliESDtrack *track = fesd->GetTrack(itrack);
1794 if (!track) continue;
1796 AliExternalTrackParam *tpcTrack =
1797 (AliExternalTrackParam *)track->GetTPCInnerParam();
1801 PropagateTrackToBxByBz(tpcTrack,kRadius,track->GetMass(),kMaxStep,kFALSE);
1804 Int_t n=trkArray.GetEntriesFast();
1805 selectedIdx[n]=track->GetID();
1806 trkArray.AddLast(tpcTrack);
1809 //Tracks refitted by ITS should already be at the SPD vertex
1810 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1813 PropagateTrackToBxByBz(track,kRadius,track->GetMass(),kMaxStep,kFALSE);
1814 Double_t x[3]; track->GetXYZ(x);
1815 Double_t b[3]; AliTracker::GetBxByBz(x,b);
1816 track->RelateToVertexBxByBz(fesd->GetPrimaryVertexSPD(), b, kVeryBig);
1821 // Improve the reconstructed primary vertex position using the tracks
1823 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
1824 if(fesd->GetPrimaryVertexSPD()) {
1825 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
1826 if(vtitle.Contains("cosmics")) {
1827 runVertexFinderTracks=kFALSE;
1831 if (runVertexFinderTracks) {
1832 // TPC + ITS primary vertex
1833 ftVertexer->SetITSMode();
1834 ftVertexer->SetConstraintOff();
1835 // get cuts for vertexer from AliGRPRecoParam
1837 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1838 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1839 grpRecoParam->GetVertexerTracksCutsITS(cutsVertexer);
1840 ftVertexer->SetCuts(cutsVertexer);
1841 delete [] cutsVertexer; cutsVertexer = NULL;
1842 if(fDiamondProfile && grpRecoParam->GetVertexerTracksConstraintITS())
1843 ftVertexer->SetVtxStart(fDiamondProfile);
1845 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1847 if (pvtx->GetStatus()) {
1848 fesd->SetPrimaryVertexTracks(pvtx);
1849 for (Int_t i=0; i<ntracks; i++) {
1850 AliESDtrack *t = fesd->GetTrack(i);
1851 Double_t x[3]; t->GetXYZ(x);
1852 Double_t b[3]; AliTracker::GetBxByBz(x,b);
1853 t->RelateToVertexBxByBz(pvtx, b, kVeryBig);
1856 delete pvtx; pvtx=NULL;
1859 // TPC-only primary vertex
1860 ftVertexer->SetTPCMode();
1861 ftVertexer->SetConstraintOff();
1862 // get cuts for vertexer from AliGRPRecoParam
1864 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1865 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1866 grpRecoParam->GetVertexerTracksCutsTPC(cutsVertexer);
1867 ftVertexer->SetCuts(cutsVertexer);
1868 delete [] cutsVertexer; cutsVertexer = NULL;
1869 if(fDiamondProfileTPC && grpRecoParam->GetVertexerTracksConstraintTPC())
1870 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1872 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1874 if (pvtx->GetStatus()) {
1875 fesd->SetPrimaryVertexTPC(pvtx);
1876 for (Int_t i=0; i<ntracks; i++) {
1877 AliESDtrack *t = fesd->GetTrack(i);
1878 Double_t x[3]; t->GetXYZ(x);
1879 Double_t b[3]; AliTracker::GetBxByBz(x,b);
1880 t->RelateToVertexTPCBxByBz(pvtx, b, kVeryBig);
1883 delete pvtx; pvtx=NULL;
1887 delete[] selectedIdx;
1889 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1894 AliV0vertexer vtxer;
1895 vtxer.Tracks2V0vertices(fesd);
1897 if (fRunCascadeFinder) {
1899 AliCascadeVertexer cvtxer;
1900 cvtxer.V0sTracks2CascadeVertices(fesd);
1905 if (fCleanESD) CleanESD(fesd);
1907 if (fRunQA && IsInTasks(AliQAv1::kESDS)) {
1908 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1909 AliQAManager::QAManager()->RunOneEvent(fesd) ;
1912 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
1913 qadm->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1914 if (qadm && IsInTasks(AliQAv1::kESDS))
1915 qadm->Exec(AliQAv1::kESDS, fesd);
1918 if (fWriteESDfriend) {
1919 // fesdf->~AliESDfriend();
1920 // new (fesdf) AliESDfriend(); // Reset...
1921 fesd->GetESDfriend(fesdf);
1925 // Auto-save the ESD tree in case of prompt reco @P2
1926 if (fRawReader && fRawReader->UseAutoSaveESD()) {
1927 ftree->AutoSave("SaveSelf");
1928 TFile *friendfile = (TFile *)(gROOT->GetListOfFiles()->FindObject("AliESDfriends.root"));
1929 if (friendfile) friendfile->Save();
1936 if (fRunAliEVE) RunAliEVE();
1940 if (fWriteESDfriend) {
1941 fesdf->~AliESDfriend();
1942 new (fesdf) AliESDfriend(); // Reset...
1945 ProcInfo_t procInfo;
1946 gSystem->GetProcInfo(&procInfo);
1947 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, procInfo.fMemResident, procInfo.fMemVirtual));
1950 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1951 if (fReconstructor[iDet]) {
1952 fReconstructor[iDet]->SetRecoParam(NULL);
1953 fReconstructor[iDet]->SetEventInfo(NULL);
1955 if (fTracker[iDet]) fTracker[iDet]->SetEventInfo(NULL);
1958 if (fRunQA || fRunGlobalQA)
1959 AliQAManager::QAManager()->Increment() ;
1964 //_____________________________________________________________________________
1965 void AliReconstruction::SlaveTerminate()
1967 // Finalize the run on the slave side
1968 // Called after the exit
1969 // from the event loop
1970 AliCodeTimerAuto("",0);
1972 if (fIsNewRunLoader) { // galice.root didn't exist
1973 fRunLoader->WriteHeader("OVERWRITE");
1974 fRunLoader->CdGAFile();
1975 fRunLoader->Write(0, TObject::kOverwrite);
1978 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1979 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1981 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1982 cdbMapCopy->SetOwner(1);
1983 cdbMapCopy->SetName("cdbMap");
1984 TIter iter(cdbMap->GetTable());
1987 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1988 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1989 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1990 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1993 TList *cdbListCopy = new TList();
1994 cdbListCopy->SetOwner(1);
1995 cdbListCopy->SetName("cdbList");
1997 TIter iter2(cdbList);
2000 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
2001 cdbListCopy->Add(new TObjString(id->ToString().Data()));
2004 ftree->GetUserInfo()->Add(cdbMapCopy);
2005 ftree->GetUserInfo()->Add(cdbListCopy);
2010 if (fWriteESDfriend)
2011 ftree->SetBranchStatus("ESDfriend*",0);
2012 // we want to have only one tree version number
2013 ftree->Write(ftree->GetName(),TObject::kOverwrite);
2014 fhlttree->Write(fhlttree->GetName(),TObject::kOverwrite);
2016 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
2017 if (fRunPlaneEff && !FinishPlaneEff()) {
2018 AliWarning("Finish PlaneEff evaluation failed");
2021 // End of cycle for the in-loop
2023 AliQAManager::QAManager()->EndOfCycle() ;
2026 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2028 if (IsInTasks(AliQAv1::kRECPOINTS))
2029 qadm->EndOfCycle(AliQAv1::kRECPOINTS);
2030 if (IsInTasks(AliQAv1::kESDS))
2031 qadm->EndOfCycle(AliQAv1::kESDS);
2036 if (fRunQA || fRunGlobalQA) {
2038 !fProofOutputLocation.IsNull() &&
2039 fProofOutputArchive.IsNull() &&
2040 !fProofOutputDataset) {
2041 TString qaOutputFile(Form("%sMerged.%s.Data.root",
2042 fProofOutputLocation.Data(),
2043 AliQAv1::GetQADataFileName()));
2044 TProofOutputFile *qaProofFile = new TProofOutputFile(Form("Merged.%s.Data.root",
2045 AliQAv1::GetQADataFileName()));
2046 qaProofFile->SetOutputFileName(qaOutputFile.Data());
2047 if (AliDebugLevel() > 0) qaProofFile->Dump();
2048 fOutput->Add(qaProofFile);
2049 MergeQA(qaProofFile->GetFileName());
2060 if (!fProofOutputFileName.IsNull() &&
2061 !fProofOutputLocation.IsNull() &&
2062 fProofOutputDataset &&
2063 !fProofOutputArchive.IsNull()) {
2064 TProofOutputFile *zipProofFile = new TProofOutputFile(fProofOutputFileName.Data(),
2066 fProofOutputLocation.Data());
2067 if (AliDebugLevel() > 0) zipProofFile->Dump();
2068 fOutput->Add(zipProofFile);
2069 TString fileList(fProofOutputArchive.Data());
2070 fileList.ReplaceAll(","," ");
2072 #if ROOT_SVN_REVISION >= 30174
2073 command.Form("zip -n root %s/%s %s",zipProofFile->GetDir(kTRUE),zipProofFile->GetFileName(),fileList.Data());
2075 command.Form("zip -n root %s/%s %s",zipProofFile->GetDir(),zipProofFile->GetFileName(),fileList.Data());
2077 AliInfo(Form("Executing: %s",command.Data()));
2078 gSystem->Exec(command.Data());
2083 //_____________________________________________________________________________
2084 void AliReconstruction::Terminate()
2086 // Create tags for the events in the ESD tree (the ESD tree is always present)
2087 // In case of empty events the tags will contain dummy values
2088 AliCodeTimerAuto("",0);
2090 // Do not call the ESD tag creator in case of PROOF-based reconstruction
2092 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
2093 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData, AliQAv1::Instance()->GetQA(), AliQAv1::Instance()->GetEventSpecies(), AliQAv1::kNDET, AliRecoParam::kNSpecies);
2094 delete esdtagCreator;
2097 // Cleanup of CDB manager: cache and active storages!
2098 AliCDBManager::Instance()->ClearCache();
2101 //_____________________________________________________________________________
2102 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
2104 // run the local reconstruction
2106 static Int_t eventNr=0;
2107 AliCodeTimerAuto("",0)
2109 TString detStr = detectors;
2110 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2111 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2112 AliReconstructor* reconstructor = GetReconstructor(iDet);
2113 if (!reconstructor) continue;
2114 AliLoader* loader = fLoader[iDet];
2115 // Matthias April 2008: temporary fix to run HLT reconstruction
2116 // although the HLT loader is missing
2117 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
2119 reconstructor->Reconstruct(fRawReader, NULL);
2122 reconstructor->Reconstruct(dummy, NULL);
2127 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
2130 // conversion of digits
2131 if (fRawReader && reconstructor->HasDigitConversion()) {
2132 AliInfo(Form("converting raw data digits into root objects for %s",
2133 fgkDetectorName[iDet]));
2134 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
2135 // fgkDetectorName[iDet]),0);
2136 loader->LoadDigits("update");
2137 loader->CleanDigits();
2138 loader->MakeDigitsContainer();
2139 TTree* digitsTree = loader->TreeD();
2140 reconstructor->ConvertDigits(fRawReader, digitsTree);
2141 loader->WriteDigits("OVERWRITE");
2142 loader->UnloadDigits();
2144 // local reconstruction
2145 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2146 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]),0);
2147 loader->LoadRecPoints("update");
2148 loader->CleanRecPoints();
2149 loader->MakeRecPointsContainer();
2150 TTree* clustersTree = loader->TreeR();
2151 if (fRawReader && !reconstructor->HasDigitConversion()) {
2152 reconstructor->Reconstruct(fRawReader, clustersTree);
2154 loader->LoadDigits("read");
2155 TTree* digitsTree = loader->TreeD();
2157 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
2158 if (fStopOnError) return kFALSE;
2160 reconstructor->Reconstruct(digitsTree, clustersTree);
2161 if (fRunQA && IsInTasks(AliQAv1::kDIGITSR)) {
2162 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2163 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, digitsTree) ;
2166 loader->UnloadDigits();
2168 if (fRunQA && IsInTasks(AliQAv1::kRECPOINTS)) {
2169 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2170 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, clustersTree) ;
2172 loader->WriteRecPoints("OVERWRITE");
2173 loader->UnloadRecPoints();
2174 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
2176 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2177 AliError(Form("the following detectors were not found: %s",
2179 if (fStopOnError) return kFALSE;
2184 //_____________________________________________________________________________
2185 Bool_t AliReconstruction::RunSPDTrackleting(AliESDEvent*& esd)
2187 // run the SPD trackleting (for SPD efficiency purpouses)
2189 AliCodeTimerAuto("",0)
2191 Double_t vtxPos[3] = {0, 0, 0};
2192 Double_t vtxErr[3] = {0.0, 0.0, 0.0};
2194 TArrayF mcVertex(3);
2196 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
2197 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
2198 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
2201 const AliESDVertex *vertex = esd->GetVertex();
2203 AliWarning("Vertex not found");
2206 vertex->GetXYZ(vtxPos);
2207 vertex->GetSigmaXYZ(vtxErr);
2208 if (fSPDTrackleter) {
2209 AliInfo("running the SPD Trackleter for Plane Efficiency Evaluation");
2212 fLoader[0]->LoadRecPoints("read");
2213 TTree* tree = fLoader[0]->TreeR();
2215 AliError("Can't get the ITS cluster tree");
2218 fSPDTrackleter->LoadClusters(tree);
2219 fSPDTrackleter->SetVertex(vtxPos, vtxErr);
2221 if (fSPDTrackleter->Clusters2Tracks(esd) != 0) {
2222 AliError("AliITSTrackleterSPDEff Clusters2Tracks failed");
2223 // fLoader[0]->UnloadRecPoints();
2226 //fSPDTrackleter->UnloadRecPoints();
2228 AliWarning("SPDTrackleter not available");
2234 //_____________________________________________________________________________
2235 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
2237 // run the barrel tracking
2239 AliCodeTimerAuto("",0)
2241 AliVertexer *vertexer = CreateVertexer();
2242 if (!vertexer) return kFALSE;
2244 AliInfo("running the ITS vertex finder");
2245 AliESDVertex* vertex = NULL;
2247 fLoader[0]->LoadRecPoints();
2248 TTree* cltree = fLoader[0]->TreeR();
2250 if(fDiamondProfileSPD) vertexer->SetVtxStart(fDiamondProfileSPD);
2251 vertex = vertexer->FindVertexForCurrentEvent(cltree);
2254 AliError("Can't get the ITS cluster tree");
2256 fLoader[0]->UnloadRecPoints();
2259 AliError("Can't get the ITS loader");
2262 AliWarning("Vertex not found");
2263 vertex = new AliESDVertex();
2264 vertex->SetName("default");
2267 vertex->SetName("reconstructed");
2272 vertex->GetXYZ(vtxPos);
2273 vertex->GetSigmaXYZ(vtxErr);
2275 esd->SetPrimaryVertexSPD(vertex);
2276 AliESDVertex *vpileup = NULL;
2277 Int_t novertices = 0;
2278 vpileup = vertexer->GetAllVertices(novertices);
2280 for (Int_t kk=1; kk<novertices; kk++)esd->AddPileupVertexSPD(&vpileup[kk]);
2282 // if SPD multiplicity has been determined, it is stored in the ESD
2283 AliMultiplicity *mult = vertexer->GetMultiplicity();
2284 if(mult)esd->SetMultiplicity(mult);
2286 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2287 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
2296 //_____________________________________________________________________________
2297 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
2299 // run the HLT barrel tracking
2301 AliCodeTimerAuto("",0)
2304 AliError("Missing runLoader!");
2308 AliInfo("running HLT tracking");
2310 // Get a pointer to the HLT reconstructor
2311 AliReconstructor *reconstructor = GetReconstructor(kNDetectors-1);
2312 if (!reconstructor) return kFALSE;
2315 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2316 TString detName = fgkDetectorName[iDet];
2317 AliDebug(1, Form("%s HLT tracking", detName.Data()));
2318 reconstructor->SetOption(detName.Data());
2319 AliTracker *tracker = reconstructor->CreateTracker();
2321 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
2322 if (fStopOnError) return kFALSE;
2326 Double_t vtxErr[3]={0.005,0.005,0.010};
2327 const AliESDVertex *vertex = esd->GetVertex();
2328 vertex->GetXYZ(vtxPos);
2329 tracker->SetVertex(vtxPos,vtxErr);
2331 fLoader[iDet]->LoadRecPoints("read");
2332 TTree* tree = fLoader[iDet]->TreeR();
2334 AliError(Form("Can't get the %s cluster tree", detName.Data()));
2337 tracker->LoadClusters(tree);
2339 if (tracker->Clusters2Tracks(esd) != 0) {
2340 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
2344 tracker->UnloadClusters();
2352 //_____________________________________________________________________________
2353 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
2355 // run the muon spectrometer tracking
2357 AliCodeTimerAuto("",0)
2360 AliError("Missing runLoader!");
2363 Int_t iDet = 7; // for MUON
2365 AliInfo("is running...");
2367 // Get a pointer to the MUON reconstructor
2368 AliReconstructor *reconstructor = GetReconstructor(iDet);
2369 if (!reconstructor) return kFALSE;
2372 TString detName = fgkDetectorName[iDet];
2373 AliDebug(1, Form("%s tracking", detName.Data()));
2374 AliTracker *tracker = reconstructor->CreateTracker();
2376 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2381 fLoader[iDet]->LoadRecPoints("read");
2383 tracker->LoadClusters(fLoader[iDet]->TreeR());
2385 Int_t rv = tracker->Clusters2Tracks(esd);
2387 fLoader[iDet]->UnloadRecPoints();
2389 tracker->UnloadClusters();
2395 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2403 //_____________________________________________________________________________
2404 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2406 // run the barrel tracking
2407 static Int_t eventNr=0;
2408 AliCodeTimerAuto("",0)
2410 AliInfo("running tracking");
2412 // Set the event info which is used
2413 // by the trackers in order to obtain
2414 // information about read-out detectors,
2416 AliDebug(1, "Setting event info");
2417 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2418 if (!fTracker[iDet]) continue;
2419 fTracker[iDet]->SetEventInfo(&fEventInfo);
2422 //Fill the ESD with the T0 info (will be used by the TOF)
2423 if (fReconstructor[11] && fLoader[11]) {
2424 fLoader[11]->LoadRecPoints("READ");
2425 TTree *treeR = fLoader[11]->TreeR();
2427 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2431 // pass 1: TPC + ITS inwards
2432 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2433 if (!fTracker[iDet]) continue;
2434 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2437 fLoader[iDet]->LoadRecPoints("read");
2438 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2439 TTree* tree = fLoader[iDet]->TreeR();
2441 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2444 fTracker[iDet]->LoadClusters(tree);
2445 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2447 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2448 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2451 // preliminary PID in TPC needed by the ITS tracker
2453 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2454 AliESDpid::MakePID(esd);
2456 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2459 // pass 2: ALL backwards
2461 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2462 if (!fTracker[iDet]) continue;
2463 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2466 if (iDet > 1) { // all except ITS, TPC
2468 fLoader[iDet]->LoadRecPoints("read");
2469 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2470 tree = fLoader[iDet]->TreeR();
2472 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2475 fTracker[iDet]->LoadClusters(tree);
2476 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2480 if (iDet>1) // start filling residuals for the "outer" detectors
2482 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2483 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2485 AliRecoParam::EventSpecie_t es=fRecoParam.GetEventSpecie();
2486 TObjArray * elem = arr[AliRecoParam::AConvert(es)];
2487 if ( elem && (! elem->At(0)) ) {
2488 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2489 if (qadm) qadm->InitRecPointsForTracker() ;
2493 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2494 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2499 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2500 fTracker[iDet]->UnloadClusters();
2501 fLoader[iDet]->UnloadRecPoints();
2503 // updated PID in TPC needed by the ITS tracker -MI
2505 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2506 AliESDpid::MakePID(esd);
2508 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2510 //stop filling residuals for the "outer" detectors
2511 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2513 // pass 3: TRD + TPC + ITS refit inwards
2515 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2516 if (!fTracker[iDet]) continue;
2517 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2520 if (iDet<2) // start filling residuals for TPC and ITS
2522 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2523 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2525 AliRecoParam::EventSpecie_t es=fRecoParam.GetEventSpecie();
2526 TObjArray * elem = arr[AliRecoParam::AConvert(es)];
2527 if ( elem && (! elem->At(0)) ) {
2528 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2529 if (qadm) qadm->InitRecPointsForTracker() ;
2534 if (fTracker[iDet]->RefitInward(esd) != 0) {
2535 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2538 // run postprocessing
2539 if (fTracker[iDet]->PostProcess(esd) != 0) {
2540 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2543 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2546 // write space-points to the ESD in case alignment data output
2548 if (fWriteAlignmentData)
2549 WriteAlignmentData(esd);
2551 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2552 if (!fTracker[iDet]) continue;
2554 fTracker[iDet]->UnloadClusters();
2555 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2556 fLoader[iDet]->UnloadRecPoints();
2557 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2559 // stop filling residuals for TPC and ITS
2560 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2566 //_____________________________________________________________________________
2567 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2569 // Remove the data which are not needed for the physics analysis.
2572 Int_t nTracks=esd->GetNumberOfTracks();
2573 Int_t nV0s=esd->GetNumberOfV0s();
2575 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2577 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2578 Bool_t rc=esd->Clean(cleanPars);
2580 nTracks=esd->GetNumberOfTracks();
2581 nV0s=esd->GetNumberOfV0s();
2583 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2588 //_____________________________________________________________________________
2589 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2591 // fill the event summary data
2593 AliCodeTimerAuto("",0)
2594 static Int_t eventNr=0;
2595 TString detStr = detectors;
2597 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2598 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2599 AliReconstructor* reconstructor = GetReconstructor(iDet);
2600 if (!reconstructor) continue;
2601 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2602 TTree* clustersTree = NULL;
2603 if (fLoader[iDet]) {
2604 fLoader[iDet]->LoadRecPoints("read");
2605 clustersTree = fLoader[iDet]->TreeR();
2606 if (!clustersTree) {
2607 AliError(Form("Can't get the %s clusters tree",
2608 fgkDetectorName[iDet]));
2609 if (fStopOnError) return kFALSE;
2612 if (fRawReader && !reconstructor->HasDigitConversion()) {
2613 reconstructor->FillESD(fRawReader, clustersTree, esd);
2615 TTree* digitsTree = NULL;
2616 if (fLoader[iDet]) {
2617 fLoader[iDet]->LoadDigits("read");
2618 digitsTree = fLoader[iDet]->TreeD();
2620 AliError(Form("Can't get the %s digits tree",
2621 fgkDetectorName[iDet]));
2622 if (fStopOnError) return kFALSE;
2625 reconstructor->FillESD(digitsTree, clustersTree, esd);
2626 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2628 if (fLoader[iDet]) {
2629 fLoader[iDet]->UnloadRecPoints();
2633 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2634 AliError(Form("the following detectors were not found: %s",
2636 if (fStopOnError) return kFALSE;
2638 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2643 //_____________________________________________________________________________
2644 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2646 // Reads the trigger decision which is
2647 // stored in Trigger.root file and fills
2648 // the corresponding esd entries
2650 AliCodeTimerAuto("",0)
2652 AliInfo("Filling trigger information into the ESD");
2655 AliCTPRawStream input(fRawReader);
2656 if (!input.Next()) {
2657 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2660 if (esd->GetTriggerMask() != input.GetClassMask())
2661 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2662 input.GetClassMask(),esd->GetTriggerMask()));
2663 if (esd->GetOrbitNumber() != input.GetOrbitID())
2664 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2665 input.GetOrbitID(),esd->GetOrbitNumber()));
2666 if (esd->GetBunchCrossNumber() != input.GetBCID())
2667 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2668 input.GetBCID(),esd->GetBunchCrossNumber()));
2669 AliESDHeader* esdheader = esd->GetHeader();
2670 esdheader->SetL0TriggerInputs(input.GetL0Inputs());
2671 esdheader->SetL1TriggerInputs(input.GetL1Inputs());
2672 esdheader->SetL2TriggerInputs(input.GetL2Inputs());
2674 UInt_t orbit=input.GetOrbitID();
2675 for(Int_t i=0 ; i<input.GetNIRs() ; i++ )
2676 if(TMath::Abs(Int_t(orbit-(input.GetIR(i))->GetOrbit()))<=1){
2677 esdheader->AddTriggerIR(input.GetIR(i));
2683 //_____________________________________________________________________________
2684 Bool_t AliReconstruction::FillTriggerScalers(AliESDEvent*& esd)
2687 //fRunScalers->Print();
2688 if(fRunScalers && fRunScalers->CheckRunScalers()){
2689 AliTimeStamp* timestamp = new AliTimeStamp(esd->GetOrbitNumber(), esd->GetPeriodNumber(), esd->GetBunchCrossNumber());
2690 //AliTimeStamp* timestamp = new AliTimeStamp(10308000, 0, (ULong64_t)486238);
2691 AliESDHeader* esdheader = fesd->GetHeader();
2692 for(Int_t i=0;i<50;i++){
2693 if((1ull<<i) & esd->GetTriggerMask()){
2694 AliTriggerScalersESD* scalesd = fRunScalers->GetScalersForEventClass( timestamp, i+1);
2695 if(scalesd)esdheader->SetTriggerScalersRecord(scalesd);
2701 //_____________________________________________________________________________
2702 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2705 // Filling information from RawReader Header
2708 if (!fRawReader) return kFALSE;
2710 AliInfo("Filling information from RawReader Header");
2712 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2713 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2714 esd->SetPeriodNumber(fRawReader->GetPeriod());
2716 esd->SetTimeStamp(fRawReader->GetTimestamp());
2717 esd->SetEventType(fRawReader->GetType());
2723 //_____________________________________________________________________________
2724 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2726 // check whether detName is contained in detectors
2727 // if yes, it is removed from detectors
2729 // check if all detectors are selected
2730 if ((detectors.CompareTo("ALL") == 0) ||
2731 detectors.BeginsWith("ALL ") ||
2732 detectors.EndsWith(" ALL") ||
2733 detectors.Contains(" ALL ")) {
2738 // search for the given detector
2739 Bool_t result = kFALSE;
2740 if ((detectors.CompareTo(detName) == 0) ||
2741 detectors.BeginsWith(detName+" ") ||
2742 detectors.EndsWith(" "+detName) ||
2743 detectors.Contains(" "+detName+" ")) {
2744 detectors.ReplaceAll(detName, "");
2748 // clean up the detectors string
2749 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2750 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2751 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2756 //_____________________________________________________________________________
2757 Bool_t AliReconstruction::InitRunLoader()
2759 // get or create the run loader
2761 if (gAlice) delete gAlice;
2764 TFile *gafile = TFile::Open(fGAliceFileName.Data());
2765 // if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2766 if (gafile) { // galice.root exists
2770 // load all base libraries to get the loader classes
2771 TString libs = gSystem->GetLibraries();
2772 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2773 TString detName = fgkDetectorName[iDet];
2774 if (detName == "HLT") continue;
2775 if (libs.Contains("lib" + detName + "base.so")) continue;
2776 gSystem->Load("lib" + detName + "base.so");
2778 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2780 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2785 fRunLoader->CdGAFile();
2786 fRunLoader->LoadgAlice();
2788 //PH This is a temporary fix to give access to the kinematics
2789 //PH that is needed for the labels of ITS clusters
2790 fRunLoader->LoadHeader();
2791 fRunLoader->LoadKinematics();
2793 } else { // galice.root does not exist
2795 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2797 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2798 AliConfig::GetDefaultEventFolderName(),
2801 AliError(Form("could not create run loader in file %s",
2802 fGAliceFileName.Data()));
2806 fIsNewRunLoader = kTRUE;
2807 fRunLoader->MakeTree("E");
2809 if (fNumberOfEventsPerFile > 0)
2810 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2812 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2818 //_____________________________________________________________________________
2819 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2821 // get the reconstructor object and the loader for a detector
2823 if (fReconstructor[iDet]) {
2824 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2825 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2826 fReconstructor[iDet]->SetRecoParam(par);
2827 fReconstructor[iDet]->SetRunInfo(fRunInfo);
2829 return fReconstructor[iDet];
2832 // load the reconstructor object
2833 TPluginManager* pluginManager = gROOT->GetPluginManager();
2834 TString detName = fgkDetectorName[iDet];
2835 TString recName = "Ali" + detName + "Reconstructor";
2837 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2839 AliReconstructor* reconstructor = NULL;
2840 // first check if a plugin is defined for the reconstructor
2841 TPluginHandler* pluginHandler =
2842 pluginManager->FindHandler("AliReconstructor", detName);
2843 // if not, add a plugin for it
2844 if (!pluginHandler) {
2845 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2846 TString libs = gSystem->GetLibraries();
2847 if (libs.Contains("lib" + detName + "base.so") ||
2848 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2849 pluginManager->AddHandler("AliReconstructor", detName,
2850 recName, detName + "rec", recName + "()");
2852 pluginManager->AddHandler("AliReconstructor", detName,
2853 recName, detName, recName + "()");
2855 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2857 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2858 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2860 if (reconstructor) {
2861 TObject* obj = fOptions.FindObject(detName.Data());
2862 if (obj) reconstructor->SetOption(obj->GetTitle());
2863 reconstructor->SetRunInfo(fRunInfo);
2864 reconstructor->Init();
2865 fReconstructor[iDet] = reconstructor;
2868 // get or create the loader
2869 if (detName != "HLT") {
2870 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2871 if (!fLoader[iDet]) {
2872 AliConfig::Instance()
2873 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2875 // first check if a plugin is defined for the loader
2877 pluginManager->FindHandler("AliLoader", detName);
2878 // if not, add a plugin for it
2879 if (!pluginHandler) {
2880 TString loaderName = "Ali" + detName + "Loader";
2881 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2882 pluginManager->AddHandler("AliLoader", detName,
2883 loaderName, detName + "base",
2884 loaderName + "(const char*, TFolder*)");
2885 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2887 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2889 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2890 fRunLoader->GetEventFolder());
2892 if (!fLoader[iDet]) { // use default loader
2893 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2895 if (!fLoader[iDet]) {
2896 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2897 if (fStopOnError) return NULL;
2899 fRunLoader->AddLoader(fLoader[iDet]);
2900 fRunLoader->CdGAFile();
2901 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2902 fRunLoader->Write(0, TObject::kOverwrite);
2907 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2908 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2909 reconstructor->SetRecoParam(par);
2910 reconstructor->SetRunInfo(fRunInfo);
2912 return reconstructor;
2915 //_____________________________________________________________________________
2916 AliVertexer* AliReconstruction::CreateVertexer()
2918 // create the vertexer
2919 // Please note that the caller is the owner of the
2922 AliVertexer* vertexer = NULL;
2923 AliReconstructor* itsReconstructor = GetReconstructor(0);
2924 if (itsReconstructor && ((fRunLocalReconstruction.Contains("ITS")) || fRunTracking.Contains("ITS"))) {
2925 vertexer = itsReconstructor->CreateVertexer();
2928 AliWarning("couldn't create a vertexer for ITS");
2934 //_____________________________________________________________________________
2935 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2937 // create the trackers
2938 AliInfo("Creating trackers");
2940 TString detStr = detectors;
2941 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2942 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2943 AliReconstructor* reconstructor = GetReconstructor(iDet);
2944 if (!reconstructor) continue;
2945 TString detName = fgkDetectorName[iDet];
2946 if (detName == "HLT") {
2947 fRunHLTTracking = kTRUE;
2950 if (detName == "MUON") {
2951 fRunMuonTracking = kTRUE;
2956 fTracker[iDet] = reconstructor->CreateTracker();
2957 if (!fTracker[iDet] && (iDet < 7)) {
2958 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2959 if (fStopOnError) return kFALSE;
2961 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2967 //_____________________________________________________________________________
2968 void AliReconstruction::CleanUp()
2970 // delete trackers and the run loader and close and delete the file
2972 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2973 delete fReconstructor[iDet];
2974 fReconstructor[iDet] = NULL;
2975 fLoader[iDet] = NULL;
2976 delete fTracker[iDet];
2977 fTracker[iDet] = NULL;
2982 delete fSPDTrackleter;
2983 fSPDTrackleter = NULL;
2992 delete fParentRawReader;
2993 fParentRawReader=NULL;
3001 if (AliQAManager::QAManager())
3002 AliQAManager::QAManager()->ShowQA() ;
3003 AliQAManager::Destroy() ;
3007 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
3009 // Write space-points which are then used in the alignment procedures
3010 // For the moment only ITS, TPC, TRD and TOF
3012 Int_t ntracks = esd->GetNumberOfTracks();
3013 for (Int_t itrack = 0; itrack < ntracks; itrack++)
3015 AliESDtrack *track = esd->GetTrack(itrack);
3018 for (Int_t i=0; i<200; ++i) idx[i] = -1; //PH avoid uninitialized values
3019 for (Int_t iDet = 5; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
3020 nsp += track->GetNcls(iDet);
3022 if (iDet==0) { // ITS "extra" clusters
3023 track->GetClusters(iDet,idx);
3024 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
3029 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
3030 track->SetTrackPointArray(sp);
3032 for (Int_t iDet = 5; iDet >= 0; iDet--) {
3033 AliTracker *tracker = fTracker[iDet];
3034 if (!tracker) continue;
3035 Int_t nspdet = track->GetClusters(iDet,idx);
3037 if (iDet==0) // ITS "extra" clusters
3038 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
3040 if (nspdet <= 0) continue;
3044 while (isp2 < nspdet) {
3045 Bool_t isvalid=kTRUE;
3047 Int_t index=idx[isp++];
3048 if (index < 0) continue;
3050 TString dets = fgkDetectorName[iDet];
3051 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
3052 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
3053 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
3054 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
3055 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
3057 isvalid = tracker->GetTrackPoint(index,p);
3060 if (!isvalid) continue;
3061 if (iDet==0 && (isp-1)>=6) p.SetExtra();
3062 sp->AddPoint(isptrack,&p); isptrack++;
3069 //_____________________________________________________________________________
3070 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
3072 // The method reads the raw-data error log
3073 // accumulated within the rawReader.
3074 // It extracts the raw-data errors related to
3075 // the current event and stores them into
3076 // a TClonesArray inside the esd object.
3078 if (!fRawReader) return;
3080 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
3082 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
3084 if (iEvent != log->GetEventNumber()) continue;
3086 esd->AddRawDataErrorLog(log);
3091 //_____________________________________________________________________________
3092 void AliReconstruction::CheckQA()
3094 // check the QA of SIM for this run and remove the detectors
3095 // with status Fatal
3097 // TString newRunLocalReconstruction ;
3098 // TString newRunTracking ;
3099 // TString newFillESD ;
3101 // for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
3102 // TString detName(AliQAv1::GetDetName(iDet)) ;
3103 // AliQAv1 * qa = AliQAv1::Instance(AliQAv1::DETECTORINDEX_t(iDet)) ;
3104 // if ( qa->IsSet(AliQAv1::DETECTORINDEX_t(iDet), AliQAv1::kSIM, specie, AliQAv1::kFATAL)) {
3105 // AliInfo(Form("QA status for %s %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed",
3106 // detName.Data(), AliRecoParam::GetEventSpecieName(es))) ;
3108 // if ( fRunLocalReconstruction.Contains(AliQAv1::GetDetName(iDet)) ||
3109 // fRunLocalReconstruction.Contains("ALL") ) {
3110 // newRunLocalReconstruction += detName ;
3111 // newRunLocalReconstruction += " " ;
3113 // if ( fRunTracking.Contains(AliQAv1::GetDetName(iDet)) ||
3114 // fRunTracking.Contains("ALL") ) {
3115 // newRunTracking += detName ;
3116 // newRunTracking += " " ;
3118 // if ( fFillESD.Contains(AliQAv1::GetDetName(iDet)) ||
3119 // fFillESD.Contains("ALL") ) {
3120 // newFillESD += detName ;
3121 // newFillESD += " " ;
3125 // fRunLocalReconstruction = newRunLocalReconstruction ;
3126 // fRunTracking = newRunTracking ;
3127 // fFillESD = newFillESD ;
3130 //_____________________________________________________________________________
3131 Int_t AliReconstruction::GetDetIndex(const char* detector)
3133 // return the detector index corresponding to detector
3135 for (index = 0; index < kNDetectors ; index++) {
3136 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
3141 //_____________________________________________________________________________
3142 Bool_t AliReconstruction::FinishPlaneEff() {
3144 // Here execute all the necessary operationis, at the end of the tracking phase,
3145 // in case that evaluation of PlaneEfficiencies was required for some detector.
3146 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
3148 // This Preliminary version works only FOR ITS !!!!!
3149 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3152 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3155 TString detStr = fLoadCDB;
3156 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
3157 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
3158 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3159 if(fTracker[iDet] && fTracker[iDet]->GetPlaneEff()) {
3160 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
3161 TString name=planeeff->GetName();
3163 TFile* pefile = TFile::Open(name, "RECREATE");
3164 ret=(Bool_t)planeeff->Write();
3166 if(planeeff->GetCreateHistos()) {
3167 TString hname=planeeff->GetName();
3168 hname+="Histo.root";
3169 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
3172 if(fSPDTrackleter) {
3173 AliPlaneEff *planeeff=fSPDTrackleter->GetPlaneEff();
3174 TString name="AliITSPlaneEffSPDtracklet.root";
3175 TFile* pefile = TFile::Open(name, "RECREATE");
3176 ret=(Bool_t)planeeff->Write();
3178 AliESDEvent *dummy=NULL;
3179 ret=(Bool_t)fSPDTrackleter->PostProcess(dummy); // take care of writing other files
3184 //_____________________________________________________________________________
3185 Bool_t AliReconstruction::InitPlaneEff() {
3187 // Here execute all the necessary operations, before of the tracking phase,
3188 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
3189 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
3190 // which should be updated/recalculated.
3192 // This Preliminary version will work only FOR ITS !!!!!
3193 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3196 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3198 AliWarning(Form("Implementation of this method not yet completed !! Method return kTRUE"));
3200 fSPDTrackleter = NULL;
3201 TString detStr = fLoadCDB;
3202 if (IsSelected(fgkDetectorName[0], detStr)) {
3203 AliReconstructor* itsReconstructor = GetReconstructor(0);
3204 if (itsReconstructor) {
3205 fSPDTrackleter = itsReconstructor->CreateTrackleter(); // this is NULL unless required in RecoParam
3207 if (fSPDTrackleter) {
3208 AliInfo("Trackleter for SPD has been created");
3214 //_____________________________________________________________________________
3215 Bool_t AliReconstruction::InitAliEVE()
3217 // This method should be called only in case
3218 // AliReconstruction is run
3219 // within the alieve environment.
3220 // It will initialize AliEVE in a way
3221 // so that it can visualize event processed
3222 // by AliReconstruction.
3223 // The return flag shows whenever the
3224 // AliEVE initialization was successful or not.
3227 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
3228 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
3229 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
3231 gROOT->ProcessLine("if (!AliEveEventManager::GetMaster()){new AliEveEventManager();AliEveEventManager::GetMaster()->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(AliEveEventManager::GetMaster());};");
3232 gROOT->ProcessLine("alieve_online_init()");
3237 //_____________________________________________________________________________
3238 void AliReconstruction::RunAliEVE()
3240 // Runs AliEVE visualisation of
3241 // the current event.
3242 // Should be executed only after
3243 // successful initialization of AliEVE.
3245 AliInfo("Running AliEVE...");
3246 gROOT->ProcessLine(Form("AliEveEventManager::GetMaster()->SetEvent((AliRunLoader*)0x%lx,(AliRawReader*)0x%lx,(AliESDEvent*)0x%lx,(AliESDfriend*)0x%lx);",fRunLoader,fRawReader,fesd,fesdf));
3250 //_____________________________________________________________________________
3251 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
3253 // Allows to run QA for a selected set of detectors
3254 // and a selected set of tasks among RAWS, DIGITSR, RECPOINTS and ESDS
3255 // all selected detectors run the same selected tasks
3257 if (!detAndAction.Contains(":")) {
3258 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
3262 Int_t colon = detAndAction.Index(":") ;
3263 fQADetectors = detAndAction(0, colon) ;
3264 if (fQADetectors.Contains("ALL") ){
3265 TString tmp = fFillESD ;
3266 Int_t minus = fQADetectors.Last('-') ;
3267 TString toKeep = fFillESD ;
3268 TString toRemove("") ;
3269 while (minus >= 0) {
3270 toRemove = fQADetectors(minus+1, fQADetectors.Length()) ;
3271 toRemove = toRemove.Strip() ;
3272 toKeep.ReplaceAll(toRemove, "") ;
3273 fQADetectors.ReplaceAll(Form("-%s", toRemove.Data()), "") ;
3274 minus = fQADetectors.Last('-') ;
3276 fQADetectors = toKeep ;
3278 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
3279 if (fQATasks.Contains("ALL") ) {
3280 fQATasks = Form("%d %d %d %d", AliQAv1::kRAWS, AliQAv1::kDIGITSR, AliQAv1::kRECPOINTS, AliQAv1::kESDS) ;
3282 fQATasks.ToUpper() ;
3284 if ( fQATasks.Contains("RAW") )
3285 tempo = Form("%d ", AliQAv1::kRAWS) ;
3286 if ( fQATasks.Contains("DIGIT") )
3287 tempo += Form("%d ", AliQAv1::kDIGITSR) ;
3288 if ( fQATasks.Contains("RECPOINT") )
3289 tempo += Form("%d ", AliQAv1::kRECPOINTS) ;
3290 if ( fQATasks.Contains("ESD") )
3291 tempo += Form("%d ", AliQAv1::kESDS) ;
3293 if (fQATasks.IsNull()) {
3294 AliInfo("No QA requested\n") ;
3299 TString tempo(fQATasks) ;
3300 tempo.ReplaceAll(Form("%d", AliQAv1::kRAWS), AliQAv1::GetTaskName(AliQAv1::kRAWS)) ;
3301 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITSR), AliQAv1::GetTaskName(AliQAv1::kDIGITSR)) ;
3302 tempo.ReplaceAll(Form("%d", AliQAv1::kRECPOINTS), AliQAv1::GetTaskName(AliQAv1::kRECPOINTS)) ;
3303 tempo.ReplaceAll(Form("%d", AliQAv1::kESDS), AliQAv1::GetTaskName(AliQAv1::kESDS)) ;
3304 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
3309 //_____________________________________________________________________________
3310 Bool_t AliReconstruction::InitRecoParams()
3312 // The method accesses OCDB and retrieves all
3313 // the available reco-param objects from there.
3315 Bool_t isOK = kTRUE;
3317 if (fRecoParam.GetDetRecoParamArray(kNDetectors)) {
3318 AliInfo("Using custom GRP reconstruction parameters");
3321 AliInfo("Loading GRP reconstruction parameter objects");
3323 AliCDBPath path("GRP","Calib","RecoParam");
3324 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3326 AliWarning("Couldn't find GRP RecoParam entry in OCDB");
3330 TObject *recoParamObj = entry->GetObject();
3331 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3332 // GRP has a normal TobjArray of AliDetectorRecoParam objects
3333 // Registering them in AliRecoParam
3334 fRecoParam.AddDetRecoParamArray(kNDetectors,dynamic_cast<TObjArray*>(recoParamObj));
3336 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3337 // GRP has only onse set of reco parameters
3338 // Registering it in AliRecoParam
3339 AliInfo("Single set of GRP reconstruction parameters found");
3340 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3341 fRecoParam.AddDetRecoParam(kNDetectors,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3344 AliError("No valid GRP RecoParam object found in the OCDB");
3351 TString detStr = fLoadCDB;
3352 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3354 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3356 if (fRecoParam.GetDetRecoParamArray(iDet)) {
3357 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
3361 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
3363 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
3364 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3366 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
3370 TObject *recoParamObj = entry->GetObject();
3371 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3372 // The detector has a normal TobjArray of AliDetectorRecoParam objects
3373 // Registering them in AliRecoParam
3374 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
3376 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3377 // The detector has only onse set of reco parameters
3378 // Registering it in AliRecoParam
3379 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
3380 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3381 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3384 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
3388 // FIX ME: We have to disable the unloading of reco-param CDB
3389 // entries because QA framework is using them. Has to be fix in
3390 // a way that the QA takes the objects already constructed in
3392 // AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
3396 if (AliDebugLevel() > 0) fRecoParam.Print();
3401 //_____________________________________________________________________________
3402 Bool_t AliReconstruction::GetEventInfo()
3404 // Fill the event info object
3406 AliCodeTimerAuto("",0)
3408 AliCentralTrigger *aCTP = NULL;
3410 fEventInfo.SetEventType(fRawReader->GetType());
3412 ULong64_t mask = fRawReader->GetClassMask();
3413 fEventInfo.SetTriggerMask(mask);
3414 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
3415 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
3417 aCTP = new AliCentralTrigger();
3418 TString configstr("");
3419 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
3420 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
3424 aCTP->SetClassMask(mask);
3425 aCTP->SetClusterMask(clmask);
3428 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
3430 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
3431 aCTP = fRunLoader->GetTrigger();
3432 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
3433 // get inputs from actp - just get
3434 AliESDHeader* esdheader = fesd->GetHeader();
3435 esdheader->SetL0TriggerInputs(aCTP->GetL0TriggerInputs());
3436 esdheader->SetL1TriggerInputs(aCTP->GetL1TriggerInputs());
3437 esdheader->SetL2TriggerInputs(aCTP->GetL2TriggerInputs());
3438 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
3441 AliWarning("No trigger can be loaded! The trigger information will not be used!");
3446 AliTriggerConfiguration *config = aCTP->GetConfiguration();
3448 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
3449 if (fRawReader) delete aCTP;
3453 UChar_t clustmask = 0;
3455 ULong64_t trmask = fEventInfo.GetTriggerMask();
3456 const TObjArray& classesArray = config->GetClasses();
3457 Int_t nclasses = classesArray.GetEntriesFast();
3458 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
3459 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
3461 Int_t trindex = TMath::Nint(TMath::Log2(trclass->GetMask()));
3462 fesd->SetTriggerClass(trclass->GetName(),trindex);
3463 if (fRawReader) fRawReader->LoadTriggerClass(trclass->GetName(),trindex);
3464 if (trmask & (1ull << trindex)) {
3466 trclasses += trclass->GetName();
3468 clustmask |= trclass->GetCluster()->GetClusterMask();
3472 fEventInfo.SetTriggerClasses(trclasses);
3474 // Set the information in ESD
3475 fesd->SetTriggerMask(trmask);
3476 fesd->SetTriggerCluster(clustmask);
3478 if (!aCTP->CheckTriggeredDetectors()) {
3479 if (fRawReader) delete aCTP;
3483 if (fRawReader) delete aCTP;
3485 // We have to fill also the HLT decision here!!
3491 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3493 // Match the detector list found in the rec.C or the default 'ALL'
3494 // to the list found in the GRP (stored there by the shuttle PP which
3495 // gets the information from ECS)
3496 static TString resultList;
3497 TString detList = detectorList;
3501 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3502 if ((detectorMask >> iDet) & 0x1) {
3503 TString det = AliDAQ::OfflineModuleName(iDet);
3504 if ((detList.CompareTo("ALL") == 0) ||
3505 ((detList.BeginsWith("ALL ") ||
3506 detList.EndsWith(" ALL") ||
3507 detList.Contains(" ALL ")) &&
3508 !(detList.BeginsWith("-"+det+" ") ||
3509 detList.EndsWith(" -"+det) ||
3510 detList.Contains(" -"+det+" "))) ||
3511 (detList.CompareTo(det) == 0) ||
3512 detList.BeginsWith(det+" ") ||
3513 detList.EndsWith(" "+det) ||
3514 detList.Contains( " "+det+" " )) {
3515 if (!resultList.EndsWith(det + " ")) {
3524 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3525 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3526 if ((detList.CompareTo("ALL") == 0) ||
3527 ((detList.BeginsWith("ALL ") ||
3528 detList.EndsWith(" ALL") ||
3529 detList.Contains(" ALL ")) &&
3530 !(detList.BeginsWith("-"+hltDet+" ") ||
3531 detList.EndsWith(" -"+hltDet) ||
3532 detList.Contains(" -"+hltDet+" "))) ||
3533 (detList.CompareTo(hltDet) == 0) ||
3534 detList.BeginsWith(hltDet+" ") ||
3535 detList.EndsWith(" "+hltDet) ||
3536 detList.Contains( " "+hltDet+" " )) {
3537 resultList += hltDet;
3541 return resultList.Data();
3545 //______________________________________________________________________________
3546 void AliReconstruction::Abort(const char *method, EAbort what)
3548 // Abort processing. If what = kAbortProcess, the Process() loop will be
3549 // aborted. If what = kAbortFile, the current file in a chain will be
3550 // aborted and the processing will continue with the next file, if there
3551 // is no next file then Process() will be aborted. Abort() can also be
3552 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3553 // the SlaveTerminate() and Terminate() are always called. The abort flag
3554 // can be checked in these methods using GetAbort().
3556 // The method is overwritten in AliReconstruction for better handling of
3557 // reco specific errors
3559 if (!fStopOnError) return;
3563 TString whyMess = method;
3564 whyMess += " failed! Aborting...";
3566 AliError(whyMess.Data());
3569 TString mess = "Abort";
3570 if (fAbort == kAbortProcess)
3571 mess = "AbortProcess";
3572 else if (fAbort == kAbortFile)
3575 Info(mess, whyMess.Data());
3578 //______________________________________________________________________________
3579 Bool_t AliReconstruction::ProcessEvent(void* event)
3581 // Method that is used in case the event loop
3582 // is steered from outside, for example by AMORE
3583 // 'event' is a pointer to the DATE event in the memory
3585 if (fRawReader) delete fRawReader;
3586 fRawReader = new AliRawReaderDate(event);
3587 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
3594 //______________________________________________________________________________
3595 Bool_t AliReconstruction::ParseOutput()
3597 // The method parses the output file
3598 // location string in order to steer
3599 // properly the selector
3601 TPMERegexp re1("(\\w+\\.zip#\\w+\\.root):([,*\\w+\\.root,*]+)@dataset://(\\w++)");
3602 TPMERegexp re2("(\\w+\\.root)?@?dataset://(\\w++)");
3604 if (re1.Match(fESDOutput) == 4) {
3605 // root archive with output files stored and regustered
3607 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",re1[1].Data()));
3608 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",re1[3].Data()));
3609 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_DATASET",""));
3610 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_ARCHIVE",re1[2].Data()));
3611 AliInfo(Form("%s files will be stored within %s in dataset %s",
3616 else if (re2.Match(fESDOutput) == 3) {
3617 // output file stored and registered
3619 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",(re2[1].IsNull()) ? "AliESDs.root" : re2[1].Data()));
3620 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",re2[2].Data()));
3621 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_DATASET",""));
3622 AliInfo(Form("%s will be stored in dataset %s",
3623 (re2[1].IsNull()) ? "AliESDs.root" : re2[1].Data(),
3627 if (fESDOutput.IsNull()) {
3628 // Output location not given.
3629 // Assuming xrootd has been already started and
3630 // the output file has to be sent back
3631 // to the client machine
3632 TString esdUrl(Form("root://%s/%s/",
3633 TUrl(gSystem->HostName()).GetHostFQDN(),
3635 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE","AliESDs.root"));
3636 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",esdUrl.Data()));
3637 AliInfo(Form("AliESDs.root will be stored in %s",
3641 // User specified an output location.
3642 // Ones has just to parse it here
3643 TUrl outputUrl(fESDOutput.Data());
3644 TString outputFile(gSystem->BaseName(outputUrl.GetFile()));
3645 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",outputFile.IsNull() ? "AliESDs.root" : outputFile.Data()));
3646 TString outputLocation(outputUrl.GetUrl());
3647 outputLocation.ReplaceAll(outputFile.Data(),"");
3648 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",outputLocation.Data()));
3649 AliInfo(Form("%s will be stored in %s",
3650 outputFile.IsNull() ? "AliESDs.root" : outputFile.Data(),
3651 outputLocation.Data()));