1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for running the reconstruction //
22 // Clusters and tracks are created for all detectors and all events by //
25 // AliReconstruction rec; //
28 // The Run method returns kTRUE in case of successful execution. //
30 // If the input to the reconstruction are not simulated digits but raw data, //
31 // this can be specified by an argument of the Run method or by the method //
33 // rec.SetInput("..."); //
35 // The input formats and the corresponding argument are: //
36 // - DDL raw data files: directory name, ends with "/" //
37 // - raw data root file: root file name, extension ".root" //
38 // - raw data DATE file: DATE file name, any other non-empty string //
39 // - MC root files : empty string, default //
41 // By default all events are reconstructed. The reconstruction can be //
42 // limited to a range of events by giving the index of the first and the //
43 // last event as an argument to the Run method or by calling //
45 // rec.SetEventRange(..., ...); //
47 // The index -1 (default) can be used for the last event to indicate no //
48 // upper limit of the event range. //
50 // In case of raw-data reconstruction the user can modify the default //
51 // number of events per digits/clusters/tracks file. In case the option //
52 // is not used the number is set 1. In case the user provides 0, than //
53 // the number of events is equal to the number of events inside the //
54 // raw-data file (i.e. one digits/clusters/tracks file): //
56 // rec.SetNumberOfEventsPerFile(...); //
59 // The name of the galice file can be changed from the default //
60 // "galice.root" by passing it as argument to the AliReconstruction //
61 // constructor or by //
63 // rec.SetGAliceFile("..."); //
65 // The local reconstruction can be switched on or off for individual //
68 // rec.SetRunLocalReconstruction("..."); //
70 // The argument is a (case sensitive) string with the names of the //
71 // detectors separated by a space. The special string "ALL" selects all //
72 // available detectors. This is the default. //
74 // The reconstruction of the primary vertex position can be switched off by //
76 // rec.SetRunVertexFinder(kFALSE); //
78 // The tracking and the creation of ESD tracks can be switched on for //
79 // selected detectors by //
81 // rec.SetRunTracking("..."); //
83 // Uniform/nonuniform field tracking switches (default: uniform field) //
85 // rec.SetUniformFieldTracking(); ( rec.SetUniformFieldTracking(kFALSE); ) //
87 // The filling of additional ESD information can be steered by //
89 // rec.SetFillESD("..."); //
91 // Again, for both methods the string specifies the list of detectors. //
92 // The default is "ALL". //
94 // The call of the shortcut method //
96 // rec.SetRunReconstruction("..."); //
98 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
99 // SetFillESD with the same detector selecting string as argument. //
101 // The reconstruction requires digits or raw data as input. For the creation //
102 // of digits and raw data have a look at the class AliSimulation. //
104 // The input data of a detector can be replaced by the corresponding HLT //
105 // data by calling (usual detector string) //
106 // SetUseHLTData("..."); //
108 // For debug purposes the method SetCheckPointLevel can be used. If the //
109 // argument is greater than 0, files with ESD events will be written after //
110 // selected steps of the reconstruction for each event: //
111 // level 1: after tracking and after filling of ESD (final) //
112 // level 2: in addition after each tracking step //
113 // level 3: in addition after the filling of ESD for each detector //
114 // If a final check point file exists for an event, this event will be //
115 // skipped in the reconstruction. The tracking and the filling of ESD for //
116 // a detector will be skipped as well, if the corresponding check point //
117 // file exists. The ESD event will then be loaded from the file instead. //
119 ///////////////////////////////////////////////////////////////////////////////
126 #include <TPluginManager.h>
127 #include <TGeoManager.h>
128 #include <TLorentzVector.h>
131 #include <TObjArray.h>
133 #include "AliReconstruction.h"
134 #include "AliCodeTimer.h"
135 #include "AliReconstructor.h"
137 #include "AliRunLoader.h"
139 #include "AliRawReaderFile.h"
140 #include "AliRawReaderDate.h"
141 #include "AliRawReaderRoot.h"
142 #include "AliRawEventHeaderBase.h"
143 #include "AliESDEvent.h"
144 #include "AliESDMuonTrack.h"
145 #include "AliESDfriend.h"
146 #include "AliESDVertex.h"
147 #include "AliESDcascade.h"
148 #include "AliESDkink.h"
149 #include "AliESDtrack.h"
150 #include "AliESDCaloCluster.h"
151 #include "AliESDCaloCells.h"
152 #include "AliMultiplicity.h"
153 #include "AliTracker.h"
154 #include "AliVertexer.h"
155 #include "AliVertexerTracks.h"
156 #include "AliV0vertexer.h"
157 #include "AliCascadeVertexer.h"
158 #include "AliHeader.h"
159 #include "AliGenEventHeader.h"
161 #include "AliESDpid.h"
162 #include "AliESDtrack.h"
163 #include "AliESDPmdTrack.h"
165 #include "AliESDTagCreator.h"
166 #include "AliAODTagCreator.h"
168 #include "AliGeomManager.h"
169 #include "AliTrackPointArray.h"
170 #include "AliCDBManager.h"
171 #include "AliCDBStorage.h"
172 #include "AliCDBEntry.h"
173 #include "AliAlignObj.h"
175 #include "AliCentralTrigger.h"
176 #include "AliTriggerConfiguration.h"
177 #include "AliTriggerClass.h"
178 #include "AliCTPRawStream.h"
180 #include "AliQADataMakerRec.h"
181 #include "AliGlobalQADataMaker.h"
183 #include "AliQADataMakerSteer.h"
185 #include "AliPlaneEff.h"
187 #include "AliSysInfo.h" // memory snapshots
188 #include "AliRawHLTManager.h"
191 ClassImp(AliReconstruction)
194 //_____________________________________________________________________________
195 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "HMPID", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "T0", "VZERO", "ACORDE", "HLT"};
197 //_____________________________________________________________________________
198 AliReconstruction::AliReconstruction(const char* gAliceFilename,
199 const char* name, const char* title) :
202 fUniformField(kTRUE),
203 fRunVertexFinder(kTRUE),
204 fRunVertexFinderTracks(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"),
240 fParentRawReader(NULL),
243 fDiamondProfile(NULL),
244 fDiamondProfileTPC(NULL),
245 fMeanVertexConstraint(kTRUE),
249 fAlignObjArray(NULL),
252 fInitCDBCalled(kFALSE),
253 fSetRunNumberFromDataCalled(kFALSE),
257 fSameQACycle(kFALSE),
259 fRunPlaneEff(kFALSE),
271 fIsNewRunLoader(kFALSE),
274 // create reconstruction object with default parameters
276 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
277 fReconstructor[iDet] = NULL;
278 fLoader[iDet] = NULL;
279 fTracker[iDet] = NULL;
280 fQADataMaker[iDet] = NULL;
281 fQACycles[iDet] = 999999;
283 fQADataMaker[fgkNDetectors]=NULL; //Global QA
287 //_____________________________________________________________________________
288 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
291 fUniformField(rec.fUniformField),
292 fRunVertexFinder(rec.fRunVertexFinder),
293 fRunVertexFinderTracks(rec.fRunVertexFinderTracks),
294 fRunHLTTracking(rec.fRunHLTTracking),
295 fRunMuonTracking(rec.fRunMuonTracking),
296 fRunV0Finder(rec.fRunV0Finder),
297 fRunCascadeFinder(rec.fRunCascadeFinder),
298 fStopOnError(rec.fStopOnError),
299 fWriteAlignmentData(rec.fWriteAlignmentData),
300 fWriteESDfriend(rec.fWriteESDfriend),
301 fWriteAOD(rec.fWriteAOD),
302 fFillTriggerESD(rec.fFillTriggerESD),
304 fCleanESD(rec.fCleanESD),
305 fV0DCAmax(rec.fV0DCAmax),
306 fV0CsPmin(rec.fV0CsPmin),
310 fRunLocalReconstruction(rec.fRunLocalReconstruction),
311 fRunTracking(rec.fRunTracking),
312 fFillESD(rec.fFillESD),
313 fUseTrackingErrorsForAlignment(rec.fUseTrackingErrorsForAlignment),
314 fGAliceFileName(rec.fGAliceFileName),
316 fEquipIdMap(rec.fEquipIdMap),
317 fFirstEvent(rec.fFirstEvent),
318 fLastEvent(rec.fLastEvent),
319 fNumberOfEventsPerFile(rec.fNumberOfEventsPerFile),
322 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
323 fLoadAlignData(rec.fLoadAlignData),
324 fESDPar(rec.fESDPar),
325 fUseHLTData(rec.fUseHLTData),
329 fParentRawReader(NULL),
332 fDiamondProfile(NULL),
333 fDiamondProfileTPC(NULL),
334 fMeanVertexConstraint(rec.fMeanVertexConstraint),
338 fAlignObjArray(rec.fAlignObjArray),
339 fCDBUri(rec.fCDBUri),
341 fInitCDBCalled(rec.fInitCDBCalled),
342 fSetRunNumberFromDataCalled(rec.fSetRunNumberFromDataCalled),
344 fRunGlobalQA(rec.fRunGlobalQA),
345 fInLoopQA(rec.fInLoopQA),
346 fSameQACycle(rec.fSameQACycle),
347 fRunPlaneEff(rec.fRunPlaneEff),
359 fIsNewRunLoader(rec.fIsNewRunLoader),
364 for (Int_t i = 0; i < rec.fOptions.GetEntriesFast(); i++) {
365 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
367 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
368 fReconstructor[iDet] = NULL;
369 fLoader[iDet] = NULL;
370 fTracker[iDet] = NULL;
371 fQADataMaker[iDet] = NULL;
372 fQACycles[iDet] = rec.fQACycles[iDet];
374 fQADataMaker[fgkNDetectors]=NULL; //Global QA
375 for (Int_t i = 0; i < rec.fSpecCDBUri.GetEntriesFast(); i++) {
376 if (rec.fSpecCDBUri[i]) fSpecCDBUri.Add(rec.fSpecCDBUri[i]->Clone());
380 //_____________________________________________________________________________
381 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
383 // assignment operator
385 this->~AliReconstruction();
386 new(this) AliReconstruction(rec);
390 //_____________________________________________________________________________
391 AliReconstruction::~AliReconstruction()
397 fSpecCDBUri.Delete();
399 AliCodeTimer::Instance()->Print();
402 //_____________________________________________________________________________
403 void AliReconstruction::InitCDB()
405 // activate a default CDB storage
406 // First check if we have any CDB storage set, because it is used
407 // to retrieve the calibration and alignment constants
409 if (fInitCDBCalled) return;
410 fInitCDBCalled = kTRUE;
412 AliCDBManager* man = AliCDBManager::Instance();
413 if (man->IsDefaultStorageSet())
415 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
416 AliWarning("Default CDB storage has been already set !");
417 AliWarning(Form("Ignoring the default storage declared in AliReconstruction: %s",fCDBUri.Data()));
418 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
419 fCDBUri = man->GetDefaultStorage()->GetURI();
422 if (fCDBUri.Length() > 0)
424 AliDebug(2,"!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
425 AliDebug(2, Form("Default CDB storage is set to: %s", fCDBUri.Data()));
426 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
428 fCDBUri="local://$ALICE_ROOT";
429 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
430 AliWarning("Default CDB storage not yet set !!!!");
431 AliWarning(Form("Setting it now to: %s", fCDBUri.Data()));
432 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
435 man->SetDefaultStorage(fCDBUri);
438 // Now activate the detector specific CDB storage locations
439 for (Int_t i = 0; i < fSpecCDBUri.GetEntriesFast(); i++) {
440 TObject* obj = fSpecCDBUri[i];
442 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
443 AliDebug(2, Form("Specific CDB storage for %s is set to: %s",obj->GetName(),obj->GetTitle()));
444 AliDebug(2, "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
445 man->SetSpecificStorage(obj->GetName(), obj->GetTitle());
450 //_____________________________________________________________________________
451 void AliReconstruction::SetDefaultStorage(const char* uri) {
452 // Store the desired default CDB storage location
453 // Activate it later within the Run() method
459 //_____________________________________________________________________________
460 void AliReconstruction::SetSpecificStorage(const char* calibType, const char* uri) {
461 // Store a detector-specific CDB storage location
462 // Activate it later within the Run() method
464 AliCDBPath aPath(calibType);
465 if(!aPath.IsValid()){
466 // if calibType is not wildcard but it is a valid detector, add "/*" to make it a valid path
467 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
468 if(!strcmp(calibType, fgkDetectorName[iDet])) {
469 aPath.SetPath(Form("%s/*", calibType));
470 AliInfo(Form("Path for specific storage set to %s", aPath.GetPath().Data()));
474 if(!aPath.IsValid()){
475 AliError(Form("Not a valid path or detector: %s", calibType));
480 // // check that calibType refers to a "valid" detector name
481 // Bool_t isDetector = kFALSE;
482 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
483 // TString detName = fgkDetectorName[iDet];
484 // if(aPath.GetLevel0() == detName) {
485 // isDetector = kTRUE;
491 // AliError(Form("Not a valid detector: %s", aPath.GetLevel0().Data()));
495 TObject* obj = fSpecCDBUri.FindObject(aPath.GetPath().Data());
496 if (obj) fSpecCDBUri.Remove(obj);
497 fSpecCDBUri.Add(new TNamed(aPath.GetPath().Data(), uri));
501 //_____________________________________________________________________________
502 Bool_t AliReconstruction::SetRunNumberFromData()
504 // The method is called in Run() in order
505 // to set a correct run number.
506 // In case of raw data reconstruction the
507 // run number is taken from the raw data header
509 if (fSetRunNumberFromDataCalled) return kTRUE;
510 fSetRunNumberFromDataCalled = kTRUE;
512 AliCDBManager* man = AliCDBManager::Instance();
514 if(man->GetRun() > 0) {
515 AliWarning("Run number is taken from raw-event header! Ignoring settings in AliCDBManager!");
519 AliError("No run loader is found !");
522 // read run number from gAlice
523 if(fRunLoader->GetAliRun())
524 AliCDBManager::Instance()->SetRun(fRunLoader->GetHeader()->GetRun());
527 if(fRawReader->NextEvent()) {
528 AliCDBManager::Instance()->SetRun(fRawReader->GetRunNumber());
529 fRawReader->RewindEvents();
532 if(man->GetRun() > 0) {
533 AliWarning("No raw events is found ! Using settings in AliCDBManager !");
538 AliWarning("Neither raw events nor settings in AliCDBManager are found !");
544 AliError("Neither gAlice nor RawReader objects are found !");
554 //_____________________________________________________________________________
555 void AliReconstruction::SetCDBLock() {
556 // Set CDB lock: from now on it is forbidden to reset the run number
557 // or the default storage or to activate any further storage!
559 AliCDBManager::Instance()->SetLock(1);
562 //_____________________________________________________________________________
563 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
565 // Read the alignment objects from CDB.
566 // Each detector is supposed to have the
567 // alignment objects in DET/Align/Data CDB path.
568 // All the detector objects are then collected,
569 // sorted by geometry level (starting from ALIC) and
570 // then applied to the TGeo geometry.
571 // Finally an overlaps check is performed.
573 // Load alignment data from CDB and fill fAlignObjArray
574 if(fLoadAlignFromCDB){
576 TString detStr = detectors;
577 TString loadAlObjsListOfDets = "";
579 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
580 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
581 loadAlObjsListOfDets += fgkDetectorName[iDet];
582 loadAlObjsListOfDets += " ";
583 } // end loop over detectors
584 loadAlObjsListOfDets.Prepend("GRP "); //add alignment objects for non-sensitive modules
585 AliGeomManager::ApplyAlignObjsFromCDB(loadAlObjsListOfDets.Data());
587 // Check if the array with alignment objects was
588 // provided by the user. If yes, apply the objects
589 // to the present TGeo geometry
590 if (fAlignObjArray) {
591 if (gGeoManager && gGeoManager->IsClosed()) {
592 if (AliGeomManager::ApplyAlignObjsToGeom(*fAlignObjArray) == kFALSE) {
593 AliError("The misalignment of one or more volumes failed!"
594 "Compare the list of simulated detectors and the list of detector alignment data!");
599 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
605 delete fAlignObjArray; fAlignObjArray=0;
610 //_____________________________________________________________________________
611 void AliReconstruction::SetGAliceFile(const char* fileName)
613 // set the name of the galice file
615 fGAliceFileName = fileName;
618 //_____________________________________________________________________________
619 void AliReconstruction::SetInput(const char* input)
621 // In case the input string starts with 'mem://', we run in an online mode
622 // and AliRawReaderDateOnline object is created. In all other cases a raw-data
623 // file is assumed. One can give as an input:
624 // mem://: - events taken from DAQ monitoring libs online
626 // mem://<filename> - emulation of the above mode (via DATE monitoring libs)
630 //_____________________________________________________________________________
631 void AliReconstruction::SetOption(const char* detector, const char* option)
633 // set options for the reconstruction of a detector
635 TObject* obj = fOptions.FindObject(detector);
636 if (obj) fOptions.Remove(obj);
637 fOptions.Add(new TNamed(detector, option));
640 //_____________________________________________________________________________
641 Bool_t AliReconstruction::Run(const char* input)
644 AliCodeTimerAuto("");
646 if (!InitRun(input)) return kFALSE;
648 //******* The loop over events
650 while ((iEvent < fRunLoader->GetNumberOfEvents()) ||
651 (fRawReader && fRawReader->NextEvent())) {
652 if (!RunEvent(iEvent)) return kFALSE;
656 if (!FinishRun()) return kFALSE;
661 //_____________________________________________________________________________
662 Bool_t AliReconstruction::InitRun(const char* input)
664 // Initialize all the stuff before
665 // going into the event loop
666 // If the second argument is given, the first one is ignored and
667 // the reconstruction works in an online mode
668 AliCodeTimerAuto("");
670 // Overwrite the previous setting
671 if (input) fInput = input;
673 // set the input in case of raw data
674 fRawReader = AliRawReader::Create(fInput.Data());
676 AliInfo("Reconstruction will run over digits");
678 if (!fEquipIdMap.IsNull() && fRawReader)
679 fRawReader->LoadEquipmentIdsMap(fEquipIdMap);
681 if (!fUseHLTData.IsNull()) {
682 // create the RawReaderHLT which performs redirection of HLT input data for
683 // the specified detectors
684 AliRawReader* pRawReader=AliRawHLTManager::CreateRawReaderHLT(fRawReader, fUseHLTData.Data());
686 fParentRawReader=fRawReader;
687 fRawReader=pRawReader;
689 AliError(Form("can not create Raw Reader for HLT input %s", fUseHLTData.Data()));
693 AliSysInfo::AddStamp("Start");
694 // get the run loader
695 if (!InitRunLoader()) return kFALSE;
696 AliSysInfo::AddStamp("LoadLoader");
698 // Initialize the CDB storage
701 AliSysInfo::AddStamp("LoadCDB");
703 // Set run number in CDBManager (if it is not already set by the user)
704 if (!SetRunNumberFromData()) if (fStopOnError) return kFALSE;
706 // Set CDB lock: from now on it is forbidden to reset the run number
707 // or the default storage or to activate any further storage!
710 // Import ideal TGeo geometry and apply misalignment
712 TString geom(gSystem->DirName(fGAliceFileName));
713 geom += "/geometry.root";
714 AliGeomManager::LoadGeometry(geom.Data());
715 if (!gGeoManager) if (fStopOnError) return kFALSE;
718 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
719 AliSysInfo::AddStamp("LoadGeom");
722 AliQADataMakerSteer qas ;
723 if (fRunQA && fRawReader) {
724 qas.Run(fRunLocalReconstruction, fRawReader) ;
725 fSameQACycle = kTRUE ;
727 // checking the QA of previous steps
731 // local reconstruction
732 if (!fRunLocalReconstruction.IsNull()) {
733 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
734 if (fStopOnError) {CleanUp(); return kFALSE;}
740 if (fRunVertexFinder && !CreateVertexer()) {
746 AliSysInfo::AddStamp("Vertexer");
749 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
755 AliSysInfo::AddStamp("LoadTrackers");
757 // get the possibly already existing ESD file and tree
758 fesd = new AliESDEvent(); fhltesd = new AliESDEvent();
759 if (!gSystem->AccessPathName("AliESDs.root")){
760 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
761 ffileOld = TFile::Open("AliESDs.old.root");
762 if (ffileOld && ffileOld->IsOpen()) {
763 ftreeOld = (TTree*) ffileOld->Get("esdTree");
764 if (ftreeOld)fesd->ReadFromTree(ftreeOld);
765 fhlttreeOld = (TTree*) ffileOld->Get("HLTesdTree");
766 if (fhlttreeOld) fhltesd->ReadFromTree(fhlttreeOld);
770 // create the ESD output file and tree
771 ffile = TFile::Open("AliESDs.root", "RECREATE");
772 ffile->SetCompressionLevel(2);
773 if (!ffile->IsOpen()) {
774 AliError("opening AliESDs.root failed");
775 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
778 ftree = new TTree("esdTree", "Tree with ESD objects");
779 fesd = new AliESDEvent();
780 fesd->CreateStdContent();
781 fesd->WriteToTree(ftree);
783 fhlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
784 fhltesd = new AliESDEvent();
785 fhltesd->CreateStdContent();
786 fhltesd->WriteToTree(fhlttree);
789 delete esd; delete hltesd;
790 esd = NULL; hltesd = NULL;
792 // create the branch with ESD additions
796 if (fWriteESDfriend) {
797 fesdf = new AliESDfriend();
798 TBranch *br=ftree->Branch("ESDfriend.","AliESDfriend", &fesdf);
799 br->SetFile("AliESDfriends.root");
800 fesd->AddObject(fesdf);
803 // Get the GRP CDB entry
804 AliCDBEntry* entryGRP = AliCDBManager::Instance()->Get("GRP/GRP/Data");
807 fGRPData = dynamic_cast<TMap*> (entryGRP->GetObject());
810 AliError("No GRP entry found in OCDB!");
812 // Get the diamond profile from OCDB
813 AliCDBEntry* entry = AliCDBManager::Instance()
814 ->Get("GRP/Calib/MeanVertex");
817 fDiamondProfile = dynamic_cast<AliESDVertex*> (entry->GetObject());
819 AliError("No diamond profile found in OCDB!");
823 entry = AliCDBManager::Instance()
824 ->Get("GRP/Calib/MeanVertexTPC");
827 fDiamondProfileTPC = dynamic_cast<AliESDVertex*> (entry->GetObject());
829 AliError("No diamond profile found in OCDB!");
832 ftVertexer = new AliVertexerTracks(AliTracker::GetBz());
833 if(fDiamondProfile && fMeanVertexConstraint) ftVertexer->SetVtxStart(fDiamondProfile);
835 if (fRawReader) fRawReader->RewindEvents();
838 gSystem->GetProcInfo(&ProcInfo);
839 AliInfo(Form("Current memory usage %d %d", ProcInfo.fMemResident, ProcInfo.fMemVirtual));
842 //Initialize the QA and start of cycle for out-of-cycle QA
844 TString detStr(fFillESD);
845 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
846 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
847 AliQADataMakerRec *qadm = GetQADataMaker(iDet);
849 AliInfo(Form("Initializing the QA data maker for %s",
850 fgkDetectorName[iDet]));
851 qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun());
852 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
854 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
855 qadm->StartOfCycle(AliQA::kESDS,"same");
859 AliQADataMakerRec *qadm = GetQADataMaker(AliQA::kGLOBAL);
860 AliInfo(Form("Initializing the global QA data maker"));
862 qadm->Init(AliQA::kRECPOINTS, AliCDBManager::Instance()->GetRun());
863 AliTracker::SetResidualsArray(arr);
864 qadm->Init(AliQA::kESDS, AliCDBManager::Instance()->GetRun());
866 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
867 qadm->StartOfCycle(AliQA::kESDS, "same");
871 fSameQACycle = kTRUE;
874 //Initialize the Plane Efficiency framework
875 if (fRunPlaneEff && !InitPlaneEff()) {
876 if(fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
879 if (strcmp(gProgName,"alieve") == 0)
880 fRunAliEVE = InitAliEVE();
885 //_____________________________________________________________________________
886 Bool_t AliReconstruction::RunEvent(Int_t iEvent)
888 // run the reconstruction over a single event
889 // The event loop is steered in Run method
891 AliCodeTimerAuto("");
893 if (iEvent >= fRunLoader->GetNumberOfEvents()) {
894 fRunLoader->SetEventNumber(iEvent);
895 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
897 //?? fRunLoader->MakeTree("H");
898 fRunLoader->TreeE()->Fill();
901 if ((iEvent < fFirstEvent) || ((fLastEvent >= 0) && (iEvent > fLastEvent))) {
902 // copy old ESD to the new one
904 fesd->ReadFromTree(ftreeOld);
905 ftreeOld->GetEntry(iEvent);
909 fhltesd->ReadFromTree(fhlttreeOld);
910 fhlttreeOld->GetEntry(iEvent);
916 AliInfo(Form("processing event %d", iEvent));
918 //Start of cycle for the in-loop QA
921 TString detStr(fFillESD);
922 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
923 if (!IsSelected(fgkDetectorName[iDet], detStr))
925 AliQADataMakerRec *qadm = GetQADataMaker(iDet);
928 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
929 qadm->StartOfCycle(AliQA::kESDS, "same") ;
932 AliQADataMakerRec *qadm = GetQADataMaker(AliQA::kGLOBAL);
933 qadm->StartOfCycle(AliQA::kRECPOINTS, fSameQACycle);
934 qadm->StartOfCycle(AliQA::kESDS, "same");
939 fRunLoader->GetEvent(iEvent);
942 sprintf(aFileName, "ESD_%d.%d_final.root",
943 fRunLoader->GetHeader()->GetRun(),
944 fRunLoader->GetHeader()->GetEventNrInRun());
945 if (!gSystem->AccessPathName(aFileName)) return kTRUE;
947 // local single event reconstruction
948 if (!fRunLocalReconstruction.IsNull()) {
949 TString detectors=fRunLocalReconstruction;
950 // run HLT event reconstruction first
951 // ;-( IsSelected changes the string
952 if (IsSelected("HLT", detectors) &&
953 !RunLocalEventReconstruction("HLT")) {
954 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
956 detectors=fRunLocalReconstruction;
957 detectors.ReplaceAll("HLT", "");
958 if (!RunLocalEventReconstruction(detectors)) {
959 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
963 fesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
964 fhltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
965 fesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
966 fhltesd->SetEventNumberInFile(fRunLoader->GetHeader()->GetEventNrInRun());
968 // Set magnetic field from the tracker
969 fesd->SetMagneticField(AliTracker::GetBz());
970 fhltesd->SetMagneticField(AliTracker::GetBz());
974 // Fill raw-data error log into the ESD
975 if (fRawReader) FillRawDataErrorLog(iEvent,fesd);
978 if (fRunVertexFinder) {
979 if (!ReadESD(fesd, "vertex")) {
980 if (!RunVertexFinder(fesd)) {
981 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
983 if (fCheckPointLevel > 0) WriteESD(fesd, "vertex");
988 if (!fRunTracking.IsNull()) {
989 if (fRunMuonTracking) {
990 if (!RunMuonTracking(fesd)) {
991 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
997 if (!fRunTracking.IsNull()) {
998 if (!ReadESD(fesd, "tracking")) {
999 if (!RunTracking(fesd)) {
1000 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
1002 if (fCheckPointLevel > 0) WriteESD(fesd, "tracking");
1007 if (!fFillESD.IsNull()) {
1008 TString detectors=fFillESD;
1009 // run HLT first and on hltesd
1010 // ;-( IsSelected changes the string
1011 if (IsSelected("HLT", detectors) &&
1012 !FillESD(fhltesd, "HLT")) {
1013 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
1016 detectors.ReplaceAll("HLT", "");
1017 if (!FillESD(fesd, detectors)) {
1018 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
1022 // fill Event header information from the RawEventHeader
1023 if (fRawReader){FillRawEventHeaderESD(fesd);}
1026 AliESDpid::MakePID(fesd);
1027 if (fCheckPointLevel > 1) WriteESD(fesd, "PID");
1029 if (fFillTriggerESD) {
1030 if (!ReadESD(fesd, "trigger")) {
1031 if (!FillTriggerESD(fesd)) {
1032 if (fStopOnError) {CleanUp(ffile, ffileOld); return kFALSE;}
1034 if (fCheckPointLevel > 1) WriteESD(fesd, "trigger");
1041 // Propagate track to the beam pipe (if not already done by ITS)
1043 const Int_t ntracks = fesd->GetNumberOfTracks();
1044 const Double_t kBz = fesd->GetMagneticField();
1045 const Double_t kRadius = 2.8; //something less than the beam pipe radius
1048 UShort_t *selectedIdx=new UShort_t[ntracks];
1050 for (Int_t itrack=0; itrack<ntracks; itrack++){
1051 const Double_t kMaxStep = 5; //max step over the material
1054 AliESDtrack *track = fesd->GetTrack(itrack);
1055 if (!track) continue;
1057 AliExternalTrackParam *tpcTrack =
1058 (AliExternalTrackParam *)track->GetTPCInnerParam();
1062 PropagateTrackTo(tpcTrack,kRadius,track->GetMass(),kMaxStep,kTRUE);
1067 Int_t n=trkArray.GetEntriesFast();
1068 selectedIdx[n]=track->GetID();
1069 trkArray.AddLast(tpcTrack);
1072 if (track->GetX() < kRadius) continue;
1075 PropagateTrackTo(track,kRadius,track->GetMass(),kMaxStep,kTRUE);
1077 track->RelateToVertex(fesd->GetPrimaryVertexSPD(), kBz, kRadius);
1082 // Improve the reconstructed primary vertex position using the tracks
1084 TObject *obj = fOptions.FindObject("ITS");
1086 TString optITS = obj->GetTitle();
1087 if (optITS.Contains("cosmics") || optITS.Contains("COSMICS"))
1088 fRunVertexFinderTracks=kFALSE;
1090 if (fRunVertexFinderTracks) {
1091 // TPC + ITS primary vertex
1092 ftVertexer->SetITSrefitRequired();
1093 if(fDiamondProfile && fMeanVertexConstraint) {
1094 ftVertexer->SetVtxStart(fDiamondProfile);
1096 ftVertexer->SetConstraintOff();
1098 AliESDVertex *pvtx=ftVertexer->FindPrimaryVertex(fesd);
1100 if (pvtx->GetStatus()) {
1101 fesd->SetPrimaryVertex(pvtx);
1102 for (Int_t i=0; i<ntracks; i++) {
1103 AliESDtrack *t = fesd->GetTrack(i);
1104 t->RelateToVertex(pvtx, kBz, kRadius);
1109 // TPC-only primary vertex
1110 ftVertexer->SetITSrefitNotRequired();
1111 if(fDiamondProfileTPC && fMeanVertexConstraint) {
1112 ftVertexer->SetVtxStart(fDiamondProfileTPC);
1114 ftVertexer->SetConstraintOff();
1116 pvtx=ftVertexer->FindPrimaryVertex(&trkArray,selectedIdx);
1118 if (pvtx->GetStatus()) {
1119 fesd->SetPrimaryVertexTPC(pvtx);
1120 Int_t nsel=trkArray.GetEntriesFast();
1121 for (Int_t i=0; i<nsel; i++) {
1122 AliExternalTrackParam *t =
1123 (AliExternalTrackParam *)trkArray.UncheckedAt(i);
1124 t->PropagateToDCA(pvtx, kBz, kRadius);
1130 delete[] selectedIdx;
1132 if(fDiamondProfile) fesd->SetDiamond(fDiamondProfile);
1137 AliV0vertexer vtxer;
1138 vtxer.Tracks2V0vertices(fesd);
1140 if (fRunCascadeFinder) {
1142 AliCascadeVertexer cvtxer;
1143 cvtxer.V0sTracks2CascadeVertices(fesd);
1148 if (fCleanESD) CleanESD(fesd);
1152 AliQADataMakerRec *qadm = GetQADataMaker(AliQA::kGLOBAL);
1153 if (qadm) qadm->Exec(AliQA::kESDS, fesd);
1157 if (fWriteESDfriend) {
1158 fesdf->~AliESDfriend();
1159 new (fesdf) AliESDfriend(); // Reset...
1160 fesd->GetESDfriend(fesdf);
1168 if (fRunAliEVE) RunAliEVE();
1170 if (fCheckPointLevel > 0) WriteESD(fesd, "final");
1173 if (fWriteESDfriend) {
1174 fesdf->~AliESDfriend();
1175 new (fesdf) AliESDfriend(); // Reset...
1178 ProcInfo_t ProcInfo;
1179 gSystem->GetProcInfo(&ProcInfo);
1180 AliInfo(Form("Event %d -> Current memory usage %d %d",iEvent, ProcInfo.fMemResident, ProcInfo.fMemVirtual));
1183 // End of cycle for the in-loop QA
1186 RunQA(fFillESD.Data(), fesd);
1187 TString detStr(fFillESD);
1188 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1189 if (!IsSelected(fgkDetectorName[iDet], detStr))
1191 AliQADataMakerRec * qadm = GetQADataMaker(iDet);
1194 qadm->EndOfCycle(AliQA::kRECPOINTS);
1195 qadm->EndOfCycle(AliQA::kESDS);
1200 AliQADataMakerRec *qadm = GetQADataMaker(AliQA::kGLOBAL);
1202 qadm->EndOfCycle(AliQA::kRECPOINTS);
1203 qadm->EndOfCycle(AliQA::kESDS);
1212 //_____________________________________________________________________________
1213 Bool_t AliReconstruction::FinishRun()
1216 // Called after the exit
1217 // from the event loop
1218 AliCodeTimerAuto("");
1220 if (fIsNewRunLoader) { // galice.root didn't exist
1221 fRunLoader->WriteHeader("OVERWRITE");
1222 fRunLoader->CdGAFile();
1223 fRunLoader->Write(0, TObject::kOverwrite);
1226 ftree->GetUserInfo()->Add(fesd);
1227 fhlttree->GetUserInfo()->Add(fhltesd);
1229 const TMap *cdbMap = AliCDBManager::Instance()->GetStorageMap();
1230 const TList *cdbList = AliCDBManager::Instance()->GetRetrievedIds();
1232 TMap *cdbMapCopy = new TMap(cdbMap->GetEntries());
1233 cdbMapCopy->SetOwner(1);
1234 cdbMapCopy->SetName("cdbMap");
1235 TIter iter(cdbMap->GetTable());
1238 while((pair = dynamic_cast<TPair*> (iter.Next()))){
1239 TObjString* keyStr = dynamic_cast<TObjString*> (pair->Key());
1240 TObjString* valStr = dynamic_cast<TObjString*> (pair->Value());
1241 cdbMapCopy->Add(new TObjString(keyStr->GetName()), new TObjString(valStr->GetName()));
1244 TList *cdbListCopy = new TList();
1245 cdbListCopy->SetOwner(1);
1246 cdbListCopy->SetName("cdbList");
1248 TIter iter2(cdbList);
1251 while((id = dynamic_cast<AliCDBId*> (iter2.Next()))){
1252 cdbListCopy->Add(new TObjString(id->ToString().Data()));
1255 ftree->GetUserInfo()->Add(cdbMapCopy);
1256 ftree->GetUserInfo()->Add(cdbListCopy);
1259 if(fESDPar.Contains("ESD.par")){
1260 AliInfo("Attaching ESD.par to Tree");
1261 TNamed *fn = CopyFileToTNamed(fESDPar.Data(),"ESD.par");
1262 ftree->GetUserInfo()->Add(fn);
1268 if (fWriteESDfriend)
1269 ftree->SetBranchStatus("ESDfriend*",0);
1270 // we want to have only one tree version number
1271 ftree->Write(ftree->GetName(),TObject::kOverwrite);
1274 // Finish with Plane Efficiency evaluation: before of CleanUp !!!
1275 if (fRunPlaneEff && !FinishPlaneEff()) {
1276 AliWarning("Finish PlaneEff evaluation failed");
1280 CleanUp(ffile, ffileOld);
1283 AliWarning("AOD creation not supported anymore during reconstruction. See ANALYSIS/AliAnalysisTaskESDfilter.cxx instead.");
1286 // Create tags for the events in the ESD tree (the ESD tree is always present)
1287 // In case of empty events the tags will contain dummy values
1288 AliESDTagCreator *esdtagCreator = new AliESDTagCreator();
1289 esdtagCreator->CreateESDTags(fFirstEvent,fLastEvent,fGRPData);
1291 AliWarning("AOD tag creation not supported anymore during reconstruction.");
1294 //Finish QA and end of cycle for out-of-loop QA
1297 AliQADataMakerSteer qas;
1298 qas.Run(fRunLocalReconstruction.Data(), AliQA::kRECPOINTS, fSameQACycle);
1300 qas.Run(fRunTracking.Data(), AliQA::kESDS, fSameQACycle);
1302 AliQADataMakerRec *qadm = GetQADataMaker(AliQA::kGLOBAL);
1304 qadm->EndOfCycle(AliQA::kRECPOINTS);
1305 qadm->EndOfCycle(AliQA::kESDS);
1312 // Cleanup of CDB manager: cache and active storages!
1313 AliCDBManager::Instance()->ClearCache();
1319 //_____________________________________________________________________________
1320 Bool_t AliReconstruction::RunLocalReconstruction(const TString& /*detectors*/)
1322 // run the local reconstruction
1323 static Int_t eventNr=0;
1324 AliCodeTimerAuto("")
1326 // AliCDBManager* man = AliCDBManager::Instance();
1327 // Bool_t origCache = man->GetCacheFlag();
1329 // TString detStr = detectors;
1330 // for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1331 // if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1332 // AliReconstructor* reconstructor = GetReconstructor(iDet);
1333 // if (!reconstructor) continue;
1334 // if (reconstructor->HasLocalReconstruction()) continue;
1336 // AliCodeTimerStart(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1337 // AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1339 // AliCodeTimerStart(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1340 // AliInfo(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1342 // man->SetCacheFlag(kTRUE);
1343 // TString calibPath = Form("%s/Calib/*", fgkDetectorName[iDet]);
1344 // man->GetAll(calibPath); // entries are cached!
1346 // AliCodeTimerStop(Form("Loading calibration data from OCDB for %s", fgkDetectorName[iDet]));
1348 // if (fRawReader) {
1349 // fRawReader->RewindEvents();
1350 // reconstructor->Reconstruct(fRunLoader, fRawReader);
1352 // reconstructor->Reconstruct(fRunLoader);
1355 // AliCodeTimerStop(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1356 // AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr));
1358 // // unload calibration data
1359 // man->UnloadFromCache(calibPath);
1360 // //man->ClearCache();
1363 // man->SetCacheFlag(origCache);
1365 // if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1366 // AliError(Form("the following detectors were not found: %s",
1368 // if (fStopOnError) return kFALSE;
1375 //_____________________________________________________________________________
1376 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
1378 // run the local reconstruction
1380 static Int_t eventNr=0;
1381 AliCodeTimerAuto("")
1383 TString detStr = detectors;
1384 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1385 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1386 AliReconstructor* reconstructor = GetReconstructor(iDet);
1387 if (!reconstructor) continue;
1388 AliLoader* loader = fLoader[iDet];
1389 // Matthias April 2008: temporary fix to run HLT reconstruction
1390 // although the HLT loader is missing
1391 if (strcmp(fgkDetectorName[iDet], "HLT")==0) {
1393 reconstructor->Reconstruct(fRawReader, NULL);
1396 reconstructor->Reconstruct(dummy, NULL);
1401 AliWarning(Form("No loader is defined for %s!",fgkDetectorName[iDet]));
1404 // conversion of digits
1405 if (fRawReader && reconstructor->HasDigitConversion()) {
1406 AliInfo(Form("converting raw data digits into root objects for %s",
1407 fgkDetectorName[iDet]));
1408 AliCodeTimerAuto(Form("converting raw data digits into root objects for %s",
1409 fgkDetectorName[iDet]));
1410 loader->LoadDigits("update");
1411 loader->CleanDigits();
1412 loader->MakeDigitsContainer();
1413 TTree* digitsTree = loader->TreeD();
1414 reconstructor->ConvertDigits(fRawReader, digitsTree);
1415 loader->WriteDigits("OVERWRITE");
1416 loader->UnloadDigits();
1418 // local reconstruction
1419 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1420 AliCodeTimerAuto(Form("running reconstruction for %s", fgkDetectorName[iDet]));
1421 loader->LoadRecPoints("update");
1422 loader->CleanRecPoints();
1423 loader->MakeRecPointsContainer();
1424 TTree* clustersTree = loader->TreeR();
1425 if (fRawReader && !reconstructor->HasDigitConversion()) {
1426 reconstructor->Reconstruct(fRawReader, clustersTree);
1428 loader->LoadDigits("read");
1429 TTree* digitsTree = loader->TreeD();
1431 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
1432 if (fStopOnError) return kFALSE;
1434 reconstructor->Reconstruct(digitsTree, clustersTree);
1436 loader->UnloadDigits();
1439 // In-loop QA for local reconstrucion
1440 if (fRunQA && fInLoopQA) {
1441 AliQADataMakerRec * qadm = GetQADataMaker(iDet);
1444 //(Form("Running QA data maker for %s", fgkDetectorName[iDet]));
1446 //(Form("Running QA data maker for %s", fgkDetectorName[iDet]));
1448 qadm->Exec(AliQA::kRECPOINTS, clustersTree) ;
1451 //(Form("Running QA data maker for %s", fgkDetectorName[iDet]));
1455 loader->WriteRecPoints("OVERWRITE");
1456 loader->UnloadRecPoints();
1457 AliSysInfo::AddStamp(Form("LRec%s_%d",fgkDetectorName[iDet],eventNr), iDet,1,eventNr);
1460 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1461 AliError(Form("the following detectors were not found: %s",
1463 if (fStopOnError) return kFALSE;
1469 //_____________________________________________________________________________
1470 Bool_t AliReconstruction::RunVertexFinder(AliESDEvent*& esd)
1472 // run the barrel tracking
1474 AliCodeTimerAuto("")
1476 AliESDVertex* vertex = NULL;
1477 Double_t vtxPos[3] = {0, 0, 0};
1478 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
1479 TArrayF mcVertex(3);
1480 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
1481 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
1482 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
1486 if(fDiamondProfile) fVertexer->SetVtxStart(fDiamondProfile);
1487 AliInfo("running the ITS vertex finder");
1488 if (fLoader[0]) fLoader[0]->LoadRecPoints();
1489 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
1490 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
1492 AliWarning("Vertex not found");
1493 vertex = new AliESDVertex();
1494 vertex->SetName("default");
1497 vertex->SetName("reconstructed");
1501 AliInfo("getting the primary vertex from MC");
1502 vertex = new AliESDVertex(vtxPos, vtxErr);
1506 vertex->GetXYZ(vtxPos);
1507 vertex->GetSigmaXYZ(vtxErr);
1509 AliWarning("no vertex reconstructed");
1510 vertex = new AliESDVertex(vtxPos, vtxErr);
1512 esd->SetPrimaryVertexSPD(vertex);
1513 // if SPD multiplicity has been determined, it is stored in the ESD
1514 AliMultiplicity *mult = fVertexer->GetMultiplicity();
1515 if(mult)esd->SetMultiplicity(mult);
1517 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1518 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
1525 //_____________________________________________________________________________
1526 Bool_t AliReconstruction::RunHLTTracking(AliESDEvent*& esd)
1528 // run the HLT barrel tracking
1530 AliCodeTimerAuto("")
1533 AliError("Missing runLoader!");
1537 AliInfo("running HLT tracking");
1539 // Get a pointer to the HLT reconstructor
1540 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
1541 if (!reconstructor) return kFALSE;
1544 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1545 TString detName = fgkDetectorName[iDet];
1546 AliDebug(1, Form("%s HLT tracking", detName.Data()));
1547 reconstructor->SetOption(detName.Data());
1548 AliTracker *tracker = reconstructor->CreateTracker();
1550 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
1551 if (fStopOnError) return kFALSE;
1555 Double_t vtxErr[3]={0.005,0.005,0.010};
1556 const AliESDVertex *vertex = esd->GetVertex();
1557 vertex->GetXYZ(vtxPos);
1558 tracker->SetVertex(vtxPos,vtxErr);
1560 fLoader[iDet]->LoadRecPoints("read");
1561 TTree* tree = fLoader[iDet]->TreeR();
1563 AliError(Form("Can't get the %s cluster tree", detName.Data()));
1566 tracker->LoadClusters(tree);
1568 if (tracker->Clusters2Tracks(esd) != 0) {
1569 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
1573 tracker->UnloadClusters();
1581 //_____________________________________________________________________________
1582 Bool_t AliReconstruction::RunMuonTracking(AliESDEvent*& esd)
1584 // run the muon spectrometer tracking
1586 AliCodeTimerAuto("")
1589 AliError("Missing runLoader!");
1592 Int_t iDet = 7; // for MUON
1594 AliInfo("is running...");
1596 // Get a pointer to the MUON reconstructor
1597 AliReconstructor *reconstructor = GetReconstructor(iDet);
1598 if (!reconstructor) return kFALSE;
1601 TString detName = fgkDetectorName[iDet];
1602 AliDebug(1, Form("%s tracking", detName.Data()));
1603 AliTracker *tracker = reconstructor->CreateTracker();
1605 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1610 fLoader[iDet]->LoadRecPoints("read");
1612 tracker->LoadClusters(fLoader[iDet]->TreeR());
1614 Int_t rv = tracker->Clusters2Tracks(esd);
1618 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1622 fLoader[iDet]->UnloadRecPoints();
1624 tracker->UnloadClusters();
1632 //_____________________________________________________________________________
1633 Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd)
1635 // run the barrel tracking
1636 static Int_t eventNr=0;
1637 AliCodeTimerAuto("")
1639 AliInfo("running tracking");
1641 //Fill the ESD with the T0 info (will be used by the TOF)
1642 if (fReconstructor[11] && fLoader[11]) {
1643 fLoader[11]->LoadRecPoints("READ");
1644 TTree *treeR = fLoader[11]->TreeR();
1645 GetReconstructor(11)->FillESD((TTree *)NULL,treeR,esd);
1648 // pass 1: TPC + ITS inwards
1649 for (Int_t iDet = 1; iDet >= 0; iDet--) {
1650 if (!fTracker[iDet]) continue;
1651 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
1654 fLoader[iDet]->LoadRecPoints("read");
1655 AliSysInfo::AddStamp(Form("RLoadCluster%s_%d",fgkDetectorName[iDet],eventNr),iDet,1, eventNr);
1656 TTree* tree = fLoader[iDet]->TreeR();
1658 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1661 fTracker[iDet]->LoadClusters(tree);
1662 AliSysInfo::AddStamp(Form("TLoadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1664 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
1665 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
1668 if (fCheckPointLevel > 1) {
1669 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
1671 // preliminary PID in TPC needed by the ITS tracker
1673 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1674 AliESDpid::MakePID(esd);
1676 AliSysInfo::AddStamp(Form("Tracking0%s_%d",fgkDetectorName[iDet],eventNr), iDet,3,eventNr);
1679 // pass 2: ALL backwards
1681 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1682 if (!fTracker[iDet]) continue;
1683 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
1686 if (iDet > 1) { // all except ITS, TPC
1688 fLoader[iDet]->LoadRecPoints("read");
1689 AliSysInfo::AddStamp(Form("RLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,1, eventNr);
1690 tree = fLoader[iDet]->TreeR();
1692 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
1695 fTracker[iDet]->LoadClusters(tree);
1696 AliSysInfo::AddStamp(Form("TLoadCluster0%s_%d",fgkDetectorName[iDet],eventNr), iDet,2, eventNr);
1700 if (iDet>1) // start filling residuals for the "outer" detectors
1701 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
1703 if (fTracker[iDet]->PropagateBack(esd) != 0) {
1704 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
1707 if (fCheckPointLevel > 1) {
1708 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
1712 if (iDet > 2) { // all except ITS, TPC, TRD
1713 fTracker[iDet]->UnloadClusters();
1714 fLoader[iDet]->UnloadRecPoints();
1716 // updated PID in TPC needed by the ITS tracker -MI
1718 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
1719 AliESDpid::MakePID(esd);
1721 AliSysInfo::AddStamp(Form("Tracking1%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1723 //stop filling residuals for the "outer" detectors
1724 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
1726 // write space-points to the ESD in case alignment data output
1728 if (fWriteAlignmentData)
1729 WriteAlignmentData(esd);
1731 // pass 3: TRD + TPC + ITS refit inwards
1733 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1734 if (!fTracker[iDet]) continue;
1735 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1738 if (iDet<2) // start filling residuals for TPC and ITS
1739 if (fRunGlobalQA) AliTracker::SetFillResiduals(kTRUE);
1741 if (fTracker[iDet]->RefitInward(esd) != 0) {
1742 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1745 // run postprocessing
1746 if (fTracker[iDet]->PostProcess(esd) != 0) {
1747 AliError(Form("%s postprocessing failed", fgkDetectorName[iDet]));
1750 if (fCheckPointLevel > 1) {
1751 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1753 AliSysInfo::AddStamp(Form("Tracking2%s_%d",fgkDetectorName[iDet],eventNr), iDet,3, eventNr);
1755 fTracker[iDet]->UnloadClusters();
1756 AliSysInfo::AddStamp(Form("TUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,4, eventNr);
1757 fLoader[iDet]->UnloadRecPoints();
1758 AliSysInfo::AddStamp(Form("RUnloadCluster%s_%d",fgkDetectorName[iDet],eventNr), iDet,5, eventNr);
1760 // stop filling residuals for TPC and ITS
1761 if (fRunGlobalQA) AliTracker::SetFillResiduals(kFALSE);
1767 //_____________________________________________________________________________
1768 Bool_t AliReconstruction::CleanESD(AliESDEvent *esd){
1770 // Remove the data which are not needed for the physics analysis.
1773 Int_t nTracks=esd->GetNumberOfTracks();
1774 Int_t nV0s=esd->GetNumberOfV0s();
1776 (Form("Number of ESD tracks and V0s before cleaning: %d %d",nTracks,nV0s));
1778 Float_t cleanPars[]={fV0DCAmax,fV0CsPmin,fDmax,fZmax};
1779 Bool_t rc=esd->Clean(cleanPars);
1781 nTracks=esd->GetNumberOfTracks();
1782 nV0s=esd->GetNumberOfV0s();
1784 (Form("Number of ESD tracks and V0s after cleaning %d %d",nTracks,nV0s));
1789 //_____________________________________________________________________________
1790 Bool_t AliReconstruction::FillESD(AliESDEvent*& esd, const TString& detectors)
1792 // fill the event summary data
1794 AliCodeTimerAuto("")
1795 static Int_t eventNr=0;
1796 TString detStr = detectors;
1798 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1799 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1800 AliReconstructor* reconstructor = GetReconstructor(iDet);
1801 if (!reconstructor) continue;
1802 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1803 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1804 TTree* clustersTree = NULL;
1805 if (fLoader[iDet]) {
1806 fLoader[iDet]->LoadRecPoints("read");
1807 clustersTree = fLoader[iDet]->TreeR();
1808 if (!clustersTree) {
1809 AliError(Form("Can't get the %s clusters tree",
1810 fgkDetectorName[iDet]));
1811 if (fStopOnError) return kFALSE;
1814 if (fRawReader && !reconstructor->HasDigitConversion()) {
1815 reconstructor->FillESD(fRawReader, clustersTree, esd);
1817 TTree* digitsTree = NULL;
1818 if (fLoader[iDet]) {
1819 fLoader[iDet]->LoadDigits("read");
1820 digitsTree = fLoader[iDet]->TreeD();
1822 AliError(Form("Can't get the %s digits tree",
1823 fgkDetectorName[iDet]));
1824 if (fStopOnError) return kFALSE;
1827 reconstructor->FillESD(digitsTree, clustersTree, esd);
1828 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1830 if (fLoader[iDet]) {
1831 fLoader[iDet]->UnloadRecPoints();
1834 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1838 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1839 AliError(Form("the following detectors were not found: %s",
1841 if (fStopOnError) return kFALSE;
1843 AliSysInfo::AddStamp(Form("FillESD%d",eventNr), 0,1, eventNr);
1848 //_____________________________________________________________________________
1849 Bool_t AliReconstruction::FillTriggerESD(AliESDEvent*& esd)
1851 // Reads the trigger decision which is
1852 // stored in Trigger.root file and fills
1853 // the corresponding esd entries
1855 AliCodeTimerAuto("")
1857 AliInfo("Filling trigger information into the ESD");
1859 AliCentralTrigger *aCTP = NULL;
1862 AliCTPRawStream input(fRawReader);
1863 if (!input.Next()) {
1864 AliError("No valid CTP (trigger) DDL raw data is found ! The trigger information is not stored in the ESD !");
1867 esd->SetTriggerMask(input.GetClassMask());
1868 esd->SetTriggerCluster(input.GetClusterMask());
1870 aCTP = new AliCentralTrigger();
1871 TString configstr("");
1872 if (!aCTP->LoadConfiguration(configstr)) { // Load CTP config from OCDB
1873 AliError("No trigger configuration found in OCDB! The trigger classes information will no be stored in ESD!");
1879 AliRunLoader *runloader = AliRunLoader::GetRunLoader();
1881 if (!runloader->LoadTrigger()) {
1882 aCTP = runloader->GetTrigger();
1883 esd->SetTriggerMask(aCTP->GetClassMask());
1884 esd->SetTriggerCluster(aCTP->GetClusterMask());
1887 AliWarning("No trigger can be loaded! The trigger information is not stored in the ESD !");
1892 AliError("No run loader is available! The trigger information is not stored in the ESD !");
1897 // Now fill the trigger class names into AliESDRun object
1898 AliTriggerConfiguration *config = aCTP->GetConfiguration();
1900 AliError("No trigger configuration has been found! The trigger classes information will no be stored in ESD!");
1901 if (fRawReader) delete aCTP;
1905 const TObjArray& classesArray = config->GetClasses();
1906 Int_t nclasses = classesArray.GetEntriesFast();
1907 for( Int_t j=0; j<nclasses; j++ ) {
1908 AliTriggerClass* trclass = (AliTriggerClass*)classesArray.At( j );
1909 Int_t trindex = (Int_t)TMath::Log2(trclass->GetMask());
1910 esd->SetTriggerClass(trclass->GetName(),trindex);
1913 if (fRawReader) delete aCTP;
1921 //_____________________________________________________________________________
1922 Bool_t AliReconstruction::FillRawEventHeaderESD(AliESDEvent*& esd)
1925 // Filling information from RawReader Header
1928 AliInfo("Filling information from RawReader Header");
1929 esd->SetBunchCrossNumber(0);
1930 esd->SetOrbitNumber(0);
1931 esd->SetPeriodNumber(0);
1932 esd->SetTimeStamp(0);
1933 esd->SetEventType(0);
1934 const AliRawEventHeaderBase * eventHeader = fRawReader->GetEventHeader();
1937 const UInt_t *id = eventHeader->GetP("Id");
1938 esd->SetBunchCrossNumber((id)[1]&0x00000fff);
1939 esd->SetOrbitNumber((((id)[0]<<20)&0xf00000)|(((id)[1]>>12)&0xfffff));
1940 esd->SetPeriodNumber(((id)[0]>>4)&0x0fffffff);
1942 esd->SetTimeStamp((eventHeader->Get("Timestamp")));
1943 esd->SetEventType((eventHeader->Get("Type")));
1950 //_____________________________________________________________________________
1951 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1953 // check whether detName is contained in detectors
1954 // if yes, it is removed from detectors
1956 // check if all detectors are selected
1957 if ((detectors.CompareTo("ALL") == 0) ||
1958 detectors.BeginsWith("ALL ") ||
1959 detectors.EndsWith(" ALL") ||
1960 detectors.Contains(" ALL ")) {
1965 // search for the given detector
1966 Bool_t result = kFALSE;
1967 if ((detectors.CompareTo(detName) == 0) ||
1968 detectors.BeginsWith(detName+" ") ||
1969 detectors.EndsWith(" "+detName) ||
1970 detectors.Contains(" "+detName+" ")) {
1971 detectors.ReplaceAll(detName, "");
1975 // clean up the detectors string
1976 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1977 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1978 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1983 //_____________________________________________________________________________
1984 Bool_t AliReconstruction::InitRunLoader()
1986 // get or create the run loader
1988 if (gAlice) delete gAlice;
1991 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1992 // load all base libraries to get the loader classes
1993 TString libs = gSystem->GetLibraries();
1994 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1995 TString detName = fgkDetectorName[iDet];
1996 if (detName == "HLT") continue;
1997 if (libs.Contains("lib" + detName + "base.so")) continue;
1998 gSystem->Load("lib" + detName + "base.so");
2000 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
2002 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
2006 fRunLoader->CdGAFile();
2007 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
2008 if (fRunLoader->LoadgAlice() == 0) {
2009 gAlice = fRunLoader->GetAliRun();
2010 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
2013 if (!gAlice && !fRawReader) {
2014 AliError(Form("no gAlice object found in file %s",
2015 fGAliceFileName.Data()));
2020 //PH This is a temporary fix to give access to the kinematics
2021 //PH that is needed for the labels of ITS clusters
2022 fRunLoader->LoadHeader();
2023 fRunLoader->LoadKinematics();
2025 } else { // galice.root does not exist
2027 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
2031 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
2032 AliConfig::GetDefaultEventFolderName(),
2035 AliError(Form("could not create run loader in file %s",
2036 fGAliceFileName.Data()));
2040 fIsNewRunLoader = kTRUE;
2041 fRunLoader->MakeTree("E");
2043 if (fNumberOfEventsPerFile > 0)
2044 fRunLoader->SetNumberOfEventsPerFile(fNumberOfEventsPerFile);
2046 fRunLoader->SetNumberOfEventsPerFile((UInt_t)-1);
2052 //_____________________________________________________________________________
2053 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
2055 // get the reconstructor object and the loader for a detector
2057 if (fReconstructor[iDet]) return fReconstructor[iDet];
2059 // load the reconstructor object
2060 TPluginManager* pluginManager = gROOT->GetPluginManager();
2061 TString detName = fgkDetectorName[iDet];
2062 TString recName = "Ali" + detName + "Reconstructor";
2064 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT")) return NULL;
2066 AliReconstructor* reconstructor = NULL;
2067 // first check if a plugin is defined for the reconstructor
2068 TPluginHandler* pluginHandler =
2069 pluginManager->FindHandler("AliReconstructor", detName);
2070 // if not, add a plugin for it
2071 if (!pluginHandler) {
2072 AliDebug(1, Form("defining plugin for %s", recName.Data()));
2073 TString libs = gSystem->GetLibraries();
2074 if (libs.Contains("lib" + detName + "base.so") ||
2075 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2076 pluginManager->AddHandler("AliReconstructor", detName,
2077 recName, detName + "rec", recName + "()");
2079 pluginManager->AddHandler("AliReconstructor", detName,
2080 recName, detName, recName + "()");
2082 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
2084 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2085 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
2087 if (reconstructor) {
2088 TObject* obj = fOptions.FindObject(detName.Data());
2089 if (obj) reconstructor->SetOption(obj->GetTitle());
2090 reconstructor->Init();
2091 fReconstructor[iDet] = reconstructor;
2094 // get or create the loader
2095 if (detName != "HLT") {
2096 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
2097 if (!fLoader[iDet]) {
2098 AliConfig::Instance()
2099 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
2101 // first check if a plugin is defined for the loader
2103 pluginManager->FindHandler("AliLoader", detName);
2104 // if not, add a plugin for it
2105 if (!pluginHandler) {
2106 TString loaderName = "Ali" + detName + "Loader";
2107 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
2108 pluginManager->AddHandler("AliLoader", detName,
2109 loaderName, detName + "base",
2110 loaderName + "(const char*, TFolder*)");
2111 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
2113 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2115 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
2116 fRunLoader->GetEventFolder());
2118 if (!fLoader[iDet]) { // use default loader
2119 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
2121 if (!fLoader[iDet]) {
2122 AliWarning(Form("couldn't get loader for %s", detName.Data()));
2123 if (fStopOnError) return NULL;
2125 fRunLoader->AddLoader(fLoader[iDet]);
2126 fRunLoader->CdGAFile();
2127 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
2128 fRunLoader->Write(0, TObject::kOverwrite);
2133 return reconstructor;
2136 //_____________________________________________________________________________
2137 Bool_t AliReconstruction::CreateVertexer()
2139 // create the vertexer
2142 AliReconstructor* itsReconstructor = GetReconstructor(0);
2143 if (itsReconstructor) {
2144 fVertexer = itsReconstructor->CreateVertexer();
2147 AliWarning("couldn't create a vertexer for ITS");
2148 if (fStopOnError) return kFALSE;
2154 //_____________________________________________________________________________
2155 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
2157 // create the trackers
2159 TString detStr = detectors;
2160 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2161 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2162 AliReconstructor* reconstructor = GetReconstructor(iDet);
2163 if (!reconstructor) continue;
2164 TString detName = fgkDetectorName[iDet];
2165 if (detName == "HLT") {
2166 fRunHLTTracking = kTRUE;
2169 if (detName == "MUON") {
2170 fRunMuonTracking = kTRUE;
2175 fTracker[iDet] = reconstructor->CreateTracker();
2176 if (!fTracker[iDet] && (iDet < 7)) {
2177 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
2178 if (fStopOnError) return kFALSE;
2180 AliSysInfo::AddStamp(Form("LTracker%s",fgkDetectorName[iDet]), iDet,0);
2186 //_____________________________________________________________________________
2187 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
2189 // delete trackers and the run loader and close and delete the file
2191 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2192 delete fReconstructor[iDet];
2193 fReconstructor[iDet] = NULL;
2194 fLoader[iDet] = NULL;
2195 delete fTracker[iDet];
2196 fTracker[iDet] = NULL;
2197 // delete fQADataMaker[iDet];
2198 // fQADataMaker[iDet] = NULL;
2203 if (ftVertexer) delete ftVertexer;
2206 if(!(AliCDBManager::Instance()->GetCacheFlag())) {
2207 delete fDiamondProfile;
2208 fDiamondProfile = NULL;
2209 delete fDiamondProfileTPC;
2210 fDiamondProfileTPC = NULL;
2220 if (fParentRawReader) delete fParentRawReader;
2221 fParentRawReader=NULL;
2231 gSystem->Unlink("AliESDs.old.root");
2236 //_____________________________________________________________________________
2238 Bool_t AliReconstruction::ReadESD(AliESDEvent*& esd, const char* recStep) const
2240 // read the ESD event from a file
2242 if (!esd) return kFALSE;
2244 sprintf(fileName, "ESD_%d.%d_%s.root",
2245 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
2246 if (gSystem->AccessPathName(fileName)) return kFALSE;
2248 AliInfo(Form("reading ESD from file %s", fileName));
2249 AliDebug(1, Form("reading ESD from file %s", fileName));
2250 TFile* file = TFile::Open(fileName);
2251 if (!file || !file->IsOpen()) {
2252 AliError(Form("opening %s failed", fileName));
2259 esd = (AliESDEvent*) file->Get("ESD");
2268 //_____________________________________________________________________________
2269 void AliReconstruction::WriteESD(AliESDEvent* esd, const char* recStep) const
2271 // write the ESD event to a file
2275 sprintf(fileName, "ESD_%d.%d_%s.root",
2276 esd->GetRunNumber(), esd->GetEventNumberInFile(), recStep);
2278 AliDebug(1, Form("writing ESD to file %s", fileName));
2279 TFile* file = TFile::Open(fileName, "recreate");
2280 if (!file || !file->IsOpen()) {
2281 AliError(Form("opening %s failed", fileName));
2290 void AliReconstruction::WriteAlignmentData(AliESDEvent* esd)
2292 // Write space-points which are then used in the alignment procedures
2293 // For the moment only ITS, TRD and TPC
2295 // Load TOF clusters
2297 fLoader[3]->LoadRecPoints("read");
2298 TTree* tree = fLoader[3]->TreeR();
2300 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
2303 fTracker[3]->LoadClusters(tree);
2305 Int_t ntracks = esd->GetNumberOfTracks();
2306 for (Int_t itrack = 0; itrack < ntracks; itrack++)
2308 AliESDtrack *track = esd->GetTrack(itrack);
2311 for (Int_t iDet = 3; iDet >= 0; iDet--)
2312 nsp += track->GetNcls(iDet);
2314 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
2315 track->SetTrackPointArray(sp);
2317 for (Int_t iDet = 3; iDet >= 0; iDet--) {
2318 AliTracker *tracker = fTracker[iDet];
2319 if (!tracker) continue;
2320 Int_t nspdet = track->GetNcls(iDet);
2321 if (nspdet <= 0) continue;
2322 track->GetClusters(iDet,idx);
2326 while (isp2 < nspdet) {
2328 TString dets = fgkDetectorName[iDet];
2329 if ((fUseTrackingErrorsForAlignment.CompareTo(dets) == 0) ||
2330 fUseTrackingErrorsForAlignment.BeginsWith(dets+" ") ||
2331 fUseTrackingErrorsForAlignment.EndsWith(" "+dets) ||
2332 fUseTrackingErrorsForAlignment.Contains(" "+dets+" ")) {
2333 isvalid = tracker->GetTrackPointTrackingError(idx[isp2],p,track);
2335 isvalid = tracker->GetTrackPoint(idx[isp2],p);
2338 const Int_t kNTPCmax = 159;
2339 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
2340 if (!isvalid) continue;
2341 sp->AddPoint(isptrack,&p); isptrack++; isp++;
2347 fTracker[3]->UnloadClusters();
2348 fLoader[3]->UnloadRecPoints();
2352 //_____________________________________________________________________________
2353 void AliReconstruction::FillRawDataErrorLog(Int_t iEvent, AliESDEvent* esd)
2355 // The method reads the raw-data error log
2356 // accumulated within the rawReader.
2357 // It extracts the raw-data errors related to
2358 // the current event and stores them into
2359 // a TClonesArray inside the esd object.
2361 if (!fRawReader) return;
2363 for(Int_t i = 0; i < fRawReader->GetNumberOfErrorLogs(); i++) {
2365 AliRawDataErrorLog *log = fRawReader->GetErrorLog(i);
2367 if (iEvent != log->GetEventNumber()) continue;
2369 esd->AddRawDataErrorLog(log);
2374 TNamed* AliReconstruction::CopyFileToTNamed(TString fPath,TString pName){
2375 // Dump a file content into a char in TNamed
2377 in.open(fPath.Data(),ios::in | ios::binary|ios::ate);
2378 Int_t kBytes = (Int_t)in.tellg();
2379 printf("Size: %d \n",kBytes);
2382 char* memblock = new char [kBytes];
2383 in.seekg (0, ios::beg);
2384 in.read (memblock, kBytes);
2386 TString fData(memblock,kBytes);
2387 fn = new TNamed(pName,fData);
2388 printf("fData Size: %d \n",fData.Sizeof());
2389 printf("pName Size: %d \n",pName.Sizeof());
2390 printf("fn Size: %d \n",fn->Sizeof());
2394 AliInfo(Form("Could not Open %s\n",fPath.Data()));
2400 void AliReconstruction::TNamedToFile(TTree* fTree, TString pName){
2401 // This is not really needed in AliReconstruction at the moment
2402 // but can serve as a template
2404 TList *fList = fTree->GetUserInfo();
2405 TNamed *fn = (TNamed*)fList->FindObject(pName.Data());
2406 printf("fn Size: %d \n",fn->Sizeof());
2408 TString fTmp(fn->GetName()); // to be 100% sure in principle pName also works
2409 const char* cdata = fn->GetTitle();
2410 printf("fTmp Size %d\n",fTmp.Sizeof());
2412 int size = fn->Sizeof()-fTmp.Sizeof()-sizeof(UChar_t)-sizeof(Int_t); // see dfinition of TString::SizeOf()...
2413 printf("calculated size %d\n",size);
2414 ofstream out(pName.Data(),ios::out | ios::binary);
2415 out.write(cdata,size);
2420 //_____________________________________________________________________________
2421 AliQADataMakerRec * AliReconstruction::GetQADataMaker(Int_t iDet)
2423 // get the quality assurance data maker object and the loader for a detector
2425 if (fQADataMaker[iDet])
2426 return fQADataMaker[iDet];
2428 AliQADataMakerRec * qadm = NULL;
2429 if (iDet == fgkNDetectors) { //Global QA
2430 qadm = new AliGlobalQADataMaker();
2431 fQADataMaker[iDet] = qadm;
2435 // load the QA data maker object
2436 TPluginManager* pluginManager = gROOT->GetPluginManager();
2437 TString detName = fgkDetectorName[iDet];
2438 TString qadmName = "Ali" + detName + "QADataMakerRec";
2439 if (!fIsNewRunLoader && !fRunLoader->GetLoader(detName+"Loader") && (detName != "HLT"))
2442 // first check if a plugin is defined for the quality assurance data maker
2443 TPluginHandler* pluginHandler = pluginManager->FindHandler("AliQADataMakerRec", detName);
2444 // if not, add a plugin for it
2445 if (!pluginHandler) {
2446 AliDebug(1, Form("defining plugin for %s", qadmName.Data()));
2447 TString libs = gSystem->GetLibraries();
2448 if (libs.Contains("lib" + detName + "base.so") ||
2449 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
2450 pluginManager->AddHandler("AliQADataMakerRec", detName,
2451 qadmName, detName + "qadm", qadmName + "()");
2453 pluginManager->AddHandler("AliQADataMakerRec", detName,
2454 qadmName, detName, qadmName + "()");
2456 pluginHandler = pluginManager->FindHandler("AliQADataMakerRec", detName);
2458 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
2459 qadm = (AliQADataMakerRec *) pluginHandler->ExecPlugin(0);
2462 fQADataMaker[iDet] = qadm;
2467 //_____________________________________________________________________________
2468 Bool_t AliReconstruction::RunQA(const char* detectors, AliESDEvent *& esd)
2470 // run the Quality Assurance data producer
2472 AliCodeTimerAuto("")
2473 TString detStr = detectors;
2474 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2475 if (!IsSelected(fgkDetectorName[iDet], detStr))
2477 AliQADataMakerRec * qadm = GetQADataMaker(iDet);
2480 AliCodeTimerStart(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2481 AliInfo(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2483 qadm->Exec(AliQA::kESDS, esd) ;
2486 AliCodeTimerStop(Form("running quality assurance data maker for %s", fgkDetectorName[iDet]));
2488 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
2489 AliError(Form("the following detectors were not found: %s",
2499 //_____________________________________________________________________________
2500 void AliReconstruction::CheckQA()
2502 // check the QA of SIM for this run and remove the detectors
2503 // with status Fatal
2505 TString newRunLocalReconstruction ;
2506 TString newRunTracking ;
2507 TString newFillESD ;
2509 for (Int_t iDet = 0; iDet < AliQA::kNDET; iDet++) {
2510 TString detName(AliQA::GetDetName(iDet)) ;
2511 AliQA * qa = AliQA::Instance(AliQA::DETECTORINDEX_t(iDet)) ;
2512 if ( qa->IsSet(AliQA::DETECTORINDEX_t(iDet), AliQA::kSIM, AliQA::kFATAL)) {
2513 AliInfo(Form("QA status for %s in Hits and/or SDIGITS and/or Digits was Fatal; No reconstruction performed", detName.Data())) ;
2515 if ( fRunLocalReconstruction.Contains(AliQA::GetDetName(iDet)) ||
2516 fRunLocalReconstruction.Contains("ALL") ) {
2517 newRunLocalReconstruction += detName ;
2518 newRunLocalReconstruction += " " ;
2520 if ( fRunTracking.Contains(AliQA::GetDetName(iDet)) ||
2521 fRunTracking.Contains("ALL") ) {
2522 newRunTracking += detName ;
2523 newRunTracking += " " ;
2525 if ( fFillESD.Contains(AliQA::GetDetName(iDet)) ||
2526 fFillESD.Contains("ALL") ) {
2527 newFillESD += detName ;
2532 fRunLocalReconstruction = newRunLocalReconstruction ;
2533 fRunTracking = newRunTracking ;
2534 fFillESD = newFillESD ;
2537 //_____________________________________________________________________________
2538 Int_t AliReconstruction::GetDetIndex(const char* detector)
2540 // return the detector index corresponding to detector
2542 for (index = 0; index < fgkNDetectors ; index++) {
2543 if ( strcmp(detector, fgkDetectorName[index]) == 0 )
2548 //_____________________________________________________________________________
2549 Bool_t AliReconstruction::FinishPlaneEff() {
2551 // Here execute all the necessary operationis, at the end of the tracking phase,
2552 // in case that evaluation of PlaneEfficiencies was required for some detector.
2553 // E.g., write into a DataBase file the PlaneEfficiency which have been evaluated.
2555 // This Preliminary version works only FOR ITS !!!!!
2556 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2559 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2562 //for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
2563 for (Int_t iDet = 0; iDet < 1; iDet++) { // for the time being only ITS
2564 //if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
2565 if(fTracker[iDet]) {
2566 AliPlaneEff *planeeff=fTracker[iDet]->GetPlaneEff();
2567 ret=planeeff->WriteIntoCDB();
2568 if(planeeff->GetCreateHistos()) {
2569 TString name="PlaneEffHisto";
2570 name+=fgkDetectorName[iDet];
2572 ret*=planeeff->WriteHistosToFile(name,"RECREATE");
2578 //_____________________________________________________________________________
2579 Bool_t AliReconstruction::InitPlaneEff() {
2581 // Here execute all the necessary operations, before of the tracking phase,
2582 // for the evaluation of PlaneEfficiencies, in case required for some detectors.
2583 // E.g., read from a DataBase file a first evaluation of the PlaneEfficiency
2584 // which should be updated/recalculated.
2586 // This Preliminary version will work only FOR ITS !!!!!
2587 // other detectors (TOF,TRD, etc. have to develop their specific codes)
2590 // Return: kTRUE if all operations have been done properly, kFALSE otherwise
2592 AliWarning(Form("Implementation of this method not yet done !! Method return kTRUE"));
2596 //_____________________________________________________________________________
2597 Bool_t AliReconstruction::InitAliEVE()
2599 // This method should be called only in case
2600 // AliReconstruction is run
2601 // within the alieve environment.
2602 // It will initialize AliEVE in a way
2603 // so that it can visualize event processed
2604 // by AliReconstruction.
2605 // The return flag shows whenever the
2606 // AliEVE initialization was successful or not.
2609 macroStr.Form("%s/EVE/macros/alieve_online.C",gSystem->ExpandPathName("$ALICE_ROOT"));
2610 AliInfo(Form("Loading AliEVE macro: %s",macroStr.Data()));
2611 if (gROOT->LoadMacro(macroStr.Data()) != 0) return kFALSE;
2613 gROOT->ProcessLine("if (!gAliEveEvent) {gAliEveEvent = new AliEveEventManager();gAliEveEvent->SetAutoLoad(kTRUE);gAliEveEvent->AddNewEventCommand(\"alieve_online()\");gEve->AddEvent(gAliEveEvent);};");
2618 //_____________________________________________________________________________
2619 void AliReconstruction::RunAliEVE()
2621 // Runs AliEVE visualisation of
2622 // the current event.
2623 // Should be executed only after
2624 // successful initialization of AliEVE.
2626 AliInfo("Running AliEVE...");
2627 gROOT->ProcessLine(Form("gAliEveEvent->SetEvent((AliRunLoader*)%p,(AliRawReader*)%p,(AliESDEvent*)%p);",fRunLoader,fRawReader,fesd));
2628 gROOT->ProcessLine("gAliEveEvent->StartStopAutoLoadTimer();");