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 ///////////////////////////////////////////////////////////////////////////////
121 #include <TPluginManager.h>
122 #include <TGeoManager.h>
123 #include <TLorentzVector.h>
127 #include "AliReconstruction.h"
128 #include "AliCodeTimer.h"
129 #include "AliReconstructor.h"
131 #include "AliRunLoader.h"
133 #include "AliRawReaderFile.h"
134 #include "AliRawReaderDate.h"
135 #include "AliRawReaderRoot.h"
136 #include "AliRawEventHeaderBase.h"
137 #include "AliESDEvent.h"
138 #include "AliESDMuonTrack.h"
139 #include "AliESDfriend.h"
140 #include "AliESDVertex.h"
141 #include "AliESDcascade.h"
142 #include "AliESDkink.h"
143 #include "AliESDtrack.h"
144 #include "AliESDCaloCluster.h"
145 #include "AliMultiplicity.h"
146 #include "AliTracker.h"
147 #include "AliVertexer.h"
148 #include "AliVertexerTracks.h"
149 #include "AliV0vertexer.h"
150 #include "AliCascadeVertexer.h"
151 #include "AliHeader.h"
152 #include "AliGenEventHeader.h"
154 #include "AliESDpid.h"
155 #include "AliESDtrack.h"
156 #include "AliESDPmdTrack.h"
158 #include "AliESDTagCreator.h"
159 #include "AliAODTagCreator.h"
161 #include "AliGeomManager.h"
162 #include "AliTrackPointArray.h"
163 #include "AliCDBManager.h"
164 #include "AliCDBEntry.h"
165 #include "AliAlignObj.h"
167 #include "AliCentralTrigger.h"
168 #include "AliCTPRawStream.h"
170 #include "AliAODEvent.h"
171 #include "AliAODHeader.h"
172 #include "AliAODTrack.h"
173 #include "AliAODVertex.h"
174 #include "AliAODv0.h"
175 #include "AliAODJet.h"
176 #include "AliAODCaloCells.h"
177 #include "AliAODCaloCluster.h"
178 #include "AliAODPmdCluster.h"
179 #include "AliAODFmdCluster.h"
180 #include "AliAODTracklets.h"
182 #include "AliQADataMaker.h"
184 #include "AliSysInfo.h" // memory snapshots
187 ClassImp(AliReconstruction)
190 //_____________________________________________________________________________
191 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
193 //_____________________________________________________________________________
194 AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
195 const char* name, const char* title) :
198 fUniformField(kTRUE),
199 fRunVertexFinder(kTRUE),
200 fRunHLTTracking(kFALSE),
201 fRunMuonTracking(kFALSE),
203 fRunCascadeFinder(kTRUE),
204 fStopOnError(kFALSE),
205 fWriteAlignmentData(kFALSE),
206 fWriteESDfriend(kFALSE),
208 fFillTriggerESD(kTRUE),
214 fRunLocalReconstruction("ALL"),
217 fUseTrackingErrorsForAlignment(""),
218 fGAliceFileName(gAliceFilename),
223 fNumberOfEventsPerFile(1),
226 fLoadAlignFromCDB(kTRUE),
227 fLoadAlignData("ALL"),
234 fDiamondProfile(NULL),
236 fAlignObjArray(NULL),
241 // create reconstruction object with default parameters
243 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
244 fReconstructor[iDet] = NULL;
245 fLoader[iDet] = NULL;
246 fTracker[iDet] = NULL;
247 fQADataMaker[iDet] = NULL;
248 fQACycles[iDet] = 999999;
253 //_____________________________________________________________________________
254 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
257 fUniformField(rec.fUniformField),
258 fRunVertexFinder(rec.fRunVertexFinder),
259 fRunHLTTracking(rec.fRunHLTTracking),
260 fRunMuonTracking(rec.fRunMuonTracking),
261 fRunV0Finder(rec.fRunV0Finder),
262 fRunCascadeFinder(rec.fRunCascadeFinder),
263 fStopOnError(rec.fStopOnError),
264 fWriteAlignmentData(rec.fWriteAlignmentData),
265 fWriteESDfriend(rec.fWriteESDfriend),
266 fWriteAOD(rec.fWriteAOD),
267 fFillTriggerESD(rec.fFillTriggerESD),
269 fCleanESD(rec.fCleanESD),
273 fRunLocalReconstruction(rec.fRunLocalReconstruction),
274 fRunTracking(rec.fRunTracking),
275 fFillESD(rec.fFillESD),
276 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
277 fGAliceFileName(rec.fGAliceFileName),
279 fEquipIdMap(rec.fEquipIdMap),
280 fFirstEvent(rec.fFirstEvent),
281 fLastEvent(rec.fLastEvent),
282 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
285 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
286 fLoadAlignData(rec.fLoadAlignData),
287 fESDPar(rec.fESDPar),
293 fDiamondProfile(NULL),
295 fAlignObjArray(rec.fAlignObjArray),
296 fCDBUri(rec.fCDBUri),
297 fRemoteCDBUri(rec.fRemoteCDBUri),
302 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
303 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
305 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
306 fReconstructor[iDet] = NULL;
307 fLoader[iDet] = NULL;
308 fTracker[iDet] = NULL;
309 fQADataMaker[iDet] = NULL;
310 fQACycles[iDet] = rec.fQACycles[iDet];
312 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
313 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
317 //_____________________________________________________________________________
318 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
320 // assignment operator
322 this->~AliReconstruction();
323 new(this) AliReconstruction(rec);
327 //_____________________________________________________________________________
328 AliReconstruction::~AliReconstruction()
334 fSpecCDBUri.Delete();
336 AliCodeTimer::Instance()->Print();
339 //_____________________________________________________________________________
340 void AliReconstruction::InitCDBStorage()
342 // activate a default CDB storage
343 // First check if we have any CDB storage set, because it is used
344 // to retrieve the calibration and alignment constants
346 AliCDBManager* man = AliCDBManager::Instance();
347 if (man->IsDefaultStorageSet())
349 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
350 AliWarning("Default CDB storage has been already set !");
351 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
352 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
356 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
357 AliDebug(2, Form("Default CDB storage is set to: %s",fCDBUri.Data()));
358 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
359 man->SetDefaultStorage(fCDBUri);
362 // Remote storage (the Grid storage) is used if it is activated
363 // and if the object is not found in the default storage
365 // if (man->IsRemoteStorageSet())
367 // AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
368 // AliWarning("Remote CDB storage has been already set !");
369 // AliWarning(Form("Ignoring the remote storage declared in AliReconstruction: %s",fRemoteCDBUri.Data()));
370 // AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
371 // fRemoteCDBUri = "";
374 // AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
375 // AliDebug(2, Form("Remote CDB storage is set to: %s",fRemoteCDBUri.Data()));
376 // AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
377 // man->SetRemoteStorage(fRemoteCDBUri);
380 // Now activate the detector specific CDB storage locations
381 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
382 TObject* obj = fSpecCDBUri[i];
384 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
385 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
386 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
387 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
392 //_____________________________________________________________________________
393 void AliReconstruction::SetDefaultStorage(const char* uri) {
394 // Store the desired default CDB storage location
395 // Activate it later within the Run() method
401 //_____________________________________________________________________________
402 void AliReconstruction::SetRemoteStorage(const char* uri) {
403 // Store the desired remote CDB storage location
404 // Activate it later within the Run() method
405 // Remote storage (the Grid storage) is used if it is activated
406 // and if the object is not found in the default storage
412 //_____________________________________________________________________________
413 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
414 // Store a detector-specific CDB storage location
415 // Activate it later within the Run() method
417 AliCDBPath aPath(calibType);
418 if(!aPath.IsValid()){
419 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
420 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
421 if(!strcmp(calibType, fgkDetectorName[iDet])) {
422 aPath.SetPath(Form("%s/*", calibType));
423 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
427 if(!aPath.IsValid()){
428 AliError(Form("Not a valid path or detector: %s", calibType));
433 // // check that calibType refers to a "valid" detector name
434 // Bool_t isDetector = kFALSE;
435 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
436 // TString detName = fgkDetectorName[iDet];
437 // if(aPath.GetLevel0() == detName) {
438 // isDetector = kTRUE;
444 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
448 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
449 if (obj) fSpecCDBUri.Remove(obj);
450 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
457 //_____________________________________________________________________________
458 Bool_t AliReconstruction::SetRunNumber()
460 // The method is called in Run() in order
461 // to set a correct run number.
462 // In case of raw data reconstruction the
463 // run number is taken from the raw data header
465 if(AliCDBManager::Instance()->GetRun() < 0) {
467 AliError("No run loader is found !");
470 // read run number from gAlice
471 if(fRunLoader->GetAliRun())
472 AliCDBManager::Instance()->SetRun(fRunLoader->GetAliRun()->GetRunNumber());
475 if(fRawReader->NextEvent()) {
476 AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
477 fRawReader->RewindEvents();
480 AliError("No raw-data events found !");
485 AliError("Neither gAlice nor RawReader objects are found !");
489 AliInfo(Form("CDB Run number: %d",AliCDBManager::Instance()->GetRun()));
494 //_____________________________________________________________________________
495 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
497 // Read the alignment objects from CDB.
498 // Each detector is supposed to have the
499 // alignment objects in DET/Align/Data CDB path.
500 // All the detector objects are then collected,
501 // sorted by geometry level (starting from ALIC) and
502 // then applied to the TGeo geometry.
503 // Finally an overlaps check is performed.
505 // Load alignment data from CDB and fill fAlignObjArray
506 if(fLoadAlignFromCDB){
508 TString detStr = detectors;
509 TString loadAlObjsListOfDets = "";
511 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
512 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
513 loadAlObjsListOfDets += fgkDetectorName[iDet];
514 loadAlObjsListOfDets += " ";
515 } // end loop over detectors
516 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
517 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
519 // Check if the array with alignment objects was
520 // provided by the user. If yes, apply the objects
521 // to the present TGeo geometry
522 if (fAlignObjArray) {
523 if (gGeoManager && gGeoManager->IsClosed()) {
524 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
525 AliError("The misalignment of one or more volumes failed!"
526 "Compare the list of simulated detectors and the list of detector alignment data!");
531 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
537 delete fAlignObjArray; fAlignObjArray=0;
542 //_____________________________________________________________________________
543 void AliReconstruction::SetGAliceFile(const char* fileName)
545 // set the name of the galice file
547 fGAliceFileName = fileName;
550 //_____________________________________________________________________________
551 void AliReconstruction::SetOption(const char* detector, const char* option)
553 // set options for the reconstruction of a detector
555 TObject* obj = fOptions.FindObject(detector);
556 if (obj) fOptions.Remove(obj);
557 fOptions.Add(new TNamed(detector, option));
561 //_____________________________________________________________________________
562 Bool_t AliReconstruction::Run(const char* input)
564 // run the reconstruction
569 if (!input) input = fInput.Data();
570 TString fileName(input);
571 if (fileName.EndsWith("/")) {
572 fRawReader = new AliRawReaderFile(fileName);
573 } else if (fileName.EndsWith(".root")) {
574 fRawReader = new AliRawReaderRoot(fileName);
575 } else if (!fileName.IsNull()) {
576 fRawReader = new AliRawReaderDate(fileName);
577 fRawReader->SelectEvents(7);
579 if (!fEquipIdMap.IsNull() && fRawReader)
580 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
582 AliSysInfo::AddStamp("Start");
583 // get the run loader
584 if (!InitRunLoader()) return kFALSE;
585 AliSysInfo::AddStamp("LoadLoader");
587 // Initialize the CDB storage
589 AliSysInfo::AddStamp("LoadCDB");
591 // Set run number in CDBManager (if it is not already set by the user)
592 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
594 // Import ideal TGeo geometry and apply misalignment
596 TString geom(gSystem->DirName(fGAliceFileName));
597 geom += "/geometry.root";
598 AliGeomManager::LoadGeometry(geom.Data());
599 if (!gGeoManager) if (fStopOnError) return kFALSE;
602 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
603 AliSysInfo::AddStamp("LoadGeom");
605 // local reconstruction
606 if (!fRunLocalReconstruction.IsNull()) {
607 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
608 if (fStopOnError) {CleanUp(); return kFALSE;}
611 // if (!fRunVertexFinder && fRunTracking.IsNull() &&
612 // fFillESD.IsNull()) return kTRUE;
615 if (fRunVertexFinder && !CreateVertexer()) {
621 AliSysInfo::AddStamp("Vertexer");
624 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
630 AliSysInfo::AddStamp("LoadTrackers");
632 // get the possibly already existing ESD file and tree
633 AliESDEvent* esd = new AliESDEvent(); AliESDEvent* hltesd = new AliESDEvent();
634 TFile* fileOld = NULL;
635 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
636 if (!gSystem->AccessPathName("AliESDs.root")){
637 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
638 fileOld = TFile::Open("AliESDs.old.root");
639 if (fileOld && fileOld->IsOpen()) {
640 treeOld = (TTree*) fileOld->Get("esdTree");
641 if (treeOld)esd->ReadFromTree(treeOld);
642 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
643 if (hlttreeOld) hltesd->ReadFromTree(hlttreeOld);
647 // create the ESD output file and tree
648 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
649 file->SetCompressionLevel(2);
650 if (!file->IsOpen()) {
651 AliError("opening AliESDs.root failed");
652 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
655 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
656 esd = new AliESDEvent();
657 esd->CreateStdContent();
658 esd->WriteToTree(tree);
660 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
661 hltesd = new AliESDEvent();
662 hltesd->CreateStdContent();
663 hltesd->WriteToTree(hlttree);
666 delete esd; delete hltesd;
667 esd = NULL; hltesd = NULL;
669 // create the branch with ESD additions
673 AliESDfriend *esdf = 0;
674 if (fWriteESDfriend) {
675 esdf = new AliESDfriend();
676 TBranch *br=tree->Branch("ESDfriend.","AliESDfriend", &esdf);
677 br->SetFile("AliESDfriends.root");
678 esd->AddObject(esdf);
682 // Get the diamond profile from OCDB
683 AliCDBEntry* entry = AliCDBManager::Instance()
684 ->Get("GRP/Calib/MeanVertex");
687 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
689 AliError("No diamond profile found in OCDB!");
692 AliVertexerTracks tVertexer(AliTracker::GetBz());
693 if(fDiamondProfile) tVertexer.SetVtxStart(fDiamondProfile);
697 if (fRawReader) fRawReader->RewindEvents();
698 TString detStr(fFillESD) ;
701 gSystem->GetProcInfo(&ProcInfo);
702 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
704 // checking the QA of previous steps
707 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
708 if (fRawReader) fRawReader->NextEvent();
709 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
710 // copy old ESD to the new one
712 esd->ReadFromTree(treeOld);
713 treeOld->GetEntry(iEvent);
717 esd->ReadFromTree(hlttreeOld);
718 hlttreeOld->GetEntry(iEvent);
724 AliInfo(Form("processing event %d", iEvent));
725 fRunLoader->GetEvent(iEvent);
728 sprintf(aFileName, "ESD_%d.%d_final.root",
729 fRunLoader->GetHeader()->GetRun(),
730 fRunLoader->GetHeader()->GetEventNrInRun());
731 if (!gSystem->AccessPathName(aFileName)) continue;
733 // local reconstruction
734 if (!fRunLocalReconstruction.IsNull()) {
735 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
736 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
741 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
742 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
743 esd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
744 hltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
746 // Set magnetic field from the tracker
747 esd->SetMagneticField(AliTracker::GetBz());
748 hltesd->SetMagneticField(AliTracker::GetBz());
752 // Fill raw-data error log into the ESD
753 if (fRawReader) FillRawDataErrorLog(iEvent,esd);
756 if (fRunVertexFinder) {
757 if (!ReadESD(esd, "vertex")) {
758 if (!RunVertexFinder(esd)) {
759 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
761 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
766 if (!fRunTracking.IsNull()) {
767 if (fRunHLTTracking) {
768 hltesd->SetVertex(esd->GetVertex());
769 if (!RunHLTTracking(hltesd)) {
770 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
776 if (!fRunTracking.IsNull()) {
777 if (fRunMuonTracking) {
778 if (!RunMuonTracking(esd)) {
779 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
785 if (!fRunTracking.IsNull()) {
786 if (!ReadESD(esd, "tracking")) {
787 if (!RunTracking(esd)) {
788 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
790 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
795 if (!fFillESD.IsNull()) {
796 if (!FillESD(esd, fFillESD)) {
797 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
801 if (!fFillESD.IsNull())
802 RunQA(fFillESD.Data(), esd);
804 // fill Event header information from the RawEventHeader
805 if (fRawReader){FillRawEventHeaderESD(esd);}
808 AliESDpid::MakePID(esd);
809 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
811 if (fFillTriggerESD) {
812 if (!ReadESD(esd, "trigger")) {
813 if (!FillTriggerESD(esd)) {
814 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
816 if (fCheckPointLevel > 1) WriteESD(esd, "trigger");
822 //Try to improve the reconstructed primary vertex position using the tracks
823 AliESDVertex *pvtx=0;
824 Bool_t dovertex=kTRUE;
825 TObject* obj = fOptions.FindObject("ITS");
827 TString optITS = obj->GetTitle();
828 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
831 if(dovertex) pvtx=tVertexer.FindPrimaryVertex(esd);
832 if(fDiamondProfile) esd->SetDiamond(fDiamondProfile);
835 if (pvtx->GetStatus()) {
836 // Store the improved primary vertex
837 esd->SetPrimaryVertex(pvtx);
838 // Propagate the tracks to the DCA to the improved primary vertex
839 Double_t somethingbig = 777.;
840 Double_t bz = esd->GetMagneticField();
841 Int_t nt=esd->GetNumberOfTracks();
843 AliESDtrack *t = esd->GetTrack(nt);
844 t->RelateToVertex(pvtx, bz, somethingbig);
851 vtxer.Tracks2V0vertices(esd);
853 if (fRunCascadeFinder) {
855 AliCascadeVertexer cvtxer;
856 cvtxer.V0sTracks2CascadeVertices(esd);
861 if (fCleanESD) CleanESD(esd);
862 if (fWriteESDfriend) {
863 esdf->~AliESDfriend();
864 new (esdf) AliESDfriend(); // Reset...
865 esd->GetESDfriend(esdf);
872 if (fCheckPointLevel > 0) WriteESD(esd, "final");
875 if (fWriteESDfriend) {
876 esdf->~AliESDfriend();
877 new (esdf) AliESDfriend(); // Reset...
880 // delete esdf; esdf = 0;
883 gSystem->GetProcInfo(&ProcInfo);
884 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
888 // write quality assurance ESDs data (one entry for all events)
889 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
890 if (!IsSelected(fgkDetectorName[iDet], detStr))
892 AliQADataMaker * qadm = GetQADataMaker(iDet);
894 qadm->EndOfCycle(AliQA::kRECPOINTS);
895 qadm->EndOfCycle(AliQA::kESDS);
896 qadm->Finish(AliQA::kRECPOINTS);
897 qadm->Finish(AliQA::kESDS) ;
900 tree->GetUserInfo()->Add(esd);
901 hlttree->GetUserInfo()->Add(hltesd);
905 if(fESDPar.Contains("ESD.par")){
906 AliInfo("Attaching ESD.par to Tree");
907 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
908 tree->GetUserInfo()->Add(fn);
914 tree->SetBranchStatus("ESDfriend*",0);
915 // we want to have only one tree version number
916 tree->Write(tree->GetName(),TObject::kOverwrite);
920 TFile *aodFile = TFile::Open("AliAOD.root", "RECREATE");
921 ESDFile2AODFile(file, aodFile);
926 CleanUp(file, fileOld);
928 // Create tags for the events in the ESD tree (the ESD tree is always present)
929 // In case of empty events the tags will contain dummy values
930 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
931 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent);
933 AliAODTagCreator *aodtagCreator = new AliAODTagCreator();
934 aodtagCreator->CreateAODTags(fFirstEvent,fLastEvent);
941 //_____________________________________________________________________________
942 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
944 // run the local reconstruction
945 static Int_t eventNr=0;
948 // AliCDBManager* man = AliCDBManager::Instance();
949 // Bool_t origCache = man->GetCacheFlag();
951 // TString detStr = detectors;
952 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
953 // if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
954 // AliReconstructor* reconstructor = GetReconstructor(iDet);
955 // if (!reconstructor) continue;
956 // if (reconstructor->HasLocalReconstruction()) continue;
958 // AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
959 // AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
961 // AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
962 // AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
964 // man->SetCacheFlag(kTRUE);
965 // TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
966 // man->GetAll(calibPath); // entries are cached!
968 // AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
971 // fRawReader->RewindEvents();
972 // reconstructor->Reconstruct(fRunLoader, fRawReader);
974 // reconstructor->Reconstruct(fRunLoader);
977 // AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
978 // AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
980 // // unload calibration data
981 // man->UnloadFromCache(calibPath);
982 // //man->ClearCache();
985 // man->SetCacheFlag(origCache);
987 // if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
988 // AliError(Form("the following detectors were not found: %s",
990 // if (fStopOnError) return kFALSE;
997 //_____________________________________________________________________________
998 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1000 // run the local reconstruction
1001 static Int_t eventNr=0;
1002 AliCodeTimerAuto("")
1004 TString detStr = detectors;
1005 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1006 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1007 AliReconstructor* reconstructor = GetReconstructor(iDet);
1008 if (!reconstructor) continue;
1009 AliLoader* loader = fLoader[iDet];
1011 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1015 // conversion of digits
1016 if (fRawReader && reconstructor->HasDigitConversion()) {
1017 AliInfo(Form("converting raw data digits into root objects for %s",
1018 fgkDetectorName[iDet]));
1019 AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1020 fgkDetectorName[iDet]));
1021 loader->LoadDigits("update");
1022 loader->CleanDigits();
1023 loader->MakeDigitsContainer();
1024 TTree* digitsTree = loader->TreeD();
1025 reconstructor->ConvertDigits(fRawReader, digitsTree);
1026 loader->WriteDigits("OVERWRITE");
1027 loader->UnloadDigits();
1030 // local reconstruction
1031 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1032 AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1033 loader->LoadRecPoints("update");
1034 loader->CleanRecPoints();
1035 loader->MakeRecPointsContainer();
1036 TTree* clustersTree = loader->TreeR();
1037 if (fRawReader && !reconstructor->HasDigitConversion()) {
1038 reconstructor->Reconstruct(fRawReader, clustersTree);
1040 loader->LoadDigits("read");
1041 TTree* digitsTree = loader->TreeD();
1043 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1044 if (fStopOnError) return kFALSE;
1046 reconstructor->Reconstruct(digitsTree, clustersTree);
1048 loader->UnloadDigits();
1051 AliQADataMaker * qadm = GetQADataMaker(iDet);
1053 AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1054 AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1056 if (qadm->IsCycleDone() ) {
1057 qadm->EndOfCycle(AliQA::kRECPOINTS) ;
1058 qadm->EndOfCycle(AliQA::kESDS) ;
1059 qadm->StartOfCycle(AliQA::kRECPOINTS) ;
1060 qadm->StartOfCycle(AliQA::kESDS, "same") ;
1062 qadm->Exec(AliQA::kRECPOINTS, clustersTree) ;
1063 AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
1066 loader->WriteRecPoints("OVERWRITE");
1067 loader->UnloadRecPoints();
1068 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1071 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1072 AliError(Form("the following detectors were not found: %s",
1074 if (fStopOnError) return kFALSE;
1080 //_____________________________________________________________________________
1081 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1083 // run the barrel tracking
1085 AliCodeTimerAuto("")
1087 AliESDVertex* vertex = NULL;
1088 Double_t vtxPos[3] = {0, 0, 0};
1089 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1090 TArrayF mcVertex(3);
1091 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1092 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1093 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1097 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1098 AliInfo("running the ITS vertex finder");
1099 if (fLoader[0]) fLoader[0]->LoadRecPoints();
1100 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
1101 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
1103 AliWarning("Vertex not found");
1104 vertex = new AliESDVertex();
1105 vertex->SetName("default");
1108 vertex->SetName("reconstructed");
1112 AliInfo("getting the primary vertex from MC");
1113 vertex = new AliESDVertex(vtxPos, vtxErr);
1117 vertex->GetXYZ(vtxPos);
1118 vertex->GetSigmaXYZ(vtxErr);
1120 AliWarning("no vertex reconstructed");
1121 vertex = new AliESDVertex(vtxPos, vtxErr);
1123 esd->SetVertex(vertex);
1124 // if SPD multiplicity has been determined, it is stored in the ESD
1125 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1126 if(mult)esd->SetMultiplicity(mult);
1128 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1129 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1136 //_____________________________________________________________________________
1137 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1139 // run the HLT barrel tracking
1141 AliCodeTimerAuto("")
1144 AliError("Missing runLoader!");
1148 AliInfo("running HLT tracking");
1150 // Get a pointer to the HLT reconstructor
1151 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1152 if (!reconstructor) return kFALSE;
1155 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1156 TString detName = fgkDetectorName[iDet];
1157 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1158 reconstructor->SetOption(detName.Data());
1159 AliTracker *tracker = reconstructor->CreateTracker();
1161 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1162 if (fStopOnError) return kFALSE;
1166 Double_t vtxErr[3]={0.005,0.005,0.010};
1167 const AliESDVertex *vertex = esd->GetVertex();
1168 vertex->GetXYZ(vtxPos);
1169 tracker->SetVertex(vtxPos,vtxErr);
1171 fLoader[iDet]->LoadRecPoints("read");
1172 TTree* tree = fLoader[iDet]->TreeR();
1174 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1177 tracker->LoadClusters(tree);
1179 if (tracker->Clusters2Tracks(esd) != 0) {
1180 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1184 tracker->UnloadClusters();
1192 //_____________________________________________________________________________
1193 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1195 // run the muon spectrometer tracking
1197 AliCodeTimerAuto("")
1200 AliError("Missing runLoader!");
1203 Int_t iDet = 7; // for MUON
1205 AliInfo("is running...");
1207 // Get a pointer to the MUON reconstructor
1208 AliReconstructor *reconstructor = GetReconstructor(iDet);
1209 if (!reconstructor) return kFALSE;
1212 TString detName = fgkDetectorName[iDet];
1213 AliDebug(1, Form("%s tracking", detName.Data()));
1214 AliTracker *tracker = reconstructor->CreateTracker();
1216 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1221 fLoader[iDet]->LoadTracks("update");
1222 fLoader[iDet]->CleanTracks();
1223 fLoader[iDet]->MakeTracksContainer();
1226 fLoader[iDet]->LoadRecPoints("read");
1227 tracker->LoadClusters(fLoader[iDet]->TreeR());
1229 Int_t rv = tracker->Clusters2Tracks(esd);
1231 fLoader[iDet]->UnloadRecPoints();
1235 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1239 tracker->UnloadClusters();
1241 fLoader[iDet]->UnloadRecPoints();
1243 fLoader[iDet]->WriteTracks("OVERWRITE");
1244 fLoader[iDet]->UnloadTracks();
1252 //_____________________________________________________________________________
1253 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1255 // run the barrel tracking
1256 static Int_t eventNr=0;
1257 AliCodeTimerAuto("")
1259 AliInfo("running tracking");
1261 //Fill the ESD with the T0 info (will be used by the TOF)
1262 if (fReconstructor[11] && fLoader[11]) {
1263 fLoader[11]->LoadRecPoints("READ");
1264 TTree *treeR = fLoader[11]->TreeR();
1265 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
1268 // pass 1: TPC + ITS inwards
1269 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1270 if (!fTracker[iDet]) continue;
1271 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1274 fLoader[iDet]->LoadRecPoints("read");
1275 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
1276 TTree* tree = fLoader[iDet]->TreeR();
1278 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1281 fTracker[iDet]->LoadClusters(tree);
1282 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1284 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1285 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1288 if (fCheckPointLevel > 1) {
1289 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1291 // preliminary PID in TPC needed by the ITS tracker
1293 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1294 AliESDpid::MakePID(esd);
1296 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
1299 // pass 2: ALL backwards
1300 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1301 if (!fTracker[iDet]) continue;
1302 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1305 if (iDet > 1) { // all except ITS, TPC
1307 fLoader[iDet]->LoadRecPoints("read");
1308 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
1309 tree = fLoader[iDet]->TreeR();
1311 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1314 fTracker[iDet]->LoadClusters(tree);
1315 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1319 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1320 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1323 if (fCheckPointLevel > 1) {
1324 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1328 if (iDet > 2) { // all except ITS, TPC, TRD
1329 fTracker[iDet]->UnloadClusters();
1330 fLoader[iDet]->UnloadRecPoints();
1332 // updated PID in TPC needed by the ITS tracker -MI
1334 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1335 AliESDpid::MakePID(esd);
1337 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1340 // write space-points to the ESD in case alignment data output
1342 if (fWriteAlignmentData)
1343 WriteAlignmentData(esd);
1345 // pass 3: TRD + TPC + ITS refit inwards
1346 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1347 if (!fTracker[iDet]) continue;
1348 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1351 if (fTracker[iDet]->RefitInward(esd) != 0) {
1352 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1355 if (fCheckPointLevel > 1) {
1356 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1358 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1360 fTracker[iDet]->UnloadClusters();
1361 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
1362 fLoader[iDet]->UnloadRecPoints();
1363 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
1366 // Propagate track to the vertex - if not done by ITS
1368 Int_t ntracks = esd->GetNumberOfTracks();
1369 for (Int_t itrack=0; itrack<ntracks; itrack++){
1370 const Double_t kRadius = 3; // beam pipe radius
1371 const Double_t kMaxStep = 5; // max step
1372 const Double_t kMaxD = 123456; // max distance to prim vertex
1373 Double_t fieldZ = AliTracker::GetBz(); //
1374 AliESDtrack * track = esd->GetTrack(itrack);
1375 if (!track) continue;
1376 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1377 AliTracker::PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1378 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1384 //_____________________________________________________________________________
1385 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1387 // Remove the data which are not needed for the physics analysis.
1390 AliInfo("Cleaning the ESD...");
1391 Int_t nTracks=esd->GetNumberOfTracks();
1392 AliInfo(Form("Number of ESD tracks before cleaning %d",nTracks));
1394 Float_t cleanPars[]={fDmax,fZmax};
1395 Bool_t rc=esd->Clean(cleanPars);
1397 nTracks=esd->GetNumberOfTracks();
1398 AliInfo(Form("Number of ESD tracks after cleaning %d",nTracks));
1403 //_____________________________________________________________________________
1404 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1406 // fill the event summary data
1408 AliCodeTimerAuto("")
1409 static Int_t eventNr=0;
1410 TString detStr = detectors;
1411 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1412 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1413 AliReconstructor* reconstructor = GetReconstructor(iDet);
1414 if (!reconstructor) continue;
1416 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1417 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1418 TTree* clustersTree = NULL;
1419 if (fLoader[iDet]) {
1420 fLoader[iDet]->LoadRecPoints("read");
1421 clustersTree = fLoader[iDet]->TreeR();
1422 if (!clustersTree) {
1423 AliError(Form("Can't get the %s clusters tree",
1424 fgkDetectorName[iDet]));
1425 if (fStopOnError) return kFALSE;
1428 if (fRawReader && !reconstructor->HasDigitConversion()) {
1429 reconstructor->FillESD(fRawReader, clustersTree, esd);
1431 TTree* digitsTree = NULL;
1432 if (fLoader[iDet]) {
1433 fLoader[iDet]->LoadDigits("read");
1434 digitsTree = fLoader[iDet]->TreeD();
1436 AliError(Form("Can't get the %s digits tree",
1437 fgkDetectorName[iDet]));
1438 if (fStopOnError) return kFALSE;
1441 reconstructor->FillESD(digitsTree, clustersTree, esd);
1442 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1444 if (fLoader[iDet]) {
1445 fLoader[iDet]->UnloadRecPoints();
1448 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1452 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1453 AliError(Form("the following detectors were not found: %s",
1455 if (fStopOnError) return kFALSE;
1457 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
1462 //_____________________________________________________________________________
1463 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
1465 // Reads the trigger decision which is
1466 // stored in Trigger.root file and fills
1467 // the corresponding esd entries
1469 AliCodeTimerAuto("")
1471 AliInfo("Filling trigger information into the ESD");
1474 AliCTPRawStream input(fRawReader);
1475 if (!input.Next()) {
1476 AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1479 esd->SetTriggerMask(input.GetClassMask());
1480 esd->SetTriggerCluster(input.GetClusterMask());
1483 AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1485 if (!runloader->LoadTrigger()) {
1486 AliCentralTrigger *aCTP = runloader->GetTrigger();
1487 esd->SetTriggerMask(aCTP->GetClassMask());
1488 esd->SetTriggerCluster(aCTP->GetClusterMask());
1491 AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1496 AliError("No run loader is available! The trigger information is not stored in the ESD !");
1508 //_____________________________________________________________________________
1509 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
1512 // Filling information from RawReader Header
1515 AliInfo("Filling information from RawReader Header");
1516 esd->SetBunchCrossNumber(0);
1517 esd->SetOrbitNumber(0);
1518 esd->SetPeriodNumber(0);
1519 esd->SetTimeStamp(0);
1520 esd->SetEventType(0);
1521 const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1524 const UInt_t *id = eventHeader->GetP("Id");
1525 esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1526 esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1527 esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1529 esd->SetTimeStamp((eventHeader->Get("Timestamp")));
1530 esd->SetEventType((eventHeader->Get("Type")));
1537 //_____________________________________________________________________________
1538 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1540 // check whether detName is contained in detectors
1541 // if yes, it is removed from detectors
1543 // check if all detectors are selected
1544 if ((detectors.CompareTo("ALL") == 0) ||
1545 detectors.BeginsWith("ALL ") ||
1546 detectors.EndsWith(" ALL") ||
1547 detectors.Contains(" ALL ")) {
1552 // search for the given detector
1553 Bool_t result = kFALSE;
1554 if ((detectors.CompareTo(detName) == 0) ||
1555 detectors.BeginsWith(detName+" ") ||
1556 detectors.EndsWith(" "+detName) ||
1557 detectors.Contains(" "+detName+" ")) {
1558 detectors.ReplaceAll(detName, "");
1562 // clean up the detectors string
1563 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1564 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1565 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1570 //_____________________________________________________________________________
1571 Bool_t AliReconstruction::InitRunLoader()
1573 // get or create the run loader
1575 if (gAlice) delete gAlice;
1578 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1579 // load all base libraries to get the loader classes
1580 TString libs = gSystem->GetLibraries();
1581 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1582 TString detName = fgkDetectorName[iDet];
1583 if (detName == "HLT") continue;
1584 if (libs.Contains("lib" + detName + "base.so")) continue;
1585 gSystem->Load("lib" + detName + "base.so");
1587 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1589 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1593 fRunLoader->CdGAFile();
1594 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1595 if (fRunLoader->LoadgAlice() == 0) {
1596 gAlice = fRunLoader->GetAliRun();
1597 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
1600 if (!gAlice && !fRawReader) {
1601 AliError(Form("no gAlice object found in file %s",
1602 fGAliceFileName.Data()));
1607 //PH This is a temporary fix to give access to the kinematics
1608 //PH that is needed for the labels of ITS clusters
1609 fRunLoader->LoadKinematics();
1611 } else { // galice.root does not exist
1613 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1617 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1618 AliConfig::GetDefaultEventFolderName(),
1621 AliError(Form("could not create run loader in file %s",
1622 fGAliceFileName.Data()));
1626 fRunLoader->MakeTree("E");
1628 while (fRawReader->NextEvent()) {
1629 fRunLoader->SetEventNumber(iEvent);
1630 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1632 fRunLoader->MakeTree("H");
1633 fRunLoader->TreeE()->Fill();
1636 fRawReader->RewindEvents();
1637 if (fNumberOfEventsPerFile > 0)
1638 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
1640 fRunLoader->SetNumberOfEventsPerFile(iEvent);
1641 fRunLoader->WriteHeader("OVERWRITE");
1642 fRunLoader->CdGAFile();
1643 fRunLoader->Write(0, TObject::kOverwrite);
1644 // AliTracker::SetFieldMap(???);
1650 //_____________________________________________________________________________
1651 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1653 // get the reconstructor object and the loader for a detector
1655 if (fReconstructor[iDet]) return fReconstructor[iDet];
1657 // load the reconstructor object
1658 TPluginManager* pluginManager = gROOT->GetPluginManager();
1659 TString detName = fgkDetectorName[iDet];
1660 TString recName = "Ali" + detName + "Reconstructor";
1661 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1663 AliReconstructor* reconstructor = NULL;
1664 // first check if a plugin is defined for the reconstructor
1665 TPluginHandler* pluginHandler =
1666 pluginManager->FindHandler("AliReconstructor", detName);
1667 // if not, add a plugin for it
1668 if (!pluginHandler) {
1669 AliDebug(1, Form("defining plugin for %s", recName.Data()));
1670 TString libs = gSystem->GetLibraries();
1671 if (libs.Contains("lib" + detName + "base.so") ||
1672 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1673 pluginManager->AddHandler("AliReconstructor", detName,
1674 recName, detName + "rec", recName + "()");
1676 pluginManager->AddHandler("AliReconstructor", detName,
1677 recName, detName, recName + "()");
1679 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1681 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1682 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1684 if (reconstructor) {
1685 TObject* obj = fOptions.FindObject(detName.Data());
1686 if (obj) reconstructor->SetOption(obj->GetTitle());
1687 reconstructor->Init();
1688 fReconstructor[iDet] = reconstructor;
1691 // get or create the loader
1692 if (detName != "HLT") {
1693 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1694 if (!fLoader[iDet]) {
1695 AliConfig::Instance()
1696 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1698 // first check if a plugin is defined for the loader
1700 pluginManager->FindHandler("AliLoader", detName);
1701 // if not, add a plugin for it
1702 if (!pluginHandler) {
1703 TString loaderName = "Ali" + detName + "Loader";
1704 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1705 pluginManager->AddHandler("AliLoader", detName,
1706 loaderName, detName + "base",
1707 loaderName + "(const char*, TFolder*)");
1708 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1710 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1712 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1713 fRunLoader->GetEventFolder());
1715 if (!fLoader[iDet]) { // use default loader
1716 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1718 if (!fLoader[iDet]) {
1719 AliWarning(Form("couldn't get loader for %s", detName.Data()));
1720 if (fStopOnError) return NULL;
1722 fRunLoader->AddLoader(fLoader[iDet]);
1723 fRunLoader->CdGAFile();
1724 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1725 fRunLoader->Write(0, TObject::kOverwrite);
1730 return reconstructor;
1733 //_____________________________________________________________________________
1734 Bool_t AliReconstruction::CreateVertexer()
1736 // create the vertexer
1739 AliReconstructor* itsReconstructor = GetReconstructor(0);
1740 if (itsReconstructor) {
1741 fVertexer = itsReconstructor->CreateVertexer();
1744 AliWarning("couldn't create a vertexer for ITS");
1745 if (fStopOnError) return kFALSE;
1751 //_____________________________________________________________________________
1752 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1754 // create the trackers
1756 TString detStr = detectors;
1757 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1758 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1759 AliReconstructor* reconstructor = GetReconstructor(iDet);
1760 if (!reconstructor) continue;
1761 TString detName = fgkDetectorName[iDet];
1762 if (detName == "HLT") {
1763 fRunHLTTracking = kTRUE;
1766 if (detName == "MUON") {
1767 fRunMuonTracking = kTRUE;
1772 fTracker[iDet] = reconstructor->CreateTracker();
1773 if (!fTracker[iDet] && (iDet < 7)) {
1774 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1775 if (fStopOnError) return kFALSE;
1777 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
1783 //_____________________________________________________________________________
1784 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1786 // delete trackers and the run loader and close and delete the file
1788 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1789 delete fReconstructor[iDet];
1790 fReconstructor[iDet] = NULL;
1791 fLoader[iDet] = NULL;
1792 delete fTracker[iDet];
1793 fTracker[iDet] = NULL;
1794 delete fQADataMaker[iDet];
1795 fQADataMaker[iDet] = NULL;
1799 delete fDiamondProfile;
1800 fDiamondProfile = NULL;
1815 gSystem->Unlink("AliESDs.old.root");
1819 //_____________________________________________________________________________
1821 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
1823 // read the ESD event from a file
1825 if (!esd) return kFALSE;
1827 sprintf(fileName, "ESD_%d.%d_%s.root",
1828 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1829 if (gSystem->AccessPathName(fileName)) return kFALSE;
1831 AliInfo(Form("reading ESD from file %s", fileName));
1832 AliDebug(1, Form("reading ESD from file %s", fileName));
1833 TFile* file = TFile::Open(fileName);
1834 if (!file || !file->IsOpen()) {
1835 AliError(Form("opening %s failed", fileName));
1842 esd = (AliESDEvent*) file->Get("ESD");
1851 //_____________________________________________________________________________
1852 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
1854 // write the ESD event to a file
1858 sprintf(fileName, "ESD_%d.%d_%s.root",
1859 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
1861 AliDebug(1, Form("writing ESD to file %s", fileName));
1862 TFile* file = TFile::Open(fileName, "recreate");
1863 if (!file || !file->IsOpen()) {
1864 AliError(Form("opening %s failed", fileName));
1876 //_____________________________________________________________________________
1877 void AliReconstruction::ESDFile2AODFile(TFile* esdFile, TFile* aodFile)
1879 // write all files from the given esd file to an aod file
1881 // create an AliAOD object
1882 AliAODEvent *aod = new AliAODEvent();
1883 aod->CreateStdContent();
1889 TTree *aodTree = new TTree("aodTree", "AliAOD tree");
1890 aodTree->Branch(aod->GetList());
1893 TTree *t = (TTree*) esdFile->Get("esdTree");
1894 AliESDEvent *esd = new AliESDEvent();
1895 esd->ReadFromTree(t);
1897 Int_t nEvents = t->GetEntries();
1899 // set arrays and pointers
1909 // loop over events and fill them
1910 for (Int_t iEvent = 0; iEvent < nEvents; ++iEvent) {
1911 //cout << "event: " << iEvent << endl;
1912 t->GetEntry(iEvent);
1914 // Multiplicity information needed by the header (to be revised!)
1915 Int_t nTracks = esd->GetNumberOfTracks();
1916 Int_t nPosTracks = 0;
1917 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack)
1918 if (esd->GetTrack(iTrack)->Charge()> 0) nPosTracks++;
1920 // Access the header
1921 AliAODHeader *header = aod->GetHeader();
1924 header->SetRunNumber (esd->GetRunNumber() );
1925 header->SetBunchCrossNumber(esd->GetBunchCrossNumber());
1926 header->SetOrbitNumber (esd->GetOrbitNumber() );
1927 header->SetPeriodNumber (esd->GetPeriodNumber() );
1928 header->SetTriggerMask (esd->GetTriggerMask() );
1929 header->SetTriggerCluster (esd->GetTriggerCluster() );
1930 header->SetEventType (esd->GetEventType() );
1931 header->SetMagneticField (esd->GetMagneticField() );
1932 header->SetZDCN1Energy (esd->GetZDCN1Energy() );
1933 header->SetZDCP1Energy (esd->GetZDCP1Energy() );
1934 header->SetZDCN2Energy (esd->GetZDCN2Energy() );
1935 header->SetZDCP2Energy (esd->GetZDCP2Energy() );
1936 header->SetZDCEMEnergy (esd->GetZDCEMEnergy() );
1937 header->SetRefMultiplicity (nTracks);
1938 header->SetRefMultiplicityPos(nPosTracks);
1939 header->SetRefMultiplicityNeg(nTracks - nPosTracks);
1940 header->SetMuonMagFieldScale(-999.); // FIXME
1941 header->SetCentrality(-999.); // FIXME
1943 Int_t nV0s = esd->GetNumberOfV0s();
1944 Int_t nCascades = esd->GetNumberOfCascades();
1945 Int_t nKinks = esd->GetNumberOfKinks();
1946 Int_t nVertices = nV0s + nCascades + nKinks + 1 /* = prim. vtx*/;
1948 Int_t nCaloClus = esd->GetNumberOfCaloClusters();
1950 Int_t nPmdClus = esd->GetNumberOfPmdTracks();
1952 aod->ResetStd(nTracks, nVertices, nV0s+nCascades, nJets, nCaloClus, nFmdClus, nPmdClus);
1954 // Array to take into account the tracks already added to the AOD
1955 Bool_t * usedTrack = NULL;
1957 usedTrack = new Bool_t[nTracks];
1958 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) usedTrack[iTrack]=kFALSE;
1960 // Array to take into account the V0s already added to the AOD
1961 Bool_t * usedV0 = NULL;
1963 usedV0 = new Bool_t[nV0s];
1964 for (Int_t iV0=0; iV0<nV0s; ++iV0) usedV0[iV0]=kFALSE;
1966 // Array to take into account the kinks already added to the AOD
1967 Bool_t * usedKink = NULL;
1969 usedKink = new Bool_t[nKinks];
1970 for (Int_t iKink=0; iKink<nKinks; ++iKink) usedKink[iKink]=kFALSE;
1973 // Access to the AOD container of vertices
1974 TClonesArray &vertices = *(aod->GetVertices());
1977 // Access to the AOD container of tracks
1978 TClonesArray &tracks = *(aod->GetTracks());
1981 // Access to the AOD container of V0s
1982 TClonesArray &V0s = *(aod->GetV0s());
1985 // Add primary vertex. The primary tracks will be defined
1986 // after the loops on the composite objects (V0, cascades, kinks)
1987 const AliESDVertex *vtx = esd->GetPrimaryVertex();
1989 vtx->GetXYZ(pos); // position
1990 vtx->GetCovMatrix(covVtx); //covariance matrix
1992 AliAODVertex * primary = new(vertices[jVertices++])
1993 AliAODVertex(pos, covVtx, vtx->GetChi2toNDF(), NULL, -1, AliAODVertex::kPrimary);
1996 AliAODTrack *aodTrack = 0x0;
1998 // Create vertices starting from the most complex objects
2001 for (Int_t nCascade = 0; nCascade < nCascades; ++nCascade) {
2002 AliESDcascade *cascade = esd->GetCascade(nCascade);
2004 cascade->GetXYZ(pos[0], pos[1], pos[2]);
2005 cascade->GetPosCovXi(covVtx);
2007 // Add the cascade vertex
2008 AliAODVertex * vcascade = new(vertices[jVertices++]) AliAODVertex(pos,
2010 cascade->GetChi2Xi(), // = chi2/NDF since NDF = 2*2-3
2013 AliAODVertex::kCascade);
2015 primary->AddDaughter(vcascade); // the cascade 'particle' (represented by a vertex) is added as a daughter to the primary vertex
2017 // Add the V0 from the cascade. The ESD class have to be optimized...
2018 // Now we have to search for the corresponding V0 in the list of V0s
2019 // using the indeces of the positive and negative tracks
2021 Int_t posFromV0 = cascade->GetPindex();
2022 Int_t negFromV0 = cascade->GetNindex();
2025 AliESDv0 * v0 = 0x0;
2028 for (Int_t iV0=0; iV0<nV0s; ++iV0) {
2030 v0 = esd->GetV0(iV0);
2031 Int_t posV0 = v0->GetPindex();
2032 Int_t negV0 = v0->GetNindex();
2034 if (posV0==posFromV0 && negV0==negFromV0) {
2040 AliAODVertex * vV0FromCascade = 0x0;
2042 if (indV0>-1 && !usedV0[indV0]) {
2044 // the V0 exists in the array of V0s and is not used
2046 usedV0[indV0] = kTRUE;
2048 v0->GetXYZ(pos[0], pos[1], pos[2]);
2049 v0->GetPosCov(covVtx);
2051 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2053 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2059 // the V0 doesn't exist in the array of V0s or was used
2060 cerr << "Error: event " << iEvent << " cascade " << nCascade
2061 << " The V0 " << indV0
2062 << " doesn't exist in the array of V0s or was used!" << endl;
2064 cascade->GetXYZ(pos[0], pos[1], pos[2]);
2065 cascade->GetPosCov(covVtx);
2067 vV0FromCascade = new(vertices[jVertices++]) AliAODVertex(pos,
2069 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2073 vcascade->AddDaughter(vV0FromCascade);
2077 // Add the positive tracks from the V0
2079 if (! usedTrack[posFromV0]) {
2081 usedTrack[posFromV0] = kTRUE;
2083 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2084 esdTrack->GetPxPyPz(p_pos);
2085 esdTrack->GetXYZ(pos);
2086 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2087 esdTrack->GetESDpid(pid);
2089 vV0FromCascade->AddDaughter(aodTrack =
2090 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2091 esdTrack->GetLabel(),
2097 (Short_t)esdTrack->Charge(),
2098 esdTrack->GetITSClusterMap(),
2101 kTRUE, // check if this is right
2102 kFALSE, // check if this is right
2103 AliAODTrack::kSecondary)
2105 aodTrack->ConvertAliPIDtoAODPID();
2108 cerr << "Error: event " << iEvent << " cascade " << nCascade
2109 << " track " << posFromV0 << " has already been used!" << endl;
2112 // Add the negative tracks from the V0
2114 if (!usedTrack[negFromV0]) {
2116 usedTrack[negFromV0] = kTRUE;
2118 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2119 esdTrack->GetPxPyPz(p_neg);
2120 esdTrack->GetXYZ(pos);
2121 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2122 esdTrack->GetESDpid(pid);
2124 vV0FromCascade->AddDaughter(aodTrack =
2125 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2126 esdTrack->GetLabel(),
2132 (Short_t)esdTrack->Charge(),
2133 esdTrack->GetITSClusterMap(),
2136 kTRUE, // check if this is right
2137 kFALSE, // check if this is right
2138 AliAODTrack::kSecondary)
2140 aodTrack->ConvertAliPIDtoAODPID();
2143 cerr << "Error: event " << iEvent << " cascade " << nCascade
2144 << " track " << negFromV0 << " has already been used!" << endl;
2147 // add it to the V0 array as well
2148 Double_t d0[2] = { -999., -99.};
2149 // counting is probably wrong
2150 new(V0s[jV0s++]) AliAODv0(vV0FromCascade, -999., -99., p_pos, p_neg, d0); // to be refined
2152 // Add the bachelor track from the cascade
2154 Int_t bachelor = cascade->GetBindex();
2156 if(!usedTrack[bachelor]) {
2158 usedTrack[bachelor] = kTRUE;
2160 AliESDtrack *esdTrack = esd->GetTrack(bachelor);
2161 esdTrack->GetPxPyPz(p);
2162 esdTrack->GetXYZ(pos);
2163 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2164 esdTrack->GetESDpid(pid);
2166 vcascade->AddDaughter(aodTrack =
2167 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2168 esdTrack->GetLabel(),
2174 (Short_t)esdTrack->Charge(),
2175 esdTrack->GetITSClusterMap(),
2178 kTRUE, // check if this is right
2179 kFALSE, // check if this is right
2180 AliAODTrack::kSecondary)
2182 aodTrack->ConvertAliPIDtoAODPID();
2185 cerr << "Error: event " << iEvent << " cascade " << nCascade
2186 << " track " << bachelor << " has already been used!" << endl;
2189 // Add the primary track of the cascade (if any)
2191 } // end of the loop on cascades
2195 for (Int_t nV0 = 0; nV0 < nV0s; ++nV0) {
2197 if (usedV0[nV0]) continue; // skip if aready added to the AOD
2199 AliESDv0 *v0 = esd->GetV0(nV0);
2201 v0->GetXYZ(pos[0], pos[1], pos[2]);
2202 v0->GetPosCov(covVtx);
2204 AliAODVertex * vV0 =
2205 new(vertices[jVertices++]) AliAODVertex(pos,
2207 v0->GetChi2V0(), // = chi2/NDF since NDF = 2*2-3
2211 primary->AddDaughter(vV0);
2213 Int_t posFromV0 = v0->GetPindex();
2214 Int_t negFromV0 = v0->GetNindex();
2216 // Add the positive tracks from the V0
2218 if (!usedTrack[posFromV0]) {
2220 usedTrack[posFromV0] = kTRUE;
2222 AliESDtrack *esdTrack = esd->GetTrack(posFromV0);
2223 esdTrack->GetPxPyPz(p_pos);
2224 esdTrack->GetXYZ(pos);
2225 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2226 esdTrack->GetESDpid(pid);
2228 vV0->AddDaughter(aodTrack =
2229 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2230 esdTrack->GetLabel(),
2236 (Short_t)esdTrack->Charge(),
2237 esdTrack->GetITSClusterMap(),
2240 kTRUE, // check if this is right
2241 kFALSE, // check if this is right
2242 AliAODTrack::kSecondary)
2244 aodTrack->ConvertAliPIDtoAODPID();
2247 cerr << "Error: event " << iEvent << " V0 " << nV0
2248 << " track " << posFromV0 << " has already been used!" << endl;
2251 // Add the negative tracks from the V0
2253 if (!usedTrack[negFromV0]) {
2255 usedTrack[negFromV0] = kTRUE;
2257 AliESDtrack *esdTrack = esd->GetTrack(negFromV0);
2258 esdTrack->GetPxPyPz(p_neg);
2259 esdTrack->GetXYZ(pos);
2260 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2261 esdTrack->GetESDpid(pid);
2263 vV0->AddDaughter(aodTrack =
2264 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2265 esdTrack->GetLabel(),
2271 (Short_t)esdTrack->Charge(),
2272 esdTrack->GetITSClusterMap(),
2275 kTRUE, // check if this is right
2276 kFALSE, // check if this is right
2277 AliAODTrack::kSecondary)
2279 aodTrack->ConvertAliPIDtoAODPID();
2282 cerr << "Error: event " << iEvent << " V0 " << nV0
2283 << " track " << negFromV0 << " has already been used!" << endl;
2286 // add it to the V0 array as well
2287 Double_t d0[2] = { 999., 99.};
2288 new(V0s[jV0s++]) AliAODv0(vV0, 999., 99., p_pos, p_neg, d0); // to be refined
2289 } // end of the loop on V0s
2291 // Kinks: it is a big mess the access to the information in the kinks
2292 // The loop is on the tracks in order to find the mother and daugther of each kink
2295 for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
2297 AliESDtrack * esdTrack = esd->GetTrack(iTrack);
2299 Int_t ikink = esdTrack->GetKinkIndex(0);
2302 // Negative kink index: mother, positive: daughter
2304 // Search for the second track of the kink
2306 for (Int_t jTrack = iTrack+1; jTrack<nTracks; ++jTrack) {
2308 AliESDtrack * esdTrack1 = esd->GetTrack(jTrack);
2310 Int_t jkink = esdTrack1->GetKinkIndex(0);
2312 if ( TMath::Abs(ikink)==TMath::Abs(jkink) ) {
2314 // The two tracks are from the same kink
2316 if (usedKink[TMath::Abs(ikink)-1]) continue; // skip used kinks
2319 Int_t idaughter = -1;
2321 if (ikink<0 && jkink>0) {
2326 else if (ikink>0 && jkink<0) {
2332 cerr << "Error: Wrong combination of kink indexes: "
2333 << ikink << " " << jkink << endl;
2337 // Add the mother track
2339 AliAODTrack * mother = NULL;
2341 if (!usedTrack[imother]) {
2343 usedTrack[imother] = kTRUE;
2345 AliESDtrack *esdTrack = esd->GetTrack(imother);
2346 esdTrack->GetPxPyPz(p);
2347 esdTrack->GetXYZ(pos);
2348 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2349 esdTrack->GetESDpid(pid);
2352 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2353 esdTrack->GetLabel(),
2359 (Short_t)esdTrack->Charge(),
2360 esdTrack->GetITSClusterMap(),
2363 kTRUE, // check if this is right
2364 kTRUE, // check if this is right
2365 AliAODTrack::kPrimary);
2366 primary->AddDaughter(mother);
2367 mother->ConvertAliPIDtoAODPID();
2370 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2371 << " track " << imother << " has already been used!" << endl;
2374 // Add the kink vertex
2375 AliESDkink * kink = esd->GetKink(TMath::Abs(ikink)-1);
2377 AliAODVertex * vkink =
2378 new(vertices[jVertices++]) AliAODVertex(kink->GetPosition(),
2382 esdTrack->GetID(), // This is the track ID of the mother's track!
2383 AliAODVertex::kKink);
2384 // Add the daughter track
2386 AliAODTrack * daughter = NULL;
2388 if (!usedTrack[idaughter]) {
2390 usedTrack[idaughter] = kTRUE;
2392 AliESDtrack *esdTrack = esd->GetTrack(idaughter);
2393 esdTrack->GetPxPyPz(p);
2394 esdTrack->GetXYZ(pos);
2395 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2396 esdTrack->GetESDpid(pid);
2399 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2400 esdTrack->GetLabel(),
2406 (Short_t)esdTrack->Charge(),
2407 esdTrack->GetITSClusterMap(),
2410 kTRUE, // check if this is right
2411 kTRUE, // check if this is right
2412 AliAODTrack::kPrimary);
2413 vkink->AddDaughter(daughter);
2414 daughter->ConvertAliPIDtoAODPID();
2417 cerr << "Error: event " << iEvent << " kink " << TMath::Abs(ikink)-1
2418 << " track " << idaughter << " has already been used!" << endl;
2425 // Tracks (primary and orphan)
2426 for (Int_t nTrack = 0; nTrack < nTracks; ++nTrack) {
2429 if (usedTrack[nTrack]) continue;
2431 AliESDtrack *esdTrack = esd->GetTrack(nTrack);
2432 esdTrack->GetPxPyPz(p);
2433 esdTrack->GetXYZ(pos);
2434 esdTrack->GetCovarianceXYZPxPyPz(covTr);
2435 esdTrack->GetESDpid(pid);
2437 Float_t impactXY, impactZ;
2439 esdTrack->GetImpactParameters(impactXY,impactZ);
2442 // track inside the beam pipe
2444 primary->AddDaughter(aodTrack =
2445 new(tracks[jTracks++]) AliAODTrack(esdTrack->GetID(),
2446 esdTrack->GetLabel(),
2452 (Short_t)esdTrack->Charge(),
2453 esdTrack->GetITSClusterMap(),
2456 kTRUE, // check if this is right
2457 kTRUE, // check if this is right
2458 AliAODTrack::kPrimary)
2460 aodTrack->ConvertAliPIDtoAODPID();
2463 // outside the beam pipe: orphan track
2464 // Don't write them anymore!
2467 } // end of loop on tracks
2470 Int_t nMuTracks = esd->GetNumberOfMuonTracks();
2471 for (Int_t nMuTrack = 0; nMuTrack < nMuTracks; ++nMuTrack) {
2473 AliESDMuonTrack *esdMuTrack = esd->GetMuonTrack(nMuTrack);
2474 p[0] = esdMuTrack->Px();
2475 p[1] = esdMuTrack->Py();
2476 p[2] = esdMuTrack->Pz();
2477 pos[0] = primary->GetX();
2478 pos[1] = primary->GetY();
2479 pos[2] = primary->GetZ();
2481 // has to be changed once the muon pid is provided by the ESD
2482 for (Int_t i = 0; i < 10; pid[i++] = 0.); pid[AliAODTrack::kMuon]=1.;
2484 primary->AddDaughter(aodTrack =
2485 new(tracks[jTracks++]) AliAODTrack(0, // no ID provided
2486 0, // no label provided
2491 NULL, // no covariance matrix provided
2492 esdMuTrack->Charge(),
2493 0, // ITSClusterMap is set below
2496 kFALSE, // muon tracks are not used to fit the primary vtx
2497 kFALSE, // not used for vertex fit
2498 AliAODTrack::kPrimary)
2501 aodTrack->SetHitsPatternInTrigCh(esdMuTrack->GetHitsPatternInTrigCh());
2502 Int_t track2Trigger = esdMuTrack->GetMatchTrigger();
2503 aodTrack->SetMatchTrigger(track2Trigger);
2505 aodTrack->SetChi2MatchTrigger(esdMuTrack->GetChi2MatchTrigger());
2507 aodTrack->SetChi2MatchTrigger(0.);
2510 // Access to the AOD container of PMD clusters
2511 TClonesArray &pmdClusters = *(aod->GetPmdClusters());
2512 Int_t jPmdClusters=0;
2514 for (Int_t iPmd = 0; iPmd < nPmdClus; ++iPmd) {
2515 // file pmd clusters, to be revised!
2516 AliESDPmdTrack *pmdTrack = esd->GetPmdTrack(iPmd);
2519 Double_t pos[3] = { pmdTrack->GetClusterX(), pmdTrack->GetClusterY(), pmdTrack->GetClusterZ() };
2520 Double_t pid[9] = { 0., 0., 0., 0., 0., 0., 0., 0., 0. }; // to be revised!
2522 // assoc cluster not set
2523 new(pmdClusters[jPmdClusters++]) AliAODPmdCluster(iPmd, nLabel, label, pmdTrack->GetClusterADC(), pos, pid);
2526 // Access to the AOD container of clusters
2527 TClonesArray &caloClusters = *(aod->GetCaloClusters());
2531 TArrayS EMCCellNumber(15000);
2532 TArrayD EMCCellAmplitude(15000);
2533 Int_t nEMCCells = 0;
2534 const Float_t fEMCAmpScale = 1./500;
2536 for (Int_t iClust=0; iClust<nCaloClus; ++iClust) {
2538 AliESDCaloCluster * cluster = esd->GetCaloCluster(iClust);
2540 Int_t id = cluster->GetID();
2543 Float_t energy = cluster->E();
2544 cluster->GetPosition(posF);
2545 Char_t ttype=AliAODCluster::kUndef;
2547 if (cluster->GetClusterType() == AliESDCaloCluster::kPHOSCluster) {
2548 ttype=AliAODCluster::kPHOSNeutral;
2550 else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALClusterv1) {
2551 ttype = AliAODCluster::kEMCALClusterv1;
2553 else if (cluster->GetClusterType() == AliESDCaloCluster::kEMCALPseudoCluster) {
2554 // Collect raw tower info
2555 for (Int_t iDig = 0; iDig < cluster->GetNumberOfDigits(); iDig++) {
2556 EMCCellNumber[nEMCCells] = cluster->GetDigitIndex()->At(iDig);
2557 EMCCellAmplitude[nEMCCells] = fEMCAmpScale*cluster->GetDigitAmplitude()->At(iDig);
2560 // don't write cluster data (it's just a pseudo cluster, holding the tower information)
2564 AliAODCaloCluster *caloCluster = new(caloClusters[jClusters++]) AliAODCaloCluster(id,
2572 caloCluster->SetCaloCluster(); // to be refined!
2574 } // end of loop on calo clusters
2576 // fill EMC cell info
2577 AliAODCaloCells &EMCCells = *(aod->GetCaloCells());
2578 EMCCells.CreateContainer(nEMCCells);
2579 EMCCells.SetType(AliAODCaloCells::kEMCAL);
2580 for (Int_t iCell = 0; iCell < nEMCCells; iCell++) {
2581 EMCCells.SetCell(iCell,EMCCellNumber[iCell],EMCCellAmplitude[iCell]);
2586 AliAODTracklets &SPDTracklets = *(aod->GetTracklets());
2587 const AliMultiplicity *mult = esd->GetMultiplicity();
2589 if (mult->GetNumberOfTracklets()>0) {
2590 SPDTracklets.CreateContainer(mult->GetNumberOfTracklets());
2592 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
2593 SPDTracklets.SetTracklet(n, mult->GetTheta(n), mult->GetPhi(n), mult->GetDeltaPhi(n), mult->GetLabel(n));
2597 Printf("ERROR: AliMultiplicity could not be retrieved from ESD");
2600 delete [] usedTrack;
2604 // fill the tree for this event
2606 } // end of event loop
2608 aodTree->GetUserInfo()->Add(aod);
2613 // write the tree to the specified file
2614 aodFile = aodTree->GetCurrentFile();
2621 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2623 // Write space-points which are then used in the alignment procedures
2624 // For the moment only ITS, TRD and TPC
2626 // Load TOF clusters
2628 fLoader[3]->LoadRecPoints("read");
2629 TTree* tree = fLoader[3]->TreeR();
2631 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2634 fTracker[3]->LoadClusters(tree);
2636 Int_t ntracks = esd->GetNumberOfTracks();
2637 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2639 AliESDtrack *track = esd->GetTrack(itrack);
2642 for (Int_t iDet = 3; iDet >= 0; iDet--)
2643 nsp += track->GetNcls(iDet);
2645 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2646 track->SetTrackPointArray(sp);
2648 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2649 AliTracker *tracker = fTracker[iDet];
2650 if (!tracker) continue;
2651 Int_t nspdet = track->GetNcls(iDet);
2652 if (nspdet <= 0) continue;
2653 track->GetClusters(iDet,idx);
2657 while (isp < nspdet) {
2659 if(IsSelected(fgkDetectorName[iDet],fUseTrackingErrorsForAlignment)) {
2660 isvalid = tracker->GetTrackPointTrackingError(idx[isp2],p,track);
2662 isvalid = tracker->GetTrackPoint(idx[isp2],p);
2665 const Int_t kNTPCmax = 159;
2666 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
2667 if (!isvalid) continue;
2668 sp->AddPoint(isptrack,&p); isptrack++; isp++;
2674 fTracker[3]->UnloadClusters();
2675 fLoader[3]->UnloadRecPoints();
2679 //_____________________________________________________________________________
2680 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2682 // The method reads the raw-data error log
2683 // accumulated within the rawReader.
2684 // It extracts the raw-data errors related to
2685 // the current event and stores them into
2686 // a TClonesArray inside the esd object.
2688 if (!fRawReader) return;
2690 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2692 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2694 if (iEvent != log->GetEventNumber()) continue;
2696 esd->AddRawDataErrorLog(log);
2701 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString fName){
2702 // Dump a file content into a char in TNamed
2704 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2705 Int_t kBytes = (Int_t)in.tellg();
2706 printf("Size: %d \n",kBytes);
2709 char* memblock = new char [kBytes];
2710 in.seekg (0, ios::beg);
2711 in.read (memblock, kBytes);
2713 TString fData(memblock,kBytes);
2714 fn = new TNamed(fName,fData);
2715 printf("fData Size: %d \n",fData.Sizeof());
2716 printf("fName Size: %d \n",fName.Sizeof());
2717 printf("fn Size: %d \n",fn->Sizeof());
2721 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2727 void AliReconstruction::TNamedToFile(TTree* fTree, TString fName){
2728 // This is not really needed in AliReconstruction at the moment
2729 // but can serve as a template
2731 TList *fList = fTree->GetUserInfo();
2732 TNamed *fn = (TNamed*)fList->FindObject(fName.Data());
2733 printf("fn Size: %d \n",fn->Sizeof());
2735 TString fTmp(fn->GetName()); // to be 100% sure in principle fName also works
2736 const char* cdata = fn->GetTitle();
2737 printf("fTmp Size %d\n",fTmp.Sizeof());
2739 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2740 printf("calculated size %d\n",size);
2741 ofstream out(fName.Data(),ios::out | ios::binary);
2742 out.write(cdata,size);
2749 //_____________________________________________________________________________
2750 AliQADataMaker * AliReconstruction::GetQADataMaker(Int_t iDet)
2752 // get the quality assurance data maker object and the loader for a detector
2754 if (fQADataMaker[iDet])
2755 return fQADataMaker[iDet];
2757 // load the QA data maker object
2758 TPluginManager* pluginManager = gROOT->GetPluginManager();
2759 TString detName = fgkDetectorName[iDet];
2760 TString qadmName = "Ali" + detName + "QADataMaker";
2761 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT"))
2764 AliQADataMaker * qadm = NULL;
2765 // first check if a plugin is defined for the quality assurance data maker
2766 TPluginHandler* pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2767 // if not, add a plugin for it
2768 if (!pluginHandler) {
2769 AliDebug(1, Form("defining plugin for %s", qadmName.Data()));
2770 TString libs = gSystem->GetLibraries();
2771 if (libs.Contains("lib" + detName + "base.so") ||
2772 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2773 pluginManager->AddHandler("AliQADataMaker", detName,
2774 qadmName, detName + "qadm", qadmName + "()");
2776 pluginManager->AddHandler("AliQADataMaker", detName,
2777 qadmName, detName, qadmName + "()");
2779 pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName);
2781 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2782 qadm = (AliQADataMaker *) pluginHandler->ExecPlugin(0);
2785 AliInfo(Form("Initializing quality assurance data maker for %s", fgkDetectorName[iDet]));
2786 qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun(), GetQACycles(fgkDetectorName[iDet]));
2787 qadm->StartOfCycle(AliQA::kRECPOINTS);
2788 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
2789 qadm->StartOfCycle(AliQA::kESDS, "same") ;
2790 fQADataMaker[iDet] = qadm;
2796 //_____________________________________________________________________________
2797 Bool_t AliReconstruction::RunQA(const char* detectors, AliESDEvent *& esd)
2799 // run the Quality Assurance data producer
2801 AliCodeTimerAuto("")
2802 TString detStr = detectors;
2803 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2804 if (!IsSelected(fgkDetectorName[iDet], detStr))
2806 AliQADataMaker * qadm = GetQADataMaker(iDet);
2809 AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2810 AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2812 qadm->Exec(AliQA::kESDS, esd) ;
2815 AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2817 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2818 AliError(Form("the following detectors were not found: %s",
2829 //_____________________________________________________________________________
2830 void AliReconstruction::CheckQA()
2832 // check the QA of SIM for this run and remove the detectors
2833 // with status Fatal
2835 TString newDetList ;
2836 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2837 TString detName(AliQA::GetDetName(iDet)) ;
2838 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2839 fRunLocalReconstruction.Contains("ALL") ) {
2840 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX(iDet)) ;
2841 if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2842 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2843 } else if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kERROR)) {
2844 AliError(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was ERROR", detName.Data())) ;
2845 newDetList += detName ;
2847 } else if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kWARNING) ) {
2848 AliWarning(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was WARNING", detName.Data())) ;
2849 newDetList += detName ;
2851 } else if ( qa->IsSet(AliQA::DETECTORINDEX(iDet), AliQA::kSIM, AliQA::kINFO) ) {
2852 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was INFO", detName.Data())) ;
2853 newDetList += detName ;
2858 fRunLocalReconstruction = newDetList ;
2861 //_____________________________________________________________________________
2862 Int_t AliReconstruction::GetDetIndex(const char* detector)
2864 // return the detector index corresponding to detector
2866 for (index = 0; index < fgkNDetectors ; index++) {
2867 if ( strcmp(detector, fgkDetectorName[index]) == 0 )