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 // For debug purposes the method SetCheckPointLevel can be used. If the //
105 // argument is greater than 0, files with ESD events will be written after //
106 // selected steps of the reconstruction for each event: //
107 // level 1: after tracking and after filling of ESD (final) //
108 // level 2: in addition after each tracking step //
109 // level 3: in addition after the filling of ESD for each detector //
110 // If a final check point file exists for an event, this event will be //
111 // skipped in the reconstruction. The tracking and the filling of ESD for //
112 // a detector will be skipped as well, if the corresponding check point //
113 // file exists. The ESD event will then be loaded from the file instead. //
115 ///////////////////////////////////////////////////////////////////////////////
122 #include <TPluginManager.h>
123 #include <TGeoManager.h>
124 #include <TLorentzVector.h>
128 #include "AliReconstruction.h"
129 #include "AliCodeTimer.h"
130 #include "AliReconstructor.h"
132 #include "AliRunLoader.h"
134 #include "AliRawReaderFile.h"
135 #include "AliRawReaderDate.h"
136 #include "AliRawReaderRoot.h"
137 #include "AliRawEventHeaderBase.h"
138 #include "AliESDEvent.h"
139 #include "AliESDMuonTrack.h"
140 #include "AliESDfriend.h"
141 #include "AliESDVertex.h"
142 #include "AliESDcascade.h"
143 #include "AliESDkink.h"
144 #include "AliESDtrack.h"
145 #include "AliESDCaloCluster.h"
146 #include "AliESDCaloCells.h"
147 #include "AliMultiplicity.h"
148 #include "AliTracker.h"
149 #include "AliVertexer.h"
150 #include "AliVertexerTracks.h"
151 #include "AliV0vertexer.h"
152 #include "AliCascadeVertexer.h"
153 #include "AliHeader.h"
154 #include "AliGenEventHeader.h"
156 #include "AliESDpid.h"
157 #include "AliESDtrack.h"
158 #include "AliESDPmdTrack.h"
160 #include "AliESDTagCreator.h"
161 #include "AliAODTagCreator.h"
163 #include "AliGeomManager.h"
164 #include "AliTrackPointArray.h"
165 #include "AliCDBManager.h"
166 #include "AliCDBStorage.h"
167 #include "AliCDBEntry.h"
168 #include "AliAlignObj.h"
170 #include "AliCentralTrigger.h"
171 #include "AliCTPRawStream.h"
173 #include "AliAODEvent.h"
174 #include "AliAODHeader.h"
175 #include "AliAODTrack.h"
176 #include "AliAODVertex.h"
177 #include "AliAODv0.h"
178 #include "AliAODJet.h"
179 #include "AliAODCaloCells.h"
180 #include "AliAODCaloCluster.h"
181 #include "AliAODPmdCluster.h"
182 #include "AliAODFmdCluster.h"
183 #include "AliAODTracklets.h"
185 #include "AliQADataMaker.h"
187 #include "AliQADataMakerSteer.h"
189 #include "AliSysInfo.h" // memory snapshots
192 ClassImp(AliReconstruction)
195 //_____________________________________________________________________________
196 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
198 //_____________________________________________________________________________
199 AliReconstruction::AliReconstruction(const char* gAliceFilename,
200 const char* name, const char* title) :
203 fUniformField(kTRUE),
204 fRunVertexFinder(kTRUE),
205 fRunHLTTracking(kFALSE),
206 fRunMuonTracking(kFALSE),
208 fRunCascadeFinder(kTRUE),
209 fStopOnError(kFALSE),
210 fWriteAlignmentData(kFALSE),
211 fWriteESDfriend(kFALSE),
213 fFillTriggerESD(kTRUE),
221 fRunLocalReconstruction("ALL"),
224 fUseTrackingErrorsForAlignment(""),
225 fGAliceFileName(gAliceFilename),
230 fNumberOfEventsPerFile(1),
233 fLoadAlignFromCDB(kTRUE),
234 fLoadAlignData("ALL"),
241 fDiamondProfile(NULL),
242 fMeanVertexConstraint(kTRUE),
246 fAlignObjArray(NULL),
249 fInitCDBCalled(kFALSE),
250 fSetRunNumberFromDataCalled(kFALSE),
255 // create reconstruction object with default parameters
257 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
258 fReconstructor[iDet] = NULL;
259 fLoader[iDet] = NULL;
260 fTracker[iDet] = NULL;
261 fQADataMaker[iDet] = NULL;
262 fQACycles[iDet] = 999999;
267 //_____________________________________________________________________________
268 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
271 fUniformField(rec.fUniformField),
272 fRunVertexFinder(rec.fRunVertexFinder),
273 fRunHLTTracking(rec.fRunHLTTracking),
274 fRunMuonTracking(rec.fRunMuonTracking),
275 fRunV0Finder(rec.fRunV0Finder),
276 fRunCascadeFinder(rec.fRunCascadeFinder),
277 fStopOnError(rec.fStopOnError),
278 fWriteAlignmentData(rec.fWriteAlignmentData),
279 fWriteESDfriend(rec.fWriteESDfriend),
280 fWriteAOD(rec.fWriteAOD),
281 fFillTriggerESD(rec.fFillTriggerESD),
283 fCleanESD(rec.fCleanESD),
284 fV0DCAmax(rec.fV0DCAmax),
285 fV0CsPmin(rec.fV0CsPmin),
289 fRunLocalReconstruction(rec.fRunLocalReconstruction),
290 fRunTracking(rec.fRunTracking),
291 fFillESD(rec.fFillESD),
292 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
293 fGAliceFileName(rec.fGAliceFileName),
295 fEquipIdMap(rec.fEquipIdMap),
296 fFirstEvent(rec.fFirstEvent),
297 fLastEvent(rec.fLastEvent),
298 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
301 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
302 fLoadAlignData(rec.fLoadAlignData),
303 fESDPar(rec.fESDPar),
309 fDiamondProfile(NULL),
310 fMeanVertexConstraint(rec.fMeanVertexConstraint),
314 fAlignObjArray(rec.fAlignObjArray),
315 fCDBUri(rec.fCDBUri),
317 fInitCDBCalled(rec.fInitCDBCalled),
318 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
325 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
326 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
328 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
329 fReconstructor[iDet] = NULL;
330 fLoader[iDet] = NULL;
331 fTracker[iDet] = NULL;
332 fQADataMaker[iDet] = NULL;
333 fQACycles[iDet] = rec.fQACycles[iDet];
335 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
336 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
340 //_____________________________________________________________________________
341 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
343 // assignment operator
345 this->~AliReconstruction();
346 new(this) AliReconstruction(rec);
350 //_____________________________________________________________________________
351 AliReconstruction::~AliReconstruction()
357 fSpecCDBUri.Delete();
359 AliCodeTimer::Instance()->Print();
362 //_____________________________________________________________________________
363 void AliReconstruction::InitCDB()
365 // activate a default CDB storage
366 // First check if we have any CDB storage set, because it is used
367 // to retrieve the calibration and alignment constants
369 if (fInitCDBCalled) return;
370 fInitCDBCalled = kTRUE;
372 AliCDBManager* man = AliCDBManager::Instance();
373 if (man->IsDefaultStorageSet())
375 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
376 AliWarning("Default CDB storage has been already set !");
377 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
378 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
379 fCDBUri = man->GetDefaultStorage()->GetURI();
382 if (fCDBUri.Length() > 0)
384 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
385 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
386 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
388 fCDBUri="local://$ALICE_ROOT";
389 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
390 AliWarning("Default CDB storage not yet set !!!!");
391 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
392 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
395 man->SetDefaultStorage(fCDBUri);
398 // Now activate the detector specific CDB storage locations
399 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
400 TObject* obj = fSpecCDBUri[i];
402 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
403 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
404 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
405 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
410 //_____________________________________________________________________________
411 void AliReconstruction::SetDefaultStorage(const char* uri) {
412 // Store the desired default CDB storage location
413 // Activate it later within the Run() method
419 //_____________________________________________________________________________
420 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
421 // Store a detector-specific CDB storage location
422 // Activate it later within the Run() method
424 AliCDBPath aPath(calibType);
425 if(!aPath.IsValid()){
426 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
427 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
428 if(!strcmp(calibType, fgkDetectorName[iDet])) {
429 aPath.SetPath(Form("%s/*", calibType));
430 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
434 if(!aPath.IsValid()){
435 AliError(Form("Not a valid path or detector: %s", calibType));
440 // // check that calibType refers to a "valid" detector name
441 // Bool_t isDetector = kFALSE;
442 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
443 // TString detName = fgkDetectorName[iDet];
444 // if(aPath.GetLevel0() == detName) {
445 // isDetector = kTRUE;
451 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
455 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
456 if (obj) fSpecCDBUri.Remove(obj);
457 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
461 //_____________________________________________________________________________
462 Bool_t AliReconstruction::SetRunNumberFromData()
464 // The method is called in Run() in order
465 // to set a correct run number.
466 // In case of raw data reconstruction the
467 // run number is taken from the raw data header
469 if (fSetRunNumberFromDataCalled) return kTRUE;
470 fSetRunNumberFromDataCalled = kTRUE;
472 AliCDBManager* man = AliCDBManager::Instance();
474 if(man->GetRun() > 0) {
475 AliWarning("Run number is taken from event header! Ignoring settings in AliCDBManager!");
479 AliError("No run loader is found !");
482 // read run number from gAlice
483 if(fRunLoader->GetAliRun())
484 AliCDBManager::Instance()->SetRun(fRunLoader->GetHeader()->GetRun());
487 if(fRawReader->NextEvent()) {
488 AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
489 fRawReader->RewindEvents();
492 AliError("No raw-data events found !");
497 AliError("Neither gAlice nor RawReader objects are found !");
507 //_____________________________________________________________________________
508 void AliReconstruction::SetCDBLock() {
509 // Set CDB lock: from now on it is forbidden to reset the run number
510 // or the default storage or to activate any further storage!
512 AliCDBManager::Instance()->SetLock(1);
515 //_____________________________________________________________________________
516 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
518 // Read the alignment objects from CDB.
519 // Each detector is supposed to have the
520 // alignment objects in DET/Align/Data CDB path.
521 // All the detector objects are then collected,
522 // sorted by geometry level (starting from ALIC) and
523 // then applied to the TGeo geometry.
524 // Finally an overlaps check is performed.
526 // Load alignment data from CDB and fill fAlignObjArray
527 if(fLoadAlignFromCDB){
529 TString detStr = detectors;
530 TString loadAlObjsListOfDets = "";
532 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
533 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
534 loadAlObjsListOfDets += fgkDetectorName[iDet];
535 loadAlObjsListOfDets += " ";
536 } // end loop over detectors
537 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
538 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
540 // Check if the array with alignment objects was
541 // provided by the user. If yes, apply the objects
542 // to the present TGeo geometry
543 if (fAlignObjArray) {
544 if (gGeoManager && gGeoManager->IsClosed()) {
545 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
546 AliError("The misalignment of one or more volumes failed!"
547 "Compare the list of simulated detectors and the list of detector alignment data!");
552 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
558 delete fAlignObjArray; fAlignObjArray=0;
563 //_____________________________________________________________________________
564 void AliReconstruction::SetGAliceFile(const char* fileName)
566 // set the name of the galice file
568 fGAliceFileName = fileName;
571 //_____________________________________________________________________________
572 void AliReconstruction::SetOption(const char* detector, const char* option)
574 // set options for the reconstruction of a detector
576 TObject* obj = fOptions.FindObject(detector);
577 if (obj) fOptions.Remove(obj);
578 fOptions.Add(new TNamed(detector, option));
582 //_____________________________________________________________________________
583 Bool_t AliReconstruction::Run(const char* input)
585 // run the reconstruction
590 if (!input) input = fInput.Data();
591 TString fileName(input);
592 if (fileName.EndsWith("/")) {
593 fRawReader = new AliRawReaderFile(fileName);
594 } else if (fileName.EndsWith(".root")) {
595 fRawReader = new AliRawReaderRoot(fileName);
596 } else if (!fileName.IsNull()) {
597 fRawReader = new AliRawReaderDate(fileName);
598 fRawReader->SelectEvents(7);
600 if (!fEquipIdMap.IsNull() && fRawReader)
601 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
603 AliSysInfo::AddStamp("Start");
604 // get the run loader
605 if (!InitRunLoader()) return kFALSE;
606 AliSysInfo::AddStamp("LoadLoader");
608 // Initialize the CDB storage
611 AliSysInfo::AddStamp("LoadCDB");
613 // Set run number in CDBManager (if it is not already set by the user)
614 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
616 // Set CDB lock: from now on it is forbidden to reset the run number
617 // or the default storage or to activate any further storage!
620 // Import ideal TGeo geometry and apply misalignment
622 TString geom(gSystem->DirName(fGAliceFileName));
623 geom += "/geometry.root";
624 AliGeomManager::LoadGeometry(geom.Data());
625 if (!gGeoManager) if (fStopOnError) return kFALSE;
628 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
629 AliSysInfo::AddStamp("LoadGeom");
632 AliQADataMakerSteer qas ;
633 if ( fRunQA && fRawReader)
634 qas.Run(fRunLocalReconstruction, fRawReader) ;
636 // checking the QA of previous steps
639 // local reconstruction
640 if (!fRunLocalReconstruction.IsNull()) {
641 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
642 if (fStopOnError) {CleanUp(); return kFALSE;}
645 // if (!fRunVertexFinder && fRunTracking.IsNull() &&
646 // fFillESD.IsNull()) return kTRUE;
649 if (fRunVertexFinder && !CreateVertexer()) {
655 AliSysInfo::AddStamp("Vertexer");
658 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
664 AliSysInfo::AddStamp("LoadTrackers");
666 // get the possibly already existing ESD file and tree
667 AliESDEvent* esd = new AliESDEvent(); AliESDEvent* hltesd = new AliESDEvent();
668 TFile* fileOld = NULL;
669 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
670 if (!gSystem->AccessPathName("AliESDs.root")){
671 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
672 fileOld = TFile::Open("AliESDs.old.root");
673 if (fileOld && fileOld->IsOpen()) {
674 treeOld = (TTree*) fileOld->Get("esdTree");
675 if (treeOld)esd->ReadFromTree(treeOld);
676 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
677 if (hlttreeOld) hltesd->ReadFromTree(hlttreeOld);
681 // create the ESD output file and tree
682 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
683 file->SetCompressionLevel(2);
684 if (!file->IsOpen()) {
685 AliError("opening AliESDs.root failed");
686 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
689 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
690 esd = new AliESDEvent();
691 esd->CreateStdContent();
692 esd->WriteToTree(tree);
694 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
695 hltesd = new AliESDEvent();
696 hltesd->CreateStdContent();
697 hltesd->WriteToTree(hlttree);
700 delete esd; delete hltesd;
701 esd = NULL; hltesd = NULL;
703 // create the branch with ESD additions
707 AliESDfriend *esdf = 0;
708 if (fWriteESDfriend) {
709 esdf = new AliESDfriend();
710 TBranch *br=tree->Branch("ESDfriend.","AliESDfriend", &esdf);
711 br->SetFile("AliESDfriends.root");
712 esd->AddObject(esdf);
716 // Get the GRP CDB entry
717 AliCDBEntry* entryGRP = AliCDBManager::Instance()->Get("GRP/GRP/Data");
720 fGRPList = dynamic_cast<TList*> (entryGRP->GetObject());
722 AliError("No GRP entry found in OCDB!");
725 // Get the diamond profile from OCDB
726 AliCDBEntry* entry = AliCDBManager::Instance()
727 ->Get("GRP/Calib/MeanVertex");
730 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
732 AliError("No diamond profile found in OCDB!");
735 AliVertexerTracks tVertexer(AliTracker::GetBz());
736 if(fDiamondProfile && fMeanVertexConstraint) tVertexer.SetVtxStart(fDiamondProfile);
740 if (fRawReader) fRawReader->RewindEvents();
741 TString detStr(fFillESD) ;
744 gSystem->GetProcInfo(&ProcInfo);
745 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
747 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
748 if (fRawReader) fRawReader->NextEvent();
749 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
750 // copy old ESD to the new one
752 esd->ReadFromTree(treeOld);
753 treeOld->GetEntry(iEvent);
757 esd->ReadFromTree(hlttreeOld);
758 hlttreeOld->GetEntry(iEvent);
764 AliInfo(Form("processing event %d", iEvent));
765 fRunLoader->GetEvent(iEvent);
768 sprintf(aFileName, "ESD_%d.%d_final.root",
769 fRunLoader->GetHeader()->GetRun(),
770 fRunLoader->GetHeader()->GetEventNrInRun());
771 if (!gSystem->AccessPathName(aFileName)) continue;
773 // local reconstruction
774 if (!fRunLocalReconstruction.IsNull()) {
775 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
776 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
781 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
782 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
783 esd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
784 hltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
786 // Set magnetic field from the tracker
787 esd->SetMagneticField(AliTracker::GetBz());
788 hltesd->SetMagneticField(AliTracker::GetBz());
792 // Fill raw-data error log into the ESD
793 if (fRawReader) FillRawDataErrorLog(iEvent,esd);
796 if (fRunVertexFinder) {
797 if (!ReadESD(esd, "vertex")) {
798 if (!RunVertexFinder(esd)) {
799 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
801 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
806 if (!fRunTracking.IsNull()) {
807 if (fRunHLTTracking) {
808 hltesd->SetVertex(esd->GetVertex());
809 if (!RunHLTTracking(hltesd)) {
810 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
816 if (!fRunTracking.IsNull()) {
817 if (fRunMuonTracking) {
818 if (!RunMuonTracking(esd)) {
819 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
825 if (!fRunTracking.IsNull()) {
826 if (!ReadESD(esd, "tracking")) {
827 if (!RunTracking(esd)) {
828 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
830 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
835 if (!fFillESD.IsNull()) {
836 if (!FillESD(esd, fFillESD)) {
837 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
841 // if (!fFillESD.IsNull())
843 // do QA in the event loop if requested
845 RunQA(fFillESD.Data(), esd);
847 // fill Event header information from the RawEventHeader
848 if (fRawReader){FillRawEventHeaderESD(esd);}
851 AliESDpid::MakePID(esd);
852 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
854 if (fFillTriggerESD) {
855 if (!ReadESD(esd, "trigger")) {
856 if (!FillTriggerESD(esd)) {
857 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
859 if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
865 //Try to improve the reconstructed primary vertex position using the tracks
866 AliESDVertex *pvtx=0;
867 Bool_t dovertex=kTRUE;
868 TObject* obj = fOptions.FindObject("ITS");
870 TString optITS = obj->GetTitle();
871 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
874 if(dovertex) pvtx=tVertexer.FindPrimaryVertex(esd);
875 if(fDiamondProfile) esd->SetDiamond(fDiamondProfile);
878 if (pvtx->GetStatus()) {
879 // Store the improved primary vertex
880 esd->SetPrimaryVertex(pvtx);
881 // Propagate the tracks to the DCA to the improved primary vertex
882 Double_t somethingbig = 777.;
883 Double_t bz = esd->GetMagneticField();
884 Int_t nt=esd->GetNumberOfTracks();
886 AliESDtrack *t = esd->GetTrack(nt);
887 t->RelateToVertex(pvtx, bz, somethingbig);
894 vtxer.Tracks2V0vertices(esd);
896 if (fRunCascadeFinder) {
898 AliCascadeVertexer cvtxer;
899 cvtxer.V0sTracks2CascadeVertices(esd);
904 if (fCleanESD) CleanESD(esd);
905 if (fWriteESDfriend) {
906 esdf->~AliESDfriend();
907 new (esdf) AliESDfriend(); // Reset...
908 esd->GetESDfriend(esdf);
915 if (fCheckPointLevel > 0) WriteESD(esd, "final");
918 if (fWriteESDfriend) {
919 esdf->~AliESDfriend();
920 new (esdf) AliESDfriend(); // Reset...
923 // delete esdf; esdf = 0;
926 gSystem->GetProcInfo(&ProcInfo);
927 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
933 // write quality assurance ESDs data (one entry for all events)
935 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
936 if (!IsSelected(fgkDetectorName[iDet], detStr))
938 AliQADataMaker * qadm = GetQADataMaker(iDet);
940 qadm->EndOfCycle(AliQA::kRECPOINTS);
941 qadm->EndOfCycle(AliQA::kESDS);
946 tree->GetUserInfo()->Add(esd);
947 hlttree->GetUserInfo()->Add(hltesd);
949 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
950 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
952 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
953 cdbMapCopy->SetOwner(1);
954 cdbMapCopy->SetName("cdbMap");
955 TIter iter(cdbMap->GetTable());
958 while((pair = dynamic_cast<TPair*> (iter.Next()))){
959 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
960 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
961 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
964 TList *cdbListCopy = new TList();
965 cdbListCopy->SetOwner(1);
966 cdbListCopy->SetName("cdbList");
968 TIter iter2(cdbList);
971 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
972 cdbListCopy->Add(new TObjString(id->ToString().Data()));
975 tree->GetUserInfo()->Add(cdbMapCopy);
976 tree->GetUserInfo()->Add(cdbListCopy);
979 if(fESDPar.Contains("ESD.par")){
980 AliInfo("Attaching ESD.par to Tree");
981 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
982 tree->GetUserInfo()->Add(fn);
988 tree->SetBranchStatus("ESDfriend*",0);
989 // we want to have only one tree version number
990 tree->Write(tree->GetName(),TObject::kOverwrite);
994 TFile *aodFile = TFile::Open("AliAOD.root", "RECREATE");
995 ESDFile2AODFile(file, aodFile);
1000 CleanUp(file, fileOld);
1002 // Create tags for the events in the ESD tree (the ESD tree is always present)
1003 // In case of empty events the tags will contain dummy values
1004 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1005 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPList);
1007 AliAODTagCreator *aodtagCreator = new AliAODTagCreator();
1008 aodtagCreator->CreateAODTags(fFirstEvent,fLastEvent,fGRPList);
1011 //QA fone outside the event loop
1014 qas.Run(fRunLocalReconstruction.Data(), AliQA::kRECPOINTS) ;
1016 qas.Run(fRunTracking.Data(), AliQA::kESDS) ;
1020 // Cleanup of CDB manager: cache and active storages!
1021 AliCDBManager::Instance()->ClearCache();
1028 //_____________________________________________________________________________
1029 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
1031 // run the local reconstruction
1032 static Int_t eventNr=0;
1033 AliCodeTimerAuto("")
1035 // AliCDBManager* man = AliCDBManager::Instance();
1036 // Bool_t origCache = man->GetCacheFlag();
1038 // TString detStr = detectors;
1039 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1040 // if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1041 // AliReconstructor* reconstructor = GetReconstructor(iDet);
1042 // if (!reconstructor) continue;
1043 // if (reconstructor->HasLocalReconstruction()) continue;
1045 // AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1046 // AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1048 // AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1049 // AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1051 // man->SetCacheFlag(kTRUE);
1052 // TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
1053 // man->GetAll(calibPath); // entries are cached!
1055 // AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1057 // if (fRawReader) {
1058 // fRawReader->RewindEvents();
1059 // reconstructor->Reconstruct(fRunLoader, fRawReader);
1061 // reconstructor->Reconstruct(fRunLoader);
1064 // AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1065 // AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
1067 // // unload calibration data
1068 // man->UnloadFromCache(calibPath);
1069 // //man->ClearCache();
1072 // man->SetCacheFlag(origCache);
1074 // if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1075 // AliError(Form("the following detectors were not found: %s",
1077 // if (fStopOnError) return kFALSE;
1084 //_____________________________________________________________________________
1085 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1087 // run the local reconstruction
1088 static Int_t eventNr=0;
1089 AliCodeTimerAuto("")
1091 TString detStr = detectors;
1092 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1093 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1094 AliReconstructor* reconstructor = GetReconstructor(iDet);
1095 if (!reconstructor) continue;
1096 AliLoader* loader = fLoader[iDet];
1098 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1102 // conversion of digits
1103 if (fRawReader && reconstructor->HasDigitConversion()) {
1104 AliInfo(Form("converting raw data digits into root objects for %s",
1105 fgkDetectorName[iDet]));
1106 AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1107 fgkDetectorName[iDet]));
1108 loader->LoadDigits("update");
1109 loader->CleanDigits();
1110 loader->MakeDigitsContainer();
1111 TTree* digitsTree = loader->TreeD();
1112 reconstructor->ConvertDigits(fRawReader, digitsTree);
1113 loader->WriteDigits("OVERWRITE");
1114 loader->UnloadDigits();
1117 // local reconstruction
1118 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1119 AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1120 loader->LoadRecPoints("update");
1121 loader->CleanRecPoints();
1122 loader->MakeRecPointsContainer();
1123 TTree* clustersTree = loader->TreeR();
1124 if (fRawReader && !reconstructor->HasDigitConversion()) {
1125 reconstructor->Reconstruct(fRawReader, clustersTree);
1127 loader->LoadDigits("read");
1128 TTree* digitsTree = loader->TreeD();
1130 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1131 if (fStopOnError) return kFALSE;
1133 reconstructor->Reconstruct(digitsTree, clustersTree);
1135 loader->UnloadDigits();
1138 // do QA in the event loop if requested
1140 AliQADataMaker * qadm = GetQADataMaker(iDet);
1142 AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1143 AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1144 if (qadm->IsCycleDone() ) {
1145 qadm->EndOfCycle(AliQA::kRECPOINTS) ;
1146 qadm->EndOfCycle(AliQA::kESDS) ;
1147 qadm->StartOfCycle(AliQA::kRECPOINTS) ;
1148 qadm->StartOfCycle(AliQA::kESDS, "same") ;
1150 qadm->Exec(AliQA::kRECPOINTS, clustersTree) ;
1151 AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1155 loader->WriteRecPoints("OVERWRITE");
1156 loader->UnloadRecPoints();
1157 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1160 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1161 AliError(Form("the following detectors were not found: %s",
1163 if (fStopOnError) return kFALSE;
1169 //_____________________________________________________________________________
1170 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1172 // run the barrel tracking
1174 AliCodeTimerAuto("")
1176 AliESDVertex* vertex = NULL;
1177 Double_t vtxPos[3] = {0, 0, 0};
1178 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1179 TArrayF mcVertex(3);
1180 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1181 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1182 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1186 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1187 AliInfo("running the ITS vertex finder");
1188 if (fLoader[0]) fLoader[0]->LoadRecPoints();
1189 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
1190 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
1192 AliWarning("Vertex not found");
1193 vertex = new AliESDVertex();
1194 vertex->SetName("default");
1197 vertex->SetName("reconstructed");
1201 AliInfo("getting the primary vertex from MC");
1202 vertex = new AliESDVertex(vtxPos, vtxErr);
1206 vertex->GetXYZ(vtxPos);
1207 vertex->GetSigmaXYZ(vtxErr);
1209 AliWarning("no vertex reconstructed");
1210 vertex = new AliESDVertex(vtxPos, vtxErr);
1212 esd->SetVertex(vertex);
1213 // if SPD multiplicity has been determined, it is stored in the ESD
1214 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1215 if(mult)esd->SetMultiplicity(mult);
1217 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1218 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1225 //_____________________________________________________________________________
1226 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1228 // run the HLT barrel tracking
1230 AliCodeTimerAuto("")
1233 AliError("Missing runLoader!");
1237 AliInfo("running HLT tracking");
1239 // Get a pointer to the HLT reconstructor
1240 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1241 if (!reconstructor) return kFALSE;
1244 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1245 TString detName = fgkDetectorName[iDet];
1246 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1247 reconstructor->SetOption(detName.Data());
1248 AliTracker *tracker = reconstructor->CreateTracker();
1250 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1251 if (fStopOnError) return kFALSE;
1255 Double_t vtxErr[3]={0.005,0.005,0.010};
1256 const AliESDVertex *vertex = esd->GetVertex();
1257 vertex->GetXYZ(vtxPos);
1258 tracker->SetVertex(vtxPos,vtxErr);
1260 fLoader[iDet]->LoadRecPoints("read");
1261 TTree* tree = fLoader[iDet]->TreeR();
1263 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1266 tracker->LoadClusters(tree);
1268 if (tracker->Clusters2Tracks(esd) != 0) {
1269 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1273 tracker->UnloadClusters();
1281 //_____________________________________________________________________________
1282 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1284 // run the muon spectrometer tracking
1286 AliCodeTimerAuto("")
1289 AliError("Missing runLoader!");
1292 Int_t iDet = 7; // for MUON
1294 AliInfo("is running...");
1296 // Get a pointer to the MUON reconstructor
1297 AliReconstructor *reconstructor = GetReconstructor(iDet);
1298 if (!reconstructor) return kFALSE;
1301 TString detName = fgkDetectorName[iDet];
1302 AliDebug(1, Form("%s tracking", detName.Data()));
1303 AliTracker *tracker = reconstructor->CreateTracker();
1305 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1310 fLoader[iDet]->LoadRecPoints("read");
1312 tracker->LoadClusters(fLoader[iDet]->TreeR());
1314 Int_t rv = tracker->Clusters2Tracks(esd);
1318 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1322 fLoader[iDet]->UnloadRecPoints();
1324 tracker->UnloadClusters();
1332 //_____________________________________________________________________________
1333 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1335 // run the barrel tracking
1336 static Int_t eventNr=0;
1337 AliCodeTimerAuto("")
1339 AliInfo("running tracking");
1341 //Fill the ESD with the T0 info (will be used by the TOF)
1342 if (fReconstructor[11] && fLoader[11]) {
1343 fLoader[11]->LoadRecPoints("READ");
1344 TTree *treeR = fLoader[11]->TreeR();
1345 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
1348 // pass 1: TPC + ITS inwards
1349 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1350 if (!fTracker[iDet]) continue;
1351 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1354 fLoader[iDet]->LoadRecPoints("read");
1355 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
1356 TTree* tree = fLoader[iDet]->TreeR();
1358 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1361 fTracker[iDet]->LoadClusters(tree);
1362 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1364 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1365 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1368 if (fCheckPointLevel > 1) {
1369 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1371 // preliminary PID in TPC needed by the ITS tracker
1373 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1374 AliESDpid::MakePID(esd);
1376 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
1379 // pass 2: ALL backwards
1380 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1381 if (!fTracker[iDet]) continue;
1382 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1385 if (iDet > 1) { // all except ITS, TPC
1387 fLoader[iDet]->LoadRecPoints("read");
1388 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
1389 tree = fLoader[iDet]->TreeR();
1391 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1394 fTracker[iDet]->LoadClusters(tree);
1395 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1399 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1400 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1403 if (fCheckPointLevel > 1) {
1404 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1408 if (iDet > 2) { // all except ITS, TPC, TRD
1409 fTracker[iDet]->UnloadClusters();
1410 fLoader[iDet]->UnloadRecPoints();
1412 // updated PID in TPC needed by the ITS tracker -MI
1414 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1415 AliESDpid::MakePID(esd);
1417 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1420 // write space-points to the ESD in case alignment data output
1422 if (fWriteAlignmentData)
1423 WriteAlignmentData(esd);
1425 // pass 3: TRD + TPC + ITS refit inwards
1426 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1427 if (!fTracker[iDet]) continue;
1428 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1431 if (fTracker[iDet]->RefitInward(esd) != 0) {
1432 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1435 if (fCheckPointLevel > 1) {
1436 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1438 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1440 fTracker[iDet]->UnloadClusters();
1441 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
1442 fLoader[iDet]->UnloadRecPoints();
1443 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
1446 // Propagate track to the vertex - if not done by ITS
1448 Int_t ntracks = esd->GetNumberOfTracks();
1449 for (Int_t itrack=0; itrack<ntracks; itrack++){
1450 const Double_t kRadius = 3; // beam pipe radius
1451 const Double_t kMaxStep = 5; // max step
1452 const Double_t kMaxD = 123456; // max distance to prim vertex
1453 Double_t fieldZ = AliTracker::GetBz(); //
1454 AliESDtrack * track = esd->GetTrack(itrack);
1455 if (!track) continue;
1456 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1457 AliTracker::PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1458 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1464 //_____________________________________________________________________________
1465 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1467 // Remove the data which are not needed for the physics analysis.
1470 Int_t nTracks=esd->GetNumberOfTracks();
1471 Int_t nV0s=esd->GetNumberOfV0s();
1473 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
1475 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
1476 Bool_t rc=esd->Clean(cleanPars);
1478 nTracks=esd->GetNumberOfTracks();
1479 nV0s=esd->GetNumberOfV0s();
1481 (Form("Number of ESD tracks and V0s after cleaning %d",nTracks,nV0s));
1486 //_____________________________________________________________________________
1487 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1489 // fill the event summary data
1491 AliCodeTimerAuto("")
1492 static Int_t eventNr=0;
1493 TString detStr = detectors;
1495 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1496 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1497 AliReconstructor* reconstructor = GetReconstructor(iDet);
1498 if (!reconstructor) continue;
1499 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1500 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1501 TTree* clustersTree = NULL;
1502 if (fLoader[iDet]) {
1503 fLoader[iDet]->LoadRecPoints("read");
1504 clustersTree = fLoader[iDet]->TreeR();
1505 if (!clustersTree) {
1506 AliError(Form("Can't get the %s clusters tree",
1507 fgkDetectorName[iDet]));
1508 if (fStopOnError) return kFALSE;
1511 if (fRawReader && !reconstructor->HasDigitConversion()) {
1512 reconstructor->FillESD(fRawReader, clustersTree, esd);
1514 TTree* digitsTree = NULL;
1515 if (fLoader[iDet]) {
1516 fLoader[iDet]->LoadDigits("read");
1517 digitsTree = fLoader[iDet]->TreeD();
1519 AliError(Form("Can't get the %s digits tree",
1520 fgkDetectorName[iDet]));
1521 if (fStopOnError) return kFALSE;
1524 reconstructor->FillESD(digitsTree, clustersTree, esd);
1525 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1527 if (fLoader[iDet]) {
1528 fLoader[iDet]->UnloadRecPoints();
1531 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1535 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1536 AliError(Form("the following detectors were not found: %s",
1538 if (fStopOnError) return kFALSE;
1540 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
1545 //_____________________________________________________________________________
1546 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
1548 // Reads the trigger decision which is
1549 // stored in Trigger.root file and fills
1550 // the corresponding esd entries
1552 AliCodeTimerAuto("")
1554 AliInfo("Filling trigger information into the ESD");
1557 AliCTPRawStream input(fRawReader);
1558 if (!input.Next()) {
1559 AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1562 esd->SetTriggerMask(input.GetClassMask());
1563 esd->SetTriggerCluster(input.GetClusterMask());
1566 AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1568 if (!runloader->LoadTrigger()) {
1569 AliCentralTrigger *aCTP = runloader->GetTrigger();
1570 esd->SetTriggerMask(aCTP->GetClassMask());
1571 esd->SetTriggerCluster(aCTP->GetClusterMask());
1574 AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1579 AliError("No run loader is available! The trigger information is not stored in the ESD !");
1591 //_____________________________________________________________________________
1592 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
1595 // Filling information from RawReader Header
1598 AliInfo("Filling information from RawReader Header");
1599 esd->SetBunchCrossNumber(0);
1600 esd->SetOrbitNumber(0);
1601 esd->SetPeriodNumber(0);
1602 esd->SetTimeStamp(0);
1603 esd->SetEventType(0);
1604 const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1607 const UInt_t *id = eventHeader->GetP("Id");
1608 esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1609 esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1610 esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1612 esd->SetTimeStamp((eventHeader->Get("Timestamp")));
1613 esd->SetEventType((eventHeader->Get("Type")));
1620 //_____________________________________________________________________________
1621 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1623 // check whether detName is contained in detectors
1624 // if yes, it is removed from detectors
1626 // check if all detectors are selected
1627 if ((detectors.CompareTo("ALL") == 0) ||
1628 detectors.BeginsWith("ALL ") ||
1629 detectors.EndsWith(" ALL") ||
1630 detectors.Contains(" ALL ")) {
1635 // search for the given detector
1636 Bool_t result = kFALSE;
1637 if ((detectors.CompareTo(detName) == 0) ||
1638 detectors.BeginsWith(detName+" ") ||
1639 detectors.EndsWith(" "+detName) ||
1640 detectors.Contains(" "+detName+" ")) {
1641 detectors.ReplaceAll(detName, "");
1645 // clean up the detectors string
1646 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1647 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1648 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1653 //_____________________________________________________________________________
1654 Bool_t AliReconstruction::InitRunLoader()
1656 // get or create the run loader
1658 if (gAlice) delete gAlice;
1661 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1662 // load all base libraries to get the loader classes
1663 TString libs = gSystem->GetLibraries();
1664 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1665 TString detName = fgkDetectorName[iDet];
1666 if (detName == "HLT") continue;
1667 if (libs.Contains("lib" + detName + "base.so")) continue;
1668 gSystem->Load("lib" + detName + "base.so");
1670 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1672 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1676 fRunLoader->CdGAFile();
1677 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1678 if (fRunLoader->LoadgAlice() == 0) {
1679 gAlice = fRunLoader->GetAliRun();
1680 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
1683 if (!gAlice && !fRawReader) {
1684 AliError(Form("no gAlice object found in file %s",
1685 fGAliceFileName.Data()));
1690 //PH This is a temporary fix to give access to the kinematics
1691 //PH that is needed for the labels of ITS clusters
1692 fRunLoader->LoadHeader();
1693 fRunLoader->LoadKinematics();
1695 } else { // galice.root does not exist
1697 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1701 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1702 AliConfig::GetDefaultEventFolderName(),
1705 AliError(Form("could not create run loader in file %s",
1706 fGAliceFileName.Data()));
1710 fRunLoader->MakeTree("E");
1712 while (fRawReader->NextEvent()) {
1713 fRunLoader->SetEventNumber(iEvent);
1714 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1716 fRunLoader->MakeTree("H");
1717 fRunLoader->TreeE()->Fill();
1720 fRawReader->RewindEvents();
1721 if (fNumberOfEventsPerFile > 0)
1722 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
1724 fRunLoader->SetNumberOfEventsPerFile(iEvent);
1725 fRunLoader->WriteHeader("OVERWRITE");
1726 fRunLoader->CdGAFile();
1727 fRunLoader->Write(0, TObject::kOverwrite);
1728 // AliTracker::SetFieldMap(???);
1734 //_____________________________________________________________________________
1735 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1737 // get the reconstructor object and the loader for a detector
1739 if (fReconstructor[iDet]) return fReconstructor[iDet];
1741 // load the reconstructor object
1742 TPluginManager* pluginManager = gROOT->GetPluginManager();
1743 TString detName = fgkDetectorName[iDet];
1744 TString recName = "Ali" + detName + "Reconstructor";
1745 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1747 AliReconstructor* reconstructor = NULL;
1748 // first check if a plugin is defined for the reconstructor
1749 TPluginHandler* pluginHandler =
1750 pluginManager->FindHandler("AliReconstructor", detName);
1751 // if not, add a plugin for it
1752 if (!pluginHandler) {
1753 AliDebug(1, Form("defining plugin for %s", recName.Data()));
1754 TString libs = gSystem->GetLibraries();
1755 if (libs.Contains("lib" + detName + "base.so") ||
1756 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1757 pluginManager->AddHandler("AliReconstructor", detName,
1758 recName, detName + "rec", recName + "()");
1760 pluginManager->AddHandler("AliReconstructor", detName,
1761 recName, detName, recName + "()");
1763 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1765 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1766 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1768 if (reconstructor) {
1769 TObject* obj = fOptions.FindObject(detName.Data());
1770 if (obj) reconstructor->SetOption(obj->GetTitle());
1771 reconstructor->Init();
1772 fReconstructor[iDet] = reconstructor;
1775 // get or create the loader
1776 if (detName != "HLT") {
1777 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1778 if (!fLoader[iDet]) {
1779 AliConfig::Instance()
1780 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1782 // first check if a plugin is defined for the loader
1784 pluginManager->FindHandler("AliLoader", detName);
1785 // if not, add a plugin for it
1786 if (!pluginHandler) {
1787 TString loaderName = "Ali" + detName + "Loader";
1788 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1789 pluginManager->AddHandler("AliLoader", detName,
1790 loaderName, detName + "base",
1791 loaderName + "(const char*, TFolder*)");
1792 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1794 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1796 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1797 fRunLoader->GetEventFolder());
1799 if (!fLoader[iDet]) { // use default loader
1800 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1802 if (!fLoader[iDet]) {
1803 AliWarning(Form("couldn't get loader for %s", detName.Data()));
1804 if (fStopOnError) return NULL;
1806 fRunLoader->AddLoader(fLoader[iDet]);
1807 fRunLoader->CdGAFile();
1808 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1809 fRunLoader->Write(0, TObject::kOverwrite);
1814 return reconstructor;
1817 //_____________________________________________________________________________
1818 Bool_t AliReconstruction::CreateVertexer()
1820 // create the vertexer
1823 AliReconstructor* itsReconstructor = GetReconstructor(0);
1824 if (itsReconstructor) {
1825 fVertexer = itsReconstructor->CreateVertexer();
1828 AliWarning("couldn't create a vertexer for ITS");
1829 if (fStopOnError) return kFALSE;
1835 //_____________________________________________________________________________
1836 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1838 // create the trackers
1840 TString detStr = detectors;
1841 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1842 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1843 AliReconstructor* reconstructor = GetReconstructor(iDet);
1844 if (!reconstructor) continue;
1845 TString detName = fgkDetectorName[iDet];
1846 if (detName == "HLT") {
1847 fRunHLTTracking = kTRUE;
1850 if (detName == "MUON") {
1851 fRunMuonTracking = kTRUE;
1856 fTracker[iDet] = reconstructor->CreateTracker();
1857 if (!fTracker[iDet] && (iDet < 7)) {
1858 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1859 if (fStopOnError) return kFALSE;
1861 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
1867 //_____________________________________________________________________________
1868 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1870 // delete trackers and the run loader and close and delete the file
1872 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1873 delete fReconstructor[iDet];
1874 fReconstructor[iDet] = NULL;
1875 fLoader[iDet] = NULL;
1876 delete fTracker[iDet];
1877 fTracker[iDet] = NULL;
1878 // delete fQADataMaker[iDet];
1879 // fQADataMaker[iDet] = NULL;
1884 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
1885 delete fDiamondProfile;
1886 fDiamondProfile = NULL;
1905 gSystem->Unlink("AliESDs.old.root");
1909 //_____________________________________________________________________________
1911 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
1913 // read the ESD event from a file
1915 if (!esd) return kFALSE;
1917 sprintf(fileName, "ESD_%d.%d_%s.root",
1918 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1919 if (gSystem->AccessPathName(fileName)) return kFALSE;
1921 AliInfo(Form("reading ESD from file %s", fileName));
1922 AliDebug(1, Form("reading ESD from file %s", fileName));
1923 TFile* file = TFile::Open(fileName);
1924 if (!file || !file->IsOpen()) {
1925 AliError(Form("opening %s failed", fileName));
1932 esd = (AliESDEvent*) file->Get("ESD");
1941 //_____________________________________________________________________________
1942 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
1944 // write the ESD event to a file
1948 sprintf(fileName, "ESD_%d.%d_%s.root",
1949 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1951 AliDebug(1, Form("writing ESD to file %s", fileName));
1952 TFile* file = TFile::Open(fileName, "recreate");
1953 if (!file || !file->IsOpen()) {
1954 AliError(Form("opening %s failed", fileName));
1966 //_____________________________________________________________________________
1967 void AliReconstruction::ESDFile2AODFile(TFile* esdFile, TFile* aodFile)
1969 // write all files from the given esd file to an aod file
1971 // create an AliAOD object
1972 AliAODEvent *aod = new AliAODEvent();
1973 aod->CreateStdContent();
1979 TTree *aodTree = new TTree("aodTree", "AliAOD tree");
1980 aodTree->Branch(aod->GetList());
1983 TTree *t = (TTree*) esdFile->Get("esdTree");
1984 AliESDEvent *esd = new AliESDEvent();
1985 esd->ReadFromTree(t);
1987 Int_t nEvents = t->GetEntries();
1989 // set arrays and pointers
1999 // loop over events and fill them
2000 for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
2001 //cout << "event: " << iEvent << endl;
2002 t->GetEntry(iEvent);
2004 // Multiplicity information needed by the header (to be revised!)
2005 Int_t nTracks = esd->GetNumberOfTracks();
2006 Int_t nPosTracks = 0;
2007 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack)
2008 if (esd->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
2010 // Access the header
2011 AliAODHeader *header = aod->GetHeader();
2014 header->SetRunNumber (esd->GetRunNumber() );
2015 header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
2016 header->SetOrbitNumber (esd->GetOrbitNumber() );
2017 header->SetPeriodNumber (esd->GetPeriodNumber() );
2018 header->SetTriggerMask (esd->GetTriggerMask() );
2019 header->SetTriggerCluster (esd->GetTriggerCluster() );
2020 header->SetEventType (esd->GetEventType() );
2021 header->SetMagneticField (esd->GetMagneticField() );
2022 header->SetZDCN1Energy (esd->GetZDCN1Energy() );
2023 header->SetZDCP1Energy (esd->GetZDCP1Energy() );
2024 header->SetZDCN2Energy (esd->GetZDCN2Energy() );
2025 header->SetZDCP2Energy (esd->GetZDCP2Energy() );
2026 header->SetZDCEMEnergy (esd->GetZDCEMEnergy(0),esd->GetZDCEMEnergy(1));
2027 header->SetRefMultiplicity (nTracks);
2028 header->SetRefMultiplicityPos(nPosTracks);
2029 header->SetRefMultiplicityNeg(nTracks - nPosTracks);
2030 header->SetMuonMagFieldScale(-999.); // FIXME
2031 header->SetCentrality(-999.); // FIXME
2033 Int_t nV0s = esd->GetNumberOfV0s();
2034 Int_t nCascades = esd->GetNumberOfCascades();
2035 Int_t nKinks = esd->GetNumberOfKinks();
2036 Int_t nVertices = nV0s + 2*nCascades /*could lead to two vertices, one V0 and the Xi */+ nKinks + 1 /* = prim. vtx*/;
2038 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
2040 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
2042 aod->ResetStd(nTracks, nVertices, nV0s+nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);
2044 // Array to take into account the tracks already added to the AOD
2045 Bool_t * usedTrack = NULL;
2047 usedTrack = new Bool_t[nTracks];
2048 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
2050 // Array to take into account the V0s already added to the AOD
2051 Bool_t * usedV0 = NULL;
2053 usedV0 = new Bool_t[nV0s];
2054 for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
2056 // Array to take into account the kinks already added to the AOD
2057 Bool_t * usedKink = NULL;
2059 usedKink = new Bool_t[nKinks];
2060 for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
2063 // Access to the AOD container of vertices
2064 TClonesArray &vertices = *(aod->GetVertices());
2067 // Access to the AOD container of tracks
2068 TClonesArray &tracks = *(aod->GetTracks());
2071 // Access to the AOD container of V0s
2072 TClonesArray &V0s = *(aod->GetV0s());
2075 // Add primary vertex. The primary tracks will be defined
2076 // after the loops on the composite objects (V0, cascades, kinks)
2077 const AliESDVertex *vtx = esd->GetPrimaryVertex();
2079 vtx->GetXYZ(pos); // position
2080 vtx->GetCovMatrix(covVtx); //covariance matrix
2082 AliAODVertex * primary = new(vertices[jVertices++])
2083 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
2086 AliAODTrack *aodTrack = 0x0;
2088 // Create vertices starting from the most complex objects
2091 for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
2092 AliESDcascade *cascade = esd->GetCascade(nCascade);
2094 cascade->GetXYZ(pos[0], pos[1], pos[2]);
2095 cascade->GetPosCovXi(covVtx);
2097 // Add the cascade vertex
2098 AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
2100 cascade->GetChi2Xi(), // = chi2/NDF since NDF = 2*2-3
2103 AliAODVertex::kCascade);
2105 primary->AddDaughter(vcascade); // the cascade 'particle' (represented by a vertex) is added as a daughter to the primary vertex
2107 // Add the V0 from the cascade. The ESD class have to be optimized...
2108 // Now we have to search for the corresponding V0 in the list of V0s
2109 // using the indeces of the positive and negative tracks
2111 Int_t posFromV0 = cascade->GetPindex();
2112 Int_t negFromV0 = cascade->GetNindex();
2115 AliESDv0 * v0 = 0x0;
2118 for (Int_t iV0=0; iV0<nV0s; ++iV0) {
2120 v0 = esd->GetV0(iV0);
2121 Int_t posV0 = v0->GetPindex();
2122 Int_t negV0 = v0->GetNindex();
2124 if (posV0==posFromV0 && negV0==negFromV0) {
2130 AliAODVertex * vV0FromCascade = 0x0;
2132 if (indV0>-1 && !usedV0[indV0]) {
2134 // the V0 exists in the array of V0s and is not used
2136 usedV0[indV0] = kTRUE;
2138 v0->GetXYZ(pos[0], pos[1], pos[2]);
2139 v0->GetPosCov(covVtx);
2141 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2143 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2149 // the V0 doesn't exist in the array of V0s or was used
2150 cerr << "Error: event " << iEvent << " cascade " << nCascade
2151 << " The V0 " << indV0
2152 << " doesn't exist in the array of V0s or was used!" << endl;
2154 cascade->GetXYZ(pos[0], pos[1], pos[2]);
2155 cascade->GetPosCov(covVtx);
2157 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2159 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2163 vcascade->AddDaughter(vV0FromCascade);
2167 // Add the positive tracks from the V0
2169 if (! usedTrack[posFromV0]) {
2171 usedTrack[posFromV0] = kTRUE;
2173 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2174 esdTrack->GetPxPyPz(p_pos);
2175 esdTrack->GetXYZ(pos);
2176 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2177 esdTrack->GetESDpid(pid);
2179 vV0FromCascade->AddDaughter(aodTrack =
2180 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2181 esdTrack->GetLabel(),
2187 (Short_t)esdTrack->Charge(),
2188 esdTrack->GetITSClusterMap(),
2191 kTRUE, // check if this is right
2192 kFALSE, // check if this is right
2193 AliAODTrack::kSecondary)
2195 aodTrack->ConvertAliPIDtoAODPID();
2198 cerr << "Error: event " << iEvent << " cascade " << nCascade
2199 << " track " << posFromV0 << " has already been used!" << endl;
2202 // Add the negative tracks from the V0
2204 if (!usedTrack[negFromV0]) {
2206 usedTrack[negFromV0] = kTRUE;
2208 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2209 esdTrack->GetPxPyPz(p_neg);
2210 esdTrack->GetXYZ(pos);
2211 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2212 esdTrack->GetESDpid(pid);
2214 vV0FromCascade->AddDaughter(aodTrack =
2215 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2216 esdTrack->GetLabel(),
2222 (Short_t)esdTrack->Charge(),
2223 esdTrack->GetITSClusterMap(),
2226 kTRUE, // check if this is right
2227 kFALSE, // check if this is right
2228 AliAODTrack::kSecondary)
2230 aodTrack->ConvertAliPIDtoAODPID();
2233 cerr << "Error: event " << iEvent << " cascade " << nCascade
2234 << " track " << negFromV0 << " has already been used!" << endl;
2237 // add it to the V0 array as well
2238 Double_t d0[2] = { -999., -99.};
2239 // counting is probably wrong
2240 new(V0s[jV0s++]) AliAODv0(vV0FromCascade, -999., -99., p_pos, p_neg, d0); // to be refined
2242 // Add the bachelor track from the cascade
2244 Int_t bachelor = cascade->GetBindex();
2246 if(!usedTrack[bachelor]) {
2248 usedTrack[bachelor] = kTRUE;
2250 AliESDtrack *esdTrack = esd->GetTrack(bachelor);
2251 esdTrack->GetPxPyPz(p);
2252 esdTrack->GetXYZ(pos);
2253 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2254 esdTrack->GetESDpid(pid);
2256 vcascade->AddDaughter(aodTrack =
2257 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2258 esdTrack->GetLabel(),
2264 (Short_t)esdTrack->Charge(),
2265 esdTrack->GetITSClusterMap(),
2268 kTRUE, // check if this is right
2269 kFALSE, // check if this is right
2270 AliAODTrack::kSecondary)
2272 aodTrack->ConvertAliPIDtoAODPID();
2275 cerr << "Error: event " << iEvent << " cascade " << nCascade
2276 << " track " << bachelor << " has already been used!" << endl;
2279 // Add the primary track of the cascade (if any)
2281 } // end of the loop on cascades
2285 for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
2287 if (usedV0[nV0]) continue; // skip if aready added to the AOD
2289 AliESDv0 *v0 = esd->GetV0(nV0);
2291 v0->GetXYZ(pos[0], pos[1], pos[2]);
2292 v0->GetPosCov(covVtx);
2294 AliAODVertex * vV0 =
2295 new(vertices[jVertices++]) AliAODVertex(pos,
2297 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2301 primary->AddDaughter(vV0);
2303 Int_t posFromV0 = v0->GetPindex();
2304 Int_t negFromV0 = v0->GetNindex();
2306 // Add the positive tracks from the V0
2308 if (!usedTrack[posFromV0]) {
2310 usedTrack[posFromV0] = kTRUE;
2312 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2313 esdTrack->GetPxPyPz(p_pos);
2314 esdTrack->GetXYZ(pos);
2315 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2316 esdTrack->GetESDpid(pid);
2318 vV0->AddDaughter(aodTrack =
2319 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2320 esdTrack->GetLabel(),
2326 (Short_t)esdTrack->Charge(),
2327 esdTrack->GetITSClusterMap(),
2330 kTRUE, // check if this is right
2331 kFALSE, // check if this is right
2332 AliAODTrack::kSecondary)
2334 aodTrack->ConvertAliPIDtoAODPID();
2337 cerr << "Error: event " << iEvent << " V0 " << nV0
2338 << " track " << posFromV0 << " has already been used!" << endl;
2341 // Add the negative tracks from the V0
2343 if (!usedTrack[negFromV0]) {
2345 usedTrack[negFromV0] = kTRUE;
2347 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2348 esdTrack->GetPxPyPz(p_neg);
2349 esdTrack->GetXYZ(pos);
2350 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2351 esdTrack->GetESDpid(pid);
2353 vV0->AddDaughter(aodTrack =
2354 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2355 esdTrack->GetLabel(),
2361 (Short_t)esdTrack->Charge(),
2362 esdTrack->GetITSClusterMap(),
2365 kTRUE, // check if this is right
2366 kFALSE, // check if this is right
2367 AliAODTrack::kSecondary)
2369 aodTrack->ConvertAliPIDtoAODPID();
2372 cerr << "Error: event " << iEvent << " V0 " << nV0
2373 << " track " << negFromV0 << " has already been used!" << endl;
2376 // add it to the V0 array as well
2377 Double_t d0[2] = { 999., 99.};
2378 new(V0s[jV0s++]) AliAODv0(vV0, 999., 99., p_pos, p_neg, d0); // to be refined
2381 // end of the loop on V0s
2383 // Kinks: it is a big mess the access to the information in the kinks
2384 // The loop is on the tracks in order to find the mother and daugther of each kink
2387 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2389 AliESDtrack * esdTrack = esd->GetTrack(iTrack);
2391 Int_t ikink = esdTrack->GetKinkIndex(0);
2394 // Negative kink index: mother, positive: daughter
2396 // Search for the second track of the kink
2398 for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
2400 AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
2402 Int_t jkink = esdTrack1->GetKinkIndex(0);
2404 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
2406 // The two tracks are from the same kink
2408 if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
2411 Int_t idaughter = -1;
2413 if (ikink<0 && jkink>0) {
2418 else if (ikink>0 && jkink<0) {
2424 cerr << "Error: Wrong combination of kink indexes: "
2425 << ikink << " " << jkink << endl;
2429 // Add the mother track
2431 AliAODTrack * mother = NULL;
2433 if (!usedTrack[imother]) {
2435 usedTrack[imother] = kTRUE;
2437 AliESDtrack *esdTrack = esd->GetTrack(imother);
2438 esdTrack->GetPxPyPz(p);
2439 esdTrack->GetXYZ(pos);
2440 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2441 esdTrack->GetESDpid(pid);
2444 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2445 esdTrack->GetLabel(),
2451 (Short_t)esdTrack->Charge(),
2452 esdTrack->GetITSClusterMap(),
2455 kTRUE, // check if this is right
2456 kTRUE, // check if this is right
2457 AliAODTrack::kPrimary);
2458 primary->AddDaughter(mother);
2459 mother->ConvertAliPIDtoAODPID();
2462 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2463 << " track " << imother << " has already been used!" << endl;
2466 // Add the kink vertex
2467 AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
2469 AliAODVertex * vkink =
2470 new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
2474 esdTrack->GetID(), // This is the track ID of the mother's track!
2475 AliAODVertex::kKink);
2476 // Add the daughter track
2478 AliAODTrack * daughter = NULL;
2480 if (!usedTrack[idaughter]) {
2482 usedTrack[idaughter] = kTRUE;
2484 AliESDtrack *esdTrack = esd->GetTrack(idaughter);
2485 esdTrack->GetPxPyPz(p);
2486 esdTrack->GetXYZ(pos);
2487 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2488 esdTrack->GetESDpid(pid);
2491 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2492 esdTrack->GetLabel(),
2498 (Short_t)esdTrack->Charge(),
2499 esdTrack->GetITSClusterMap(),
2502 kTRUE, // check if this is right
2503 kTRUE, // check if this is right
2504 AliAODTrack::kPrimary);
2505 vkink->AddDaughter(daughter);
2506 daughter->ConvertAliPIDtoAODPID();
2509 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2510 << " track " << idaughter << " has already been used!" << endl;
2516 vertices.Expand(jVertices);
2518 // Tracks (primary and orphan)
2519 for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
2521 if (usedTrack[nTrack]) continue;
2523 AliESDtrack *esdTrack = esd->GetTrack(nTrack);
2524 esdTrack->GetPxPyPz(p);
2525 esdTrack->GetXYZ(pos);
2526 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2527 esdTrack->GetESDpid(pid);
2529 Float_t impactXY, impactZ;
2531 esdTrack->GetImpactParameters(impactXY,impactZ);
2534 // track inside the beam pipe
2536 primary->AddDaughter(aodTrack =
2537 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2538 esdTrack->GetLabel(),
2544 (Short_t)esdTrack->Charge(),
2545 esdTrack->GetITSClusterMap(),
2548 kTRUE, // check if this is right
2549 kTRUE, // check if this is right
2550 AliAODTrack::kPrimary)
2552 aodTrack->ConvertAliPIDtoAODPID();
2555 // outside the beam pipe: orphan track
2556 // Don't write them anymore!
2559 } // end of loop on tracks
2562 Int_t nMuTracks = esd->GetNumberOfMuonTracks();
2563 for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
2565 AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack);
2566 p[0] = esdMuTrack->Px();
2567 p[1] = esdMuTrack->Py();
2568 p[2] = esdMuTrack->Pz();
2569 pos[0] = primary->GetX();
2570 pos[1] = primary->GetY();
2571 pos[2] = primary->GetZ();
2573 // has to be changed once the muon pid is provided by the ESD
2574 for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
2576 primary->AddDaughter(aodTrack =
2577 new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
2578 0, // no label provided
2583 NULL, // no covariance matrix provided
2584 esdMuTrack->Charge(),
2585 0, // ITSClusterMap is set below
2588 kFALSE, // muon tracks are not used to fit the primary vtx
2589 kFALSE, // not used for vertex fit
2590 AliAODTrack::kPrimary)
2593 aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
2594 Int_t track2Trigger = esdMuTrack->GetMatchTrigger();
2595 aodTrack->SetMatchTrigger(track2Trigger);
2597 aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
2599 aodTrack->SetChi2MatchTrigger(0.);
2601 tracks.Expand(jTracks); // remove 'empty slots' due to unwritten tracks
2603 // Access to the AOD container of PMD clusters
2604 TClonesArray &pmdClusters = *(aod->GetPmdClusters());
2605 Int_t jPmdClusters=0;
2607 for (Int_t iPmd = 0; iPmd < nPmdClus; ++iPmd) {
2608 // file pmd clusters, to be revised!
2609 AliESDPmdTrack *pmdTrack = esd->GetPmdTrack(iPmd);
2612 Double_t pos[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ() };
2613 Double_t pid[9] = { 0., 0., 0., 0., 0., 0., 0., 0., 0. }; // to be revised!
2615 // assoc cluster not set
2616 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), pos, pid);
2619 // Access to the AOD container of clusters
2620 TClonesArray &caloClusters = *(aod->GetCaloClusters());
2623 for (Int_t iClust=0; iClust<nCaloClus; ++iClust) {
2625 AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
2627 Int_t id = cluster->GetID();
2630 Float_t energy = cluster->E();
2631 cluster->GetPosition(posF);
2632 Char_t ttype=AliAODCluster::kUndef;
2634 if (cluster->GetClusterType() == AliESDCaloCluster::kPHOSCluster) {
2635 ttype=AliAODCluster::kPHOSNeutral;
2637 else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALClusterv1) {
2638 ttype = AliAODCluster::kEMCALClusterv1;
2642 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
2650 caloCluster->SetCaloCluster(); // to be refined!
2653 caloClusters.Expand(jClusters); // resize TObjArray to 'remove' slots for pseudo clusters
2654 // end of loop on calo clusters
2656 // fill EMCAL cell info
2657 if (esd->GetEMCALCells()) { // protection against missing ESD information
2658 AliESDCaloCells &esdEMcells = *(esd->GetEMCALCells());
2659 Int_t nEMcell = esdEMcells.GetNumberOfCells() ;
2661 AliAODCaloCells &aodEMcells = *(aod->GetEMCALCells());
2662 aodEMcells.CreateContainer(nEMcell);
2663 aodEMcells.SetType(AliAODCaloCells::kEMCAL);
2664 for (Int_t iCell = 0; iCell < nEMcell; iCell++) {
2665 aodEMcells.SetCell(iCell,esdEMcells.GetCellNumber(iCell),esdEMcells.GetAmplitude(iCell));
2670 // fill PHOS cell info
2671 if (esd->GetPHOSCells()) { // protection against missing ESD information
2672 AliESDCaloCells &esdPHcells = *(esd->GetPHOSCells());
2673 Int_t nPHcell = esdPHcells.GetNumberOfCells() ;
2675 AliAODCaloCells &aodPHcells = *(aod->GetPHOSCells());
2676 aodPHcells.CreateContainer(nPHcell);
2677 aodPHcells.SetType(AliAODCaloCells::kPHOS);
2678 for (Int_t iCell = 0; iCell < nPHcell; iCell++) {
2679 aodPHcells.SetCell(iCell,esdPHcells.GetCellNumber(iCell),esdPHcells.GetAmplitude(iCell));
2685 AliAODTracklets &SPDTracklets = *(aod->GetTracklets());
2686 const AliMultiplicity *mult = esd->GetMultiplicity();
2688 if (mult->GetNumberOfTracklets()>0) {
2689 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
2691 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
2692 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n));
2696 Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
2699 delete [] usedTrack;
2703 // fill the tree for this event
2705 } // end of event loop
2707 aodTree->GetUserInfo()->Add(aod);
2709 // write the tree to the specified file
2710 aodFile = aodTree->GetCurrentFile();
2717 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2719 // Write space-points which are then used in the alignment procedures
2720 // For the moment only ITS, TRD and TPC
2722 // Load TOF clusters
2724 fLoader[3]->LoadRecPoints("read");
2725 TTree* tree = fLoader[3]->TreeR();
2727 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2730 fTracker[3]->LoadClusters(tree);
2732 Int_t ntracks = esd->GetNumberOfTracks();
2733 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2735 AliESDtrack *track = esd->GetTrack(itrack);
2738 for (Int_t iDet = 3; iDet >= 0; iDet--)
2739 nsp += track->GetNcls(iDet);
2741 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2742 track->SetTrackPointArray(sp);
2744 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2745 AliTracker *tracker = fTracker[iDet];
2746 if (!tracker) continue;
2747 Int_t nspdet = track->GetNcls(iDet);
2748 if (nspdet <= 0) continue;
2749 track->GetClusters(iDet,idx);
2753 while (isp < nspdet) {
2755 TString dets = fgkDetectorName[iDet];
2756 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2757 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2758 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2759 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2760 isvalid = tracker->GetTrackPointTrackingError(idx[isp2],p,track);
2762 isvalid = tracker->GetTrackPoint(idx[isp2],p);
2765 const Int_t kNTPCmax = 159;
2766 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
2767 if (!isvalid) continue;
2768 sp->AddPoint(isptrack,&p); isptrack++; isp++;
2774 fTracker[3]->UnloadClusters();
2775 fLoader[3]->UnloadRecPoints();
2779 //_____________________________________________________________________________
2780 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2782 // The method reads the raw-data error log
2783 // accumulated within the rawReader.
2784 // It extracts the raw-data errors related to
2785 // the current event and stores them into
2786 // a TClonesArray inside the esd object.
2788 if (!fRawReader) return;
2790 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2792 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2794 if (iEvent != log->GetEventNumber()) continue;
2796 esd->AddRawDataErrorLog(log);
2801 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString fName){
2802 // Dump a file content into a char in TNamed
2804 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2805 Int_t kBytes = (Int_t)in.tellg();
2806 printf("Size: %d \n",kBytes);
2809 char* memblock = new char [kBytes];
2810 in.seekg (0, ios::beg);
2811 in.read (memblock, kBytes);
2813 TString fData(memblock,kBytes);
2814 fn = new TNamed(fName,fData);
2815 printf("fData Size: %d \n",fData.Sizeof());
2816 printf("fName Size: %d \n",fName.Sizeof());
2817 printf("fn Size: %d \n",fn->Sizeof());
2821 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2827 void AliReconstruction::TNamedToFile(TTree* fTree, TString fName){
2828 // This is not really needed in AliReconstruction at the moment
2829 // but can serve as a template
2831 TList *fList = fTree->GetUserInfo();
2832 TNamed *fn = (TNamed*)fList->FindObject(fName.Data());
2833 printf("fn Size: %d \n",fn->Sizeof());
2835 TString fTmp(fn->GetName()); // to be 100% sure in principle fName also works
2836 const char* cdata = fn->GetTitle();
2837 printf("fTmp Size %d\n",fTmp.Sizeof());
2839 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2840 printf("calculated size %d\n",size);
2841 ofstream out(fName.Data(),ios::out | ios::binary);
2842 out.write(cdata,size);
2847 //_____________________________________________________________________________
2848 AliQADataMaker * AliReconstruction::GetQADataMaker(Int_t iDet)
2850 // get the quality assurance data maker object and the loader for a detector
2852 if (fQADataMaker[iDet])
2853 return fQADataMaker[iDet];
2855 // load the QA data maker object
2856 TPluginManager* pluginManager = gROOT->GetPluginManager();
2857 TString detName = fgkDetectorName[iDet];
2858 TString qadmName = "Ali" + detName + "QADataMaker";
2859 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT"))
2862 AliQADataMaker * qadm = NULL;
2863 // first check if a plugin is defined for the quality assurance data maker
2864 TPluginHandler* pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2865 // if not, add a plugin for it
2866 if (!pluginHandler) {
2867 AliDebug(1, Form("defining plugin for %s", qadmName.Data()));
2868 TString libs = gSystem->GetLibraries();
2869 if (libs.Contains("lib" + detName + "base.so") ||
2870 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2871 pluginManager->AddHandler("AliQADataMaker", detName,
2872 qadmName, detName + "qadm", qadmName + "()");
2874 pluginManager->AddHandler("AliQADataMaker", detName,
2875 qadmName, detName, qadmName + "()");
2877 pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2879 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2880 qadm = (AliQADataMaker *) pluginHandler->ExecPlugin(0);
2883 AliInfo(Form("Initializing quality assurance data maker for %s", fgkDetectorName[iDet]));
2884 qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun(), GetQACycles(fgkDetectorName[iDet]));
2885 qadm->StartOfCycle(AliQA::kRECPOINTS);
2886 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
2887 qadm->StartOfCycle(AliQA::kESDS, "same") ;
2888 fQADataMaker[iDet] = qadm;
2894 //_____________________________________________________________________________
2895 Bool_t AliReconstruction::RunQA(const char* detectors, AliESDEvent *& esd)
2897 // run the Quality Assurance data producer
2899 AliCodeTimerAuto("")
2900 TString detStr = detectors;
2901 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2902 if (!IsSelected(fgkDetectorName[iDet], detStr))
2904 AliQADataMaker * qadm = GetQADataMaker(iDet);
2907 AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2908 AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2910 qadm->Exec(AliQA::kESDS, esd) ;
2913 AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2915 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2916 AliError(Form("the following detectors were not found: %s",
2926 //_____________________________________________________________________________
2927 void AliReconstruction::CheckQA()
2929 // check the QA of SIM for this run and remove the detectors
2930 // with status Fatal
2932 TString newRunLocalReconstruction ;
2933 TString newRunTracking ;
2934 TString newFillESD ;
2936 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2937 TString detName(AliQA::GetDetName(iDet)) ;
2938 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX(iDet)) ;
2939 if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2940 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2942 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2943 fRunLocalReconstruction.Contains("ALL") ) {
2944 newRunLocalReconstruction += detName ;
2945 newRunLocalReconstruction += " " ;
2947 if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) ||
2948 fRunTracking.Contains("ALL") ) {
2949 newRunTracking += detName ;
2950 newRunTracking += " " ;
2952 if ( fFillESD.Contains(AliQA::GetDetName(iDet)) ||
2953 fFillESD.Contains("ALL") ) {
2954 newFillESD += detName ;
2959 fRunLocalReconstruction = newRunLocalReconstruction ;
2960 fRunTracking = newRunTracking ;
2961 fFillESD = newFillESD ;
2964 //_____________________________________________________________________________
2965 Int_t AliReconstruction::GetDetIndex(const char* detector)
2967 // return the detector index corresponding to detector
2969 for (index = 0; index < fgkNDetectors ; index++) {
2970 if ( strcmp(detector, fgkDetectorName[index]) == 0 )