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"
193 ClassImp(AliReconstruction)
195 //_____________________________________________________________________________
196 const char* AliReconstruction::fgkDetectorName[AliReconstruction::kNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
198 //_____________________________________________________________________________
199 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
201 fRunVertexFinder(kTRUE),
202 fRunVertexFinderTracks(kTRUE),
203 fRunHLTTracking(kFALSE),
204 fRunMuonTracking(kFALSE),
206 fRunCascadeFinder(kTRUE),
207 fStopOnError(kFALSE),
208 fWriteAlignmentData(kFALSE),
209 fWriteESDfriend(kFALSE),
210 fFillTriggerESD(kTRUE),
218 fRunLocalReconstruction("ALL"),
222 fUseTrackingErrorsForAlignment(""),
223 fGAliceFileName(gAliceFilename),
226 fProofOutputFileName(""),
227 fProofOutputLocation(""),
228 fProofOutputDataset(kFALSE),
229 fProofOutputArchive(""),
233 fNumberOfEventsPerFile((UInt_t)-1),
235 fLoadAlignFromCDB(kTRUE),
236 fLoadAlignData("ALL"),
244 fParentRawReader(NULL),
248 fSPDTrackleter(NULL),
250 fDiamondProfileSPD(NULL),
251 fDiamondProfile(NULL),
252 fDiamondProfileTPC(NULL),
253 fListOfCosmicTriggers(NULL),
257 fAlignObjArray(NULL),
261 fInitCDBCalled(kFALSE),
262 fSetRunNumberFromDataCalled(kFALSE),
267 fSameQACycle(kFALSE),
268 fInitQACalled(kFALSE),
269 fWriteQAExpertData(kTRUE),
270 fRunPlaneEff(kFALSE),
279 fIsNewRunLoader(kFALSE),
283 // create reconstruction object with default parameters
286 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
287 fReconstructor[iDet] = NULL;
288 fLoader[iDet] = NULL;
289 fTracker[iDet] = NULL;
291 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
292 fQACycles[iDet] = 999999 ;
293 fQAWriteExpert[iDet] = kFALSE ;
299 //_____________________________________________________________________________
300 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
302 fRunVertexFinder(rec.fRunVertexFinder),
303 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
304 fRunHLTTracking(rec.fRunHLTTracking),
305 fRunMuonTracking(rec.fRunMuonTracking),
306 fRunV0Finder(rec.fRunV0Finder),
307 fRunCascadeFinder(rec.fRunCascadeFinder),
308 fStopOnError(rec.fStopOnError),
309 fWriteAlignmentData(rec.fWriteAlignmentData),
310 fWriteESDfriend(rec.fWriteESDfriend),
311 fFillTriggerESD(rec.fFillTriggerESD),
313 fCleanESD(rec.fCleanESD),
314 fV0DCAmax(rec.fV0DCAmax),
315 fV0CsPmin(rec.fV0CsPmin),
319 fRunLocalReconstruction(rec.fRunLocalReconstruction),
320 fRunTracking(rec.fRunTracking),
321 fFillESD(rec.fFillESD),
322 fLoadCDB(rec.fLoadCDB),
323 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
324 fGAliceFileName(rec.fGAliceFileName),
325 fRawInput(rec.fRawInput),
326 fESDOutput(rec.fESDOutput),
327 fProofOutputFileName(rec.fProofOutputFileName),
328 fProofOutputLocation(rec.fProofOutputLocation),
329 fProofOutputDataset(rec.fProofOutputDataset),
330 fProofOutputArchive(rec.fProofOutputArchive),
331 fEquipIdMap(rec.fEquipIdMap),
332 fFirstEvent(rec.fFirstEvent),
333 fLastEvent(rec.fLastEvent),
334 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
336 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
337 fLoadAlignData(rec.fLoadAlignData),
338 fUseHLTData(rec.fUseHLTData),
345 fParentRawReader(NULL),
347 fRecoParam(rec.fRecoParam),
349 fSPDTrackleter(NULL),
351 fDiamondProfileSPD(rec.fDiamondProfileSPD),
352 fDiamondProfile(rec.fDiamondProfile),
353 fDiamondProfileTPC(rec.fDiamondProfileTPC),
354 fListOfCosmicTriggers(NULL),
358 fAlignObjArray(rec.fAlignObjArray),
359 fCDBUri(rec.fCDBUri),
360 fQARefUri(rec.fQARefUri),
362 fInitCDBCalled(rec.fInitCDBCalled),
363 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
364 fQADetectors(rec.fQADetectors),
365 fQATasks(rec.fQATasks),
367 fRunGlobalQA(rec.fRunGlobalQA),
368 fSameQACycle(rec.fSameQACycle),
369 fInitQACalled(rec.fInitQACalled),
370 fWriteQAExpertData(rec.fWriteQAExpertData),
371 fRunPlaneEff(rec.fRunPlaneEff),
380 fIsNewRunLoader(rec.fIsNewRunLoader),
386 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
387 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
389 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
390 fReconstructor[iDet] = NULL;
391 fLoader[iDet] = NULL;
392 fTracker[iDet] = NULL;
395 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
396 fQACycles[iDet] = rec.fQACycles[iDet];
397 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
400 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
401 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
405 //_____________________________________________________________________________
406 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
408 // assignment operator
409 // Used in PROOF mode
410 // Be very careful while modifing it!
411 // Simple rules to follow:
412 // for persistent data members - use their assignment operators
413 // for non-persistent ones - do nothing or take the default values from constructor
414 // TSelector members should not be touched
415 if(&rec == this) return *this;
417 fRunVertexFinder = rec.fRunVertexFinder;
418 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
419 fRunHLTTracking = rec.fRunHLTTracking;
420 fRunMuonTracking = rec.fRunMuonTracking;
421 fRunV0Finder = rec.fRunV0Finder;
422 fRunCascadeFinder = rec.fRunCascadeFinder;
423 fStopOnError = rec.fStopOnError;
424 fWriteAlignmentData = rec.fWriteAlignmentData;
425 fWriteESDfriend = rec.fWriteESDfriend;
426 fFillTriggerESD = rec.fFillTriggerESD;
428 fCleanESD = rec.fCleanESD;
429 fV0DCAmax = rec.fV0DCAmax;
430 fV0CsPmin = rec.fV0CsPmin;
434 fRunLocalReconstruction = rec.fRunLocalReconstruction;
435 fRunTracking = rec.fRunTracking;
436 fFillESD = rec.fFillESD;
437 fLoadCDB = rec.fLoadCDB;
438 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
439 fGAliceFileName = rec.fGAliceFileName;
440 fRawInput = rec.fRawInput;
441 fESDOutput = rec.fESDOutput;
442 fProofOutputFileName = rec.fProofOutputFileName;
443 fProofOutputLocation = rec.fProofOutputLocation;
444 fProofOutputDataset = rec.fProofOutputDataset;
445 fProofOutputArchive = rec.fProofOutputArchive;
446 fEquipIdMap = rec.fEquipIdMap;
447 fFirstEvent = rec.fFirstEvent;
448 fLastEvent = rec.fLastEvent;
449 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
451 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
452 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
455 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
456 fLoadAlignData = rec.fLoadAlignData;
457 fUseHLTData = rec.fUseHLTData;
459 delete fRunInfo; fRunInfo = NULL;
460 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
462 fEventInfo = rec.fEventInfo;
464 delete fRunScalers; fRunScalers = NULL;
465 if (rec.fRunScalers) fRunScalers = new AliTriggerRunScalers(*rec.fRunScalers);
470 fParentRawReader = NULL;
472 fRecoParam = rec.fRecoParam;
474 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
475 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
476 delete fLoader[iDet]; fLoader[iDet] = NULL;
477 delete fTracker[iDet]; fTracker[iDet] = NULL;
480 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
481 fQACycles[iDet] = rec.fQACycles[iDet];
482 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
485 delete fSPDTrackleter; fSPDTrackleter = NULL;
487 delete fDiamondProfileSPD; fDiamondProfileSPD = NULL;
488 if (rec.fDiamondProfileSPD) fDiamondProfileSPD = new AliESDVertex(*rec.fDiamondProfileSPD);
489 delete fDiamondProfile; fDiamondProfile = NULL;
490 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
491 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
492 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
494 delete fListOfCosmicTriggers; fListOfCosmicTriggers = NULL;
495 if (rec.fListOfCosmicTriggers) fListOfCosmicTriggers = (THashTable*)((rec.fListOfCosmicTriggers)->Clone());
497 delete fGRPData; fGRPData = NULL;
498 // if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
499 if (rec.fGRPData) fGRPData = (AliGRPObject*)((rec.fGRPData)->Clone());
501 delete fAlignObjArray; fAlignObjArray = NULL;
504 fQARefUri = rec.fQARefUri;
505 fSpecCDBUri.Delete();
506 fInitCDBCalled = rec.fInitCDBCalled;
507 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
508 fQADetectors = rec.fQADetectors;
509 fQATasks = rec.fQATasks;
511 fRunGlobalQA = rec.fRunGlobalQA;
512 fSameQACycle = rec.fSameQACycle;
513 fInitQACalled = rec.fInitQACalled;
514 fWriteQAExpertData = rec.fWriteQAExpertData;
515 fRunPlaneEff = rec.fRunPlaneEff;
524 fIsNewRunLoader = rec.fIsNewRunLoader;
531 //_____________________________________________________________________________
532 AliReconstruction::~AliReconstruction()
537 if (fListOfCosmicTriggers) {
538 fListOfCosmicTriggers->Delete();
539 delete fListOfCosmicTriggers;
544 if (fAlignObjArray) {
545 fAlignObjArray->Delete();
546 delete fAlignObjArray;
548 fSpecCDBUri.Delete();
550 AliCodeTimer::Instance()->Print();
553 //_____________________________________________________________________________
554 void AliReconstruction::InitQA()
556 //Initialize the QA and start of cycle
557 AliCodeTimerAuto("");
559 if (fInitQACalled) return;
560 fInitQACalled = kTRUE;
562 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kRECMODE) ;
563 if (fWriteQAExpertData)
564 qam->SetWriteExpert() ;
566 if (qam->IsDefaultStorageSet()) {
567 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
568 AliWarning("Default QA reference storage has been already set !");
569 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fQARefUri.Data()));
570 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
571 fQARefUri = qam->GetDefaultStorage()->GetURI();
573 if (fQARefUri.Length() > 0) {
574 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
575 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
576 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
578 fQARefUri="local://$ALICE_ROOT/QAref";
579 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
580 AliWarning("Default QA refeference storage not yet set !!!!");
581 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
582 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
585 qam->SetDefaultStorage(fQARefUri);
589 qam->SetActiveDetectors(fQADetectors) ;
590 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++) {
591 qam->SetCycleLength(AliQAv1::DETECTORINDEX_t(det), fQACycles[det]) ;
592 qam->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
594 if (!fRawReader && !fInput && IsInTasks(AliQAv1::kRAWS))
595 fQATasks.ReplaceAll(Form("%d",AliQAv1::kRAWS), "") ;
596 qam->SetTasks(fQATasks) ;
597 qam->InitQADataMaker(AliCDBManager::Instance()->GetRun()) ;
600 Bool_t sameCycle = kFALSE ;
601 AliQADataMaker *qadm = qam->GetQADataMaker(AliQAv1::kGLOBAL);
602 AliInfo(Form("Initializing the global QA data maker"));
603 if (IsInTasks(AliQAv1::kRECPOINTS)) {
604 qadm->StartOfCycle(AliQAv1::kRECPOINTS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
605 TObjArray **arr=qadm->Init(AliQAv1::kRECPOINTS);
606 AliTracker::SetResidualsArray(arr);
609 if (IsInTasks(AliQAv1::kESDS)) {
610 qadm->StartOfCycle(AliQAv1::kESDS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
611 qadm->Init(AliQAv1::kESDS);
614 AliSysInfo::AddStamp("InitQA") ;
617 //_____________________________________________________________________________
618 void AliReconstruction::MergeQA(const char *fileName)
620 //Initialize the QA and start of cycle
621 AliCodeTimerAuto("") ;
622 AliQAManager::QAManager()->Merge(AliCDBManager::Instance()->GetRun(),fileName) ;
623 AliSysInfo::AddStamp("MergeQA") ;
626 //_____________________________________________________________________________
627 void AliReconstruction::InitCDB()
629 // activate a default CDB storage
630 // First check if we have any CDB storage set, because it is used
631 // to retrieve the calibration and alignment constants
632 AliCodeTimerAuto("");
634 if (fInitCDBCalled) return;
635 fInitCDBCalled = kTRUE;
637 AliCDBManager* man = AliCDBManager::Instance();
638 if (man->IsDefaultStorageSet())
640 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
641 AliWarning("Default CDB storage has been already set !");
642 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
643 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
644 fCDBUri = man->GetDefaultStorage()->GetURI();
647 if (fCDBUri.Length() > 0)
649 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
650 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
651 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
653 fCDBUri="local://$ALICE_ROOT/OCDB";
654 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
655 AliWarning("Default CDB storage not yet set !!!!");
656 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
657 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
660 man->SetDefaultStorage(fCDBUri);
663 // Now activate the detector specific CDB storage locations
664 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
665 TObject* obj = fSpecCDBUri[i];
667 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
668 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
669 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
670 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
672 AliSysInfo::AddStamp("InitCDB");
675 //_____________________________________________________________________________
676 void AliReconstruction::SetDefaultStorage(const char* uri) {
677 // Store the desired default CDB storage location
678 // Activate it later within the Run() method
684 //_____________________________________________________________________________
685 void AliReconstruction::SetQARefDefaultStorage(const char* uri) {
686 // Store the desired default CDB storage location
687 // Activate it later within the Run() method
690 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
693 //_____________________________________________________________________________
694 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
695 // Store a detector-specific CDB storage location
696 // Activate it later within the Run() method
698 AliCDBPath aPath(calibType);
699 if(!aPath.IsValid()){
700 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
701 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
702 if(!strcmp(calibType, fgkDetectorName[iDet])) {
703 aPath.SetPath(Form("%s/*", calibType));
704 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
708 if(!aPath.IsValid()){
709 AliError(Form("Not a valid path or detector: %s", calibType));
714 // // check that calibType refers to a "valid" detector name
715 // Bool_t isDetector = kFALSE;
716 // for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
717 // TString detName = fgkDetectorName[iDet];
718 // if(aPath.GetLevel0() == detName) {
719 // isDetector = kTRUE;
725 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
729 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
730 if (obj) fSpecCDBUri.Remove(obj);
731 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
735 //_____________________________________________________________________________
736 Bool_t AliReconstruction::SetRunNumberFromData()
738 // The method is called in Run() in order
739 // to set a correct run number.
740 // In case of raw data reconstruction the
741 // run number is taken from the raw data header
743 if (fSetRunNumberFromDataCalled) return kTRUE;
744 fSetRunNumberFromDataCalled = kTRUE;
746 AliCDBManager* man = AliCDBManager::Instance();
749 if(fRawReader->NextEvent()) {
750 if(man->GetRun() > 0) {
751 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
753 man->SetRun(fRawReader->GetRunNumber());
754 fRawReader->RewindEvents();
757 if(man->GetRun() > 0) {
758 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
761 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
767 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
769 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
774 // read run number from gAlice
775 if(rl->GetHeader()) {
776 man->SetRun(rl->GetHeader()->GetRun());
781 AliError("Neither run-loader header nor RawReader objects are found !");
793 //_____________________________________________________________________________
794 void AliReconstruction::SetCDBLock() {
795 // Set CDB lock: from now on it is forbidden to reset the run number
796 // or the default storage or to activate any further storage!
798 AliCDBManager::Instance()->SetLock(1);
801 //_____________________________________________________________________________
802 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
804 // Read the alignment objects from CDB.
805 // Each detector is supposed to have the
806 // alignment objects in DET/Align/Data CDB path.
807 // All the detector objects are then collected,
808 // sorted by geometry level (starting from ALIC) and
809 // then applied to the TGeo geometry.
810 // Finally an overlaps check is performed.
812 // Load alignment data from CDB and fill fAlignObjArray
813 if(fLoadAlignFromCDB){
815 TString detStr = detectors;
816 TString loadAlObjsListOfDets = "";
818 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
819 if(!IsSelected(fgkDetectorName[iDet], detStr)) continue;
820 if(!strcmp(fgkDetectorName[iDet],"HLT")) continue;
822 if(AliGeomManager::GetNalignable(fgkDetectorName[iDet]) != 0)
824 loadAlObjsListOfDets += fgkDetectorName[iDet];
825 loadAlObjsListOfDets += " ";
827 } // end loop over detectors
829 if(AliGeomManager::GetNalignable("GRP") != 0)
830 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
831 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
832 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
834 // Check if the array with alignment objects was
835 // provided by the user. If yes, apply the objects
836 // to the present TGeo geometry
837 if (fAlignObjArray) {
838 if (gGeoManager && gGeoManager->IsClosed()) {
839 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
840 AliError("The misalignment of one or more volumes failed!"
841 "Compare the list of simulated detectors and the list of detector alignment data!");
846 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
852 if (fAlignObjArray) {
853 fAlignObjArray->Delete();
854 delete fAlignObjArray; fAlignObjArray=NULL;
860 //_____________________________________________________________________________
861 void AliReconstruction::SetGAliceFile(const char* fileName)
863 // set the name of the galice file
865 fGAliceFileName = fileName;
868 //_____________________________________________________________________________
869 void AliReconstruction::SetInput(const char* input)
871 // In case the input string starts with 'mem://', we run in an online mode
872 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
873 // file is assumed. One can give as an input:
874 // mem://: - events taken from DAQ monitoring libs online
876 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
877 if (input) fRawInput = input;
880 //_____________________________________________________________________________
881 void AliReconstruction::SetOutput(const char* output)
883 // Set the output ESD filename
884 // 'output' is a normalt ROOT url
885 // The method is used in case of raw-data reco with PROOF
886 if (output) fESDOutput = output;
889 //_____________________________________________________________________________
890 void AliReconstruction::SetOption(const char* detector, const char* option)
892 // set options for the reconstruction of a detector
894 TObject* obj = fOptions.FindObject(detector);
895 if (obj) fOptions.Remove(obj);
896 fOptions.Add(new TNamed(detector, option));
899 //_____________________________________________________________________________
900 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
902 // Set custom reconstruction parameters for a given detector
903 // Single set of parameters for all the events
905 // First check if the reco-params are global
906 if(!strcmp(detector, "GRP")) {
908 fRecoParam.AddDetRecoParam(kNDetectors,par);
912 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
913 if(!strcmp(detector, fgkDetectorName[iDet])) {
915 fRecoParam.AddDetRecoParam(iDet,par);
922 //_____________________________________________________________________________
923 Bool_t AliReconstruction::InitGRP() {
924 //------------------------------------
925 // Initialization of the GRP entry
926 //------------------------------------
927 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
931 TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
934 AliInfo("Found a TMap in GRP/GRP/Data, converting it into an AliGRPObject");
936 fGRPData = new AliGRPObject();
937 fGRPData->ReadValuesFromMap(m);
941 AliInfo("Found an AliGRPObject in GRP/GRP/Data, reading it");
942 fGRPData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
946 // FIX ME: The unloading of GRP entry is temporarily disabled
947 // because ZDC and VZERO are using it in order to initialize
948 // their reconstructor objects. In the future one has to think
949 // of propagating AliRunInfo to the reconstructors.
950 // AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
954 AliError("No GRP entry found in OCDB!");
958 TString lhcState = fGRPData->GetLHCState();
959 if (lhcState==AliGRPObject::GetInvalidString()) {
960 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
961 lhcState = "UNKNOWN";
964 TString beamType = fGRPData->GetBeamType();
965 if (beamType==AliGRPObject::GetInvalidString()) {
966 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
967 beamType = "UNKNOWN";
970 Float_t beamEnergy = fGRPData->GetBeamEnergy();
971 if (beamEnergy==AliGRPObject::GetInvalidFloat()) {
972 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
975 // LHC: "multiply by 120 to get the energy in MeV"
978 TString runType = fGRPData->GetRunType();
979 if (runType==AliGRPObject::GetInvalidString()) {
980 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
984 Int_t activeDetectors = fGRPData->GetDetectorMask();
985 if (activeDetectors==AliGRPObject::GetInvalidUInt()) {
986 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
987 activeDetectors = 1074790399;
990 fRunInfo = new AliRunInfo(lhcState, beamType, beamEnergy, runType, activeDetectors);
994 // Process the list of active detectors
995 if (activeDetectors) {
996 UInt_t detMask = activeDetectors;
997 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
998 fRunTracking = MatchDetectorList(fRunTracking,detMask);
999 fFillESD = MatchDetectorList(fFillESD,detMask);
1000 fQADetectors = MatchDetectorList(fQADetectors,detMask);
1001 fLoadCDB.Form("%s %s %s %s",
1002 fRunLocalReconstruction.Data(),
1003 fRunTracking.Data(),
1005 fQADetectors.Data());
1006 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
1007 if (!((detMask >> AliDAQ::DetectorID("ITSSPD")) & 0x1) &&
1008 !((detMask >> AliDAQ::DetectorID("ITSSDD")) & 0x1) &&
1009 !((detMask >> AliDAQ::DetectorID("ITSSSD")) & 0x1) ) {
1010 // switch off the vertexer
1011 AliInfo("SPD,SDD,SSD is not in the list of active detectors. Vertexer switched off.");
1012 fRunVertexFinder = kFALSE;
1014 if (!((detMask >> AliDAQ::DetectorID("TRG")) & 0x1)) {
1015 // switch off the reading of CTP raw-data payload
1016 if (fFillTriggerESD) {
1017 AliInfo("CTP is not in the list of active detectors. CTP data reading switched off.");
1018 fFillTriggerESD = kFALSE;
1023 AliInfo("===================================================================================");
1024 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
1025 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
1026 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
1027 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
1028 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
1029 AliInfo("===================================================================================");
1031 //*** Dealing with the magnetic field map
1032 if ( TGeoGlobalMagField::Instance()->IsLocked() ) {
1033 if (TGeoGlobalMagField::Instance()->GetField()->TestBit(AliMagF::kOverrideGRP)) {
1034 AliInfo("ExpertMode!!! GRP information will be ignored !");
1035 AliInfo("ExpertMode!!! Running with the externally locked B field !");
1038 AliInfo("Destroying existing B field instance!");
1039 delete TGeoGlobalMagField::Instance();
1042 if ( !TGeoGlobalMagField::Instance()->IsLocked() ) {
1043 // Construct the field map out of the information retrieved from GRP.
1046 Float_t l3Current = fGRPData->GetL3Current((AliGRPObject::Stats)0);
1047 if (l3Current == AliGRPObject::GetInvalidFloat()) {
1048 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1052 Char_t l3Polarity = fGRPData->GetL3Polarity();
1053 if (l3Polarity == AliGRPObject::GetInvalidChar()) {
1054 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1059 Float_t diCurrent = fGRPData->GetDipoleCurrent((AliGRPObject::Stats)0);
1060 if (diCurrent == AliGRPObject::GetInvalidFloat()) {
1061 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1065 Char_t diPolarity = fGRPData->GetDipolePolarity();
1066 if (diPolarity == AliGRPObject::GetInvalidChar()) {
1067 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1071 // read special bits for the polarity convention and map type
1072 Int_t polConvention = fGRPData->IsPolarityConventionLHC() ? AliMagF::kConvLHC : AliMagF::kConvDCS2008;
1073 Bool_t uniformB = fGRPData->IsUniformBMap();
1076 AliMagF* fld = AliMagF::CreateFieldMap(TMath::Abs(l3Current) * (l3Polarity ? -1:1),
1077 TMath::Abs(diCurrent) * (diPolarity ? -1:1),
1078 polConvention,uniformB,beamEnergy, beamType.Data());
1080 TGeoGlobalMagField::Instance()->SetField( fld );
1081 TGeoGlobalMagField::Instance()->Lock();
1082 AliInfo("Running with the B field constructed out of GRP !");
1084 else AliFatal("Failed to create a B field map !");
1086 else AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
1089 //*** Get the diamond profiles from OCDB
1090 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1092 fDiamondProfileSPD = dynamic_cast<AliESDVertex*> (entry->GetObject());
1094 AliError("No SPD diamond profile found in OCDB!");
1097 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1099 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
1101 AliError("No diamond profile found in OCDB!");
1104 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
1106 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
1108 AliError("No TPC diamond profile found in OCDB!");
1111 entry = AliCDBManager::Instance()->Get("GRP/Calib/CosmicTriggers");
1113 fListOfCosmicTriggers = dynamic_cast<THashTable*>(entry->GetObject());
1115 AliCDBManager::Instance()->UnloadFromCache("GRP/Calib/CosmicTriggers");
1118 if (!fListOfCosmicTriggers) {
1119 AliWarning("Can not get list of cosmic triggers from OCDB! Cosmic event specie will be effectively disabled!");
1125 //_____________________________________________________________________________
1126 Bool_t AliReconstruction::LoadCDB()
1128 AliCodeTimerAuto("");
1130 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1132 TString detStr = fLoadCDB;
1133 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1134 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1135 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1139 //_____________________________________________________________________________
1140 Bool_t AliReconstruction::LoadTriggerScalersCDB()
1142 AliCodeTimerAuto("");
1144 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/CTP/Scalers");
1148 AliInfo("Found an AliTriggerRunScalers in GRP/CTP/Scalers, reading it");
1149 fRunScalers = dynamic_cast<AliTriggerRunScalers*> (entry->GetObject());
1151 if (fRunScalers->CorrectScalersOverflow() == 0) AliInfo("32bit Trigger counters corrected for overflow");
1156 //_____________________________________________________________________________
1157 Bool_t AliReconstruction::Run(const char* input)
1160 AliCodeTimerAuto("");
1163 if (GetAbort() != TSelector::kContinue) return kFALSE;
1165 TChain *chain = NULL;
1166 if (fRawReader && (chain = fRawReader->GetChain())) {
1167 Long64_t nEntries = (fLastEvent < 0) ? (TChain::kBigNumber) : (fLastEvent - fFirstEvent + 1);
1172 gProof->Exec("TGrid::Connect(\"alien://\")",kTRUE);
1174 TMessage::EnableSchemaEvolutionForAll(kTRUE);
1175 gProof->Exec("TMessage::EnableSchemaEvolutionForAll(kTRUE)",kTRUE);
1177 gProof->AddInput(this);
1179 if (!ParseOutput()) return kFALSE;
1181 gProof->SetParameter("PROOF_MaxSlavesPerNode", 9999);
1183 chain->Process("AliReconstruction","",nEntries,fFirstEvent);
1186 chain->Process(this,"",nEntries,fFirstEvent);
1191 if (GetAbort() != TSelector::kContinue) return kFALSE;
1193 if (GetAbort() != TSelector::kContinue) return kFALSE;
1194 //******* The loop over events
1195 AliInfo("Starting looping over events");
1197 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1198 (fRawReader && fRawReader->NextEvent())) {
1199 if (!ProcessEvent(iEvent)) {
1200 Abort("ProcessEvent",TSelector::kAbortFile);
1206 if (GetAbort() != TSelector::kContinue) return kFALSE;
1208 if (GetAbort() != TSelector::kContinue) return kFALSE;
1214 //_____________________________________________________________________________
1215 void AliReconstruction::InitRawReader(const char* input)
1217 AliCodeTimerAuto("");
1219 // Init raw-reader and
1220 // set the input in case of raw data
1221 if (input) fRawInput = input;
1222 fRawReader = AliRawReader::Create(fRawInput.Data());
1224 AliInfo("Reconstruction will run over digits");
1226 if (!fEquipIdMap.IsNull() && fRawReader)
1227 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1229 if (!fUseHLTData.IsNull()) {
1230 // create the RawReaderHLT which performs redirection of HLT input data for
1231 // the specified detectors
1232 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1234 fParentRawReader=fRawReader;
1235 fRawReader=pRawReader;
1237 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1240 AliSysInfo::AddStamp("CreateRawReader");
1243 //_____________________________________________________________________________
1244 void AliReconstruction::InitRun(const char* input)
1246 // Initialization of raw-reader,
1247 // run number, CDB etc.
1248 AliCodeTimerAuto("");
1249 AliSysInfo::AddStamp("Start");
1251 // Initialize raw-reader if any
1252 InitRawReader(input);
1254 // Initialize the CDB storage
1257 // Set run number in CDBManager (if it is not already set by the user)
1258 if (!SetRunNumberFromData()) {
1259 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1263 // Set CDB lock: from now on it is forbidden to reset the run number
1264 // or the default storage or to activate any further storage!
1269 //_____________________________________________________________________________
1270 void AliReconstruction::Begin(TTree *)
1272 // Initialize AlReconstruction before
1273 // going into the event loop
1274 // Should follow the TSelector convention
1275 // i.e. initialize only the object on the client side
1276 AliCodeTimerAuto("");
1278 AliReconstruction *reco = NULL;
1280 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1283 AliSysInfo::AddStamp("ReadInputInBegin");
1286 // Import ideal TGeo geometry and apply misalignment
1288 TString geom(gSystem->DirName(fGAliceFileName));
1289 geom += "/geometry.root";
1290 AliGeomManager::LoadGeometry(geom.Data());
1292 Abort("LoadGeometry", TSelector::kAbortProcess);
1295 AliSysInfo::AddStamp("LoadGeom");
1296 TString detsToCheck=fRunLocalReconstruction;
1297 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1298 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1301 AliSysInfo::AddStamp("CheckGeom");
1304 if (!MisalignGeometry(fLoadAlignData)) {
1305 Abort("MisalignGeometry", TSelector::kAbortProcess);
1308 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1309 AliSysInfo::AddStamp("MisalignGeom");
1312 Abort("InitGRP", TSelector::kAbortProcess);
1315 AliSysInfo::AddStamp("InitGRP");
1318 Abort("LoadCDB", TSelector::kAbortProcess);
1321 AliSysInfo::AddStamp("LoadCDB");
1323 if (!LoadTriggerScalersCDB()) {
1324 Abort("LoadTriggerScalersCDB", TSelector::kAbortProcess);
1327 AliSysInfo::AddStamp("LoadTriggerScalersCDB");
1330 // Read the reconstruction parameters from OCDB
1331 if (!InitRecoParams()) {
1332 AliWarning("Not all detectors have correct RecoParam objects initialized");
1334 AliSysInfo::AddStamp("InitRecoParams");
1336 if (fInput && gProof) {
1337 if (reco) *reco = *this;
1339 gGeoManager->SetName("Geometry");
1340 gProof->AddInputData(gGeoManager,kTRUE);
1342 gProof->AddInputData(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()),kTRUE);
1343 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1344 AliMagF *magFieldMap = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1345 magFieldMap->SetName("MagneticFieldMap");
1346 gProof->AddInputData(magFieldMap,kTRUE);
1351 //_____________________________________________________________________________
1352 void AliReconstruction::SlaveBegin(TTree*)
1354 // Initialization related to run-loader,
1355 // vertexer, trackers, recontructors
1356 // In proof mode it is executed on the slave
1357 AliCodeTimerAuto("");
1359 TProofOutputFile *outProofFile = NULL;
1361 if (AliDebugLevel() > 0) fInput->Print();
1362 if (AliDebugLevel() > 10) fInput->Dump();
1363 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1366 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1368 AliGeomManager::SetGeometry(tgeo);
1370 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1371 Int_t runNumber = -1;
1372 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1373 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1374 man->SetCacheFlag(kTRUE);
1375 man->SetLock(kTRUE);
1379 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1380 AliMagF *newMap = new AliMagF(*map);
1381 if (!newMap->LoadParameterization()) {
1382 Abort("AliMagF::LoadParameterization", TSelector::kAbortProcess);
1385 TGeoGlobalMagField::Instance()->SetField(newMap);
1386 TGeoGlobalMagField::Instance()->Lock();
1388 if (TNamed *outputFileName = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE"))
1389 fProofOutputFileName = outputFileName->GetTitle();
1390 if (TNamed *outputLocation = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE_LOCATION"))
1391 fProofOutputLocation = outputLocation->GetTitle();
1392 if (fInput->FindObject("PROOF_OUTPUTFILE_DATASET"))
1393 fProofOutputDataset = kTRUE;
1394 if (TNamed *archiveList = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE_ARCHIVE"))
1395 fProofOutputArchive = archiveList->GetTitle();
1396 if (!fProofOutputFileName.IsNull() &&
1397 !fProofOutputLocation.IsNull() &&
1398 fProofOutputArchive.IsNull()) {
1399 if (!fProofOutputDataset) {
1400 outProofFile = new TProofOutputFile(fProofOutputFileName.Data(),"M");
1401 outProofFile->SetOutputFileName(Form("%s%s",fProofOutputLocation.Data(),fProofOutputFileName.Data()));
1404 outProofFile = new TProofOutputFile(fProofOutputFileName.Data(),"DROV",fProofOutputLocation.Data());
1406 if (AliDebugLevel() > 0) outProofFile->Dump();
1407 fOutput->Add(outProofFile);
1409 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1412 // get the run loader
1413 if (!InitRunLoader()) {
1414 Abort("InitRunLoader", TSelector::kAbortProcess);
1417 AliSysInfo::AddStamp("LoadLoader");
1419 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1422 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1423 Abort("CreateTrackers", TSelector::kAbortProcess);
1426 AliSysInfo::AddStamp("CreateTrackers");
1428 // create the ESD output file and tree
1429 if (!outProofFile) {
1430 ffile = TFile::Open("AliESDs.root", "RECREATE");
1431 ffile->SetCompressionLevel(2);
1432 if (!ffile->IsOpen()) {
1433 Abort("OpenESDFile", TSelector::kAbortProcess);
1438 AliInfo(Form("Opening output PROOF file: %s/%s",
1439 outProofFile->GetDir(), outProofFile->GetFileName()));
1440 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1441 Abort(Form("Problems opening output PROOF file: %s/%s",
1442 outProofFile->GetDir(), outProofFile->GetFileName()),
1443 TSelector::kAbortProcess);
1448 ftree = new TTree("esdTree", "Tree with ESD objects");
1449 fesd = new AliESDEvent();
1450 fesd->CreateStdContent();
1452 fesd->WriteToTree(ftree);
1453 if (fWriteESDfriend) {
1455 // Since we add the branch manually we must
1456 // book and add it after WriteToTree
1457 // otherwise it is created twice,
1458 // once via writetotree and once here.
1459 // The case for AliESDfriend is now
1460 // caught also in AlIESDEvent::WriteToTree but
1461 // be careful when changing the name (AliESDfriend is not
1462 // a TNamed so we had to hardwire it)
1463 fesdf = new AliESDfriend();
1464 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
1465 br->SetFile("AliESDfriends.root");
1466 fesd->AddObject(fesdf);
1468 ftree->GetUserInfo()->Add(fesd);
1470 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1471 fhltesd = new AliESDEvent();
1472 fhltesd->CreateStdContent();
1474 // read the ESD template from CDB
1475 // HLT is allowed to put non-std content to its ESD, the non-std
1476 // objects need to be created before invocation of WriteToTree in
1477 // order to create all branches. Initialization is done from an
1478 // ESD layout template in CDB
1479 AliCDBManager* man = AliCDBManager::Instance();
1480 AliCDBPath hltESDConfigPath("HLT/ConfigHLT/esdLayout");
1481 AliCDBEntry* hltESDConfig=NULL;
1482 if (man->GetId(hltESDConfigPath)!=NULL &&
1483 (hltESDConfig=man->Get(hltESDConfigPath))!=NULL) {
1484 AliESDEvent* pESDLayout=dynamic_cast<AliESDEvent*>(hltESDConfig->GetObject());
1486 // init all internal variables from the list of objects
1487 pESDLayout->GetStdContent();
1489 // copy content and create non-std objects
1490 *fhltesd=*pESDLayout;
1493 AliError(Form("error setting hltEsd layout from %s: invalid object type",
1494 hltESDConfigPath.GetPath().Data()));
1498 fhltesd->WriteToTree(fhlttree);
1499 fhlttree->GetUserInfo()->Add(fhltesd);
1501 ProcInfo_t procInfo;
1502 gSystem->GetProcInfo(&procInfo);
1503 AliInfo(Form("Current memory usage %d %d", procInfo.fMemResident, procInfo.fMemVirtual));
1506 //Initialize the QA and start of cycle
1507 if (fRunQA || fRunGlobalQA)
1510 //Initialize the Plane Efficiency framework
1511 if (fRunPlaneEff && !InitPlaneEff()) {
1512 Abort("InitPlaneEff", TSelector::kAbortProcess);
1516 if (strcmp(gProgName,"alieve") == 0)
1517 fRunAliEVE = InitAliEVE();
1522 //_____________________________________________________________________________
1523 Bool_t AliReconstruction::Process(Long64_t entry)
1525 // run the reconstruction over a single entry
1526 // from the chain with raw data
1527 AliCodeTimerAuto("");
1529 TTree *currTree = fChain->GetTree();
1530 AliRawVEvent *event = NULL;
1531 currTree->SetBranchAddress("rawevent",&event);
1532 currTree->GetEntry(entry);
1533 fRawReader = new AliRawReaderRoot(event);
1534 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1542 //_____________________________________________________________________________
1543 void AliReconstruction::Init(TTree *tree)
1546 AliError("The input tree is not found!");
1552 //_____________________________________________________________________________
1553 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1555 // run the reconstruction over a single event
1556 // The event loop is steered in Run method
1558 AliCodeTimerAuto("");
1560 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1561 fRunLoader->SetEventNumber(iEvent);
1562 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1564 fRunLoader->TreeE()->Fill();
1565 if (fRawReader && fRawReader->UseAutoSaveESD())
1566 fRunLoader->TreeE()->AutoSave("SaveSelf");
1569 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1573 AliInfo(Form("processing event %d", iEvent));
1575 fRunLoader->GetEvent(iEvent);
1577 // Fill Event-info object
1579 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo,fListOfCosmicTriggers);
1580 AliInfo(Form("Current event specie: %s",fRecoParam.PrintEventSpecie()));
1582 // Set the reco-params
1584 TString detStr = fLoadCDB;
1585 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1586 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1587 AliReconstructor *reconstructor = GetReconstructor(iDet);
1588 if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
1589 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
1590 reconstructor->SetRecoParam(par);
1591 reconstructor->SetEventInfo(&fEventInfo);
1593 AliQAManager::QAManager()->SetRecoParam(iDet, par) ;
1594 AliQAManager::QAManager()->SetEventSpecie(AliRecoParam::Convert(par->GetEventSpecie())) ;
1601 if (fRunQA && IsInTasks(AliQAv1::kRAWS)) {
1602 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1603 AliQAManager::QAManager()->RunOneEvent(fRawReader) ;
1605 // local single event reconstruction
1606 if (!fRunLocalReconstruction.IsNull()) {
1607 TString detectors=fRunLocalReconstruction;
1608 // run HLT event reconstruction first
1609 // ;-( IsSelected changes the string
1610 if (IsSelected("HLT", detectors) &&
1611 !RunLocalEventReconstruction("HLT")) {
1612 if (fStopOnError) {CleanUp(); return kFALSE;}
1614 detectors=fRunLocalReconstruction;
1615 detectors.ReplaceAll("HLT", "");
1616 if (!RunLocalEventReconstruction(detectors)) {
1617 if (fStopOnError) {CleanUp(); return kFALSE;}
1621 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1622 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1623 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1624 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1626 // Set magnetic field from the tracker
1627 fesd->SetMagneticField(AliTracker::GetBz());
1628 fhltesd->SetMagneticField(AliTracker::GetBz());
1630 AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1631 if (fld) { // set info needed for field initialization
1632 fesd->SetCurrentL3(fld->GetCurrentSol());
1633 fesd->SetCurrentDip(fld->GetCurrentDip());
1634 fesd->SetBeamEnergy(fld->GetBeamEnergy());
1635 fesd->SetBeamType(fld->GetBeamTypeText());
1636 fesd->SetUniformBMap(fld->IsUniform());
1637 fesd->SetBInfoStored();
1639 fhltesd->SetCurrentL3(fld->GetCurrentSol());
1640 fhltesd->SetCurrentDip(fld->GetCurrentDip());
1641 fhltesd->SetBeamEnergy(fld->GetBeamEnergy());
1642 fhltesd->SetBeamType(fld->GetBeamTypeText());
1643 fhltesd->SetUniformBMap(fld->IsUniform());
1644 fhltesd->SetBInfoStored();
1647 // Set most probable pt, for B=0 tracking
1648 // Get the global reco-params. They are atposition 16 inside the array of detectors in fRecoParam
1649 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(kNDetectors));
1650 if (grpRecoParam) AliExternalTrackParam::SetMostProbablePt(grpRecoParam->GetMostProbablePt());
1652 // Fill raw-data error log into the ESD
1653 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1656 if (fRunVertexFinder) {
1657 if (!RunVertexFinder(fesd)) {
1658 if (fStopOnError) {CleanUp(); return kFALSE;}
1662 // For Plane Efficiency: run the SPD trackleter
1663 if (fRunPlaneEff && fSPDTrackleter) {
1664 if (!RunSPDTrackleting(fesd)) {
1665 if (fStopOnError) {CleanUp(); return kFALSE;}
1670 if (!fRunTracking.IsNull()) {
1671 if (fRunMuonTracking) {
1672 if (!RunMuonTracking(fesd)) {
1673 if (fStopOnError) {CleanUp(); return kFALSE;}
1679 if (!fRunTracking.IsNull()) {
1680 if (!RunTracking(fesd)) {
1681 if (fStopOnError) {CleanUp(); return kFALSE;}
1686 if (!fFillESD.IsNull()) {
1687 TString detectors=fFillESD;
1688 // run HLT first and on hltesd
1689 // ;-( IsSelected changes the string
1690 if (IsSelected("HLT", detectors) &&
1691 !FillESD(fhltesd, "HLT")) {
1692 if (fStopOnError) {CleanUp(); return kFALSE;}
1695 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1696 if (detectors.Contains("ALL")) {
1698 for (Int_t idet=0; idet<kNDetectors; ++idet){
1699 detectors += fgkDetectorName[idet];
1703 detectors.ReplaceAll("HLT", "");
1704 if (!FillESD(fesd, detectors)) {
1705 if (fStopOnError) {CleanUp(); return kFALSE;}
1709 // fill Event header information from the RawEventHeader
1710 if (fRawReader){FillRawEventHeaderESD(fesd);}
1711 if (fRawReader){FillRawEventHeaderESD(fhltesd);}
1714 AliESDpid::MakePID(fesd);
1716 if (fFillTriggerESD) {
1717 if (!FillTriggerESD(fesd)) {
1718 if (fStopOnError) {CleanUp(); return kFALSE;}
1721 // Always fill scalers
1722 if (!FillTriggerScalers(fesd)) {
1723 if (fStopOnError) {CleanUp(); return kFALSE;}
1730 // Propagate track to the beam pipe (if not already done by ITS)
1732 const Int_t ntracks = fesd->GetNumberOfTracks();
1733 const Double_t kBz = fesd->GetMagneticField();
1734 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1737 UShort_t *selectedIdx=new UShort_t[ntracks];
1739 for (Int_t itrack=0; itrack<ntracks; itrack++){
1740 const Double_t kMaxStep = 1; //max step over the material
1743 AliESDtrack *track = fesd->GetTrack(itrack);
1744 if (!track) continue;
1746 AliExternalTrackParam *tpcTrack =
1747 (AliExternalTrackParam *)track->GetTPCInnerParam();
1751 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kFALSE);
1754 Int_t n=trkArray.GetEntriesFast();
1755 selectedIdx[n]=track->GetID();
1756 trkArray.AddLast(tpcTrack);
1759 //Tracks refitted by ITS should already be at the SPD vertex
1760 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1763 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kFALSE);
1764 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kVeryBig);
1769 // Improve the reconstructed primary vertex position using the tracks
1771 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
1772 if(fesd->GetPrimaryVertexSPD()) {
1773 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
1774 if(vtitle.Contains("cosmics")) {
1775 runVertexFinderTracks=kFALSE;
1779 if (runVertexFinderTracks) {
1780 // TPC + ITS primary vertex
1781 ftVertexer->SetITSMode();
1782 ftVertexer->SetConstraintOff();
1783 // get cuts for vertexer from AliGRPRecoParam
1785 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1786 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1787 grpRecoParam->GetVertexerTracksCutsITS(cutsVertexer);
1788 ftVertexer->SetCuts(cutsVertexer);
1789 delete [] cutsVertexer; cutsVertexer = NULL;
1790 if(fDiamondProfile && grpRecoParam->GetVertexerTracksConstraintITS())
1791 ftVertexer->SetVtxStart(fDiamondProfile);
1793 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1795 if (pvtx->GetStatus()) {
1796 fesd->SetPrimaryVertexTracks(pvtx);
1797 for (Int_t i=0; i<ntracks; i++) {
1798 AliESDtrack *t = fesd->GetTrack(i);
1799 t->RelateToVertex(pvtx, kBz, kVeryBig);
1804 // TPC-only primary vertex
1805 ftVertexer->SetTPCMode();
1806 ftVertexer->SetConstraintOff();
1807 // get cuts for vertexer from AliGRPRecoParam
1809 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
1810 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
1811 grpRecoParam->GetVertexerTracksCutsTPC(cutsVertexer);
1812 ftVertexer->SetCuts(cutsVertexer);
1813 delete [] cutsVertexer; cutsVertexer = NULL;
1814 if(fDiamondProfileTPC && grpRecoParam->GetVertexerTracksConstraintTPC())
1815 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1817 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1819 if (pvtx->GetStatus()) {
1820 fesd->SetPrimaryVertexTPC(pvtx);
1821 for (Int_t i=0; i<ntracks; i++) {
1822 AliESDtrack *t = fesd->GetTrack(i);
1823 t->RelateToVertexTPC(pvtx, kBz, kVeryBig);
1829 delete[] selectedIdx;
1831 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1836 AliV0vertexer vtxer;
1837 vtxer.Tracks2V0vertices(fesd);
1839 if (fRunCascadeFinder) {
1841 AliCascadeVertexer cvtxer;
1842 cvtxer.V0sTracks2CascadeVertices(fesd);
1847 if (fCleanESD) CleanESD(fesd);
1849 if (fRunQA && IsInTasks(AliQAv1::kESDS)) {
1850 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1851 AliQAManager::QAManager()->RunOneEvent(fesd) ;
1854 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
1855 qadm->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1856 if (qadm && IsInTasks(AliQAv1::kESDS))
1857 qadm->Exec(AliQAv1::kESDS, fesd);
1860 if (fWriteESDfriend) {
1861 // fesdf->~AliESDfriend();
1862 // new (fesdf) AliESDfriend(); // Reset...
1863 fesd->GetESDfriend(fesdf);
1867 // Auto-save the ESD tree in case of prompt reco @P2
1868 if (fRawReader && fRawReader->UseAutoSaveESD()) {
1869 ftree->AutoSave("SaveSelf");
1870 TFile *friendfile = (TFile *)(gROOT->GetListOfFiles()->FindObject("AliESDfriends.root"));
1871 if (friendfile) friendfile->Save();
1878 if (fRunAliEVE) RunAliEVE();
1882 if (fWriteESDfriend) {
1883 fesdf->~AliESDfriend();
1884 new (fesdf) AliESDfriend(); // Reset...
1887 ProcInfo_t procInfo;
1888 gSystem->GetProcInfo(&procInfo);
1889 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, procInfo.fMemResident, procInfo.fMemVirtual));
1892 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1893 if (fReconstructor[iDet]) {
1894 fReconstructor[iDet]->SetRecoParam(NULL);
1895 fReconstructor[iDet]->SetEventInfo(NULL);
1897 if (fTracker[iDet]) fTracker[iDet]->SetEventInfo(NULL);
1900 if (fRunQA || fRunGlobalQA)
1901 AliQAManager::QAManager()->Increment() ;
1906 //_____________________________________________________________________________
1907 void AliReconstruction::SlaveTerminate()
1909 // Finalize the run on the slave side
1910 // Called after the exit
1911 // from the event loop
1912 AliCodeTimerAuto("");
1914 if (fIsNewRunLoader) { // galice.root didn't exist
1915 fRunLoader->WriteHeader("OVERWRITE");
1916 fRunLoader->CdGAFile();
1917 fRunLoader->Write(0, TObject::kOverwrite);
1920 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1921 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1923 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1924 cdbMapCopy->SetOwner(1);
1925 cdbMapCopy->SetName("cdbMap");
1926 TIter iter(cdbMap->GetTable());
1929 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1930 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1931 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1932 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1935 TList *cdbListCopy = new TList();
1936 cdbListCopy->SetOwner(1);
1937 cdbListCopy->SetName("cdbList");
1939 TIter iter2(cdbList);
1942 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1943 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1946 ftree->GetUserInfo()->Add(cdbMapCopy);
1947 ftree->GetUserInfo()->Add(cdbListCopy);
1952 if (fWriteESDfriend)
1953 ftree->SetBranchStatus("ESDfriend*",0);
1954 // we want to have only one tree version number
1955 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1956 fhlttree->Write(fhlttree->GetName(),TObject::kOverwrite);
1958 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1959 if (fRunPlaneEff && !FinishPlaneEff()) {
1960 AliWarning("Finish PlaneEff evaluation failed");
1963 // End of cycle for the in-loop
1965 AliQAManager::QAManager()->EndOfCycle() ;
1968 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
1970 if (IsInTasks(AliQAv1::kRECPOINTS))
1971 qadm->EndOfCycle(AliQAv1::kRECPOINTS);
1972 if (IsInTasks(AliQAv1::kESDS))
1973 qadm->EndOfCycle(AliQAv1::kESDS);
1978 if (fRunQA || fRunGlobalQA) {
1980 !fProofOutputLocation.IsNull() &&
1981 fProofOutputArchive.IsNull() &&
1982 !fProofOutputDataset) {
1983 TString qaOutputFile(Form("%sMerged.%s.Data.root",
1984 fProofOutputLocation.Data(),
1985 AliQAv1::GetQADataFileName()));
1986 TProofOutputFile *qaProofFile = new TProofOutputFile(Form("Merged.%s.Data.root",
1987 AliQAv1::GetQADataFileName()));
1988 qaProofFile->SetOutputFileName(qaOutputFile.Data());
1989 if (AliDebugLevel() > 0) qaProofFile->Dump();
1990 fOutput->Add(qaProofFile);
1991 MergeQA(qaProofFile->GetFileName());
2002 if (!fProofOutputFileName.IsNull() &&
2003 !fProofOutputLocation.IsNull() &&
2004 fProofOutputDataset &&
2005 !fProofOutputArchive.IsNull()) {
2006 TProofOutputFile *zipProofFile = new TProofOutputFile(fProofOutputFileName.Data(),
2008 fProofOutputLocation.Data());
2009 if (AliDebugLevel() > 0) zipProofFile->Dump();
2010 fOutput->Add(zipProofFile);
2011 TString fileList(fProofOutputArchive.Data());
2012 fileList.ReplaceAll(","," ");
2014 #if ROOT_SVN_REVISION >= 30174
2015 command.Form("zip -n root %s/%s %s",zipProofFile->GetDir(kTRUE),zipProofFile->GetFileName(),fileList.Data());
2017 command.Form("zip -n root %s/%s %s",zipProofFile->GetDir(),zipProofFile->GetFileName(),fileList.Data());
2019 AliInfo(Form("Executing: %s",command.Data()));
2020 gSystem->Exec(command.Data());
2025 //_____________________________________________________________________________
2026 void AliReconstruction::Terminate()
2028 // Create tags for the events in the ESD tree (the ESD tree is always present)
2029 // In case of empty events the tags will contain dummy values
2030 AliCodeTimerAuto("");
2032 // Do not call the ESD tag creator in case of PROOF-based reconstruction
2034 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
2035 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData, AliQAv1::Instance()->GetQA(), AliQAv1::Instance()->GetEventSpecies(), AliQAv1::kNDET, AliRecoParam::kNSpecies);
2038 // Cleanup of CDB manager: cache and active storages!
2039 AliCDBManager::Instance()->ClearCache();
2042 //_____________________________________________________________________________
2043 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
2045 // run the local reconstruction
2047 static Int_t eventNr=0;
2048 AliCodeTimerAuto("")
2050 TString detStr = detectors;
2051 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2052 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2053 AliReconstructor* reconstructor = GetReconstructor(iDet);
2054 if (!reconstructor) continue;
2055 AliLoader* loader = fLoader[iDet];
2056 // Matthias April 2008: temporary fix to run HLT reconstruction
2057 // although the HLT loader is missing
2058 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
2060 reconstructor->Reconstruct(fRawReader, NULL);
2063 reconstructor->Reconstruct(dummy, NULL);
2068 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
2071 // conversion of digits
2072 if (fRawReader && reconstructor->HasDigitConversion()) {
2073 AliInfo(Form("converting raw data digits into root objects for %s",
2074 fgkDetectorName[iDet]));
2075 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
2076 // fgkDetectorName[iDet]));
2077 loader->LoadDigits("update");
2078 loader->CleanDigits();
2079 loader->MakeDigitsContainer();
2080 TTree* digitsTree = loader->TreeD();
2081 reconstructor->ConvertDigits(fRawReader, digitsTree);
2082 loader->WriteDigits("OVERWRITE");
2083 loader->UnloadDigits();
2085 // local reconstruction
2086 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2087 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2088 loader->LoadRecPoints("update");
2089 loader->CleanRecPoints();
2090 loader->MakeRecPointsContainer();
2091 TTree* clustersTree = loader->TreeR();
2092 if (fRawReader && !reconstructor->HasDigitConversion()) {
2093 reconstructor->Reconstruct(fRawReader, clustersTree);
2095 loader->LoadDigits("read");
2096 TTree* digitsTree = loader->TreeD();
2098 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
2099 if (fStopOnError) return kFALSE;
2101 reconstructor->Reconstruct(digitsTree, clustersTree);
2102 if (fRunQA && IsInTasks(AliQAv1::kDIGITSR)) {
2103 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2104 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, digitsTree) ;
2107 loader->UnloadDigits();
2109 if (fRunQA && IsInTasks(AliQAv1::kRECPOINTS)) {
2110 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2111 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, clustersTree) ;
2113 loader->WriteRecPoints("OVERWRITE");
2114 loader->UnloadRecPoints();
2115 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
2117 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2118 AliError(Form("the following detectors were not found: %s",
2120 if (fStopOnError) return kFALSE;
2125 //_____________________________________________________________________________
2126 Bool_t AliReconstruction::RunSPDTrackleting(AliESDEvent*& esd)
2128 // run the SPD trackleting (for SPD efficiency purpouses)
2130 AliCodeTimerAuto("")
2132 Double_t vtxPos[3] = {0, 0, 0};
2133 Double_t vtxErr[3] = {0.0, 0.0, 0.0};
2135 TArrayF mcVertex(3);
2137 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
2138 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
2139 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
2142 const AliESDVertex *vertex = esd->GetVertex();
2144 AliWarning("Vertex not found");
2147 vertex->GetXYZ(vtxPos);
2148 vertex->GetSigmaXYZ(vtxErr);
2149 if (fSPDTrackleter) {
2150 AliInfo("running the SPD Trackleter for Plane Efficiency Evaluation");
2153 fLoader[0]->LoadRecPoints("read");
2154 TTree* tree = fLoader[0]->TreeR();
2156 AliError("Can't get the ITS cluster tree");
2159 fSPDTrackleter->LoadClusters(tree);
2160 fSPDTrackleter->SetVertex(vtxPos, vtxErr);
2162 if (fSPDTrackleter->Clusters2Tracks(esd) != 0) {
2163 AliError("AliITSTrackleterSPDEff Clusters2Tracks failed");
2164 // fLoader[0]->UnloadRecPoints();
2167 //fSPDTrackleter->UnloadRecPoints();
2169 AliWarning("SPDTrackleter not available");
2175 //_____________________________________________________________________________
2176 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
2178 // run the barrel tracking
2180 AliCodeTimerAuto("")
2182 AliVertexer *vertexer = CreateVertexer();
2183 if (!vertexer) return kFALSE;
2185 AliInfo("running the ITS vertex finder");
2186 AliESDVertex* vertex = NULL;
2188 fLoader[0]->LoadRecPoints();
2189 TTree* cltree = fLoader[0]->TreeR();
2191 if(fDiamondProfileSPD) vertexer->SetVtxStart(fDiamondProfileSPD);
2192 vertex = vertexer->FindVertexForCurrentEvent(cltree);
2195 AliError("Can't get the ITS cluster tree");
2197 fLoader[0]->UnloadRecPoints();
2200 AliError("Can't get the ITS loader");
2203 AliWarning("Vertex not found");
2204 vertex = new AliESDVertex();
2205 vertex->SetName("default");
2208 vertex->SetName("reconstructed");
2213 vertex->GetXYZ(vtxPos);
2214 vertex->GetSigmaXYZ(vtxErr);
2216 esd->SetPrimaryVertexSPD(vertex);
2217 AliESDVertex *vpileup = NULL;
2218 Int_t novertices = 0;
2219 vpileup = vertexer->GetAllVertices(novertices);
2221 for (Int_t kk=1; kk<novertices; kk++)esd->AddPileupVertexSPD(&vpileup[kk]);
2223 // if SPD multiplicity has been determined, it is stored in the ESD
2224 AliMultiplicity *mult = vertexer->GetMultiplicity();
2225 if(mult)esd->SetMultiplicity(mult);
2227 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2228 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
2237 //_____________________________________________________________________________
2238 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
2240 // run the HLT barrel tracking
2242 AliCodeTimerAuto("")
2245 AliError("Missing runLoader!");
2249 AliInfo("running HLT tracking");
2251 // Get a pointer to the HLT reconstructor
2252 AliReconstructor *reconstructor = GetReconstructor(kNDetectors-1);
2253 if (!reconstructor) return kFALSE;
2256 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2257 TString detName = fgkDetectorName[iDet];
2258 AliDebug(1, Form("%s HLT tracking", detName.Data()));
2259 reconstructor->SetOption(detName.Data());
2260 AliTracker *tracker = reconstructor->CreateTracker();
2262 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
2263 if (fStopOnError) return kFALSE;
2267 Double_t vtxErr[3]={0.005,0.005,0.010};
2268 const AliESDVertex *vertex = esd->GetVertex();
2269 vertex->GetXYZ(vtxPos);
2270 tracker->SetVertex(vtxPos,vtxErr);
2272 fLoader[iDet]->LoadRecPoints("read");
2273 TTree* tree = fLoader[iDet]->TreeR();
2275 AliError(Form("Can't get the %s cluster tree", detName.Data()));
2278 tracker->LoadClusters(tree);
2280 if (tracker->Clusters2Tracks(esd) != 0) {
2281 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
2285 tracker->UnloadClusters();
2293 //_____________________________________________________________________________
2294 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
2296 // run the muon spectrometer tracking
2298 AliCodeTimerAuto("")
2301 AliError("Missing runLoader!");
2304 Int_t iDet = 7; // for MUON
2306 AliInfo("is running...");
2308 // Get a pointer to the MUON reconstructor
2309 AliReconstructor *reconstructor = GetReconstructor(iDet);
2310 if (!reconstructor) return kFALSE;
2313 TString detName = fgkDetectorName[iDet];
2314 AliDebug(1, Form("%s tracking", detName.Data()));
2315 AliTracker *tracker = reconstructor->CreateTracker();
2317 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2322 fLoader[iDet]->LoadRecPoints("read");
2324 tracker->LoadClusters(fLoader[iDet]->TreeR());
2326 Int_t rv = tracker->Clusters2Tracks(esd);
2330 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2334 fLoader[iDet]->UnloadRecPoints();
2336 tracker->UnloadClusters();
2344 //_____________________________________________________________________________
2345 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
2347 // run the barrel tracking
2348 static Int_t eventNr=0;
2349 AliCodeTimerAuto("")
2351 AliInfo("running tracking");
2353 // Set the event info which is used
2354 // by the trackers in order to obtain
2355 // information about read-out detectors,
2357 AliDebug(1, "Setting event info");
2358 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2359 if (!fTracker[iDet]) continue;
2360 fTracker[iDet]->SetEventInfo(&fEventInfo);
2363 //Fill the ESD with the T0 info (will be used by the TOF)
2364 if (fReconstructor[11] && fLoader[11]) {
2365 fLoader[11]->LoadRecPoints("READ");
2366 TTree *treeR = fLoader[11]->TreeR();
2368 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2372 // pass 1: TPC + ITS inwards
2373 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2374 if (!fTracker[iDet]) continue;
2375 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2378 fLoader[iDet]->LoadRecPoints("read");
2379 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2380 TTree* tree = fLoader[iDet]->TreeR();
2382 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2385 fTracker[iDet]->LoadClusters(tree);
2386 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2388 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2389 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2392 // preliminary PID in TPC needed by the ITS tracker
2394 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2395 AliESDpid::MakePID(esd);
2397 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2400 // pass 2: ALL backwards
2402 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2403 if (!fTracker[iDet]) continue;
2404 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2407 if (iDet > 1) { // all except ITS, TPC
2409 fLoader[iDet]->LoadRecPoints("read");
2410 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2411 tree = fLoader[iDet]->TreeR();
2413 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2416 fTracker[iDet]->LoadClusters(tree);
2417 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2421 if (iDet>1) // start filling residuals for the "outer" detectors
2423 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2424 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2426 AliRecoParam::EventSpecie_t es=fRecoParam.GetEventSpecie();
2427 TObjArray * elem = arr[AliRecoParam::AConvert(es)];
2428 if ( elem && (! elem->At(0)) ) {
2429 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2430 if (qadm) qadm->InitRecPointsForTracker() ;
2434 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2435 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2440 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2441 fTracker[iDet]->UnloadClusters();
2442 fLoader[iDet]->UnloadRecPoints();
2444 // updated PID in TPC needed by the ITS tracker -MI
2446 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2447 AliESDpid::MakePID(esd);
2449 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2451 //stop filling residuals for the "outer" detectors
2452 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2454 // pass 3: TRD + TPC + ITS refit inwards
2456 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2457 if (!fTracker[iDet]) continue;
2458 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2461 if (iDet<2) // start filling residuals for TPC and ITS
2463 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2464 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2466 AliRecoParam::EventSpecie_t es=fRecoParam.GetEventSpecie();
2467 TObjArray * elem = arr[AliRecoParam::AConvert(es)];
2468 if ( elem && (! elem->At(0)) ) {
2469 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2470 if (qadm) qadm->InitRecPointsForTracker() ;
2475 if (fTracker[iDet]->RefitInward(esd) != 0) {
2476 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2479 // run postprocessing
2480 if (fTracker[iDet]->PostProcess(esd) != 0) {
2481 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2484 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2487 // write space-points to the ESD in case alignment data output
2489 if (fWriteAlignmentData)
2490 WriteAlignmentData(esd);
2492 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2493 if (!fTracker[iDet]) continue;
2495 fTracker[iDet]->UnloadClusters();
2496 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2497 fLoader[iDet]->UnloadRecPoints();
2498 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2500 // stop filling residuals for TPC and ITS
2501 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2507 //_____________________________________________________________________________
2508 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2510 // Remove the data which are not needed for the physics analysis.
2513 Int_t nTracks=esd->GetNumberOfTracks();
2514 Int_t nV0s=esd->GetNumberOfV0s();
2516 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2518 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2519 Bool_t rc=esd->Clean(cleanPars);
2521 nTracks=esd->GetNumberOfTracks();
2522 nV0s=esd->GetNumberOfV0s();
2524 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2529 //_____________________________________________________________________________
2530 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2532 // fill the event summary data
2534 AliCodeTimerAuto("")
2535 static Int_t eventNr=0;
2536 TString detStr = detectors;
2538 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2539 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2540 AliReconstructor* reconstructor = GetReconstructor(iDet);
2541 if (!reconstructor) continue;
2542 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2543 TTree* clustersTree = NULL;
2544 if (fLoader[iDet]) {
2545 fLoader[iDet]->LoadRecPoints("read");
2546 clustersTree = fLoader[iDet]->TreeR();
2547 if (!clustersTree) {
2548 AliError(Form("Can't get the %s clusters tree",
2549 fgkDetectorName[iDet]));
2550 if (fStopOnError) return kFALSE;
2553 if (fRawReader && !reconstructor->HasDigitConversion()) {
2554 reconstructor->FillESD(fRawReader, clustersTree, esd);
2556 TTree* digitsTree = NULL;
2557 if (fLoader[iDet]) {
2558 fLoader[iDet]->LoadDigits("read");
2559 digitsTree = fLoader[iDet]->TreeD();
2561 AliError(Form("Can't get the %s digits tree",
2562 fgkDetectorName[iDet]));
2563 if (fStopOnError) return kFALSE;
2566 reconstructor->FillESD(digitsTree, clustersTree, esd);
2567 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2569 if (fLoader[iDet]) {
2570 fLoader[iDet]->UnloadRecPoints();
2574 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2575 AliError(Form("the following detectors were not found: %s",
2577 if (fStopOnError) return kFALSE;
2579 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
2584 //_____________________________________________________________________________
2585 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2587 // Reads the trigger decision which is
2588 // stored in Trigger.root file and fills
2589 // the corresponding esd entries
2591 AliCodeTimerAuto("")
2593 AliInfo("Filling trigger information into the ESD");
2596 AliCTPRawStream input(fRawReader);
2597 if (!input.Next()) {
2598 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2601 if (esd->GetTriggerMask() != input.GetClassMask())
2602 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2603 input.GetClassMask(),esd->GetTriggerMask()));
2604 if (esd->GetOrbitNumber() != input.GetOrbitID())
2605 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2606 input.GetOrbitID(),esd->GetOrbitNumber()));
2607 if (esd->GetBunchCrossNumber() != input.GetBCID())
2608 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2609 input.GetBCID(),esd->GetBunchCrossNumber()));
2610 AliESDHeader* esdheader = esd->GetHeader();
2611 esdheader->SetL0TriggerInputs(input.GetL0Inputs());
2612 esdheader->SetL1TriggerInputs(input.GetL1Inputs());
2613 esdheader->SetL2TriggerInputs(input.GetL2Inputs());
2615 UInt_t orbit=input.GetOrbitID();
2616 for(Int_t i=0 ; i<input.GetNIRs() ; i++ )
2617 if(TMath::Abs(Int_t(orbit-(input.GetIR(i))->GetOrbit()))<=1){
2618 esdheader->AddTriggerIR(input.GetIR(i));
2624 //_____________________________________________________________________________
2625 Bool_t AliReconstruction::FillTriggerScalers(AliESDEvent*& esd)
2628 //fRunScalers->Print();
2629 if(fRunScalers && fRunScalers->CheckRunScalers()){
2630 AliTimeStamp* timestamp = new AliTimeStamp(esd->GetOrbitNumber(), esd->GetPeriodNumber(), esd->GetBunchCrossNumber());
2631 //AliTimeStamp* timestamp = new AliTimeStamp(10308000, 0, (ULong64_t)486238);
2632 AliESDHeader* esdheader = fesd->GetHeader();
2633 for(Int_t i=0;i<50;i++){
2634 if((1<<i) & esd->GetTriggerMask()){
2635 AliTriggerScalersESD* scalesd = fRunScalers->GetScalersForEventClass( timestamp, i+1);
2636 if(scalesd)esdheader->SetTriggerScalersRecord(scalesd);
2642 //_____________________________________________________________________________
2643 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
2646 // Filling information from RawReader Header
2649 if (!fRawReader) return kFALSE;
2651 AliInfo("Filling information from RawReader Header");
2653 esd->SetBunchCrossNumber(fRawReader->GetBCID());
2654 esd->SetOrbitNumber(fRawReader->GetOrbitID());
2655 esd->SetPeriodNumber(fRawReader->GetPeriod());
2657 esd->SetTimeStamp(fRawReader->GetTimestamp());
2658 esd->SetEventType(fRawReader->GetType());
2664 //_____________________________________________________________________________
2665 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
2667 // check whether detName is contained in detectors
2668 // if yes, it is removed from detectors
2670 // check if all detectors are selected
2671 if ((detectors.CompareTo("ALL") == 0) ||
2672 detectors.BeginsWith("ALL ") ||
2673 detectors.EndsWith(" ALL") ||
2674 detectors.Contains(" ALL ")) {
2679 // search for the given detector
2680 Bool_t result = kFALSE;
2681 if ((detectors.CompareTo(detName) == 0) ||
2682 detectors.BeginsWith(detName+" ") ||
2683 detectors.EndsWith(" "+detName) ||
2684 detectors.Contains(" "+detName+" ")) {
2685 detectors.ReplaceAll(detName, "");
2689 // clean up the detectors string
2690 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
2691 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
2692 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
2697 //_____________________________________________________________________________
2698 Bool_t AliReconstruction::InitRunLoader()
2700 // get or create the run loader
2702 if (gAlice) delete gAlice;
2705 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
2706 // load all base libraries to get the loader classes
2707 TString libs = gSystem->GetLibraries();
2708 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2709 TString detName = fgkDetectorName[iDet];
2710 if (detName == "HLT") continue;
2711 if (libs.Contains("lib" + detName + "base.so")) continue;
2712 gSystem->Load("lib" + detName + "base.so");
2714 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2716 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2721 fRunLoader->CdGAFile();
2722 fRunLoader->LoadgAlice();
2724 //PH This is a temporary fix to give access to the kinematics
2725 //PH that is needed for the labels of ITS clusters
2726 fRunLoader->LoadHeader();
2727 fRunLoader->LoadKinematics();
2729 } else { // galice.root does not exist
2731 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2733 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2734 AliConfig::GetDefaultEventFolderName(),
2737 AliError(Form("could not create run loader in file %s",
2738 fGAliceFileName.Data()));
2742 fIsNewRunLoader = kTRUE;
2743 fRunLoader->MakeTree("E");
2745 if (fNumberOfEventsPerFile > 0)
2746 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2748 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2754 //_____________________________________________________________________________
2755 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2757 // get the reconstructor object and the loader for a detector
2759 if (fReconstructor[iDet]) {
2760 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2761 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2762 fReconstructor[iDet]->SetRecoParam(par);
2763 fReconstructor[iDet]->SetRunInfo(fRunInfo);
2765 return fReconstructor[iDet];
2768 // load the reconstructor object
2769 TPluginManager* pluginManager = gROOT->GetPluginManager();
2770 TString detName = fgkDetectorName[iDet];
2771 TString recName = "Ali" + detName + "Reconstructor";
2773 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2775 AliReconstructor* reconstructor = NULL;
2776 // first check if a plugin is defined for the reconstructor
2777 TPluginHandler* pluginHandler =
2778 pluginManager->FindHandler("AliReconstructor", detName);
2779 // if not, add a plugin for it
2780 if (!pluginHandler) {
2781 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2782 TString libs = gSystem->GetLibraries();
2783 if (libs.Contains("lib" + detName + "base.so") ||
2784 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2785 pluginManager->AddHandler("AliReconstructor", detName,
2786 recName, detName + "rec", recName + "()");
2788 pluginManager->AddHandler("AliReconstructor", detName,
2789 recName, detName, recName + "()");
2791 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2793 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2794 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2796 if (reconstructor) {
2797 TObject* obj = fOptions.FindObject(detName.Data());
2798 if (obj) reconstructor->SetOption(obj->GetTitle());
2799 reconstructor->SetRunInfo(fRunInfo);
2800 reconstructor->Init();
2801 fReconstructor[iDet] = reconstructor;
2804 // get or create the loader
2805 if (detName != "HLT") {
2806 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2807 if (!fLoader[iDet]) {
2808 AliConfig::Instance()
2809 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2811 // first check if a plugin is defined for the loader
2813 pluginManager->FindHandler("AliLoader", detName);
2814 // if not, add a plugin for it
2815 if (!pluginHandler) {
2816 TString loaderName = "Ali" + detName + "Loader";
2817 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2818 pluginManager->AddHandler("AliLoader", detName,
2819 loaderName, detName + "base",
2820 loaderName + "(const char*, TFolder*)");
2821 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2823 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2825 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2826 fRunLoader->GetEventFolder());
2828 if (!fLoader[iDet]) { // use default loader
2829 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2831 if (!fLoader[iDet]) {
2832 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2833 if (fStopOnError) return NULL;
2835 fRunLoader->AddLoader(fLoader[iDet]);
2836 fRunLoader->CdGAFile();
2837 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2838 fRunLoader->Write(0, TObject::kOverwrite);
2843 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
2844 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
2845 reconstructor->SetRecoParam(par);
2846 reconstructor->SetRunInfo(fRunInfo);
2848 return reconstructor;
2851 //_____________________________________________________________________________
2852 AliVertexer* AliReconstruction::CreateVertexer()
2854 // create the vertexer
2855 // Please note that the caller is the owner of the
2858 AliVertexer* vertexer = NULL;
2859 AliReconstructor* itsReconstructor = GetReconstructor(0);
2860 if (itsReconstructor && ((fRunLocalReconstruction.Contains("ITS")) || fRunTracking.Contains("ITS"))) {
2861 vertexer = itsReconstructor->CreateVertexer();
2864 AliWarning("couldn't create a vertexer for ITS");
2870 //_____________________________________________________________________________
2871 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2873 // create the trackers
2874 AliInfo("Creating trackers");
2876 TString detStr = detectors;
2877 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2878 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2879 AliReconstructor* reconstructor = GetReconstructor(iDet);
2880 if (!reconstructor) continue;
2881 TString detName = fgkDetectorName[iDet];
2882 if (detName == "HLT") {
2883 fRunHLTTracking = kTRUE;
2886 if (detName == "MUON") {
2887 fRunMuonTracking = kTRUE;
2892 fTracker[iDet] = reconstructor->CreateTracker();
2893 if (!fTracker[iDet] && (iDet < 7)) {
2894 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2895 if (fStopOnError) return kFALSE;
2897 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2903 //_____________________________________________________________________________
2904 void AliReconstruction::CleanUp()
2906 // delete trackers and the run loader and close and delete the file
2908 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2909 delete fReconstructor[iDet];
2910 fReconstructor[iDet] = NULL;
2911 fLoader[iDet] = NULL;
2912 delete fTracker[iDet];
2913 fTracker[iDet] = NULL;
2918 delete fSPDTrackleter;
2919 fSPDTrackleter = NULL;
2928 delete fParentRawReader;
2929 fParentRawReader=NULL;
2937 if (AliQAManager::QAManager())
2938 AliQAManager::QAManager()->ShowQA() ;
2939 AliQAManager::Destroy() ;
2943 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2945 // Write space-points which are then used in the alignment procedures
2946 // For the moment only ITS, TPC, TRD and TOF
2948 Int_t ntracks = esd->GetNumberOfTracks();
2949 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2951 AliESDtrack *track = esd->GetTrack(itrack);
2954 for (Int_t i=0; i<200; ++i) idx[i] = -1; //PH avoid uninitialized values
2955 for (Int_t iDet = 5; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
2956 nsp += track->GetNcls(iDet);
2958 if (iDet==0) { // ITS "extra" clusters
2959 track->GetClusters(iDet,idx);
2960 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
2965 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2966 track->SetTrackPointArray(sp);
2968 for (Int_t iDet = 5; iDet >= 0; iDet--) {
2969 AliTracker *tracker = fTracker[iDet];
2970 if (!tracker) continue;
2971 Int_t nspdet = track->GetClusters(iDet,idx);
2973 if (iDet==0) // ITS "extra" clusters
2974 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
2976 if (nspdet <= 0) continue;
2980 while (isp2 < nspdet) {
2981 Bool_t isvalid=kTRUE;
2983 Int_t index=idx[isp++];
2984 if (index < 0) continue;
2986 TString dets = fgkDetectorName[iDet];
2987 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2988 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2989 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2990 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2991 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
2993 isvalid = tracker->GetTrackPoint(index,p);
2996 if (!isvalid) continue;
2997 if (iDet==0 && (isp-1)>=6) p.SetExtra();
2998 sp->AddPoint(isptrack,&p); isptrack++;
3005 //_____________________________________________________________________________
3006 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
3008 // The method reads the raw-data error log
3009 // accumulated within the rawReader.
3010 // It extracts the raw-data errors related to
3011 // the current event and stores them into
3012 // a TClonesArray inside the esd object.
3014 if (!fRawReader) return;
3016 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
3018 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
3020 if (iEvent != log->GetEventNumber()) continue;
3022 esd->AddRawDataErrorLog(log);
3027 //_____________________________________________________________________________
3028 void AliReconstruction::CheckQA()
3030 // check the QA of SIM for this run and remove the detectors
3031 // with status Fatal
3033 // TString newRunLocalReconstruction ;
3034 // TString newRunTracking ;
3035 // TString newFillESD ;
3037 // for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
3038 // TString detName(AliQAv1::GetDetName(iDet)) ;
3039 // AliQAv1 * qa = AliQAv1::Instance(AliQAv1::DETECTORINDEX_t(iDet)) ;
3040 // if ( qa->IsSet(AliQAv1::DETECTORINDEX_t(iDet), AliQAv1::kSIM, specie, AliQAv1::kFATAL)) {
3041 // AliInfo(Form("QA status for %s %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed",
3042 // detName.Data(), AliRecoParam::GetEventSpecieName(es))) ;
3044 // if ( fRunLocalReconstruction.Contains(AliQAv1::GetDetName(iDet)) ||
3045 // fRunLocalReconstruction.Contains("ALL") ) {
3046 // newRunLocalReconstruction += detName ;
3047 // newRunLocalReconstruction += " " ;
3049 // if ( fRunTracking.Contains(AliQAv1::GetDetName(iDet)) ||
3050 // fRunTracking.Contains("ALL") ) {
3051 // newRunTracking += detName ;
3052 // newRunTracking += " " ;
3054 // if ( fFillESD.Contains(AliQAv1::GetDetName(iDet)) ||
3055 // fFillESD.Contains("ALL") ) {
3056 // newFillESD += detName ;
3057 // newFillESD += " " ;
3061 // fRunLocalReconstruction = newRunLocalReconstruction ;
3062 // fRunTracking = newRunTracking ;
3063 // fFillESD = newFillESD ;
3066 //_____________________________________________________________________________
3067 Int_t AliReconstruction::GetDetIndex(const char* detector)
3069 // return the detector index corresponding to detector
3071 for (index = 0; index < kNDetectors ; index++) {
3072 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
3077 //_____________________________________________________________________________
3078 Bool_t AliReconstruction::FinishPlaneEff() {
3080 // Here execute all the necessary operationis, at the end of the tracking phase,
3081 // in case that evaluation of PlaneEfficiencies was required for some detector.
3082 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
3084 // This Preliminary version works only FOR ITS !!!!!
3085 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3088 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3091 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
3092 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
3093 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3094 if(fTracker[iDet] && fTracker[iDet]->GetPlaneEff()) {
3095 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
3096 TString name=planeeff->GetName();
3098 TFile* pefile = TFile::Open(name, "RECREATE");
3099 ret=(Bool_t)planeeff->Write();
3101 if(planeeff->GetCreateHistos()) {
3102 TString hname=planeeff->GetName();
3103 hname+="Histo.root";
3104 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
3107 if(fSPDTrackleter) {
3108 AliPlaneEff *planeeff=fSPDTrackleter->GetPlaneEff();
3109 TString name="AliITSPlaneEffSPDtracklet.root";
3110 TFile* pefile = TFile::Open(name, "RECREATE");
3111 ret=(Bool_t)planeeff->Write();
3113 AliESDEvent *dummy=NULL;
3114 ret=(Bool_t)fSPDTrackleter->PostProcess(dummy); // take care of writing other files
3119 //_____________________________________________________________________________
3120 Bool_t AliReconstruction::InitPlaneEff() {
3122 // Here execute all the necessary operations, before of the tracking phase,
3123 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
3124 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
3125 // which should be updated/recalculated.
3127 // This Preliminary version will work only FOR ITS !!!!!
3128 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3131 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3133 AliWarning(Form("Implementation of this method not yet completed !! Method return kTRUE"));
3135 fSPDTrackleter = NULL;
3136 AliReconstructor* itsReconstructor = GetReconstructor(0);
3137 if (itsReconstructor) {
3138 fSPDTrackleter = itsReconstructor->CreateTrackleter(); // this is NULL unless required in RecoParam
3140 if (fSPDTrackleter) {
3141 AliInfo("Trackleter for SPD has been created");
3147 //_____________________________________________________________________________
3148 Bool_t AliReconstruction::InitAliEVE()
3150 // This method should be called only in case
3151 // AliReconstruction is run
3152 // within the alieve environment.
3153 // It will initialize AliEVE in a way
3154 // so that it can visualize event processed
3155 // by AliReconstruction.
3156 // The return flag shows whenever the
3157 // AliEVE initialization was successful or not.
3160 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
3161 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
3162 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
3164 gROOT->ProcessLine("if (!AliEveEventManager::GetMaster()){new AliEveEventManager();AliEveEventManager::GetMaster()->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(AliEveEventManager::GetMaster());};");
3165 gROOT->ProcessLine("alieve_online_init()");
3170 //_____________________________________________________________________________
3171 void AliReconstruction::RunAliEVE()
3173 // Runs AliEVE visualisation of
3174 // the current event.
3175 // Should be executed only after
3176 // successful initialization of AliEVE.
3178 AliInfo("Running AliEVE...");
3179 gROOT->ProcessLine(Form("AliEveEventManager::GetMaster()->SetEvent((AliRunLoader*)0x%lx,(AliRawReader*)0x%lx,(AliESDEvent*)0x%lx,(AliESDfriend*)0x%lx);",fRunLoader,fRawReader,fesd,fesdf));
3183 //_____________________________________________________________________________
3184 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
3186 // Allows to run QA for a selected set of detectors
3187 // and a selected set of tasks among RAWS, DIGITSR, RECPOINTS and ESDS
3188 // all selected detectors run the same selected tasks
3190 if (!detAndAction.Contains(":")) {
3191 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
3195 Int_t colon = detAndAction.Index(":") ;
3196 fQADetectors = detAndAction(0, colon) ;
3197 if (fQADetectors.Contains("ALL") )
3198 fQADetectors = fFillESD ;
3199 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
3200 if (fQATasks.Contains("ALL") ) {
3201 fQATasks = Form("%d %d %d %d", AliQAv1::kRAWS, AliQAv1::kDIGITSR, AliQAv1::kRECPOINTS, AliQAv1::kESDS) ;
3203 fQATasks.ToUpper() ;
3205 if ( fQATasks.Contains("RAW") )
3206 tempo = Form("%d ", AliQAv1::kRAWS) ;
3207 if ( fQATasks.Contains("DIGIT") )
3208 tempo += Form("%d ", AliQAv1::kDIGITSR) ;
3209 if ( fQATasks.Contains("RECPOINT") )
3210 tempo += Form("%d ", AliQAv1::kRECPOINTS) ;
3211 if ( fQATasks.Contains("ESD") )
3212 tempo += Form("%d ", AliQAv1::kESDS) ;
3214 if (fQATasks.IsNull()) {
3215 AliInfo("No QA requested\n") ;
3220 TString tempo(fQATasks) ;
3221 tempo.ReplaceAll(Form("%d", AliQAv1::kRAWS), AliQAv1::GetTaskName(AliQAv1::kRAWS)) ;
3222 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITSR), AliQAv1::GetTaskName(AliQAv1::kDIGITSR)) ;
3223 tempo.ReplaceAll(Form("%d", AliQAv1::kRECPOINTS), AliQAv1::GetTaskName(AliQAv1::kRECPOINTS)) ;
3224 tempo.ReplaceAll(Form("%d", AliQAv1::kESDS), AliQAv1::GetTaskName(AliQAv1::kESDS)) ;
3225 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
3230 //_____________________________________________________________________________
3231 Bool_t AliReconstruction::InitRecoParams()
3233 // The method accesses OCDB and retrieves all
3234 // the available reco-param objects from there.
3236 Bool_t isOK = kTRUE;
3238 if (fRecoParam.GetDetRecoParamArray(kNDetectors)) {
3239 AliInfo("Using custom GRP reconstruction parameters");
3242 AliInfo("Loading GRP reconstruction parameter objects");
3244 AliCDBPath path("GRP","Calib","RecoParam");
3245 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3247 AliWarning("Couldn't find GRP RecoParam entry in OCDB");
3251 TObject *recoParamObj = entry->GetObject();
3252 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3253 // GRP has a normal TobjArray of AliDetectorRecoParam objects
3254 // Registering them in AliRecoParam
3255 fRecoParam.AddDetRecoParamArray(kNDetectors,dynamic_cast<TObjArray*>(recoParamObj));
3257 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3258 // GRP has only onse set of reco parameters
3259 // Registering it in AliRecoParam
3260 AliInfo("Single set of GRP reconstruction parameters found");
3261 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3262 fRecoParam.AddDetRecoParam(kNDetectors,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3265 AliError("No valid GRP RecoParam object found in the OCDB");
3272 TString detStr = fLoadCDB;
3273 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3275 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3277 if (fRecoParam.GetDetRecoParamArray(iDet)) {
3278 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
3282 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
3284 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
3285 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3287 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
3291 TObject *recoParamObj = entry->GetObject();
3292 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3293 // The detector has a normal TobjArray of AliDetectorRecoParam objects
3294 // Registering them in AliRecoParam
3295 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
3297 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3298 // The detector has only onse set of reco parameters
3299 // Registering it in AliRecoParam
3300 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
3301 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3302 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3305 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
3309 // FIX ME: We have to disable the unloading of reco-param CDB
3310 // entries because QA framework is using them. Has to be fix in
3311 // a way that the QA takes the objects already constructed in
3313 // AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
3317 if (AliDebugLevel() > 0) fRecoParam.Print();
3322 //_____________________________________________________________________________
3323 Bool_t AliReconstruction::GetEventInfo()
3325 // Fill the event info object
3327 AliCodeTimerAuto("")
3329 AliCentralTrigger *aCTP = NULL;
3331 fEventInfo.SetEventType(fRawReader->GetType());
3333 ULong64_t mask = fRawReader->GetClassMask();
3334 fEventInfo.SetTriggerMask(mask);
3335 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
3336 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
3338 aCTP = new AliCentralTrigger();
3339 TString configstr("");
3340 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
3341 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
3345 aCTP->SetClassMask(mask);
3346 aCTP->SetClusterMask(clmask);
3349 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
3351 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
3352 aCTP = fRunLoader->GetTrigger();
3353 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
3354 // get inputs from actp - just get
3355 AliESDHeader* esdheader = fesd->GetHeader();
3356 esdheader->SetL0TriggerInputs(aCTP->GetL0TriggerInputs());
3357 esdheader->SetL1TriggerInputs(aCTP->GetL1TriggerInputs());
3358 esdheader->SetL2TriggerInputs(aCTP->GetL2TriggerInputs());
3359 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
3362 AliWarning("No trigger can be loaded! The trigger information will not be used!");
3367 AliTriggerConfiguration *config = aCTP->GetConfiguration();
3369 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
3370 if (fRawReader) delete aCTP;
3374 UChar_t clustmask = 0;
3376 ULong64_t trmask = fEventInfo.GetTriggerMask();
3377 const TObjArray& classesArray = config->GetClasses();
3378 Int_t nclasses = classesArray.GetEntriesFast();
3379 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
3380 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
3382 Int_t trindex = TMath::Nint(TMath::Log2(trclass->GetMask()));
3383 fesd->SetTriggerClass(trclass->GetName(),trindex);
3384 if (fRawReader) fRawReader->LoadTriggerClass(trclass->GetName(),trindex);
3385 if (trmask & (1ull << trindex)) {
3387 trclasses += trclass->GetName();
3389 clustmask |= trclass->GetCluster()->GetClusterMask();
3393 fEventInfo.SetTriggerClasses(trclasses);
3395 // Set the information in ESD
3396 fesd->SetTriggerMask(trmask);
3397 fesd->SetTriggerCluster(clustmask);
3399 if (!aCTP->CheckTriggeredDetectors()) {
3400 if (fRawReader) delete aCTP;
3404 if (fRawReader) delete aCTP;
3406 // We have to fill also the HLT decision here!!
3412 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3414 // Match the detector list found in the rec.C or the default 'ALL'
3415 // to the list found in the GRP (stored there by the shuttle PP which
3416 // gets the information from ECS)
3417 static TString resultList;
3418 TString detList = detectorList;
3422 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3423 if ((detectorMask >> iDet) & 0x1) {
3424 TString det = AliDAQ::OfflineModuleName(iDet);
3425 if ((detList.CompareTo("ALL") == 0) ||
3426 ((detList.BeginsWith("ALL ") ||
3427 detList.EndsWith(" ALL") ||
3428 detList.Contains(" ALL ")) &&
3429 !(detList.BeginsWith("-"+det+" ") ||
3430 detList.EndsWith(" -"+det) ||
3431 detList.Contains(" -"+det+" "))) ||
3432 (detList.CompareTo(det) == 0) ||
3433 detList.BeginsWith(det+" ") ||
3434 detList.EndsWith(" "+det) ||
3435 detList.Contains( " "+det+" " )) {
3436 if (!resultList.EndsWith(det + " ")) {
3445 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3446 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3447 if ((detList.CompareTo("ALL") == 0) ||
3448 ((detList.BeginsWith("ALL ") ||
3449 detList.EndsWith(" ALL") ||
3450 detList.Contains(" ALL ")) &&
3451 !(detList.BeginsWith("-"+hltDet+" ") ||
3452 detList.EndsWith(" -"+hltDet) ||
3453 detList.Contains(" -"+hltDet+" "))) ||
3454 (detList.CompareTo(hltDet) == 0) ||
3455 detList.BeginsWith(hltDet+" ") ||
3456 detList.EndsWith(" "+hltDet) ||
3457 detList.Contains( " "+hltDet+" " )) {
3458 resultList += hltDet;
3462 return resultList.Data();
3466 //______________________________________________________________________________
3467 void AliReconstruction::Abort(const char *method, EAbort what)
3469 // Abort processing. If what = kAbortProcess, the Process() loop will be
3470 // aborted. If what = kAbortFile, the current file in a chain will be
3471 // aborted and the processing will continue with the next file, if there
3472 // is no next file then Process() will be aborted. Abort() can also be
3473 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3474 // the SlaveTerminate() and Terminate() are always called. The abort flag
3475 // can be checked in these methods using GetAbort().
3477 // The method is overwritten in AliReconstruction for better handling of
3478 // reco specific errors
3480 if (!fStopOnError) return;
3484 TString whyMess = method;
3485 whyMess += " failed! Aborting...";
3487 AliError(whyMess.Data());
3490 TString mess = "Abort";
3491 if (fAbort == kAbortProcess)
3492 mess = "AbortProcess";
3493 else if (fAbort == kAbortFile)
3496 Info(mess, whyMess.Data());
3499 //______________________________________________________________________________
3500 Bool_t AliReconstruction::ProcessEvent(void* event)
3502 // Method that is used in case the event loop
3503 // is steered from outside, for example by AMORE
3504 // 'event' is a pointer to the DATE event in the memory
3506 if (fRawReader) delete fRawReader;
3507 fRawReader = new AliRawReaderDate(event);
3508 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
3515 //______________________________________________________________________________
3516 Bool_t AliReconstruction::ParseOutput()
3518 // The method parses the output file
3519 // location string in order to steer
3520 // properly the selector
3522 TPMERegexp re1("(\\w+\\.zip#\\w+\\.root):([,*\\w+\\.root,*]+)@dataset://(\\w++)");
3523 TPMERegexp re2("(\\w+\\.root)?@?dataset://(\\w++)");
3525 if (re1.Match(fESDOutput) == 4) {
3526 // root archive with output files stored and regustered
3528 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",re1[1].Data()));
3529 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",re1[3].Data()));
3530 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_DATASET",""));
3531 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_ARCHIVE",re1[2].Data()));
3532 AliInfo(Form("%s files will be stored within %s in dataset %s",
3537 else if (re2.Match(fESDOutput) == 3) {
3538 // output file stored and registered
3540 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",(re2[1].IsNull()) ? "AliESDs.root" : re2[1].Data()));
3541 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",re2[2].Data()));
3542 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_DATASET",""));
3543 AliInfo(Form("%s will be stored in dataset %s",
3544 (re2[1].IsNull()) ? "AliESDs.root" : re2[1].Data(),
3548 if (fESDOutput.IsNull()) {
3549 // Output location not given.
3550 // Assuming xrootd has been already started and
3551 // the output file has to be sent back
3552 // to the client machine
3553 TString esdUrl(Form("root://%s/%s/",
3554 TUrl(gSystem->HostName()).GetHostFQDN(),
3556 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE","AliESDs.root"));
3557 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",esdUrl.Data()));
3558 AliInfo(Form("AliESDs.root will be stored in %s",
3562 // User specified an output location.
3563 // Ones has just to parse it here
3564 TUrl outputUrl(fESDOutput.Data());
3565 TString outputFile(gSystem->BaseName(outputUrl.GetFile()));
3566 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",outputFile.IsNull() ? "AliESDs.root" : outputFile.Data()));
3567 TString outputLocation(outputUrl.GetUrl());
3568 outputLocation.ReplaceAll(outputFile.Data(),"");
3569 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",outputLocation.Data()));
3570 AliInfo(Form("%s will be stored in %s",
3571 outputFile.IsNull() ? "AliESDs.root" : outputFile.Data(),
3572 outputLocation.Data()));