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>
135 #include "AliAlignObj.h"
136 #include "AliCDBEntry.h"
137 #include "AliCDBManager.h"
138 #include "AliCDBStorage.h"
139 #include "AliCTPRawStream.h"
140 #include "AliCascadeVertexer.h"
141 #include "AliCentralTrigger.h"
142 #include "AliCodeTimer.h"
144 #include "AliDetectorRecoParam.h"
145 #include "AliESDCaloCells.h"
146 #include "AliESDCaloCluster.h"
147 #include "AliESDEvent.h"
148 #include "AliESDMuonTrack.h"
149 #include "AliESDPmdTrack.h"
150 #include "AliESDTagCreator.h"
151 #include "AliESDVertex.h"
152 #include "AliESDcascade.h"
153 #include "AliESDfriend.h"
154 #include "AliESDkink.h"
155 #include "AliESDpid.h"
156 #include "AliESDtrack.h"
157 #include "AliESDtrack.h"
158 #include "AliEventInfo.h"
159 #include "AliGRPObject.h"
160 #include "AliGRPRecoParam.h"
161 #include "AliGenEventHeader.h"
162 #include "AliGeomManager.h"
163 #include "AliGlobalQADataMaker.h"
164 #include "AliHeader.h"
167 #include "AliMultiplicity.h"
169 #include "AliPlaneEff.h"
171 #include "AliQADataMakerRec.h"
172 #include "AliQAManager.h"
173 #include "AliRawVEvent.h"
174 #include "AliRawEventHeaderBase.h"
175 #include "AliRawHLTManager.h"
176 #include "AliRawReaderDate.h"
177 #include "AliRawReaderFile.h"
178 #include "AliRawReaderRoot.h"
179 #include "AliReconstruction.h"
180 #include "AliReconstructor.h"
182 #include "AliRunInfo.h"
183 #include "AliRunLoader.h"
184 #include "AliSysInfo.h" // memory snapshots
185 #include "AliTrackPointArray.h"
186 #include "AliTracker.h"
187 #include "AliTriggerClass.h"
188 #include "AliTriggerCluster.h"
189 #include "AliTriggerIR.h"
190 #include "AliTriggerConfiguration.h"
191 #include "AliV0vertexer.h"
192 #include "AliVertexer.h"
193 #include "AliTrackleter.h"
194 #include "AliVertexerTracks.h"
195 #include "AliTriggerRunScalers.h"
196 #include "AliCTPTimeParams.h"
197 #include "AliESDHLTDecision.h"
198 #include "AliTriggerInput.h"
199 #include "AliLHCData.h"
200 ClassImp(AliReconstruction)
202 //_____________________________________________________________________________
203 const char* AliReconstruction::fgkDetectorName[AliReconstruction::kNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
205 //_____________________________________________________________________________
206 AliReconstruction::AliReconstruction(const char* gAliceFilename) :
208 fRunVertexFinder(kTRUE),
209 fRunVertexFinderTracks(kTRUE),
210 fRunHLTTracking(kFALSE),
211 fRunMuonTracking(kFALSE),
213 fRunCascadeFinder(kTRUE),
214 fRunMultFinder(kTRUE),
216 fWriteAlignmentData(kFALSE),
217 fWriteESDfriend(kFALSE),
218 fFillTriggerESD(kTRUE),
226 fRunLocalReconstruction("ALL"),
230 fUseTrackingErrorsForAlignment(""),
231 fGAliceFileName(gAliceFilename),
234 fProofOutputFileName(""),
235 fProofOutputLocation(""),
236 fProofOutputDataset(kFALSE),
237 fProofOutputArchive(""),
241 fNumberOfEventsPerFile((UInt_t)-1),
242 fFractionFriends(0.04),
244 fLoadAlignFromCDB(kTRUE),
245 fLoadAlignData("ALL"),
250 fCTPTimeParams(NULL),
255 fParentRawReader(NULL),
259 fSPDTrackleter(NULL),
261 fDiamondProfileSPD(NULL),
262 fDiamondProfile(NULL),
263 fDiamondProfileTPC(NULL),
264 fListOfCosmicTriggers(NULL),
268 fAlignObjArray(NULL),
272 fInitCDBCalled(kFALSE),
273 fSetRunNumberFromDataCalled(kFALSE),
278 fSameQACycle(kFALSE),
279 fInitQACalled(kFALSE),
280 fWriteQAExpertData(kTRUE),
281 fRunPlaneEff(kFALSE),
292 fIsNewRunLoader(kFALSE),
302 // create reconstruction object with default parameters
305 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
306 fReconstructor[iDet] = NULL;
307 fUpgradeMask[iDet]=kFALSE;
308 fLoader[iDet] = NULL;
309 fTracker[iDet] = NULL;
311 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
312 fQACycles[iDet] = 999999 ;
313 fQAWriteExpert[iDet] = kFALSE ;
315 fBeamInt[0][0]=fBeamInt[0][1]=fBeamInt[1][0]=fBeamInt[1][1] = -1;
320 //_____________________________________________________________________________
321 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
323 fRunVertexFinder(rec.fRunVertexFinder),
324 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
325 fRunHLTTracking(rec.fRunHLTTracking),
326 fRunMuonTracking(rec.fRunMuonTracking),
327 fRunV0Finder(rec.fRunV0Finder),
328 fRunCascadeFinder(rec.fRunCascadeFinder),
329 fRunMultFinder(rec.fRunMultFinder),
330 fStopOnError(rec.fStopOnError),
331 fWriteAlignmentData(rec.fWriteAlignmentData),
332 fWriteESDfriend(rec.fWriteESDfriend),
333 fFillTriggerESD(rec.fFillTriggerESD),
335 fCleanESD(rec.fCleanESD),
336 fV0DCAmax(rec.fV0DCAmax),
337 fV0CsPmin(rec.fV0CsPmin),
341 fRunLocalReconstruction(rec.fRunLocalReconstruction),
342 fRunTracking(rec.fRunTracking),
343 fFillESD(rec.fFillESD),
344 fLoadCDB(rec.fLoadCDB),
345 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
346 fGAliceFileName(rec.fGAliceFileName),
347 fRawInput(rec.fRawInput),
348 fESDOutput(rec.fESDOutput),
349 fProofOutputFileName(rec.fProofOutputFileName),
350 fProofOutputLocation(rec.fProofOutputLocation),
351 fProofOutputDataset(rec.fProofOutputDataset),
352 fProofOutputArchive(rec.fProofOutputArchive),
353 fEquipIdMap(rec.fEquipIdMap),
354 fFirstEvent(rec.fFirstEvent),
355 fLastEvent(rec.fLastEvent),
356 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
357 fFractionFriends(rec.fFractionFriends),
359 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
360 fLoadAlignData(rec.fLoadAlignData),
361 fUseHLTData(rec.fUseHLTData),
365 fCTPTimeParams(NULL),
370 fParentRawReader(NULL),
372 fRecoParam(rec.fRecoParam),
374 fSPDTrackleter(NULL),
376 fDiamondProfileSPD(rec.fDiamondProfileSPD),
377 fDiamondProfile(rec.fDiamondProfile),
378 fDiamondProfileTPC(rec.fDiamondProfileTPC),
379 fListOfCosmicTriggers(NULL),
383 fAlignObjArray(rec.fAlignObjArray),
384 fCDBUri(rec.fCDBUri),
385 fQARefUri(rec.fQARefUri),
387 fInitCDBCalled(rec.fInitCDBCalled),
388 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
389 fQADetectors(rec.fQADetectors),
390 fQATasks(rec.fQATasks),
392 fRunGlobalQA(rec.fRunGlobalQA),
393 fSameQACycle(rec.fSameQACycle),
394 fInitQACalled(rec.fInitQACalled),
395 fWriteQAExpertData(rec.fWriteQAExpertData),
396 fRunPlaneEff(rec.fRunPlaneEff),
407 fIsNewRunLoader(rec.fIsNewRunLoader),
419 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
420 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
422 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
423 fReconstructor[iDet] = NULL;
424 fUpgradeMask[iDet] = kFALSE;
425 fLoader[iDet] = NULL;
426 fTracker[iDet] = NULL;
429 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
430 fQACycles[iDet] = rec.fQACycles[iDet];
431 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
434 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
435 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
438 for (int i=2;i--;) for (int j=2;j--;) fBeamInt[i][j] = rec.fBeamInt[i][j];
442 //_____________________________________________________________________________
443 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
445 // assignment operator
446 // Used in PROOF mode
447 // Be very careful while modifing it!
448 // Simple rules to follow:
449 // for persistent data members - use their assignment operators
450 // for non-persistent ones - do nothing or take the default values from constructor
451 // TSelector members should not be touched
452 if(&rec == this) return *this;
454 fRunVertexFinder = rec.fRunVertexFinder;
455 fRunVertexFinderTracks = rec.fRunVertexFinderTracks;
456 fRunHLTTracking = rec.fRunHLTTracking;
457 fRunMuonTracking = rec.fRunMuonTracking;
458 fRunV0Finder = rec.fRunV0Finder;
459 fRunCascadeFinder = rec.fRunCascadeFinder;
460 fRunMultFinder = rec.fRunMultFinder;
461 fStopOnError = rec.fStopOnError;
462 fWriteAlignmentData = rec.fWriteAlignmentData;
463 fWriteESDfriend = rec.fWriteESDfriend;
464 fFillTriggerESD = rec.fFillTriggerESD;
466 fCleanESD = rec.fCleanESD;
467 fV0DCAmax = rec.fV0DCAmax;
468 fV0CsPmin = rec.fV0CsPmin;
472 fRunLocalReconstruction = rec.fRunLocalReconstruction;
473 fRunTracking = rec.fRunTracking;
474 fFillESD = rec.fFillESD;
475 fLoadCDB = rec.fLoadCDB;
476 fUseTrackingErrorsForAlignment = rec.fUseTrackingErrorsForAlignment;
477 fGAliceFileName = rec.fGAliceFileName;
478 fRawInput = rec.fRawInput;
479 fESDOutput = rec.fESDOutput;
480 fProofOutputFileName = rec.fProofOutputFileName;
481 fProofOutputLocation = rec.fProofOutputLocation;
482 fProofOutputDataset = rec.fProofOutputDataset;
483 fProofOutputArchive = rec.fProofOutputArchive;
484 fEquipIdMap = rec.fEquipIdMap;
485 fFirstEvent = rec.fFirstEvent;
486 fLastEvent = rec.fLastEvent;
487 fNumberOfEventsPerFile = rec.fNumberOfEventsPerFile;
488 fFractionFriends = rec.fFractionFriends;
490 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
491 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
494 fLoadAlignFromCDB = rec.fLoadAlignFromCDB;
495 fLoadAlignData = rec.fLoadAlignData;
496 fUseHLTData = rec.fUseHLTData;
498 delete fRunInfo; fRunInfo = NULL;
499 if (rec.fRunInfo) fRunInfo = new AliRunInfo(*rec.fRunInfo);
501 fEventInfo = rec.fEventInfo;
503 delete fRunScalers; fRunScalers = NULL;
504 if (rec.fRunScalers) fRunScalers = new AliTriggerRunScalers(*rec.fRunScalers);
506 delete fCTPTimeParams; fCTPTimeParams = NULL;
507 if (rec.fCTPTimeParams) fCTPTimeParams = new AliCTPTimeParams(*rec.fCTPTimeParams);
508 delete fCTPTimeAlign; fCTPTimeAlign = NULL;
509 if (rec.fCTPTimeAlign) fCTPTimeAlign = new AliCTPTimeParams(*rec.fCTPTimeAlign);
513 fParentRawReader = NULL;
515 fRecoParam = rec.fRecoParam;
517 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
518 fUpgradeMask[iDet] = kFALSE;
519 delete fReconstructor[iDet]; fReconstructor[iDet] = NULL;
520 delete fLoader[iDet]; fLoader[iDet] = NULL;
521 delete fTracker[iDet]; fTracker[iDet] = NULL;
524 for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
525 fQACycles[iDet] = rec.fQACycles[iDet];
526 fQAWriteExpert[iDet] = rec.fQAWriteExpert[iDet] ;
529 delete fSPDTrackleter; fSPDTrackleter = NULL;
531 delete fDiamondProfileSPD; fDiamondProfileSPD = NULL;
532 if (rec.fDiamondProfileSPD) fDiamondProfileSPD = new AliESDVertex(*rec.fDiamondProfileSPD);
533 delete fDiamondProfile; fDiamondProfile = NULL;
534 if (rec.fDiamondProfile) fDiamondProfile = new AliESDVertex(*rec.fDiamondProfile);
535 delete fDiamondProfileTPC; fDiamondProfileTPC = NULL;
536 if (rec.fDiamondProfileTPC) fDiamondProfileTPC = new AliESDVertex(*rec.fDiamondProfileTPC);
538 delete fListOfCosmicTriggers; fListOfCosmicTriggers = NULL;
539 if (rec.fListOfCosmicTriggers) fListOfCosmicTriggers = (THashTable*)((rec.fListOfCosmicTriggers)->Clone());
541 delete fGRPData; fGRPData = NULL;
542 // if (rec.fGRPData) fGRPData = (TMap*)((rec.fGRPData)->Clone());
543 if (rec.fGRPData) fGRPData = (AliGRPObject*)((rec.fGRPData)->Clone());
545 delete fAlignObjArray; fAlignObjArray = NULL;
548 fQARefUri = rec.fQARefUri;
549 fSpecCDBUri.Delete();
550 fInitCDBCalled = rec.fInitCDBCalled;
551 fSetRunNumberFromDataCalled = rec.fSetRunNumberFromDataCalled;
552 fQADetectors = rec.fQADetectors;
553 fQATasks = rec.fQATasks;
555 fRunGlobalQA = rec.fRunGlobalQA;
556 fSameQACycle = rec.fSameQACycle;
557 fInitQACalled = rec.fInitQACalled;
558 fWriteQAExpertData = rec.fWriteQAExpertData;
559 fRunPlaneEff = rec.fRunPlaneEff;
560 for (int i=2;i--;) for (int j=2;j--;) fBeamInt[i][j] = rec.fBeamInt[i][j];
570 fIsNewRunLoader = rec.fIsNewRunLoader;
583 //_____________________________________________________________________________
584 AliReconstruction::~AliReconstruction()
589 if (fListOfCosmicTriggers) {
590 fListOfCosmicTriggers->Delete();
591 delete fListOfCosmicTriggers;
595 delete fCTPTimeParams;
596 delete fCTPTimeAlign;
598 if (fAlignObjArray) {
599 fAlignObjArray->Delete();
600 delete fAlignObjArray;
602 fSpecCDBUri.Delete();
604 AliCodeTimer::Instance()->Print();
607 //_____________________________________________________________________________
608 void AliReconstruction::InitQA()
610 //Initialize the QA and start of cycle
611 AliCodeTimerAuto("",0);
613 if (fInitQACalled) return;
614 fInitQACalled = kTRUE;
616 if (fGRPData) AliQADataMaker::SetCloningRequest( fGRPData->GetQATrigClasses(), fGRPData->GetQACloningRequest());
619 AliQAManager * qam = AliQAManager::QAManager(AliQAv1::kRECMODE) ;
620 if (fWriteQAExpertData)
621 qam->SetWriteExpert() ;
623 if (qam->IsDefaultStorageSet()) {
624 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
625 AliWarning("Default QA reference storage has been already set !");
626 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fQARefUri.Data()));
627 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
628 fQARefUri = qam->GetDefaultStorage()->GetURI();
630 if (fQARefUri.Length() > 0) {
631 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
632 AliDebug(2, Form("Default QA reference storage is set to: %s", fQARefUri.Data()));
633 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
635 fQARefUri="local://$ALICE_ROOT/QAref";
636 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
637 AliWarning("Default QA refeference storage not yet set !!!!");
638 AliWarning(Form("Setting it now to: %s", fQARefUri.Data()));
639 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
642 qam->SetDefaultStorage(fQARefUri);
646 qam->SetActiveDetectors(fQADetectors) ;
647 for (Int_t det = 0 ; det < AliQAv1::kNDET ; det++) {
648 qam->SetCycleLength(AliQAv1::DETECTORINDEX_t(det), fQACycles[det]) ;
649 qam->SetWriteExpert(AliQAv1::DETECTORINDEX_t(det)) ;
651 if (!fRawReader && !fInput && IsInTasks(AliQAv1::kRAWS))
652 fQATasks.ReplaceAll(Form("%d",AliQAv1::kRAWS), "") ;
653 qam->SetTasks(fQATasks) ;
654 qam->InitQADataMaker(AliCDBManager::Instance()->GetRun()) ;
657 Bool_t sameCycle = kFALSE ;
658 AliQADataMaker *qadm = qam->GetQADataMaker(AliQAv1::kGLOBAL);
659 AliInfo(Form("Initializing the global QA data maker"));
660 if (IsInTasks(AliQAv1::kRECPOINTS)) {
661 qadm->StartOfCycle(AliQAv1::kRECPOINTS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
662 TObjArray **arr=qadm->Init(AliQAv1::kRECPOINTS);
663 AliTracker::SetResidualsArray(arr);
666 if (IsInTasks(AliQAv1::kESDS)) {
667 qadm->StartOfCycle(AliQAv1::kESDS, AliCDBManager::Instance()->GetRun(), sameCycle) ;
668 qadm->Init(AliQAv1::kESDS);
671 AliSysInfo::AddStamp("InitQA") ;
674 //_____________________________________________________________________________
675 void AliReconstruction::MergeQA(const char *fileName)
677 //Initialize the QA and start of cycle
678 AliCodeTimerAuto("",0) ;
679 AliQAManager::QAManager()->Merge(AliCDBManager::Instance()->GetRun(),fileName) ;
680 AliSysInfo::AddStamp("MergeQA") ;
683 //_____________________________________________________________________________
684 void AliReconstruction::InitCDB()
686 // activate a default CDB storage
687 // First check if we have any CDB storage set, because it is used
688 // to retrieve the calibration and alignment constants
689 AliCodeTimerAuto("",0);
691 if (fInitCDBCalled) return;
692 fInitCDBCalled = kTRUE;
694 AliCDBManager* man = AliCDBManager::Instance();
695 if (man->IsDefaultStorageSet())
697 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
698 AliWarning("Default CDB storage has been already set !");
699 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
700 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
701 fCDBUri = man->GetDefaultStorage()->GetURI();
704 if (fCDBUri.Length() > 0)
706 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
707 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
708 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
709 man->SetDefaultStorage(fCDBUri);
711 else if (!man->GetRaw()){
712 fCDBUri="local://$ALICE_ROOT/OCDB";
713 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
714 AliWarning("Default CDB storage not yet set !!!!");
715 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
716 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
717 man->SetDefaultStorage(fCDBUri);
720 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
721 AliWarning("Default storage will be set after setting the Run Number!!!");
722 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
726 // Now activate the detector specific CDB storage locations
727 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
728 TObject* obj = fSpecCDBUri[i];
730 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
731 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
732 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
733 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
735 AliSysInfo::AddStamp("InitCDB");
738 //_____________________________________________________________________________
739 void AliReconstruction::SetDefaultStorage(const char* uri) {
740 // Store the desired default CDB storage location
741 // Activate it later within the Run() method
747 //_____________________________________________________________________________
748 void AliReconstruction::SetQARefDefaultStorage(const char* uri) {
749 // Store the desired default CDB storage location
750 // Activate it later within the Run() method
753 AliQAv1::SetQARefStorage(fQARefUri.Data()) ;
756 //_____________________________________________________________________________
757 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
758 // Store a detector-specific CDB storage location
759 // Activate it later within the Run() method
761 AliCDBPath aPath(calibType);
762 if(!aPath.IsValid()){
763 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
764 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
765 if(!strcmp(calibType, fgkDetectorName[iDet])) {
766 aPath.SetPath(Form("%s/*", calibType));
767 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
771 if(!aPath.IsValid()){
772 AliError(Form("Not a valid path or detector: %s", calibType));
777 // // check that calibType refers to a "valid" detector name
778 // Bool_t isDetector = kFALSE;
779 // for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
780 // TString detName = fgkDetectorName[iDet];
781 // if(aPath.GetLevel0() == detName) {
782 // isDetector = kTRUE;
788 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
792 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
793 if (obj) fSpecCDBUri.Remove(obj);
794 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
798 //_____________________________________________________________________________
799 Bool_t AliReconstruction::SetRunNumberFromData()
801 // The method is called in Run() in order
802 // to set a correct run number.
803 // In case of raw data reconstruction the
804 // run number is taken from the raw data header
806 if (fSetRunNumberFromDataCalled) return kTRUE;
807 fSetRunNumberFromDataCalled = kTRUE;
809 AliCDBManager* man = AliCDBManager::Instance();
812 if(fRawReader->NextEvent()) {
813 if(man->GetRun() > 0) {
814 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
816 man->SetRun(fRawReader->GetRunNumber());
817 fRawReader->RewindEvents();
820 if(man->GetRun() > 0) {
821 AliWarning("No raw-data events are found ! Using settings in AliCDBManager !");
824 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
830 AliRunLoader *rl = AliRunLoader::Open(fGAliceFileName.Data());
832 AliError(Form("No run loader found in file %s", fGAliceFileName.Data()));
837 // read run number from gAlice
838 if(rl->GetHeader()) {
839 man->SetRun(rl->GetHeader()->GetRun());
844 AliError("Neither run-loader header nor RawReader objects are found !");
856 //_____________________________________________________________________________
857 void AliReconstruction::SetCDBLock() {
858 // Set CDB lock: from now on it is forbidden to reset the run number
859 // or the default storage or to activate any further storage!
861 AliCDBManager::Instance()->SetLock(1);
864 //_____________________________________________________________________________
865 void AliReconstruction::MatchUpgradeDetector() {
866 // Translates detector name in a boolean.
867 // The boolean is used in GetReconstructor to load the
868 // upgrade reconstructor instead of the standard one.
869 for(Int_t iDet = 0; iDet < kNDetectors; iDet++) {
870 if(fUpgradeModule.Contains(fgkDetectorName[iDet])) fUpgradeMask[iDet]=kTRUE;
873 //_____________________________________________________________________________
874 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
876 // Read the alignment objects from CDB.
877 // Each detector is supposed to have the
878 // alignment objects in DET/Align/Data CDB path.
879 // All the detector objects are then collected,
880 // sorted by geometry level (starting from ALIC) and
881 // then applied to the TGeo geometry.
882 // Finally an overlaps check is performed.
884 // Load alignment data from CDB and fill fAlignObjArray
885 if(fLoadAlignFromCDB){
887 TString detStr = detectors;
888 TString loadAlObjsListOfDets = "";
890 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
891 if(!IsSelected(fgkDetectorName[iDet], detStr)) continue;
892 if(!strcmp(fgkDetectorName[iDet],"HLT")) continue;
894 if(AliGeomManager::GetNalignable(fgkDetectorName[iDet]) != 0)
896 loadAlObjsListOfDets += fgkDetectorName[iDet];
897 loadAlObjsListOfDets += " ";
899 } // end loop over detectors
901 if(AliGeomManager::GetNalignable("GRP") != 0)
902 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
903 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
904 AliCDBManager::Instance()->UnloadFromCache("*/Align/*");
906 // Check if the array with alignment objects was
907 // provided by the user. If yes, apply the objects
908 // to the present TGeo geometry
909 if (fAlignObjArray) {
910 if (gGeoManager && gGeoManager->IsClosed()) {
911 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
912 AliError("The misalignment of one or more volumes failed!"
913 "Compare the list of simulated detectors and the list of detector alignment data!");
918 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
924 if (fAlignObjArray) {
925 fAlignObjArray->Delete();
926 delete fAlignObjArray; fAlignObjArray=NULL;
932 //_____________________________________________________________________________
933 void AliReconstruction::SetGAliceFile(const char* fileName)
935 // set the name of the galice file
937 fGAliceFileName = fileName;
940 //_____________________________________________________________________________
941 void AliReconstruction::SetInput(const char* input)
943 // In case the input string starts with 'mem://', we run in an online mode
944 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
945 // file is assumed. One can give as an input:
946 // mem://: - events taken from DAQ monitoring libs online
948 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
949 if (input) fRawInput = input;
952 //_____________________________________________________________________________
953 void AliReconstruction::SetOutput(const char* output)
955 // Set the output ESD filename
956 // 'output' is a normalt ROOT url
957 // The method is used in case of raw-data reco with PROOF
958 if (output) fESDOutput = output;
961 //_____________________________________________________________________________
962 void AliReconstruction::SetOption(const char* detector, const char* option)
964 // set options for the reconstruction of a detector
966 TObject* obj = fOptions.FindObject(detector);
967 if (obj) fOptions.Remove(obj);
968 fOptions.Add(new TNamed(detector, option));
971 //_____________________________________________________________________________
972 void AliReconstruction::SetRecoParam(const char* detector, AliDetectorRecoParam *par)
974 // Set custom reconstruction parameters for a given detector
975 // Single set of parameters for all the events
977 // First check if the reco-params are global
978 if(!strcmp(detector, "GRP")) {
980 fRecoParam.AddDetRecoParam(kNDetectors,par);
984 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
985 if(!strcmp(detector, fgkDetectorName[iDet])) {
987 fRecoParam.AddDetRecoParam(iDet,par);
994 //_____________________________________________________________________________
995 Bool_t AliReconstruction::InitGRP() {
996 //------------------------------------
997 // Initialization of the GRP entry
998 //------------------------------------
999 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/Data");
1003 TMap* m = dynamic_cast<TMap*>(entry->GetObject()); // old GRP entry
1006 AliInfo("Found a TMap in GRP/GRP/Data, converting it into an AliGRPObject");
1008 fGRPData = new AliGRPObject();
1009 fGRPData->ReadValuesFromMap(m);
1013 AliInfo("Found an AliGRPObject in GRP/GRP/Data, reading it");
1014 fGRPData = dynamic_cast<AliGRPObject*>(entry->GetObject()); // new GRP entry
1018 // FIX ME: The unloading of GRP entry is temporarily disabled
1019 // because ZDC and VZERO are using it in order to initialize
1020 // their reconstructor objects. In the future one has to think
1021 // of propagating AliRunInfo to the reconstructors.
1022 // AliCDBManager::Instance()->UnloadFromCache("GRP/GRP/Data");
1026 AliError("No GRP entry found in OCDB!");
1030 TString lhcState = fGRPData->GetLHCState();
1031 if (lhcState==AliGRPObject::GetInvalidString()) {
1032 AliError("GRP/GRP/Data entry: missing value for the LHC state ! Using UNKNOWN");
1033 lhcState = "UNKNOWN";
1036 TString beamType = fGRPData->GetBeamType();
1037 if (beamType==AliGRPObject::GetInvalidString()) {
1038 AliError("GRP/GRP/Data entry: missing value for the beam type ! Using UNKNOWN");
1039 beamType = "UNKNOWN";
1042 Float_t beamEnergy = fGRPData->GetBeamEnergy();
1043 if (beamEnergy==AliGRPObject::GetInvalidFloat()) {
1044 AliError("GRP/GRP/Data entry: missing value for the beam energy ! Using 0");
1048 TString runType = fGRPData->GetRunType();
1049 if (runType==AliGRPObject::GetInvalidString()) {
1050 AliError("GRP/GRP/Data entry: missing value for the run type ! Using UNKNOWN");
1051 runType = "UNKNOWN";
1054 Int_t activeDetectors = fGRPData->GetDetectorMask();
1055 if (activeDetectors==AliGRPObject::GetInvalidUInt()) {
1056 AliError("GRP/GRP/Data entry: missing value for the detector mask ! Using 1074790399");
1057 activeDetectors = 1074790399;
1060 fRunInfo = new AliRunInfo(lhcState, beamType, beamEnergy, runType, activeDetectors);
1064 // Process the list of active detectors
1065 if (activeDetectors) {
1066 UInt_t detMask = activeDetectors;
1067 fRunLocalReconstruction = MatchDetectorList(fRunLocalReconstruction,detMask);
1068 fRunTracking = MatchDetectorList(fRunTracking,detMask);
1069 fFillESD = MatchDetectorList(fFillESD,detMask);
1070 fQADetectors = MatchDetectorList(fQADetectors,detMask);
1071 fLoadCDB.Form("%s %s %s %s",
1072 fRunLocalReconstruction.Data(),
1073 fRunTracking.Data(),
1075 fQADetectors.Data());
1076 fLoadCDB = MatchDetectorList(fLoadCDB,detMask);
1077 if (!((detMask >> AliDAQ::DetectorID("ITSSPD")) & 0x1) &&
1078 !((detMask >> AliDAQ::DetectorID("ITSSDD")) & 0x1) &&
1079 !((detMask >> AliDAQ::DetectorID("ITSSSD")) & 0x1) ) {
1080 // switch off the vertexer
1081 AliInfo("SPD,SDD,SSD is not in the list of active detectors. Vertexer and Trackleter are switched off.");
1082 fRunVertexFinder = kFALSE;
1083 fRunMultFinder = kFALSE;
1085 if (!((detMask >> AliDAQ::DetectorID("TRG")) & 0x1)) {
1086 // switch off the reading of CTP raw-data payload
1087 if (fFillTriggerESD) {
1088 AliInfo("CTP is not in the list of active detectors. CTP data reading switched off.");
1089 fFillTriggerESD = kFALSE;
1094 AliInfo("===================================================================================");
1095 AliInfo(Form("Running local reconstruction for detectors: %s",fRunLocalReconstruction.Data()));
1096 AliInfo(Form("Running tracking for detectors: %s",fRunTracking.Data()));
1097 AliInfo(Form("Filling ESD for detectors: %s",fFillESD.Data()));
1098 AliInfo(Form("Quality assurance is active for detectors: %s",fQADetectors.Data()));
1099 AliInfo(Form("CDB and reconstruction parameters are loaded for detectors: %s",fLoadCDB.Data()));
1100 AliInfo("===================================================================================");
1102 //*** Dealing with the magnetic field map
1103 if ( TGeoGlobalMagField::Instance()->IsLocked() ) {
1104 if (TGeoGlobalMagField::Instance()->GetField()->TestBit(AliMagF::kOverrideGRP)) {
1105 AliInfo("ExpertMode!!! GRP information will be ignored !");
1106 AliInfo("ExpertMode!!! Running with the externally locked B field !");
1109 AliInfo("Destroying existing B field instance!");
1110 delete TGeoGlobalMagField::Instance();
1113 if ( !TGeoGlobalMagField::Instance()->IsLocked() ) {
1114 // Construct the field map out of the information retrieved from GRP.
1117 Float_t l3Current = fGRPData->GetL3Current((AliGRPObject::Stats)0);
1118 if (l3Current == AliGRPObject::GetInvalidFloat()) {
1119 AliError("GRP/GRP/Data entry: missing value for the L3 current !");
1123 Char_t l3Polarity = fGRPData->GetL3Polarity();
1124 if (l3Polarity == AliGRPObject::GetInvalidChar()) {
1125 AliError("GRP/GRP/Data entry: missing value for the L3 polarity !");
1130 Float_t diCurrent = fGRPData->GetDipoleCurrent((AliGRPObject::Stats)0);
1131 if (diCurrent == AliGRPObject::GetInvalidFloat()) {
1132 AliError("GRP/GRP/Data entry: missing value for the dipole current !");
1136 Char_t diPolarity = fGRPData->GetDipolePolarity();
1137 if (diPolarity == AliGRPObject::GetInvalidChar()) {
1138 AliError("GRP/GRP/Data entry: missing value for the dipole polarity !");
1142 // read special bits for the polarity convention and map type
1143 Int_t polConvention = fGRPData->IsPolarityConventionLHC() ? AliMagF::kConvLHC : AliMagF::kConvDCS2008;
1144 Bool_t uniformB = fGRPData->IsUniformBMap();
1147 AliMagF* fld = AliMagF::CreateFieldMap(TMath::Abs(l3Current) * (l3Polarity ? -1:1),
1148 TMath::Abs(diCurrent) * (diPolarity ? -1:1),
1149 polConvention,uniformB,beamEnergy, beamType.Data());
1151 TGeoGlobalMagField::Instance()->SetField( fld );
1152 TGeoGlobalMagField::Instance()->Lock();
1153 AliInfo("Running with the B field constructed out of GRP !");
1155 else AliFatal("Failed to create a B field map !");
1157 else AliFatal("B field is neither set nor constructed from GRP ! Exitig...");
1160 //*** Get the diamond profiles from OCDB
1161 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
1163 fDiamondProfileSPD = dynamic_cast<AliESDVertex*> (entry->GetObject());
1165 AliError("No SPD diamond profile found in OCDB!");
1168 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
1170 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
1172 AliError("No diamond profile found in OCDB!");
1175 entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexTPC");
1177 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
1179 AliError("No TPC diamond profile found in OCDB!");
1182 entry = AliCDBManager::Instance()->Get("GRP/Calib/CosmicTriggers");
1184 fListOfCosmicTriggers = dynamic_cast<THashTable*>(entry->GetObject());
1186 AliCDBManager::Instance()->UnloadFromCache("GRP/Calib/CosmicTriggers");
1189 if (!fListOfCosmicTriggers) {
1190 AliWarning("Can not get list of cosmic triggers from OCDB! Cosmic event specie will be effectively disabled!");
1196 //_____________________________________________________________________________
1197 Bool_t AliReconstruction::LoadCDB()
1199 // Load CDB entries for all active detectors.
1200 // By default we load all the entries in <det>/Calib
1203 AliCodeTimerAuto("",0);
1205 AliCDBManager::Instance()->Get("GRP/CTP/Config");
1207 AliCDBManager::Instance()->Get("GRP/Calib/LHCClockPhase");
1209 TString detStr = fLoadCDB;
1210 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1211 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1212 AliCDBManager::Instance()->GetAll(Form("%s/Calib/*",fgkDetectorName[iDet]));
1215 // Temporary fix - one has to define the correct policy in order
1216 // to load the trigger OCDB entries only for the detectors that
1217 // in the trigger or that are needed in order to put correct
1218 // information in ESD
1219 AliCDBManager::Instance()->GetAll("TRIGGER/*/*");
1223 //_____________________________________________________________________________
1224 Bool_t AliReconstruction::LoadTriggerScalersCDB()
1226 // Load CTP scalers from OCDB.
1227 // The scalers are checked for consistency.
1229 AliCodeTimerAuto("",0);
1231 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/CTP/Scalers");
1235 AliInfo("Found an AliTriggerRunScalers in GRP/CTP/Scalers, reading it");
1236 fRunScalers = dynamic_cast<AliTriggerRunScalers*> (entry->GetObject());
1238 if (fRunScalers && (fRunScalers->CorrectScalersOverflow() == 0)) AliInfo("32bit Trigger counters corrected for overflow");
1243 //_____________________________________________________________________________
1244 Bool_t AliReconstruction::LoadCTPTimeParamsCDB()
1246 // Load CTP timing information (alignment)
1249 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/CTP/CTPtiming");
1250 if (!entry) return kFALSE;
1252 AliInfo("Found an AliCTPTimeParams in GRP/CTP/CTPtiming, reading it");
1253 fCTPTimeParams = dynamic_cast<AliCTPTimeParams*> (entry->GetObject());
1256 AliCDBEntry* entry2 = AliCDBManager::Instance()->Get("GRP/CTP/TimeAlign");
1257 if (!entry2) return kFALSE;
1259 AliInfo("Found an AliCTPTimeParams in GRP/CTP/TimeAlign, reading it");
1260 fCTPTimeAlign = dynamic_cast<AliCTPTimeParams*> (entry2->GetObject());
1261 entry2->SetOwner(0);
1266 //_____________________________________________________________________________
1267 Bool_t AliReconstruction::ReadIntensityInfoCDB()
1269 // Load LHC DIP data
1270 AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/GRP/LHCData");
1271 AliCDBEntry* entryCTP = AliCDBManager::Instance()->Get("GRP/CTP/Config");
1274 enum {kA,kB,kC,kE,kNMasks};
1275 AliTriggerConfiguration* conf = (AliTriggerConfiguration*)entryCTP->GetObject();
1276 const TObjArray& clArr = conf->GetClasses();
1277 TObjArray masks(kNMasks);
1279 AliTriggerClass* trClass = 0;
1281 masks.SetOwner(kFALSE);
1283 while ( (trClass=(AliTriggerClass*)next()) ) {
1284 TString trName = trClass->GetName();
1285 int ind = trName.Index("-"); // prefix in front of A,B,C,E
1286 if (ind<1) continue; // anomaly
1288 trName = trName.Data() + ind;
1289 AliTriggerBCMask* bcMask = trClass->GetBCMask();
1290 if (!bcMask) continue;
1292 if (trName.BeginsWith("-A-")) which |= 0x1<<kA;
1293 else if (trName.BeginsWith("-B-")) which |= 0x1<<kB;
1294 else if (trName.BeginsWith("-C-")) which |= 0x1<<kC;
1295 else if (trName.BeginsWith("-E-")) which |= 0x1<<kE;
1296 else if (trName.BeginsWith("-AC-")) which |= (0x1<<kA) | (0x1<<kC);
1297 else if (trName.BeginsWith("-CA-")) which |= (0x1<<kA) | (0x1<<kC);
1298 else { AliWarning(Form("Unknown trigger type %s\n",trClass->GetName())); continue;}
1300 for (int ip=kNMasks;ip--;) {
1301 if ( !(which&(0x1<<ip)) || masks[ip] ) continue; // does not match or already done
1302 masks[ip] = (TObject*)bcMask;
1305 if (nFound==kNMasks) break;
1308 AliInfo("Reading mean bunch intensities from GRP/GRP/LHCData");
1309 AliLHCData* dipData = dynamic_cast<AliLHCData*> (entry->GetObject());
1311 for (int ib=2;ib--;) {
1313 if (dipData && (dipData->GetMeanIntensity(ib,intI,intNI,&masks)>=0)) {
1314 fBeamInt[ib][0] = intI;
1315 fBeamInt[ib][1] = intNI;
1316 AliInfo(Form("Mean intensity for beam %d: Interacting:%.2e Non-Interacting:%.2e",ib,intI,intNI));
1324 //_____________________________________________________________________________
1325 Bool_t AliReconstruction::Run(const char* input)
1328 AliCodeTimerAuto("",0);
1331 if (GetAbort() != TSelector::kContinue) return kFALSE;
1333 TChain *chain = NULL;
1334 if (fRawReader && (chain = fRawReader->GetChain())) {
1335 Long64_t nEntries = (fLastEvent < 0) ? (TChain::kBigNumber) : (fLastEvent - fFirstEvent + 1);
1338 // Temporary fix for long raw-data runs (until socket timeout handling in PROOF is revised)
1339 gProof->Exec("gEnv->SetValue(\"Proof.SocketActivityTimeout\",-1)", kTRUE);
1342 gProof->Exec("TGrid::Connect(\"alien://\")",kTRUE);
1344 TMessage::EnableSchemaEvolutionForAll(kTRUE);
1345 gProof->Exec("TMessage::EnableSchemaEvolutionForAll(kTRUE)",kTRUE);
1347 gProof->AddInput(this);
1349 if (!ParseOutput()) return kFALSE;
1351 gProof->SetParameter("PROOF_MaxSlavesPerNode", 9999);
1353 chain->Process("AliReconstruction","",nEntries,fFirstEvent);
1356 chain->Process(this,"",nEntries,fFirstEvent);
1361 if (GetAbort() != TSelector::kContinue) return kFALSE;
1363 if (GetAbort() != TSelector::kContinue) return kFALSE;
1364 //******* The loop over events
1365 AliInfo("Starting looping over events");
1367 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
1368 (fRawReader && fRawReader->NextEvent())) {
1369 if (!ProcessEvent(iEvent)) {
1370 Abort("ProcessEvent",TSelector::kAbortFile);
1376 if (GetAbort() != TSelector::kContinue) return kFALSE;
1378 if (GetAbort() != TSelector::kContinue) return kFALSE;
1384 //_____________________________________________________________________________
1385 void AliReconstruction::InitRawReader(const char* input)
1387 // Init raw-reader and
1388 // set the input in case of raw data
1390 AliCodeTimerAuto("",0);
1392 if (input) fRawInput = input;
1393 fRawReader = AliRawReader::Create(fRawInput.Data());
1395 if (fRawInput.IsNull()) {
1396 AliInfo("Reconstruction will run over digits");
1399 AliFatal("Can not create raw-data reader ! Exiting...");
1403 if (!fEquipIdMap.IsNull() && fRawReader)
1404 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
1406 if (!fUseHLTData.IsNull()) {
1407 // create the RawReaderHLT which performs redirection of HLT input data for
1408 // the specified detectors
1409 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
1411 fParentRawReader=fRawReader;
1412 fRawReader=pRawReader;
1414 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
1417 AliSysInfo::AddStamp("CreateRawReader");
1420 //_____________________________________________________________________________
1421 void AliReconstruction::InitRun(const char* input)
1423 // Initialization of raw-reader,
1424 // run number, CDB etc.
1425 AliCodeTimerAuto("",0);
1426 AliSysInfo::AddStamp("Start");
1428 // Initialize raw-reader if any
1429 InitRawReader(input);
1431 // Initialize the CDB storage
1434 // Set run number in CDBManager (if it is not already set by the user)
1435 if (!SetRunNumberFromData()) {
1436 Abort("SetRunNumberFromData", TSelector::kAbortProcess);
1440 // Set CDB lock: from now on it is forbidden to reset the run number
1441 // or the default storage or to activate any further storage!
1446 //_____________________________________________________________________________
1447 void AliReconstruction::Begin(TTree *)
1449 // Initialize AlReconstruction before
1450 // going into the event loop
1451 // Should follow the TSelector convention
1452 // i.e. initialize only the object on the client side
1453 AliCodeTimerAuto("",0);
1455 AliReconstruction *reco = NULL;
1457 if ((reco = (AliReconstruction*)fInput->FindObject("AliReconstruction"))) {
1460 AliSysInfo::AddStamp("ReadInputInBegin");
1463 // Import ideal TGeo geometry and apply misalignment
1465 TString geom(gSystem->DirName(fGAliceFileName));
1466 geom += "/geometry.root";
1467 AliGeomManager::LoadGeometry(geom.Data());
1469 Abort("LoadGeometry", TSelector::kAbortProcess);
1472 AliSysInfo::AddStamp("LoadGeom");
1473 TString detsToCheck=fRunLocalReconstruction;
1474 if(!AliGeomManager::CheckSymNamesLUT(detsToCheck.Data())) {
1475 Abort("CheckSymNamesLUT", TSelector::kAbortProcess);
1478 AliSysInfo::AddStamp("CheckGeom");
1481 if (!MisalignGeometry(fLoadAlignData)) {
1482 Abort("MisalignGeometry", TSelector::kAbortProcess);
1485 AliCDBManager::Instance()->UnloadFromCache("GRP/Geometry/Data");
1486 AliSysInfo::AddStamp("MisalignGeom");
1489 Abort("InitGRP", TSelector::kAbortProcess);
1492 AliSysInfo::AddStamp("InitGRP");
1495 Abort("LoadCDB", TSelector::kAbortProcess);
1498 AliSysInfo::AddStamp("LoadCDB");
1500 if (!LoadTriggerScalersCDB()) {
1501 Abort("LoadTriggerScalersCDB", TSelector::kAbortProcess);
1504 AliSysInfo::AddStamp("LoadTriggerScalersCDB");
1506 if (!LoadCTPTimeParamsCDB()) {
1507 Abort("LoadCTPTimeParamsCDB", TSelector::kAbortProcess);
1510 AliSysInfo::AddStamp("LoadCTPTimeParamsCDB");
1512 if (!ReadIntensityInfoCDB()) {
1513 Abort("ReadIntensityInfoCDB", TSelector::kAbortProcess);
1516 AliSysInfo::AddStamp("ReadIntensityInfoCDB");
1518 // Read the reconstruction parameters from OCDB
1519 if (!InitRecoParams()) {
1520 AliWarning("Not all detectors have correct RecoParam objects initialized");
1522 AliSysInfo::AddStamp("InitRecoParams");
1524 if (fInput && gProof) {
1525 if (reco) *reco = *this;
1527 gGeoManager->SetName("Geometry");
1528 gProof->AddInputData(gGeoManager,kTRUE);
1530 gProof->AddInputData(const_cast<TMap*>(AliCDBManager::Instance()->GetEntryCache()),kTRUE);
1531 fInput->Add(new TParameter<Int_t>("RunNumber",AliCDBManager::Instance()->GetRun()));
1532 AliMagF *magFieldMap = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1533 magFieldMap->SetName("MagneticFieldMap");
1534 gProof->AddInputData(magFieldMap,kTRUE);
1539 //_____________________________________________________________________________
1540 void AliReconstruction::SlaveBegin(TTree*)
1542 // Initialization related to run-loader,
1543 // vertexer, trackers, recontructors
1544 // In proof mode it is executed on the slave
1545 AliCodeTimerAuto("",0);
1547 TProofOutputFile *outProofFile = NULL;
1549 if (AliDebugLevel() > 0) fInput->Print();
1550 if (AliDebugLevel() > 10) fInput->Dump();
1551 if (AliReconstruction *reco = (AliReconstruction*)fInput->FindObject("AliReconstruction")) {
1554 if (TGeoManager *tgeo = (TGeoManager*)fInput->FindObject("Geometry")) {
1556 AliGeomManager::SetGeometry(tgeo);
1558 if (TMap *entryCache = (TMap*)fInput->FindObject("CDBEntryCache")) {
1559 Int_t runNumber = -1;
1560 if (TProof::GetParameter(fInput,"RunNumber",runNumber) == 0) {
1561 AliCDBManager *man = AliCDBManager::Instance(entryCache,runNumber);
1562 man->SetCacheFlag(kTRUE);
1563 man->SetLock(kTRUE);
1567 if (AliMagF *map = (AliMagF*)fInput->FindObject("MagneticFieldMap")) {
1568 AliMagF *newMap = new AliMagF(*map);
1569 if (!newMap->LoadParameterization()) {
1570 Abort("AliMagF::LoadParameterization", TSelector::kAbortProcess);
1573 TGeoGlobalMagField::Instance()->SetField(newMap);
1574 TGeoGlobalMagField::Instance()->Lock();
1576 if (TNamed *outputFileName = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE"))
1577 fProofOutputFileName = outputFileName->GetTitle();
1578 if (TNamed *outputLocation = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE_LOCATION"))
1579 fProofOutputLocation = outputLocation->GetTitle();
1580 if (fInput->FindObject("PROOF_OUTPUTFILE_DATASET"))
1581 fProofOutputDataset = kTRUE;
1582 if (TNamed *archiveList = (TNamed*)fInput->FindObject("PROOF_OUTPUTFILE_ARCHIVE"))
1583 fProofOutputArchive = archiveList->GetTitle();
1584 if (!fProofOutputFileName.IsNull() &&
1585 !fProofOutputLocation.IsNull() &&
1586 fProofOutputArchive.IsNull()) {
1587 if (!fProofOutputDataset) {
1588 outProofFile = new TProofOutputFile(fProofOutputFileName.Data(),"M");
1589 outProofFile->SetOutputFileName(Form("%s%s",fProofOutputLocation.Data(),fProofOutputFileName.Data()));
1592 outProofFile = new TProofOutputFile(fProofOutputFileName.Data(),"DROV",fProofOutputLocation.Data());
1594 if (AliDebugLevel() > 0) outProofFile->Dump();
1595 fOutput->Add(outProofFile);
1597 AliSysInfo::AddStamp("ReadInputInSlaveBegin");
1600 // get the run loader
1601 if (!InitRunLoader()) {
1602 Abort("InitRunLoader", TSelector::kAbortProcess);
1605 AliSysInfo::AddStamp("LoadLoader");
1607 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
1610 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
1611 Abort("CreateTrackers", TSelector::kAbortProcess);
1614 AliSysInfo::AddStamp("CreateTrackers");
1616 // create the ESD output file and tree
1617 if (!outProofFile) {
1618 ffile = TFile::Open("AliESDs.root", "RECREATE");
1619 ffile->SetCompressionLevel(2);
1620 if (!ffile->IsOpen()) {
1621 Abort("OpenESDFile", TSelector::kAbortProcess);
1626 AliInfo(Form("Opening output PROOF file: %s/%s",
1627 outProofFile->GetDir(), outProofFile->GetFileName()));
1628 if (!(ffile = outProofFile->OpenFile("RECREATE"))) {
1629 Abort(Form("Problems opening output PROOF file: %s/%s",
1630 outProofFile->GetDir(), outProofFile->GetFileName()),
1631 TSelector::kAbortProcess);
1636 ftree = new TTree("esdTree", "Tree with ESD objects");
1637 fesd = new AliESDEvent();
1638 fesd->CreateStdContent();
1639 // add a so far non-std object to the ESD, this will
1640 // become part of the std content
1641 fesd->AddObject(new AliESDHLTDecision);
1643 fesd->WriteToTree(ftree);
1644 if (fWriteESDfriend) {
1645 ffileF = TFile::Open("AliESDfriends.root", "RECREATE");
1646 ftreeF = new TTree("esdFriendTree", "Tree with ESD Friend objects");
1647 fesdf = new AliESDfriend();
1648 ftreeF->Branch("ESDfriend.","AliESDfriend", &fesdf);
1649 fesd->AddObject(fesdf);
1652 ftree->GetUserInfo()->Add(fesd);
1654 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
1655 fhltesd = new AliESDEvent();
1656 fhltesd->CreateStdContent();
1657 // read the ESD template from CDB
1658 // HLT is allowed to put non-std content to its ESD, the non-std
1659 // objects need to be created before invocation of WriteToTree in
1660 // order to create all branches. Initialization is done from an
1661 // ESD layout template in CDB
1662 AliCDBManager* man = AliCDBManager::Instance();
1663 AliCDBPath hltESDConfigPath("HLT/ConfigHLT/esdLayout");
1664 AliCDBEntry* hltESDConfig=NULL;
1665 if (man->GetId(hltESDConfigPath)!=NULL &&
1666 (hltESDConfig=man->Get(hltESDConfigPath))!=NULL) {
1667 AliESDEvent* pESDLayout=dynamic_cast<AliESDEvent*>(hltESDConfig->GetObject());
1669 // init all internal variables from the list of objects
1670 pESDLayout->GetStdContent();
1672 // copy content and create non-std objects
1673 *fhltesd=*pESDLayout;
1676 AliError(Form("error setting hltEsd layout from %s: invalid object type",
1677 hltESDConfigPath.GetPath().Data()));
1681 fhltesd->WriteToTree(fhlttree);
1682 fhlttree->GetUserInfo()->Add(fhltesd);
1684 ProcInfo_t procInfo;
1685 gSystem->GetProcInfo(&procInfo);
1686 AliInfo(Form("Current memory usage %ld %ld", procInfo.fMemResident, procInfo.fMemVirtual));
1689 //Initialize the QA and start of cycle
1690 if (fRunQA || fRunGlobalQA)
1693 //Initialize the Plane Efficiency framework
1694 if (fRunPlaneEff && !InitPlaneEff()) {
1695 Abort("InitPlaneEff", TSelector::kAbortProcess);
1699 if (strcmp(gProgName,"alieve") == 0)
1700 fRunAliEVE = InitAliEVE();
1705 //_____________________________________________________________________________
1706 Bool_t AliReconstruction::Process(Long64_t entry)
1708 // run the reconstruction over a single entry
1709 // from the chain with raw data
1710 AliCodeTimerAuto("",0);
1712 TTree *currTree = fChain->GetTree();
1713 AliRawVEvent *event = NULL;
1714 currTree->SetBranchAddress("rawevent",&event);
1715 currTree->GetEntry(entry);
1716 fRawReader = new AliRawReaderRoot(event);
1717 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
1725 //_____________________________________________________________________________
1726 void AliReconstruction::Init(TTree *tree)
1728 // Implementation of TSelector::Init()
1731 AliError("The input tree is not found!");
1737 //_____________________________________________________________________________
1738 Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
1740 // run the reconstruction over a single event
1741 // The event loop is steered in Run method
1744 static Long_t oldMres=0;
1745 static Long_t oldMvir=0;
1746 static Float_t oldCPU=0;
1747 static Long_t aveDMres=0;
1748 static Long_t aveDMvir=0;
1749 static Float_t aveDCPU=0;
1751 AliCodeTimerAuto("",0);
1755 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
1756 fRunLoader->SetEventNumber(iEvent);
1758 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1760 fRunLoader->TreeE()->Fill();
1762 if (fRawReader && fRawReader->UseAutoSaveESD())
1763 fRunLoader->TreeE()->AutoSave("SaveSelf");
1766 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
1771 fRunLoader->GetEvent(iEvent);
1773 // Fill Event-info object
1775 fRecoParam.SetEventSpecie(fRunInfo,fEventInfo,fListOfCosmicTriggers);
1777 ProcInfo_t procInfo;
1778 if(iEvent==fFirstEvent) {
1779 gSystem->GetProcInfo(&procInfo);
1780 oldMres=procInfo.fMemResident;
1781 oldMvir=procInfo.fMemVirtual;
1782 oldCPU=procInfo.fCpuUser+procInfo.fCpuSys;
1784 AliInfo(Form("================================= Processing event %d of type %-10s ==================================", iEvent,fRecoParam.PrintEventSpecie()));
1786 // Set the reco-params
1788 TString detStr = fLoadCDB;
1789 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
1790 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1791 AliReconstructor *reconstructor = GetReconstructor(iDet);
1792 if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
1793 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
1794 reconstructor->SetRecoParam(par);
1795 reconstructor->GetPidSettings(&pid);
1796 reconstructor->SetEventInfo(&fEventInfo);
1798 AliQAManager::QAManager()->SetRecoParam(iDet, par) ;
1799 if (par) AliQAManager::QAManager()->SetEventSpecie(AliRecoParam::Convert(par->GetEventSpecie())) ;
1804 if (fRunQA || fRunGlobalQA) AliQADataMaker::SetEventTrigClasses(fEventInfo.GetTriggerClasses()); // RS: select which histo clones are to be filled
1807 const AliDetectorRecoParam *grppar = fRecoParam.GetDetRecoParam(kNDetectors);
1808 AliQAManager::QAManager()->SetRecoParam(AliQAv1::kGLOBAL, grppar) ;
1809 AliQAManager::QAManager()->SetEventSpecie(AliRecoParam::Convert(grppar->GetEventSpecie())) ;
1814 if (fRunQA && IsInTasks(AliQAv1::kRAWS)) {
1815 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
1816 AliQAManager::QAManager()->RunOneEvent(fRawReader) ;
1818 // local single event reconstruction
1819 if (!fRunLocalReconstruction.IsNull()) {
1820 TString detectors=fRunLocalReconstruction;
1821 // run HLT event reconstruction first
1822 // ;-( IsSelected changes the string
1823 if (IsSelected("HLT", detectors) &&
1824 !RunLocalEventReconstruction("HLT")) {
1825 if (fStopOnError) {CleanUp(); return kFALSE;}
1827 detectors=fRunLocalReconstruction;
1828 detectors.ReplaceAll("HLT", "");
1829 if (!RunLocalEventReconstruction(detectors)) {
1838 // fill Event header information from the RawEventHeader
1839 if (fRawReader){FillRawEventHeaderESD(fesd);}
1840 if (fRawReader){FillRawEventHeaderESD(fhltesd);}
1842 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1843 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
1845 ((AliESDRun*)fesd->GetESDRun())->SetDetectorsInDAQ(fRunInfo->GetDetectorMask());
1846 ((AliESDRun*)fhltesd->GetESDRun())->SetDetectorsInDAQ(fRunInfo->GetDetectorMask());
1847 ((AliESDRun*)fesd->GetESDRun())->SetDetectorsInReco(AliDAQ::DetectorPatternOffline(fFillESD.Data()));
1848 ((AliESDRun*)fhltesd->GetESDRun())->SetDetectorsInReco(AliDAQ::DetectorPatternOffline(fFillESD.Data()));
1850 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1851 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
1853 fesd->SetEventSpecie(fRecoParam.GetEventSpecie());
1854 fhltesd->SetEventSpecie(fRecoParam.GetEventSpecie());
1856 // Set magnetic field from the tracker
1857 fesd->SetMagneticField(AliTracker::GetBz());
1858 fhltesd->SetMagneticField(AliTracker::GetBz());
1860 AliESDRun *esdRun,*esdRunH;
1861 esdRun = (AliESDRun*)fesd->GetESDRun();
1862 esdRunH = (AliESDRun*)fhltesd->GetESDRun();
1863 esdRun->SetBeamEnergyIsSqrtSHalfGeV();
1864 esdRunH->SetBeamEnergyIsSqrtSHalfGeV();
1866 for (int ib=2;ib--;) for (int it=2;it--;) {
1867 esdRun->SetMeanIntensity(ib,it, fBeamInt[ib][it]);
1868 esdRunH->SetMeanIntensity(ib,it, fBeamInt[ib][it]);
1871 AliMagF* fld = (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
1872 if (fld) { // set info needed for field initialization
1873 fesd->SetCurrentL3(fld->GetCurrentSol());
1874 fesd->SetCurrentDip(fld->GetCurrentDip());
1875 fesd->SetBeamEnergy(fld->GetBeamEnergy());
1876 fesd->SetBeamType(fld->GetBeamTypeText());
1877 fesd->SetUniformBMap(fld->IsUniform());
1878 fesd->SetBInfoStored();
1880 fhltesd->SetCurrentL3(fld->GetCurrentSol());
1881 fhltesd->SetCurrentDip(fld->GetCurrentDip());
1882 fhltesd->SetBeamEnergy(fld->GetBeamEnergy());
1883 fhltesd->SetBeamType(fld->GetBeamTypeText());
1884 fhltesd->SetUniformBMap(fld->IsUniform());
1885 fhltesd->SetBInfoStored();
1888 // Set most probable pt, for B=0 tracking
1889 // Get the global reco-params. They are atposition 16 inside the array of detectors in fRecoParam
1890 const AliGRPRecoParam *grpRecoParam = dynamic_cast<const AliGRPRecoParam*>(fRecoParam.GetDetRecoParam(kNDetectors));
1891 if (grpRecoParam) AliExternalTrackParam::SetMostProbablePt(grpRecoParam->GetMostProbablePt());
1893 // Fill raw-data error log into the ESD
1894 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
1897 if (fRunVertexFinder) {
1898 if (!RunVertexFinder(fesd)) {
1899 if (fStopOnError) {CleanUp(); return kFALSE;}
1903 // For Plane Efficiency: run the SPD trackleter
1904 if (fRunPlaneEff && fSPDTrackleter) {
1905 if (!RunSPDTrackleting(fesd)) {
1906 if (fStopOnError) {CleanUp(); return kFALSE;}
1911 if (!fRunTracking.IsNull()) {
1912 if (fRunMuonTracking) {
1913 if (!RunMuonTracking(fesd)) {
1914 if (fStopOnError) {CleanUp(); return kFALSE;}
1920 if (!fRunTracking.IsNull()) {
1921 if (!RunTracking(fesd,pid)) {
1922 if (fStopOnError) {CleanUp(); return kFALSE;}
1927 if (!fFillESD.IsNull()) {
1928 TString detectors=fFillESD;
1929 // run HLT first and on hltesd
1930 // ;-( IsSelected changes the string
1931 if (IsSelected("HLT", detectors) &&
1932 !FillESD(fhltesd, "HLT")) {
1933 if (fStopOnError) {CleanUp(); return kFALSE;}
1936 // Temporary fix to avoid problems with HLT that overwrites the offline ESDs
1937 if (detectors.Contains("ALL")) {
1939 for (Int_t idet=0; idet<kNDetectors; ++idet){
1940 detectors += fgkDetectorName[idet];
1944 detectors.ReplaceAll("HLT", "");
1945 if (!FillESD(fesd, detectors)) {
1946 if (fStopOnError) {CleanUp(); return kFALSE;}
1951 if (fReconstructor[3])
1952 GetReconstructor(3)->FillEventTimeWithTOF(fesd,&pid);
1957 if (fFillTriggerESD) {
1958 if (!FillTriggerESD(fesd)) {
1959 if (fStopOnError) {CleanUp(); return kFALSE;}
1962 // Always fill scalers
1963 if (!FillTriggerScalers(fesd)) {
1964 if (fStopOnError) {CleanUp(); return kFALSE;}
1971 // Propagate track to the beam pipe (if not already done by ITS)
1973 const Int_t ntracks = fesd->GetNumberOfTracks();
1974 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1977 UShort_t *selectedIdx=new UShort_t[ntracks];
1979 for (Int_t itrack=0; itrack<ntracks; itrack++){
1980 const Double_t kMaxStep = 1; //max step over the material
1983 AliESDtrack *track = fesd->GetTrack(itrack);
1984 if (!track) continue;
1986 AliExternalTrackParam *tpcTrack =
1987 (AliExternalTrackParam *)track->GetTPCInnerParam();
1991 PropagateTrackToBxByBz(tpcTrack,kRadius,track->GetMass(),kMaxStep,kFALSE);
1994 Int_t n=trkArray.GetEntriesFast();
1995 selectedIdx[n]=track->GetID();
1996 trkArray.AddLast(tpcTrack);
1999 //Tracks refitted by ITS should already be at the SPD vertex
2000 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
2003 PropagateTrackToBxByBz(track,kRadius,track->GetMass(),kMaxStep,kFALSE);
2004 Double_t x[3]; track->GetXYZ(x);
2005 Double_t b[3]; AliTracker::GetBxByBz(x,b);
2006 track->RelateToVertexBxByBz(fesd->GetPrimaryVertexSPD(), b, kVeryBig);
2011 // Improve the reconstructed primary vertex position using the tracks
2013 Bool_t runVertexFinderTracks = fRunVertexFinderTracks;
2014 if(fesd->GetPrimaryVertexSPD()) {
2015 TString vtitle = fesd->GetPrimaryVertexSPD()->GetTitle();
2016 if(vtitle.Contains("cosmics")) {
2017 runVertexFinderTracks=kFALSE;
2021 if (runVertexFinderTracks) {
2022 // TPC + ITS primary vertex
2023 ftVertexer->SetITSMode();
2024 ftVertexer->SetConstraintOff();
2025 // get cuts for vertexer from AliGRPRecoParam
2026 Bool_t constrSPD=kFALSE;
2028 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
2029 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
2030 grpRecoParam->GetVertexerTracksCutsITS(cutsVertexer);
2031 ftVertexer->SetCuts(cutsVertexer);
2032 delete [] cutsVertexer; cutsVertexer = NULL;
2033 if(grpRecoParam->GetVertexerTracksConstraintITS()) {
2034 if(fDiamondProfile && fDiamondProfile->GetXRes()<kRadius){
2035 ftVertexer->SetVtxStart(fDiamondProfile); // apply constraint only if sigmax is smaller than the beam pipe radius
2037 if(fDiamondProfileSPD && fDiamondProfileSPD->GetXRes()<kRadius){
2038 ftVertexer->SetVtxStart(fDiamondProfileSPD);
2044 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
2047 TString title=pvtx->GetTitle();
2048 title.Append("SPD");
2049 pvtx->SetTitle(title);
2051 if (pvtx->GetStatus()) {
2052 fesd->SetPrimaryVertexTracks(pvtx);
2053 for (Int_t i=0; i<ntracks; i++) {
2054 AliESDtrack *t = fesd->GetTrack(i);
2055 Double_t x[3]; t->GetXYZ(x);
2056 Double_t b[3]; AliTracker::GetBxByBz(x,b);
2057 t->RelateToVertexBxByBz(pvtx, b, kVeryBig);
2060 delete pvtx; pvtx=NULL;
2063 // TPC-only primary vertex
2064 ftVertexer->SetTPCMode();
2065 ftVertexer->SetConstraintOff();
2066 // get cuts for vertexer from AliGRPRecoParam
2068 Int_t nCutsVertexer = grpRecoParam->GetVertexerTracksNCuts();
2069 Double_t *cutsVertexer = new Double_t[nCutsVertexer];
2070 grpRecoParam->GetVertexerTracksCutsTPC(cutsVertexer);
2071 ftVertexer->SetCuts(cutsVertexer);
2072 delete [] cutsVertexer; cutsVertexer = NULL;
2073 if(fDiamondProfileTPC && grpRecoParam->GetVertexerTracksConstraintTPC()) {
2074 if(fDiamondProfileTPC->GetXRes()<kRadius) ftVertexer->SetVtxStart(fDiamondProfileTPC); // apply constraint only if sigmax is smaller than the beam pipe radius
2077 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
2079 if (pvtx->GetStatus()) {
2080 fesd->SetPrimaryVertexTPC(pvtx);
2081 for (Int_t i=0; i<ntracks; i++) {
2082 AliESDtrack *t = fesd->GetTrack(i);
2083 Double_t x[3]; t->GetXYZ(x);
2084 Double_t b[3]; AliTracker::GetBxByBz(x,b);
2085 t->RelateToVertexTPCBxByBz(pvtx, b, kVeryBig);
2088 delete pvtx; pvtx=NULL;
2092 delete[] selectedIdx;
2094 if(fDiamondProfile && fDiamondProfile->GetXRes()<kRadius) fesd->SetDiamond(fDiamondProfile);
2095 else fesd->SetDiamond(fDiamondProfileSPD);
2099 AliV0vertexer vtxer;
2100 // get cuts for V0vertexer from AliGRPRecoParam
2102 Int_t nCutsV0vertexer = grpRecoParam->GetVertexerV0NCuts();
2103 Double_t *cutsV0vertexer = new Double_t[nCutsV0vertexer];
2104 grpRecoParam->GetVertexerV0Cuts(cutsV0vertexer);
2105 vtxer.SetCuts(cutsV0vertexer);
2106 delete [] cutsV0vertexer; cutsV0vertexer = NULL;
2108 vtxer.Tracks2V0vertices(fesd);
2110 if (fRunCascadeFinder) {
2112 AliCascadeVertexer cvtxer;
2113 // get cuts for CascadeVertexer from AliGRPRecoParam
2115 Int_t nCutsCascadeVertexer = grpRecoParam->GetVertexerCascadeNCuts();
2116 Double_t *cutsCascadeVertexer = new Double_t[nCutsCascadeVertexer];
2117 grpRecoParam->GetVertexerCascadeCuts(cutsCascadeVertexer);
2118 cvtxer.SetCuts(cutsCascadeVertexer);
2119 delete [] cutsCascadeVertexer; cutsCascadeVertexer = NULL;
2121 cvtxer.V0sTracks2CascadeVertices(fesd);
2126 if (fCleanESD) CleanESD(fesd);
2128 // RS run updated trackleter: since we want to mark the clusters used by tracks and also mark the
2129 // tracks interpreted as primary, this step should be done in the very end, when full
2130 // ESD info is available (particulalry, V0s)
2132 if (fRunMultFinder) {
2133 if (!RunMultFinder(fesd)) {
2134 if (fStopOnError) {CleanUp(); return kFALSE;}
2138 if (fRunQA && IsInTasks(AliQAv1::kESDS)) {
2139 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2140 AliQAManager::QAManager()->RunOneEvent(fesd, fhltesd) ;
2143 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2145 qadm->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2146 if (qadm && IsInTasks(AliQAv1::kESDS))
2147 qadm->Exec(AliQAv1::kESDS, fesd);
2150 // copy HLT decision from HLTesd to esd
2151 // the most relevant information is stored in a reduced container in the esd,
2152 // while the full information can be found in the HLTesd
2153 TObject* pHLTSrc=fhltesd->FindListObject(AliESDHLTDecision::Name());
2154 TObject* pHLTTgt=fesd->FindListObject(AliESDHLTDecision::Name());
2155 if (pHLTSrc && pHLTTgt) {
2156 pHLTSrc->Copy(*pHLTTgt);
2159 if (fWriteESDfriend)
2160 fesd->GetESDfriend(fesdf);
2163 if (fWriteESDfriend) {
2167 // Auto-save the ESD tree in case of prompt reco @P2
2168 if (fRawReader && fRawReader->UseAutoSaveESD()) {
2169 ftree->AutoSave("SaveSelf");
2170 if (fWriteESDfriend) ftreeF->AutoSave("SaveSelf");
2177 if (fRunAliEVE) RunAliEVE();
2181 if (fWriteESDfriend) {
2182 fesdf->~AliESDfriend();
2183 new (fesdf) AliESDfriend(); // Reset...
2186 gSystem->GetProcInfo(&procInfo);
2187 Long_t dMres=(procInfo.fMemResident-oldMres)/1024;
2188 Long_t dMvir=(procInfo.fMemVirtual-oldMvir)/1024;
2189 Float_t dCPU=procInfo.fCpuUser+procInfo.fCpuSys-oldCPU;
2190 aveDMres+=(dMres-aveDMres)/(iEvent-fFirstEvent+1);
2191 aveDMvir+=(dMvir-aveDMvir)/(iEvent-fFirstEvent+1);
2192 aveDCPU+=(dCPU-aveDCPU)/(iEvent-fFirstEvent+1);
2193 AliInfo(Form("======================= End Event %d: Res %ld(%3ld <%3ld>) Vir %ld(%3ld <%3ld>) CPU %5.2f <%5.2f> ===================",
2194 iEvent, procInfo.fMemResident/1024, dMres, aveDMres, procInfo.fMemVirtual/1024, dMvir, aveDMvir, dCPU, aveDCPU));
2195 oldMres=procInfo.fMemResident;
2196 oldMvir=procInfo.fMemVirtual;
2197 oldCPU=procInfo.fCpuUser+procInfo.fCpuSys;
2200 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2201 if (fReconstructor[iDet]) {
2202 fReconstructor[iDet]->SetRecoParam(NULL);
2203 fReconstructor[iDet]->SetEventInfo(NULL);
2205 if (fTracker[iDet]) fTracker[iDet]->SetEventInfo(NULL);
2208 if (fRunQA || fRunGlobalQA)
2209 AliQAManager::QAManager()->Increment() ;
2214 //_____________________________________________________________________________
2215 void AliReconstruction::SlaveTerminate()
2217 // Finalize the run on the slave side
2218 // Called after the exit
2219 // from the event loop
2220 AliCodeTimerAuto("",0);
2222 if (fIsNewRunLoader) { // galice.root didn't exist
2223 fRunLoader->WriteHeader("OVERWRITE");
2224 fRunLoader->WriteTrigger("OVERWRITE");
2225 fRunLoader->CdGAFile();
2226 fRunLoader->Write(0, TObject::kOverwrite);
2229 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
2230 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
2232 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
2233 cdbMapCopy->SetOwner(1);
2234 cdbMapCopy->SetName("cdbMap");
2235 TIter iter(cdbMap->GetTable());
2238 while((pair = dynamic_cast<TPair*> (iter.Next()))){
2239 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
2240 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
2241 if (keyStr && valStr)
2242 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
2245 TList *cdbListCopy = new TList();
2246 cdbListCopy->SetOwner(1);
2247 cdbListCopy->SetName("cdbList");
2249 TIter iter2(cdbList);
2252 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
2253 cdbListCopy->Add(new TObjString(id->ToString().Data()));
2256 ftree->GetUserInfo()->Add(cdbMapCopy);
2257 ftree->GetUserInfo()->Add(cdbListCopy);
2262 // we want to have only one tree version number
2263 ftree->Write(ftree->GetName(),TObject::kOverwrite);
2264 fhlttree->Write(fhlttree->GetName(),TObject::kOverwrite);
2266 if (fWriteESDfriend) {
2268 ftreeF->Write(ftreeF->GetName(),TObject::kOverwrite);
2271 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
2272 if (fRunPlaneEff && !FinishPlaneEff()) {
2273 AliWarning("Finish PlaneEff evaluation failed");
2276 // End of cycle for the in-loop
2278 if (fRunQA || fRunGlobalQA) {
2279 AliQAManager::QAManager()->EndOfCycle() ;
2281 !fProofOutputLocation.IsNull() &&
2282 fProofOutputArchive.IsNull() &&
2283 !fProofOutputDataset) {
2284 TString qaOutputFile(Form("%sMerged.%s.Data.root",
2285 fProofOutputLocation.Data(),
2286 AliQAv1::GetQADataFileName()));
2287 TProofOutputFile *qaProofFile = new TProofOutputFile(Form("Merged.%s.Data.root",
2288 AliQAv1::GetQADataFileName()));
2289 qaProofFile->SetOutputFileName(qaOutputFile.Data());
2290 if (AliDebugLevel() > 0) qaProofFile->Dump();
2291 fOutput->Add(qaProofFile);
2292 MergeQA(qaProofFile->GetFileName());
2303 if (!fProofOutputFileName.IsNull() &&
2304 !fProofOutputLocation.IsNull() &&
2305 fProofOutputDataset &&
2306 !fProofOutputArchive.IsNull()) {
2307 TProofOutputFile *zipProofFile = new TProofOutputFile(fProofOutputFileName.Data(),
2309 fProofOutputLocation.Data());
2310 if (AliDebugLevel() > 0) zipProofFile->Dump();
2311 fOutput->Add(zipProofFile);
2312 TString fileList(fProofOutputArchive.Data());
2313 fileList.ReplaceAll(","," ");
2315 #if ROOT_SVN_REVISION >= 30174
2316 command.Form("zip -n root %s/%s %s",zipProofFile->GetDir(kTRUE),zipProofFile->GetFileName(),fileList.Data());
2318 command.Form("zip -n root %s/%s %s",zipProofFile->GetDir(),zipProofFile->GetFileName(),fileList.Data());
2320 AliInfo(Form("Executing: %s",command.Data()));
2321 gSystem->Exec(command.Data());
2326 //_____________________________________________________________________________
2327 void AliReconstruction::Terminate()
2329 // Create tags for the events in the ESD tree (the ESD tree is always present)
2330 // In case of empty events the tags will contain dummy values
2331 AliCodeTimerAuto("",0);
2333 // Do not call the ESD tag creator in case of PROOF-based reconstruction
2335 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
2336 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData, AliQAv1::Instance()->GetQA(), AliQAv1::Instance()->GetEventSpecies(), AliQAv1::kNDET, AliRecoParam::kNSpecies);
2337 delete esdtagCreator;
2340 // Cleanup of CDB manager: cache and active storages!
2341 AliCDBManager::Instance()->ClearCache();
2344 //_____________________________________________________________________________
2345 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
2347 // run the local reconstruction
2349 static Int_t eventNr=0;
2350 AliCodeTimerAuto("",0)
2352 TString detStr = detectors;
2353 // execute HLT reconstruction first since other detector reconstruction
2354 // might depend on HLT data
2355 // key 'HLT' is removed from detStr by IsSelected
2356 if (!IsSelected("HLT", detStr)) {
2357 AliReconstructor* reconstructor = GetReconstructor(kNDetectors-1);
2358 if (reconstructor) {
2359 // there is no AliLoader for HLT, see
2360 // https://savannah.cern.ch/bugs/?35473
2361 AliInfo("running reconstruction for HLT");
2363 reconstructor->Reconstruct(fRawReader, NULL);
2366 reconstructor->Reconstruct(dummy, NULL);
2370 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2371 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2372 AliReconstructor* reconstructor = GetReconstructor(iDet);
2373 if (!reconstructor) continue;
2374 AliLoader* loader = fLoader[iDet];
2376 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
2379 // conversion of digits
2380 if (fRawReader && reconstructor->HasDigitConversion()) {
2381 AliInfo(Form("converting raw data digits into root objects for %s",
2382 fgkDetectorName[iDet]));
2383 // AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
2384 // fgkDetectorName[iDet]),0);
2385 loader->LoadDigits("update");
2386 loader->CleanDigits();
2387 loader->MakeDigitsContainer();
2388 TTree* digitsTree = loader->TreeD();
2389 reconstructor->ConvertDigits(fRawReader, digitsTree);
2390 loader->WriteDigits("OVERWRITE");
2391 loader->UnloadDigits();
2393 // local reconstruction
2394 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
2395 //AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]),0);
2396 loader->LoadRecPoints("update");
2397 loader->CleanRecPoints();
2398 loader->MakeRecPointsContainer();
2399 TTree* clustersTree = loader->TreeR();
2400 if (fRawReader && !reconstructor->HasDigitConversion()) {
2401 reconstructor->Reconstruct(fRawReader, clustersTree);
2403 loader->LoadDigits("read");
2404 TTree* digitsTree = loader->TreeD();
2406 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
2410 reconstructor->Reconstruct(digitsTree, clustersTree);
2411 if (fRunQA && IsInTasks(AliQAv1::kDIGITSR)) {
2412 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2413 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, digitsTree) ;
2416 loader->UnloadDigits();
2418 if (fRunQA && IsInTasks(AliQAv1::kRECPOINTS)) {
2419 AliQAManager::QAManager()->SetEventSpecie(fRecoParam.GetEventSpecie()) ;
2420 AliQAManager::QAManager()->RunOneEventInOneDetector(iDet, clustersTree) ;
2422 loader->WriteRecPoints("OVERWRITE");
2423 loader->UnloadRecPoints();
2424 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
2426 if (!IsSelected("CTP", detStr)) AliDebug(10,"No CTP");
2427 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2428 AliError(Form("the following detectors were not found: %s",
2436 //_____________________________________________________________________________
2437 Bool_t AliReconstruction::RunSPDTrackleting(AliESDEvent*& esd)
2439 // run the SPD trackleting (for SPD efficiency purpouses)
2441 AliCodeTimerAuto("",0)
2443 Double_t vtxPos[3] = {0, 0, 0};
2444 Double_t vtxErr[3] = {0.0, 0.0, 0.0};
2446 TArrayF mcVertex(3);
2448 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
2449 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
2450 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
2453 const AliESDVertex *vertex = esd->GetVertex();
2455 AliWarning("Vertex not found");
2458 vertex->GetXYZ(vtxPos);
2459 vertex->GetSigmaXYZ(vtxErr);
2460 if (fSPDTrackleter) {
2461 AliInfo("running the SPD Trackleter for Plane Efficiency Evaluation");
2464 fLoader[0]->LoadRecPoints("read");
2465 TTree* tree = fLoader[0]->TreeR();
2467 AliError("Can't get the ITS cluster tree");
2470 fSPDTrackleter->LoadClusters(tree);
2471 fSPDTrackleter->SetVertex(vtxPos, vtxErr);
2473 if (fSPDTrackleter->Clusters2Tracks(esd) != 0) {
2474 AliWarning("AliITSTrackleterSPDEff Clusters2Tracks failed");
2475 // fLoader[0]->UnloadRecPoints();
2478 //fSPDTrackleter->UnloadRecPoints();
2480 AliWarning("SPDTrackleter not available");
2486 //_____________________________________________________________________________
2487 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
2489 // run the barrel tracking
2491 AliCodeTimerAuto("",0)
2493 AliVertexer *vertexer = CreateVertexer();
2494 if (!vertexer) return kFALSE;
2496 AliInfo(Form("running the ITS vertex finder: %s",vertexer->ClassName()));
2497 AliESDVertex* vertex = NULL;
2499 fLoader[0]->LoadRecPoints();
2500 TTree* cltree = fLoader[0]->TreeR();
2502 if(fDiamondProfileSPD) vertexer->SetVtxStart(fDiamondProfileSPD);
2503 vertex = vertexer->FindVertexForCurrentEvent(cltree);
2506 AliError("Can't get the ITS cluster tree");
2508 fLoader[0]->UnloadRecPoints();
2511 AliError("Can't get the ITS loader");
2514 AliWarning("Vertex not found");
2515 vertex = new AliESDVertex();
2516 vertex->SetName("default");
2519 vertex->SetName("reconstructed");
2524 vertex->GetXYZ(vtxPos);
2525 vertex->GetSigmaXYZ(vtxErr);
2527 esd->SetPrimaryVertexSPD(vertex);
2528 AliESDVertex *vpileup = NULL;
2529 Int_t novertices = 0;
2530 vpileup = vertexer->GetAllVertices(novertices);
2532 for (Int_t kk=1; kk<novertices; kk++)esd->AddPileupVertexSPD(&vpileup[kk]);
2535 // if SPD multiplicity has been determined, it is stored in the ESD
2536 AliMultiplicity *mult = vertexer->GetMultiplicity();
2537 if(mult)esd->SetMultiplicity(mult);
2539 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2540 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
2549 //_____________________________________________________________________________
2550 Bool_t AliReconstruction::RunMultFinder(AliESDEvent*& esd)
2552 // run the trackleter for multiplicity study
2554 AliCodeTimerAuto("",0)
2556 AliTrackleter *trackleter = CreateMultFinder();
2557 if (!trackleter) return kFALSE;
2559 AliInfo(Form("running the ITS multiplicity finder: %s",trackleter->ClassName()));
2562 fLoader[0]->LoadRecPoints();
2563 TTree* cltree = fLoader[0]->TreeR();
2565 trackleter->Reconstruct(esd,cltree);
2566 AliMultiplicity *mult = trackleter->GetMultiplicity();
2567 if(mult) esd->SetMultiplicity(mult);
2570 AliError("Can't get the ITS cluster tree");
2572 fLoader[0]->UnloadRecPoints();
2575 AliError("Can't get the ITS loader");
2583 //_____________________________________________________________________________
2584 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
2586 // run the HLT barrel tracking
2588 AliCodeTimerAuto("",0)
2591 AliError("Missing runLoader!");
2595 AliInfo("running HLT tracking");
2597 // Get a pointer to the HLT reconstructor
2598 AliReconstructor *reconstructor = GetReconstructor(kNDetectors-1);
2599 if (!reconstructor) return kFALSE;
2602 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2603 TString detName = fgkDetectorName[iDet];
2604 AliDebug(1, Form("%s HLT tracking", detName.Data()));
2605 reconstructor->SetOption(detName.Data());
2606 AliTracker *tracker = reconstructor->CreateTracker();
2608 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
2609 if (fStopOnError) return kFALSE;
2613 Double_t vtxErr[3]={0.005,0.005,0.010};
2614 const AliESDVertex *vertex = esd->GetVertex();
2615 vertex->GetXYZ(vtxPos);
2616 tracker->SetVertex(vtxPos,vtxErr);
2618 fLoader[iDet]->LoadRecPoints("read");
2619 TTree* tree = fLoader[iDet]->TreeR();
2621 AliError(Form("Can't get the %s cluster tree", detName.Data()));
2624 tracker->LoadClusters(tree);
2626 if (tracker->Clusters2Tracks(esd) != 0) {
2627 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
2631 tracker->UnloadClusters();
2639 //_____________________________________________________________________________
2640 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
2642 // run the muon spectrometer tracking
2644 AliCodeTimerAuto("",0)
2647 AliError("Missing runLoader!");
2650 Int_t iDet = 7; // for MUON
2652 AliInfo("is running...");
2654 // Get a pointer to the MUON reconstructor
2655 AliReconstructor *reconstructor = GetReconstructor(iDet);
2656 if (!reconstructor) return kFALSE;
2659 TString detName = fgkDetectorName[iDet];
2660 AliDebug(1, Form("%s tracking", detName.Data()));
2661 AliTracker *tracker = reconstructor->CreateTracker();
2663 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2668 fLoader[iDet]->LoadRecPoints("read");
2670 tracker->LoadClusters(fLoader[iDet]->TreeR());
2672 Int_t rv = tracker->Clusters2Tracks(esd);
2674 fLoader[iDet]->UnloadRecPoints();
2676 tracker->UnloadClusters();
2682 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2690 //_____________________________________________________________________________
2691 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd,AliESDpid &PID)
2693 // run the barrel tracking
2694 static Int_t eventNr=0;
2695 AliCodeTimerAuto("",0)
2697 AliInfo("running tracking");
2699 // Set the event info which is used
2700 // by the trackers in order to obtain
2701 // information about read-out detectors,
2703 AliDebug(1, "Setting event info");
2704 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2705 if (!fTracker[iDet]) continue;
2706 fTracker[iDet]->SetEventInfo(&fEventInfo);
2709 //Fill the ESD with the T0 info (will be used by the TOF)
2710 if (fReconstructor[11] && fLoader[11]) {
2711 fLoader[11]->LoadRecPoints("READ");
2712 TTree *treeR = fLoader[11]->TreeR();
2714 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
2718 // pass 1: TPC + ITS inwards
2719 for (Int_t iDet = 1; iDet >= 0; iDet--) {
2720 if (!fTracker[iDet]) continue;
2721 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
2724 fLoader[iDet]->LoadRecPoints("read");
2725 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
2726 TTree* tree = fLoader[iDet]->TreeR();
2728 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2731 fTracker[iDet]->LoadClusters(tree);
2732 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2734 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
2735 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
2738 // preliminary PID in TPC needed by the ITS tracker
2740 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2741 PID.MakePID(esd,kTRUE);
2743 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
2746 // pass 2: ALL backwards
2748 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2749 if (!fTracker[iDet]) continue;
2750 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
2753 if (iDet > 1) { // all except ITS, TPC
2755 fLoader[iDet]->LoadRecPoints("read");
2756 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
2757 tree = fLoader[iDet]->TreeR();
2759 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
2762 fTracker[iDet]->LoadClusters(tree);
2763 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
2767 if (iDet>1) // start filling residuals for the "outer" detectors
2769 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2770 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2772 AliRecoParam::EventSpecie_t es=fRecoParam.GetEventSpecie();
2773 TObjArray * elem = arr[AliRecoParam::AConvert(es)];
2774 if ( elem && (! elem->At(0)) ) {
2775 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2776 if (qadm) qadm->InitRecPointsForTracker() ;
2780 if (fTracker[iDet]->PropagateBack(esd) != 0) {
2781 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
2786 if (iDet > 3) { // all except ITS, TPC, TRD and TOF
2787 fTracker[iDet]->UnloadClusters();
2788 fLoader[iDet]->UnloadRecPoints();
2790 // updated PID in TPC needed by the ITS tracker -MI
2792 //GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
2793 //AliESDpid::MakePID(esd);
2794 PID.MakePID(esd,kTRUE);
2796 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2798 //stop filling residuals for the "outer" detectors
2799 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2801 // pass 3: TRD + TPC + ITS refit inwards
2803 for (Int_t iDet = 2; iDet >= 0; iDet--) {
2804 if (!fTracker[iDet]) continue;
2805 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
2808 if (iDet<2) // start filling residuals for TPC and ITS
2810 AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kTRUE);
2811 TObjArray ** arr = AliTracker::GetResidualsArray() ;
2813 AliRecoParam::EventSpecie_t es=fRecoParam.GetEventSpecie();
2814 TObjArray * elem = arr[AliRecoParam::AConvert(es)];
2815 if ( elem && (! elem->At(0)) ) {
2816 AliQADataMaker *qadm = AliQAManager::QAManager()->GetQADataMaker(AliQAv1::kGLOBAL);
2817 if (qadm) qadm->InitRecPointsForTracker() ;
2822 if (fTracker[iDet]->RefitInward(esd) != 0) {
2823 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
2826 // run postprocessing
2827 if (fTracker[iDet]->PostProcess(esd) != 0) {
2828 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
2831 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
2834 // write space-points to the ESD in case alignment data output
2836 if (fWriteAlignmentData)
2837 WriteAlignmentData(esd);
2839 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2840 if (!fTracker[iDet]) continue;
2842 fTracker[iDet]->UnloadClusters();
2843 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
2844 fLoader[iDet]->UnloadRecPoints();
2845 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
2847 // stop filling residuals for TPC and ITS
2848 if (fRunGlobalQA) AliTracker::SetFillResiduals(fRecoParam.GetEventSpecie(), kFALSE);
2854 //_____________________________________________________________________________
2855 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
2857 // Remove the data which are not needed for the physics analysis.
2860 Int_t nTracks=esd->GetNumberOfTracks();
2861 Int_t nV0s=esd->GetNumberOfV0s();
2863 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
2865 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
2866 Bool_t rc=esd->Clean(cleanPars);
2868 nTracks=esd->GetNumberOfTracks();
2869 nV0s=esd->GetNumberOfV0s();
2871 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
2876 //_____________________________________________________________________________
2877 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
2879 // fill the event summary data
2881 AliCodeTimerAuto("",0)
2882 static Int_t eventNr=0;
2883 TString detStr = detectors;
2885 AliSysInfo::AddStamp(Form("FillESDb%d",eventNr), -19,-19, eventNr);
2886 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
2887 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2888 AliReconstructor* reconstructor = GetReconstructor(iDet);
2889 if (!reconstructor) continue;
2890 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
2891 TTree* clustersTree = NULL;
2892 if (fLoader[iDet]) {
2893 fLoader[iDet]->LoadRecPoints("read");
2894 clustersTree = fLoader[iDet]->TreeR();
2895 if (!clustersTree) {
2896 AliError(Form("Can't get the %s clusters tree",
2897 fgkDetectorName[iDet]));
2898 if (fStopOnError) return kFALSE;
2901 if (fRawReader && !reconstructor->HasDigitConversion()) {
2902 reconstructor->FillESD(fRawReader, clustersTree, esd);
2904 TTree* digitsTree = NULL;
2905 if (fLoader[iDet]) {
2906 fLoader[iDet]->LoadDigits("read");
2907 digitsTree = fLoader[iDet]->TreeD();
2909 AliError(Form("Can't get the %s digits tree",
2910 fgkDetectorName[iDet]));
2911 if (fStopOnError) return kFALSE;
2914 reconstructor->FillESD(digitsTree, clustersTree, esd);
2915 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
2917 if (fLoader[iDet]) {
2918 fLoader[iDet]->UnloadRecPoints();
2922 if (!IsSelected("CTP", detStr)) AliDebug(10,"No CTP");
2923 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2924 AliError(Form("the following detectors were not found: %s",
2926 if (fStopOnError) return kFALSE;
2928 AliSysInfo::AddStamp(Form("FillESDe%d",eventNr), -20,-20, eventNr);
2933 //_____________________________________________________________________________
2934 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
2936 // Reads the trigger decision which is
2937 // stored in Trigger.root file and fills
2938 // the corresponding esd entries
2940 AliCodeTimerAuto("",0)
2942 AliInfo("Filling trigger information into the ESD");
2945 AliCTPRawStream input(fRawReader);
2946 if (!input.Next()) {
2947 AliWarning("No valid CTP (trigger) DDL raw data is found ! The trigger info is taken from the event header!");
2950 if (esd->GetTriggerMask() != input.GetClassMask())
2951 AliError(Form("Invalid trigger pattern found in CTP raw-data: %llx %llx",
2952 input.GetClassMask(),esd->GetTriggerMask()));
2953 if (esd->GetOrbitNumber() != input.GetOrbitID())
2954 AliError(Form("Invalid orbit id found in CTP raw-data: %x %x",
2955 input.GetOrbitID(),esd->GetOrbitNumber()));
2956 if (esd->GetBunchCrossNumber() != input.GetBCID())
2957 AliError(Form("Invalid bunch-crossing id found in CTP raw-data: %x %x",
2958 input.GetBCID(),esd->GetBunchCrossNumber()));
2959 AliESDHeader* esdheader = esd->GetHeader();
2960 esdheader->SetL0TriggerInputs(input.GetL0Inputs());
2961 esdheader->SetL1TriggerInputs(input.GetL1Inputs());
2962 esdheader->SetL2TriggerInputs(input.GetL2Inputs());
2964 UInt_t orbit=input.GetOrbitID();
2965 for(Int_t i=0 ; i<input.GetNIRs() ; i++ )
2966 if(TMath::Abs(Int_t(orbit-(input.GetIR(i))->GetOrbit()))<=1){
2967 esdheader->AddTriggerIR(input.GetIR(i));
2969 AliCentralTrigger* rlCTP = fRunLoader->GetTrigger();
2970 rlCTP->SetL0TriggerInputs(input.GetL0Inputs());
2971 rlCTP->SetL1TriggerInputs(input.GetL1Inputs());
2972 rlCTP->SetL2TriggerInputs(input.GetL2Inputs());
2974 if (fIsNewRunLoader) fRunLoader->TreeCT()->Fill();
2978 //_____________________________________________________________________________
2979 Bool_t AliReconstruction::FillTriggerScalers(AliESDEvent*& esd)
2982 //fRunScalers->Print();
2983 if(fRunScalers && fRunScalers->CheckRunScalers()){
2984 AliTimeStamp* timestamp = new AliTimeStamp(esd->GetOrbitNumber(), esd->GetPeriodNumber(), esd->GetBunchCrossNumber());
2985 //AliTimeStamp* timestamp = new AliTimeStamp(10308000, 0, (ULong64_t)486238);
2986 AliESDHeader* esdheader = fesd->GetHeader();
2987 for(Int_t i=0;i<50;i++){
2988 if((1ull<<i) & esd->GetTriggerMask()){
2989 AliTriggerScalersESD* scalesd = fRunScalers->GetScalersForEventClass( timestamp, i+1);
2990 if(scalesd)esdheader->SetTriggerScalersRecord(scalesd);
2996 //_____________________________________________________________________________
2997 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
3000 // Filling information from RawReader Header
3003 if (!fRawReader) return kFALSE;
3005 AliInfo("Filling information from RawReader Header");
3007 esd->SetBunchCrossNumber(fRawReader->GetBCID());
3008 esd->SetOrbitNumber(fRawReader->GetOrbitID());
3009 esd->SetPeriodNumber(fRawReader->GetPeriod());
3011 esd->SetTimeStamp(fRawReader->GetTimestamp());
3012 esd->SetEventType(fRawReader->GetType());
3018 //_____________________________________________________________________________
3019 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
3021 // check whether detName is contained in detectors
3022 // if yes, it is removed from detectors
3024 // check if all detectors are selected
3025 if ((detectors.CompareTo("ALL") == 0) ||
3026 detectors.BeginsWith("ALL ") ||
3027 detectors.EndsWith(" ALL") ||
3028 detectors.Contains(" ALL ")) {
3033 // search for the given detector
3034 Bool_t result = kFALSE;
3035 if ((detectors.CompareTo(detName) == 0) ||
3036 detectors.BeginsWith(detName+" ") ||
3037 detectors.EndsWith(" "+detName) ||
3038 detectors.Contains(" "+detName+" ")) {
3039 detectors.ReplaceAll(detName, "");
3043 // clean up the detectors string
3044 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
3045 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
3046 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
3051 //_____________________________________________________________________________
3052 Bool_t AliReconstruction::InitRunLoader()
3054 // get or create the run loader
3056 if (gAlice) delete gAlice;
3059 TFile *gafile = TFile::Open(fGAliceFileName.Data());
3060 // if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
3061 if (gafile) { // galice.root exists
3065 // load all base libraries to get the loader classes
3066 TString libs = gSystem->GetLibraries();
3067 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3068 TString detName = fgkDetectorName[iDet];
3069 if (detName == "HLT") continue;
3070 if (libs.Contains("lib" + detName + "base.so")) continue;
3071 gSystem->Load("lib" + detName + "base.so");
3073 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
3075 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
3080 fRunLoader->CdGAFile();
3081 fRunLoader->LoadgAlice();
3083 //PH This is a temporary fix to give access to the kinematics
3084 //PH that is needed for the labels of ITS clusters
3085 fRunLoader->LoadHeader();
3086 fRunLoader->LoadKinematics();
3088 } else { // galice.root does not exist
3090 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
3092 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
3093 AliConfig::GetDefaultEventFolderName(),
3096 AliError(Form("could not create run loader in file %s",
3097 fGAliceFileName.Data()));
3101 fIsNewRunLoader = kTRUE;
3102 fRunLoader->MakeTree("E");
3103 fRunLoader->MakeTree("GG");
3105 if (fNumberOfEventsPerFile > 0)
3106 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
3108 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
3114 //_____________________________________________________________________________
3115 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
3117 // get the reconstructor object and the loader for a detector
3119 if (fReconstructor[iDet]) {
3120 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
3121 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
3122 fReconstructor[iDet]->SetRecoParam(par);
3123 fReconstructor[iDet]->SetRunInfo(fRunInfo);
3125 return fReconstructor[iDet];
3128 // load the reconstructor object
3129 TPluginManager* pluginManager = gROOT->GetPluginManager();
3130 TString detName = fgkDetectorName[iDet];
3131 TString recName = "Ali" + detName + "Reconstructor";
3133 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
3135 AliReconstructor* reconstructor = NULL;
3136 // first check if a plugin is defined for the reconstructor
3137 TPluginHandler* pluginHandler =
3138 pluginManager->FindHandler("AliReconstructor", detName);
3139 // if not, add a plugin for it
3140 if (!pluginHandler) {
3141 AliDebug(1, Form("defining plugin for %s", recName.Data()));
3142 TString libs = gSystem->GetLibraries();
3143 if (libs.Contains("lib" + detName + "base.so") ||
3144 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
3145 pluginManager->AddHandler("AliReconstructor", detName,
3146 recName, detName + "rec", recName + "()");
3148 pluginManager->AddHandler("AliReconstructor", detName,
3149 recName, detName, recName + "()");
3151 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
3153 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
3154 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
3157 // check if the upgrade reconstructor should be used instead of the standard one
3158 if(fUpgradeMask[iDet]) {
3159 if(reconstructor) delete reconstructor;
3160 TClass *cl = new TClass(Form("Ali%sUpgradeReconstructor",fgkDetectorName[iDet]));
3161 reconstructor = (AliReconstructor*)(cl->New());
3164 if (reconstructor) {
3165 TObject* obj = fOptions.FindObject(detName.Data());
3166 if (obj) reconstructor->SetOption(obj->GetTitle());
3167 reconstructor->SetRunInfo(fRunInfo);
3168 reconstructor->Init();
3169 fReconstructor[iDet] = reconstructor;
3172 // get or create the loader
3173 if (detName != "HLT") {
3174 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
3175 if (!fLoader[iDet]) {
3176 AliConfig::Instance()
3177 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
3179 // first check if a plugin is defined for the loader
3181 pluginManager->FindHandler("AliLoader", detName);
3182 // if not, add a plugin for it
3183 if (!pluginHandler) {
3184 TString loaderName = "Ali" + detName + "Loader";
3185 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
3186 pluginManager->AddHandler("AliLoader", detName,
3187 loaderName, detName + "base",
3188 loaderName + "(const char*, TFolder*)");
3189 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
3191 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
3193 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
3194 fRunLoader->GetEventFolder());
3196 if (!fLoader[iDet]) { // use default loader
3197 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
3199 if (!fLoader[iDet]) {
3200 AliWarning(Form("couldn't get loader for %s", detName.Data()));
3201 if (fStopOnError) return NULL;
3203 fRunLoader->AddLoader(fLoader[iDet]);
3204 fRunLoader->CdGAFile();
3205 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
3206 fRunLoader->Write(0, TObject::kOverwrite);
3211 if (fRecoParam.GetDetRecoParamArray(iDet) && !AliReconstructor::GetRecoParam(iDet)) {
3212 const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
3213 if (reconstructor) {
3214 reconstructor->SetRecoParam(par);
3215 reconstructor->SetRunInfo(fRunInfo);
3218 return reconstructor;
3221 //_____________________________________________________________________________
3222 AliVertexer* AliReconstruction::CreateVertexer()
3224 // create the vertexer
3225 // Please note that the caller is the owner of the
3228 AliVertexer* vertexer = NULL;
3229 AliReconstructor* itsReconstructor = GetReconstructor(0);
3230 if (itsReconstructor && ((fRunLocalReconstruction.Contains("ITS")) ||
3231 fRunTracking.Contains("ITS") || fFillESD.Contains("ITS") )) {
3232 vertexer = itsReconstructor->CreateVertexer();
3235 AliWarning("couldn't create a vertexer for ITS");
3241 //_____________________________________________________________________________
3242 AliTrackleter* AliReconstruction::CreateMultFinder()
3244 // create the ITS trackleter for mult. estimation
3245 // Please note that the caller is the owner of the
3248 AliTrackleter* trackleter = NULL;
3249 AliReconstructor* itsReconstructor = GetReconstructor(0);
3250 if (itsReconstructor && ((fRunLocalReconstruction.Contains("ITS")) ||
3251 fRunTracking.Contains("ITS") || fFillESD.Contains("ITS") )) {
3252 trackleter = itsReconstructor->CreateMultFinder();
3255 AliWarning("ITS is not in reconstruction, switching off RunMultFinder");
3256 fRunMultFinder = kFALSE;
3262 //_____________________________________________________________________________
3263 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
3265 // create the trackers
3266 AliInfo("Creating trackers");
3268 TString detStr = detectors;
3269 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3270 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3271 AliReconstructor* reconstructor = GetReconstructor(iDet);
3272 if (!reconstructor) continue;
3273 TString detName = fgkDetectorName[iDet];
3274 if (detName == "HLT") {
3275 fRunHLTTracking = kTRUE;
3278 if (detName == "MUON") {
3279 fRunMuonTracking = kTRUE;
3283 fTracker[iDet] = reconstructor->CreateTracker();
3284 if (!fTracker[iDet] && (iDet < 7)) {
3285 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
3286 if (fStopOnError) return kFALSE;
3288 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
3294 //_____________________________________________________________________________
3295 void AliReconstruction::CleanUp()
3297 // delete trackers and the run loader and close and delete the file
3299 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3300 delete fReconstructor[iDet];
3301 fReconstructor[iDet] = NULL;
3302 fLoader[iDet] = NULL;
3303 delete fTracker[iDet];
3304 fTracker[iDet] = NULL;
3311 delete fSPDTrackleter;
3312 fSPDTrackleter = NULL;
3321 delete fParentRawReader;
3322 fParentRawReader=NULL;
3330 if (AliQAManager::QAManager())
3331 AliQAManager::QAManager()->ShowQA() ;
3332 // AliQAManager::Destroy() ;
3336 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
3338 // Write space-points which are then used in the alignment procedures
3339 // For the moment only ITS, TPC, TRD and TOF
3341 Int_t ntracks = esd->GetNumberOfTracks();
3342 for (Int_t itrack = 0; itrack < ntracks; itrack++)
3344 AliESDtrack *track = esd->GetTrack(itrack);
3347 for (Int_t i=0; i<200; ++i) idx[i] = -1; //PH avoid uninitialized values
3348 for (Int_t iDet = 5; iDet >= 0; iDet--) {// TOF, TRD, TPC, ITS clusters
3349 nsp += (iDet==GetDetIndex("TRD")) ? track->GetTRDntracklets():track->GetNcls(iDet);
3351 if (iDet==GetDetIndex("ITS")) { // ITS "extra" clusters
3352 track->GetClusters(iDet,idx);
3353 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nsp++;
3358 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
3359 track->SetTrackPointArray(sp);
3361 for (Int_t iDet = 5; iDet >= 0; iDet--) {
3362 AliTracker *tracker = fTracker[iDet];
3363 if (!tracker) continue;
3364 Int_t nspdet = (iDet==GetDetIndex("TRD")) ? track->GetTRDtracklets(idx):track->GetClusters(iDet,idx);
3366 if (iDet==GetDetIndex("ITS")) // ITS "extra" clusters
3367 for (Int_t i=6; i<12; i++) if(idx[i] >= 0) nspdet++;
3369 if (nspdet <= 0) continue;
3373 while (isp2 < nspdet) {
3374 Bool_t isvalid=kTRUE;
3376 Int_t index=idx[isp++];
3377 if (index < 0) continue;
3379 TString dets = fgkDetectorName[iDet];
3380 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
3381 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
3382 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
3383 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
3384 isvalid = tracker->GetTrackPointTrackingError(index,p,track);
3386 isvalid = tracker->GetTrackPoint(index,p);
3389 if (!isvalid) continue;
3390 if (iDet==GetDetIndex("ITS") && (isp-1)>=6) p.SetExtra();
3391 sp->AddPoint(isptrack,&p); isptrack++;
3398 //_____________________________________________________________________________
3399 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
3401 // The method reads the raw-data error log
3402 // accumulated within the rawReader.
3403 // It extracts the raw-data errors related to
3404 // the current event and stores them into
3405 // a TClonesArray inside the esd object.
3407 if (!fRawReader) return;
3409 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
3411 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
3413 if (iEvent != log->GetEventNumber()) continue;
3415 esd->AddRawDataErrorLog(log);
3420 //_____________________________________________________________________________
3421 // void AliReconstruction::CheckQA()
3423 // check the QA of SIM for this run and remove the detectors
3424 // with status Fatal
3426 // TString newRunLocalReconstruction ;
3427 // TString newRunTracking ;
3428 // TString newFillESD ;
3430 // for (Int_t iDet = 0; iDet < AliQAv1::kNDET; iDet++) {
3431 // TString detName(AliQAv1::GetDetName(iDet)) ;
3432 // AliQAv1 * qa = AliQAv1::Instance(AliQAv1::DETECTORINDEX_t(iDet)) ;
3433 // if ( qa->IsSet(AliQAv1::DETECTORINDEX_t(iDet), AliQAv1::kSIM, specie, AliQAv1::kFATAL)) {
3434 // AliInfo(Form("QA status for %s %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed",
3435 // detName.Data(), AliRecoParam::GetEventSpecieName(es))) ;
3437 // if ( fRunLocalReconstruction.Contains(AliQAv1::GetDetName(iDet)) ||
3438 // fRunLocalReconstruction.Contains("ALL") ) {
3439 // newRunLocalReconstruction += detName ;
3440 // newRunLocalReconstruction += " " ;
3442 // if ( fRunTracking.Contains(AliQAv1::GetDetName(iDet)) ||
3443 // fRunTracking.Contains("ALL") ) {
3444 // newRunTracking += detName ;
3445 // newRunTracking += " " ;
3447 // if ( fFillESD.Contains(AliQAv1::GetDetName(iDet)) ||
3448 // fFillESD.Contains("ALL") ) {
3449 // newFillESD += detName ;
3450 // newFillESD += " " ;
3454 // fRunLocalReconstruction = newRunLocalReconstruction ;
3455 // fRunTracking = newRunTracking ;
3456 // fFillESD = newFillESD ;
3459 //_____________________________________________________________________________
3460 Int_t AliReconstruction::GetDetIndex(const char* detector)
3462 // return the detector index corresponding to detector
3464 for (index = 0; index < kNDetectors ; index++) {
3465 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
3470 //_____________________________________________________________________________
3471 Bool_t AliReconstruction::FinishPlaneEff() {
3473 // Here execute all the necessary operationis, at the end of the tracking phase,
3474 // in case that evaluation of PlaneEfficiencies was required for some detector.
3475 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
3477 // This Preliminary version works only FOR ITS !!!!!
3478 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3481 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3484 TString detStr = fLoadCDB;
3485 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
3486 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
3487 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3488 if(fTracker[iDet] && fTracker[iDet]->GetPlaneEff()) {
3489 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
3490 TString name=planeeff->GetName();
3492 TFile* pefile = TFile::Open(name, "RECREATE");
3493 ret=(Bool_t)planeeff->Write();
3495 if(planeeff->GetCreateHistos()) {
3496 TString hname=planeeff->GetName();
3497 hname+="Histo.root";
3498 ret*=planeeff->WriteHistosToFile(hname,"RECREATE");
3501 if(fSPDTrackleter) {
3502 AliPlaneEff *planeeff=fSPDTrackleter->GetPlaneEff();
3503 TString name="AliITSPlaneEffSPDtracklet.root";
3504 TFile* pefile = TFile::Open(name, "RECREATE");
3505 ret=(Bool_t)planeeff->Write();
3507 AliESDEvent *dummy=NULL;
3508 ret=(Bool_t)fSPDTrackleter->PostProcess(dummy); // take care of writing other files
3513 //_____________________________________________________________________________
3514 Bool_t AliReconstruction::InitPlaneEff() {
3516 // Here execute all the necessary operations, before of the tracking phase,
3517 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
3518 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
3519 // which should be updated/recalculated.
3521 // This Preliminary version will work only FOR ITS !!!!!
3522 // other detectors (TOF,TRD, etc. have to develop their specific codes)
3525 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
3528 fSPDTrackleter = NULL;
3529 TString detStr = fLoadCDB;
3530 if (IsSelected(fgkDetectorName[0], detStr)) {
3531 AliReconstructor* itsReconstructor = GetReconstructor(0);
3532 if (itsReconstructor) {
3533 fSPDTrackleter = itsReconstructor->CreateTrackleter(); // this is NULL unless required in RecoParam
3535 if (fSPDTrackleter) {
3536 AliInfo("Trackleter for SPD has been created");
3542 //_____________________________________________________________________________
3543 Bool_t AliReconstruction::InitAliEVE()
3545 // This method should be called only in case
3546 // AliReconstruction is run
3547 // within the alieve environment.
3548 // It will initialize AliEVE in a way
3549 // so that it can visualize event processed
3550 // by AliReconstruction.
3551 // The return flag shows whenever the
3552 // AliEVE initialization was successful or not.
3554 TString macroStr(gSystem->Getenv("ALIEVE_ONLINE_MACRO"));
3556 if (macroStr.IsNull())
3557 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
3559 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
3561 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
3563 gROOT->ProcessLine("if (!AliEveEventManager::GetMaster()){new AliEveEventManager();AliEveEventManager::GetMaster()->AddNewEventCommand(\"alieve_online_on_new_event()\");gEve->AddEvent(AliEveEventManager::GetMaster());};");
3564 gROOT->ProcessLine("alieve_online_init()");
3569 //_____________________________________________________________________________
3570 void AliReconstruction::RunAliEVE()
3572 // Runs AliEVE visualisation of
3573 // the current event.
3574 // Should be executed only after
3575 // successful initialization of AliEVE.
3577 AliInfo("Running AliEVE...");
3578 gROOT->ProcessLine(Form("AliEveEventManager::GetMaster()->SetEvent((AliRunLoader*)%p,(AliRawReader*)%p,(AliESDEvent*)%p,(AliESDfriend*)%p);",fRunLoader,fRawReader,fesd,fesdf));
3582 //_____________________________________________________________________________
3583 Bool_t AliReconstruction::SetRunQA(TString detAndAction)
3585 // Allows to run QA for a selected set of detectors
3586 // and a selected set of tasks among RAWS, DIGITSR, RECPOINTS and ESDS
3587 // all selected detectors run the same selected tasks
3589 if (!detAndAction.Contains(":")) {
3590 AliError( Form("%s is a wrong syntax, use \"DetectorList:ActionList\" \n", detAndAction.Data()) ) ;
3594 Int_t colon = detAndAction.Index(":") ;
3595 fQADetectors = detAndAction(0, colon) ;
3596 fQATasks = detAndAction(colon+1, detAndAction.Sizeof() ) ;
3597 if (fQATasks.Contains("ALL") ) {
3598 fQATasks = Form("%d %d %d %d", AliQAv1::kRAWS, AliQAv1::kDIGITSR, AliQAv1::kRECPOINTS, AliQAv1::kESDS) ;
3600 fQATasks.ToUpper() ;
3602 if ( fQATasks.Contains("RAW") )
3603 tempo = Form("%d ", AliQAv1::kRAWS) ;
3604 if ( fQATasks.Contains("DIGIT") )
3605 tempo += Form("%d ", AliQAv1::kDIGITSR) ;
3606 if ( fQATasks.Contains("RECPOINT") )
3607 tempo += Form("%d ", AliQAv1::kRECPOINTS) ;
3608 if ( fQATasks.Contains("ESD") )
3609 tempo += Form("%d ", AliQAv1::kESDS) ;
3611 if (fQATasks.IsNull()) {
3612 AliInfo("No QA requested\n") ;
3617 TString tempo(fQATasks) ;
3618 tempo.ReplaceAll(Form("%d", AliQAv1::kRAWS), AliQAv1::GetTaskName(AliQAv1::kRAWS)) ;
3619 tempo.ReplaceAll(Form("%d", AliQAv1::kDIGITSR), AliQAv1::GetTaskName(AliQAv1::kDIGITSR)) ;
3620 tempo.ReplaceAll(Form("%d", AliQAv1::kRECPOINTS), AliQAv1::GetTaskName(AliQAv1::kRECPOINTS)) ;
3621 tempo.ReplaceAll(Form("%d", AliQAv1::kESDS), AliQAv1::GetTaskName(AliQAv1::kESDS)) ;
3622 AliInfo( Form("QA will be done on \"%s\" for \"%s\"\n", fQADetectors.Data(), tempo.Data()) ) ;
3627 //_____________________________________________________________________________
3628 Bool_t AliReconstruction::InitRecoParams()
3630 // The method accesses OCDB and retrieves all
3631 // the available reco-param objects from there.
3633 Bool_t isOK = kTRUE;
3635 if (fRecoParam.GetDetRecoParamArray(kNDetectors)) {
3636 AliInfo("Using custom GRP reconstruction parameters");
3639 AliInfo("Loading GRP reconstruction parameter objects");
3641 AliCDBPath path("GRP","Calib","RecoParam");
3642 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3644 AliWarning("Couldn't find GRP RecoParam entry in OCDB");
3648 TObject *recoParamObj = entry->GetObject();
3649 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3650 // GRP has a normal TobjArray of AliDetectorRecoParam objects
3651 // Registering them in AliRecoParam
3652 fRecoParam.AddDetRecoParamArray(kNDetectors,dynamic_cast<TObjArray*>(recoParamObj));
3654 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3655 // GRP has only onse set of reco parameters
3656 // Registering it in AliRecoParam
3657 AliInfo("Single set of GRP reconstruction parameters found");
3658 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3659 fRecoParam.AddDetRecoParam(kNDetectors,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3662 AliError("No valid GRP RecoParam object found in the OCDB");
3669 TString detStr = fLoadCDB;
3670 for (Int_t iDet = 0; iDet < kNDetectors; iDet++) {
3672 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
3674 if (fRecoParam.GetDetRecoParamArray(iDet)) {
3675 AliInfo(Form("Using custom reconstruction parameters for detector %s",fgkDetectorName[iDet]));
3679 AliInfo(Form("Loading reconstruction parameter objects for detector %s",fgkDetectorName[iDet]));
3681 AliCDBPath path(fgkDetectorName[iDet],"Calib","RecoParam");
3682 AliCDBEntry *entry=AliCDBManager::Instance()->Get(path.GetPath());
3684 AliWarning(Form("Couldn't find RecoParam entry in OCDB for detector %s",fgkDetectorName[iDet]));
3688 TObject *recoParamObj = entry->GetObject();
3689 if (dynamic_cast<TObjArray*>(recoParamObj)) {
3690 // The detector has a normal TobjArray of AliDetectorRecoParam objects
3691 // Registering them in AliRecoParam
3692 fRecoParam.AddDetRecoParamArray(iDet,dynamic_cast<TObjArray*>(recoParamObj));
3694 else if (dynamic_cast<AliDetectorRecoParam*>(recoParamObj)) {
3695 // The detector has only onse set of reco parameters
3696 // Registering it in AliRecoParam
3697 AliInfo(Form("Single set of reconstruction parameters found for detector %s",fgkDetectorName[iDet]));
3698 dynamic_cast<AliDetectorRecoParam*>(recoParamObj)->SetAsDefault();
3699 fRecoParam.AddDetRecoParam(iDet,dynamic_cast<AliDetectorRecoParam*>(recoParamObj));
3702 AliError(Form("No valid RecoParam object found in the OCDB for detector %s",fgkDetectorName[iDet]));
3706 // FIX ME: We have to disable the unloading of reco-param CDB
3707 // entries because QA framework is using them. Has to be fix in
3708 // a way that the QA takes the objects already constructed in
3710 // AliCDBManager::Instance()->UnloadFromCache(path.GetPath());
3714 if (AliDebugLevel() > 0) fRecoParam.Print();
3719 //_____________________________________________________________________________
3720 Bool_t AliReconstruction::GetEventInfo()
3722 // Fill the event info object
3724 AliCodeTimerAuto("",0)
3726 AliCentralTrigger *aCTP = NULL;
3728 fEventInfo.SetEventType(fRawReader->GetType());
3730 ULong64_t mask = fRawReader->GetClassMask();
3731 fEventInfo.SetTriggerMask(mask);
3732 UInt_t clmask = fRawReader->GetDetectorPattern()[0];
3733 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(clmask));
3735 aCTP = new AliCentralTrigger();
3736 TString configstr("");
3737 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
3738 AliError("No trigger configuration found in OCDB! The trigger configuration information will not be used!");
3742 aCTP->SetClassMask(mask);
3743 aCTP->SetClusterMask(clmask);
3745 AliCentralTrigger* rlCTP = fRunLoader->GetTrigger();
3746 rlCTP->SetClassMask(mask);
3747 rlCTP->SetClusterMask(clmask);
3750 fEventInfo.SetEventType(AliRawEventHeaderBase::kPhysicsEvent);
3752 if (fRunLoader && (!fRunLoader->LoadTrigger())) {
3753 aCTP = fRunLoader->GetTrigger();
3754 fEventInfo.SetTriggerMask(aCTP->GetClassMask());
3755 // get inputs from actp - just get
3756 AliESDHeader* esdheader = fesd->GetHeader();
3757 esdheader->SetL0TriggerInputs(aCTP->GetL0TriggerInputs());
3758 esdheader->SetL1TriggerInputs(aCTP->GetL1TriggerInputs());
3759 esdheader->SetL2TriggerInputs(aCTP->GetL2TriggerInputs());
3760 fEventInfo.SetTriggerCluster(AliDAQ::ListOfTriggeredDetectors(aCTP->GetClusterMask()));
3763 AliWarning("No trigger can be loaded! The trigger information will not be used!");
3768 AliTriggerConfiguration *config = aCTP->GetConfiguration();
3770 AliError("No trigger configuration has been found! The trigger configuration information will not be used!");
3771 if (fRawReader) delete aCTP;
3775 UChar_t clustmask = 0;
3777 ULong64_t trmask = fEventInfo.GetTriggerMask();
3778 const TObjArray& classesArray = config->GetClasses();
3779 Int_t nclasses = classesArray.GetEntriesFast();
3780 for( Int_t iclass=0; iclass < nclasses; iclass++ ) {
3781 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At(iclass);
3782 if (trclass && trclass->GetMask()>0) {
3783 Int_t trindex = TMath::Nint(TMath::Log2(trclass->GetMask()));
3784 fesd->SetTriggerClass(trclass->GetName(),trindex);
3785 if (fRawReader) fRawReader->LoadTriggerClass(trclass->GetName(),trindex);
3786 if (trmask & (1ull << trindex)) {
3788 trclasses += trclass->GetName();
3790 clustmask |= trclass->GetCluster()->GetClusterMask();
3794 fEventInfo.SetTriggerClasses(trclasses);
3796 // Write names of active trigger inputs in ESD Header
3797 const TObjArray& inputsArray = config->GetInputs();
3798 Int_t ninputs = inputsArray.GetEntriesFast();
3799 for( Int_t iinput=0; iinput < ninputs; iinput++ ) {
3800 AliTriggerInput* trginput = (AliTriggerInput*)inputsArray.At(iinput);
3801 if (trginput && trginput->GetMask()>0) {
3802 Int_t inputIndex = (Int_t)TMath::Nint(TMath::Log2(trginput->GetMask()));
3803 AliESDHeader* headeresd = fesd->GetHeader();
3804 Int_t trglevel = (Int_t)trginput->GetLevel();
3805 if (trglevel == 0) headeresd->SetActiveTriggerInputs(trginput->GetInputName(), inputIndex);
3806 if (trglevel == 1) headeresd->SetActiveTriggerInputs(trginput->GetInputName(), inputIndex+24);
3807 if (trglevel == 2) headeresd->SetActiveTriggerInputs(trginput->GetInputName(), inputIndex+48);
3811 // Set the information in ESD
3812 fesd->SetTriggerMask(trmask);
3813 fesd->SetTriggerCluster(clustmask);
3815 if (!aCTP->CheckTriggeredDetectors()) {
3816 if (fRawReader) delete aCTP;
3820 if (fRawReader) delete aCTP;
3822 // We have to fill also the HLT decision here!!
3828 const char *AliReconstruction::MatchDetectorList(const char *detectorList, UInt_t detectorMask)
3830 // Match the detector list found in the rec.C or the default 'ALL'
3831 // to the list found in the GRP (stored there by the shuttle PP which
3832 // gets the information from ECS)
3833 static TString resultList;
3834 TString detList = detectorList;
3838 for(Int_t iDet = 0; iDet < (AliDAQ::kNDetectors-1); iDet++) {
3839 if ((detectorMask >> iDet) & 0x1) {
3840 TString det = AliDAQ::OfflineModuleName(iDet);
3841 if ((detList.CompareTo("ALL") == 0) ||
3842 ((detList.BeginsWith("ALL ") ||
3843 detList.EndsWith(" ALL") ||
3844 detList.Contains(" ALL ")) &&
3845 !(detList.BeginsWith("-"+det+" ") ||
3846 detList.EndsWith(" -"+det) ||
3847 detList.Contains(" -"+det+" "))) ||
3848 (detList.CompareTo(det) == 0) ||
3849 detList.BeginsWith(det+" ") ||
3850 detList.EndsWith(" "+det) ||
3851 detList.Contains( " "+det+" " )) {
3852 if (!resultList.EndsWith(det + " ")) {
3861 if ((detectorMask >> AliDAQ::kHLTId) & 0x1) {
3862 TString hltDet = AliDAQ::OfflineModuleName(AliDAQ::kNDetectors-1);
3863 if ((detList.CompareTo("ALL") == 0) ||
3864 ((detList.BeginsWith("ALL ") ||
3865 detList.EndsWith(" ALL") ||
3866 detList.Contains(" ALL ")) &&
3867 !(detList.BeginsWith("-"+hltDet+" ") ||
3868 detList.EndsWith(" -"+hltDet) ||
3869 detList.Contains(" -"+hltDet+" "))) ||
3870 (detList.CompareTo(hltDet) == 0) ||
3871 detList.BeginsWith(hltDet+" ") ||
3872 detList.EndsWith(" "+hltDet) ||
3873 detList.Contains( " "+hltDet+" " )) {
3874 resultList += hltDet;
3878 return resultList.Data();
3882 //______________________________________________________________________________
3883 void AliReconstruction::Abort(const char *method, EAbort what)
3885 // Abort processing. If what = kAbortProcess, the Process() loop will be
3886 // aborted. If what = kAbortFile, the current file in a chain will be
3887 // aborted and the processing will continue with the next file, if there
3888 // is no next file then Process() will be aborted. Abort() can also be
3889 // called from Begin(), SlaveBegin(), Init() and Notify(). After abort
3890 // the SlaveTerminate() and Terminate() are always called. The abort flag
3891 // can be checked in these methods using GetAbort().
3893 // The method is overwritten in AliReconstruction for better handling of
3894 // reco specific errors
3896 if (!fStopOnError) return;
3900 TString whyMess = method;
3901 whyMess += " failed! Aborting...";
3903 AliError(whyMess.Data());
3906 TString mess = "Abort";
3907 if (fAbort == kAbortProcess)
3908 mess = "AbortProcess";
3909 else if (fAbort == kAbortFile)
3912 Info(mess, whyMess.Data());
3915 //______________________________________________________________________________
3916 Bool_t AliReconstruction::ProcessEvent(void* event)
3918 // Method that is used in case the event loop
3919 // is steered from outside, for example by AMORE
3920 // 'event' is a pointer to the DATE event in the memory
3922 if (fRawReader) delete fRawReader;
3923 fRawReader = new AliRawReaderDate(event);
3924 fStatus = ProcessEvent(fRunLoader->GetNumberOfEvents());
3931 //______________________________________________________________________________
3932 Bool_t AliReconstruction::ParseOutput()
3934 // The method parses the output file
3935 // location string in order to steer
3936 // properly the selector
3938 TPMERegexp re1("(\\w+\\.zip#\\w+\\.root):([,*\\w+\\.root,*]+)@dataset://(\\w++)");
3939 TPMERegexp re2("(\\w+\\.root)?@?dataset://(\\w++)");
3941 if (re1.Match(fESDOutput) == 4) {
3942 // root archive with output files stored and regustered
3944 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",re1[1].Data()));
3945 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",re1[3].Data()));
3946 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_DATASET",""));
3947 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_ARCHIVE",re1[2].Data()));
3948 AliInfo(Form("%s files will be stored within %s in dataset %s",
3953 else if (re2.Match(fESDOutput) == 3) {
3954 // output file stored and registered
3956 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",(re2[1].IsNull()) ? "AliESDs.root" : re2[1].Data()));
3957 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",re2[2].Data()));
3958 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_DATASET",""));
3959 AliInfo(Form("%s will be stored in dataset %s",
3960 (re2[1].IsNull()) ? "AliESDs.root" : re2[1].Data(),
3964 if (fESDOutput.IsNull()) {
3965 // Output location not given.
3966 // Assuming xrootd has been already started and
3967 // the output file has to be sent back
3968 // to the client machine
3969 TString esdUrl(Form("root://%s/%s/",
3970 TUrl(gSystem->HostName()).GetHostFQDN(),
3972 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE","AliESDs.root"));
3973 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",esdUrl.Data()));
3974 AliInfo(Form("AliESDs.root will be stored in %s",
3978 // User specified an output location.
3979 // Ones has just to parse it here
3980 TUrl outputUrl(fESDOutput.Data());
3981 TString outputFile(gSystem->BaseName(outputUrl.GetFile()));
3982 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE",outputFile.IsNull() ? "AliESDs.root" : outputFile.Data()));
3983 TString outputLocation(outputUrl.GetUrl());
3984 outputLocation.ReplaceAll(outputFile.Data(),"");
3985 gProof->AddInput(new TNamed("PROOF_OUTPUTFILE_LOCATION",outputLocation.Data()));
3986 AliInfo(Form("%s will be stored in %s",
3987 outputFile.IsNull() ? "AliESDs.root" : outputFile.Data(),
3988 outputLocation.Data()));
3995 //______________________________________________________________________________
3996 Bool_t AliReconstruction::IsHighPt() const {
3997 // Selection of events containing "high" pT tracks
3998 // If at least one track is found within 1.5 and 100 GeV (pT)
3999 // that was reconstructed by both ITS and TPC, the event is accepted
4002 const Double_t pTmin = 1.5;
4003 const Double_t pTmax = 100;
4005 mask |= (AliESDtrack::kITSrefit);
4006 mask |= (AliESDtrack::kTPCrefit);
4007 const Double_t pTminCosmic = 5.;
4008 const Double_t pTmaxCosmic = 100;
4009 ULong_t maskCosmic = 0;
4010 Int_t cosmicCount=0;
4011 maskCosmic |= (AliESDtrack::kTPCrefit);
4013 Bool_t isOK = kFALSE;
4015 if (fesd && fesd->GetEventType()==AliRawEventHeaderBase::kPhysicsEvent) {
4016 // Check if this ia a physics event (code 7)
4017 Int_t ntrk = fesd->GetNumberOfTracks();
4018 for (Int_t itrk=0; itrk<ntrk; ++itrk) {
4020 AliESDtrack * trk = fesd->GetTrack(itrk);
4022 && trk->Pt() > pTmin
4023 && trk->Pt() < pTmax
4024 && (trk->GetStatus() & mask) == mask ) {
4030 && trk->GetInnerParam()
4031 && trk->GetInnerParam()->Pt() > pTminCosmic
4032 && trk->GetInnerParam()->Pt() < pTmaxCosmic
4033 && (trk->GetStatus() & maskCosmic) == maskCosmic ) {
4039 if (cosmicCount>1) isOK=kTRUE;
4044 //______________________________________________________________________________
4045 Bool_t AliReconstruction::IsCosmicOrCalibSpecie() const {
4046 // Select cosmic or calibration events
4048 Bool_t isOK = kFALSE;
4050 if (fesd && fesd->GetEventType()==AliRawEventHeaderBase::kPhysicsEvent) {
4051 // Check if this ia a physics event (code 7)
4053 UInt_t specie = fesd->GetEventSpecie();
4054 if (specie==AliRecoParam::kCosmic || specie==AliRecoParam::kCalib) {
4061 //______________________________________________________________________________
4062 void AliReconstruction::WriteESDfriend() {
4063 // Fill the ESD friend in the tree. The required fraction of ESD friends is stored
4064 // in fFractionFriends. We select events where we store the ESD friends according
4065 // to the following algorithm:
4066 // 1. Store all Cosmic or Calibration events within the required fraction
4067 // 2. Sample "high Pt" events within the remaining fraction after step 1.
4068 // 3. Sample randomly events if we still have remaining slot
4071 Bool_t isSelected = kFALSE;
4073 // Store all friends for B field OFF
4074 if (TMath::Abs(AliTrackerBase::GetBz())<0.5) isSelected=kTRUE;
4076 if (IsCosmicOrCalibSpecie()) { // Selection of calib or cosmic events
4078 Double_t curentSpecieFraction = ((Double_t)(fNspecie+1))/((Double_t)(fNall+1));
4079 // "Bayesian" estimate supposing that without events all the events are of the required type
4081 Double_t rnd = gRandom->Rndm()*curentSpecieFraction;
4082 if (rnd<fFractionFriends) {
4088 Double_t remainingFraction = fFractionFriends;
4089 remainingFraction -= ((Double_t)(fSspecie)/(Double_t)(fNall));
4091 if (IsHighPt()) { // Selection of "high Pt" events
4093 Double_t curentHighPtFraction = ((Double_t)(fNhighPt+1))/((Double_t)(fNall+1));
4094 // "Bayesian" estimate supposing that without events all the events are of the required type
4097 Double_t rnd = gRandom->Rndm()*curentHighPtFraction;
4098 if (rnd<remainingFraction) {
4104 remainingFraction -= ((Double_t)(fShighPt)/(Double_t)(fNall));
4106 // Random selection to fill the remaining fraction (if any)
4108 Double_t rnd = gRandom->Rndm();
4109 if (rnd<remainingFraction) {
4115 fesdf->~AliESDfriend();
4116 new (fesdf) AliESDfriend(); // Reset...
4117 fesdf->SetSkipBit(kTRUE);