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 // The name of the galice file can be changed from the default //
51 // "galice.root" by passing it as argument to the AliReconstruction //
52 // constructor or by //
54 // rec.SetGAliceFile("..."); //
56 // The local reconstruction can be switched on or off for individual //
59 // rec.SetRunLocalReconstruction("..."); //
61 // The argument is a (case sensitive) string with the names of the //
62 // detectors separated by a space. The special string "ALL" selects all //
63 // available detectors. This is the default. //
65 // The reconstruction of the primary vertex position can be switched off by //
67 // rec.SetRunVertexFinder(kFALSE); //
69 // The tracking and the creation of ESD tracks can be switched on for //
70 // selected detectors by //
72 // rec.SetRunTracking("..."); //
74 // The filling of additional ESD information can be steered by //
76 // rec.SetFillESD("..."); //
78 // Again, for both methods the string specifies the list of detectors. //
79 // The default is "ALL". //
81 // The call of the shortcut method //
83 // rec.SetRunReconstruction("..."); //
85 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
86 // SetFillESD with the same detector selecting string as argument. //
88 // The reconstruction requires digits or raw data as input. For the creation //
89 // of digits and raw data have a look at the class AliSimulation. //
91 // For debug purposes the method SetCheckPointLevel can be used. If the //
92 // argument is greater than 0, files with ESD events will be written after //
93 // selected steps of the reconstruction for each event: //
94 // level 1: after tracking and after filling of ESD (final) //
95 // level 2: in addition after each tracking step //
96 // level 3: in addition after the filling of ESD for each detector //
97 // If a final check point file exists for an event, this event will be //
98 // skipped in the reconstruction. The tracking and the filling of ESD for //
99 // a detector will be skipped as well, if the corresponding check point //
100 // file exists. The ESD event will then be loaded from the file instead. //
102 ///////////////////////////////////////////////////////////////////////////////
108 #include <TPluginManager.h>
109 #include <TStopwatch.h>
111 #include "AliReconstruction.h"
112 #include "AliReconstructor.h"
114 #include "AliRunLoader.h"
116 #include "AliRawReaderFile.h"
117 #include "AliRawReaderDate.h"
118 #include "AliRawReaderRoot.h"
119 #include "AliTracker.h"
121 #include "AliESDVertex.h"
122 #include "AliVertexer.h"
123 #include "AliHeader.h"
124 #include "AliGenEventHeader.h"
126 #include "AliESDpid.h"
129 ClassImp(AliReconstruction)
132 //_____________________________________________________________________________
133 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
135 //_____________________________________________________________________________
136 AliReconstruction::AliReconstruction(const char* gAliceFilename,
137 const char* name, const char* title) :
140 fRunLocalReconstruction("ALL"),
141 fRunVertexFinder(kTRUE),
144 fGAliceFileName(gAliceFilename),
148 fStopOnError(kFALSE),
157 // create reconstruction object with default parameters
159 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
160 fReconstructor[iDet] = NULL;
161 fLoader[iDet] = NULL;
162 fTracker[iDet] = NULL;
167 //_____________________________________________________________________________
168 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
171 fRunLocalReconstruction(rec.fRunLocalReconstruction),
172 fRunVertexFinder(rec.fRunVertexFinder),
173 fRunTracking(rec.fRunTracking),
174 fFillESD(rec.fFillESD),
175 fGAliceFileName(rec.fGAliceFileName),
177 fFirstEvent(rec.fFirstEvent),
178 fLastEvent(rec.fLastEvent),
179 fStopOnError(rec.fStopOnError),
190 for (Int_t i = 0; i < fOptions.GetEntriesFast(); i++) {
191 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
193 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
194 fReconstructor[iDet] = NULL;
195 fLoader[iDet] = NULL;
196 fTracker[iDet] = NULL;
200 //_____________________________________________________________________________
201 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
203 // assignment operator
205 this->~AliReconstruction();
206 new(this) AliReconstruction(rec);
210 //_____________________________________________________________________________
211 AliReconstruction::~AliReconstruction()
220 //_____________________________________________________________________________
221 void AliReconstruction::SetGAliceFile(const char* fileName)
223 // set the name of the galice file
225 fGAliceFileName = fileName;
228 //_____________________________________________________________________________
229 void AliReconstruction::SetOption(const char* detector, const char* option)
231 // set options for the reconstruction of a detector
233 TObject* obj = fOptions.FindObject(detector);
234 if (obj) fOptions.Remove(obj);
235 fOptions.Add(new TNamed(detector, option));
239 //_____________________________________________________________________________
240 Bool_t AliReconstruction::Run(const char* input,
241 Int_t firstEvent, Int_t lastEvent)
243 // run the reconstruction
246 if (!input) input = fInput.Data();
247 TString fileName(input);
248 if (fileName.EndsWith("/")) {
249 fRawReader = new AliRawReaderFile(fileName);
250 } else if (fileName.EndsWith(".root")) {
251 fRawReader = new AliRawReaderRoot(fileName);
252 } else if (!fileName.IsNull()) {
253 fRawReader = new AliRawReaderDate(fileName);
254 fRawReader->SelectEvents(7);
257 // get the run loader
258 if (!InitRunLoader()) return kFALSE;
260 // local reconstruction
261 if (!fRunLocalReconstruction.IsNull()) {
262 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
263 if (fStopOnError) {CleanUp(); return kFALSE;}
266 // if (!fRunVertexFinder && fRunTracking.IsNull() &&
267 // fFillESD.IsNull()) return kTRUE;
270 if (fRunVertexFinder && !CreateVertexer()) {
278 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
285 // get the possibly already existing ESD file and tree
286 AliESD* esd = new AliESD;
287 TFile* fileOld = NULL;
288 TTree* treeOld = NULL;
289 if (!gSystem->AccessPathName("AliESDs.root")){
290 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
291 fileOld = TFile::Open("AliESDs.old.root");
292 if (fileOld && fileOld->IsOpen()) {
293 treeOld = (TTree*) fileOld->Get("esdTree");
294 if (treeOld) treeOld->SetBranchAddress("ESD", &esd);
298 // create the ESD output file and tree
299 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
300 if (!file->IsOpen()) {
301 AliError("opening AliESDs.root failed");
302 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
304 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
305 tree->Branch("ESD", "AliESD", &esd);
311 if (fRawReader) fRawReader->RewindEvents();
313 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
314 if (fRawReader) fRawReader->NextEvent();
315 if ((iEvent < firstEvent) || ((lastEvent >= 0) && (iEvent > lastEvent))) {
316 // copy old ESD to the new one
318 treeOld->SetBranchAddress("ESD", &esd);
319 treeOld->GetEntry(iEvent);
325 AliInfo(Form("processing event %d", iEvent));
326 fRunLoader->GetEvent(iEvent);
329 sprintf(fileName, "ESD_%d.%d_final.root",
330 fRunLoader->GetHeader()->GetRun(),
331 fRunLoader->GetHeader()->GetEventNrInRun());
332 if (!gSystem->AccessPathName(fileName)) continue;
334 // local reconstruction
335 if (!fRunLocalReconstruction.IsNull()) {
336 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
337 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
342 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
343 esd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
345 esd->SetMagneticField(gAlice->Field()->SolenoidField());
351 if (fRunVertexFinder) {
352 if (!ReadESD(esd, "vertex")) {
353 if (!RunVertexFinder(esd)) {
354 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
356 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
361 if (!fRunTracking.IsNull()) {
362 if (!ReadESD(esd, "tracking")) {
363 if (!RunTracking(esd)) {
364 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
366 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
371 if (!fFillESD.IsNull()) {
372 if (!FillESD(esd, fFillESD)) {
373 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
378 AliESDpid::MakePID(esd);
379 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
384 if (fCheckPointLevel > 0) WriteESD(esd, "final");
391 CleanUp(file, fileOld);
397 //_____________________________________________________________________________
398 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
400 // run the local reconstruction
402 TStopwatch stopwatch;
405 TString detStr = detectors;
406 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
407 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
408 AliReconstructor* reconstructor = GetReconstructor(iDet);
409 if (!reconstructor) continue;
410 if (reconstructor->HasLocalReconstruction()) continue;
412 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
413 TStopwatch stopwatchDet;
414 stopwatchDet.Start();
416 fRawReader->RewindEvents();
417 reconstructor->Reconstruct(fRunLoader, fRawReader);
419 reconstructor->Reconstruct(fRunLoader);
421 AliInfo(Form("execution time for %s:", fgkDetectorName[iDet]));
422 ToAliInfo(stopwatchDet.Print());
425 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
426 AliError(Form("the following detectors were not found: %s",
428 if (fStopOnError) return kFALSE;
431 AliInfo("execution time:");
432 ToAliInfo(stopwatch.Print());
437 //_____________________________________________________________________________
438 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
440 // run the local reconstruction
442 TStopwatch stopwatch;
445 TString detStr = detectors;
446 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
447 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
448 AliReconstructor* reconstructor = GetReconstructor(iDet);
449 if (!reconstructor) continue;
450 AliLoader* loader = fLoader[iDet];
452 // conversion of digits
453 if (fRawReader && reconstructor->HasDigitConversion()) {
454 AliInfo(Form("converting raw data digits into root objects for %s",
455 fgkDetectorName[iDet]));
456 TStopwatch stopwatchDet;
457 stopwatchDet.Start();
458 loader->LoadDigits("update");
459 loader->CleanDigits();
460 loader->MakeDigitsContainer();
461 TTree* digitsTree = loader->TreeD();
462 reconstructor->ConvertDigits(fRawReader, digitsTree);
463 loader->WriteDigits("OVERWRITE");
464 loader->UnloadDigits();
465 AliDebug(1, Form("execution time for %s:", fgkDetectorName[iDet]));
466 ToAliDebug(1, stopwatchDet.Print());
469 // local reconstruction
470 if (!reconstructor->HasLocalReconstruction()) continue;
471 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
472 TStopwatch stopwatchDet;
473 stopwatchDet.Start();
474 loader->LoadRecPoints("update");
475 loader->CleanRecPoints();
476 loader->MakeRecPointsContainer();
477 TTree* clustersTree = loader->TreeR();
478 if (fRawReader && !reconstructor->HasDigitConversion()) {
479 reconstructor->Reconstruct(fRawReader, clustersTree);
481 loader->LoadDigits("read");
482 TTree* digitsTree = loader->TreeD();
484 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
485 if (fStopOnError) return kFALSE;
487 reconstructor->Reconstruct(digitsTree, clustersTree);
489 loader->UnloadDigits();
491 loader->WriteRecPoints("OVERWRITE");
492 loader->UnloadRecPoints();
493 AliDebug(1, Form("execution time for %s:", fgkDetectorName[iDet]));
494 ToAliDebug(1, stopwatchDet.Print());
497 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
498 AliError(Form("the following detectors were not found: %s",
500 if (fStopOnError) return kFALSE;
503 AliInfo("execution time:");
504 ToAliInfo(stopwatch.Print());
509 //_____________________________________________________________________________
510 Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
512 // run the barrel tracking
514 TStopwatch stopwatch;
517 AliESDVertex* vertex = NULL;
518 Double_t vtxPos[3] = {0, 0, 0};
519 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
521 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
522 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
523 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
527 AliInfo("running the ITS vertex finder");
528 if (fLoader[0]) fLoader[0]->LoadRecPoints();
529 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
530 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
532 AliWarning("Vertex not found");
533 vertex = new AliESDVertex();
536 vertex->SetTruePos(vtxPos); // store also the vertex from MC
540 AliInfo("getting the primary vertex from MC");
541 vertex = new AliESDVertex(vtxPos, vtxErr);
545 vertex->GetXYZ(vtxPos);
546 vertex->GetSigmaXYZ(vtxErr);
548 AliWarning("no vertex reconstructed");
549 vertex = new AliESDVertex(vtxPos, vtxErr);
551 esd->SetVertex(vertex);
552 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
553 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
557 AliInfo("execution time:");
558 ToAliInfo(stopwatch.Print());
563 //_____________________________________________________________________________
564 Bool_t AliReconstruction::RunTracking(AliESD*& esd)
566 // run the barrel tracking
568 TStopwatch stopwatch;
571 AliInfo("running tracking");
573 // pass 1: TPC + ITS inwards
574 for (Int_t iDet = 1; iDet >= 0; iDet--) {
575 if (!fTracker[iDet]) continue;
576 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
579 fLoader[iDet]->LoadRecPoints("read");
580 TTree* tree = fLoader[iDet]->TreeR();
582 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
585 fTracker[iDet]->LoadClusters(tree);
588 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
589 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
592 if (fCheckPointLevel > 1) {
593 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
595 // preliminary PID in TPC needed by the ITS tracker
597 GetReconstructor(1)->FillESD(fRunLoader, esd);
598 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
599 AliESDpid::MakePID(esd);
603 // pass 2: ALL backwards
604 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
605 if (!fTracker[iDet]) continue;
606 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
609 if (iDet > 1) { // all except ITS, TPC
611 if (iDet == 3) { // TOF
612 fLoader[iDet]->LoadDigits("read");
613 tree = fLoader[iDet]->TreeD();
615 fLoader[iDet]->LoadRecPoints("read");
616 tree = fLoader[iDet]->TreeR();
619 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
622 fTracker[iDet]->LoadClusters(tree);
626 if (fTracker[iDet]->PropagateBack(esd) != 0) {
627 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
630 if (fCheckPointLevel > 1) {
631 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
635 if (iDet > 2) { // all except ITS, TPC, TRD
636 fTracker[iDet]->UnloadClusters();
637 if (iDet == 3) { // TOF
638 fLoader[iDet]->UnloadDigits();
640 fLoader[iDet]->UnloadRecPoints();
645 // pass 3: TRD + TPC + ITS refit inwards
646 for (Int_t iDet = 2; iDet >= 0; iDet--) {
647 if (!fTracker[iDet]) continue;
648 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
651 if (fTracker[iDet]->RefitInward(esd) != 0) {
652 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
655 if (fCheckPointLevel > 1) {
656 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
660 fTracker[iDet]->UnloadClusters();
661 fLoader[iDet]->UnloadRecPoints();
664 AliInfo("execution time:");
665 ToAliInfo(stopwatch.Print());
670 //_____________________________________________________________________________
671 Bool_t AliReconstruction::FillESD(AliESD*& esd, const TString& detectors)
673 // fill the event summary data
675 TStopwatch stopwatch;
677 AliInfo("filling ESD");
679 TString detStr = detectors;
680 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
681 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
682 AliReconstructor* reconstructor = GetReconstructor(iDet);
683 if (!reconstructor) continue;
685 if (!ReadESD(esd, fgkDetectorName[iDet])) {
686 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
687 TTree* clustersTree = NULL;
688 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
689 fLoader[iDet]->LoadRecPoints("read");
690 clustersTree = fLoader[iDet]->TreeR();
692 AliError(Form("Can't get the %s clusters tree",
693 fgkDetectorName[iDet]));
694 if (fStopOnError) return kFALSE;
697 if (fRawReader && !reconstructor->HasDigitConversion()) {
698 reconstructor->FillESD(fRawReader, clustersTree, esd);
700 TTree* digitsTree = NULL;
702 fLoader[iDet]->LoadDigits("read");
703 digitsTree = fLoader[iDet]->TreeD();
705 AliError(Form("Can't get the %s digits tree",
706 fgkDetectorName[iDet]));
707 if (fStopOnError) return kFALSE;
710 reconstructor->FillESD(digitsTree, clustersTree, esd);
711 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
713 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
714 fLoader[iDet]->UnloadRecPoints();
718 reconstructor->FillESD(fRunLoader, fRawReader, esd);
720 reconstructor->FillESD(fRunLoader, esd);
722 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
726 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
727 AliError(Form("the following detectors were not found: %s",
729 if (fStopOnError) return kFALSE;
732 AliInfo("execution time:");
733 ToAliInfo(stopwatch.Print());
739 //_____________________________________________________________________________
740 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
742 // check whether detName is contained in detectors
743 // if yes, it is removed from detectors
745 // check if all detectors are selected
746 if ((detectors.CompareTo("ALL") == 0) ||
747 detectors.BeginsWith("ALL ") ||
748 detectors.EndsWith(" ALL") ||
749 detectors.Contains(" ALL ")) {
754 // search for the given detector
755 Bool_t result = kFALSE;
756 if ((detectors.CompareTo(detName) == 0) ||
757 detectors.BeginsWith(detName+" ") ||
758 detectors.EndsWith(" "+detName) ||
759 detectors.Contains(" "+detName+" ")) {
760 detectors.ReplaceAll(detName, "");
764 // clean up the detectors string
765 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
766 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
767 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
772 //_____________________________________________________________________________
773 Bool_t AliReconstruction::InitRunLoader()
775 // get or create the run loader
777 if (gAlice) delete gAlice;
780 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
781 // load all base libraries to get the loader classes
782 TString libs = gSystem->GetLibraries();
783 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
784 TString detName = fgkDetectorName[iDet];
785 if (detName == "HLT") continue;
786 if (libs.Contains("lib" + detName + "base.so")) continue;
787 gSystem->Load("lib" + detName + "base.so");
789 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
791 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
795 fRunLoader->CdGAFile();
796 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
797 if (fRunLoader->LoadgAlice() == 0) {
798 gAlice = fRunLoader->GetAliRun();
799 AliTracker::SetFieldMap(gAlice->Field());
802 if (!gAlice && !fRawReader) {
803 AliError(Form("no gAlice object found in file %s",
804 fGAliceFileName.Data()));
809 } else { // galice.root does not exist
811 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
815 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
816 AliConfig::GetDefaultEventFolderName(),
819 AliError(Form("could not create run loader in file %s",
820 fGAliceFileName.Data()));
824 fRunLoader->MakeTree("E");
826 while (fRawReader->NextEvent()) {
827 fRunLoader->SetEventNumber(iEvent);
828 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
830 fRunLoader->MakeTree("H");
831 fRunLoader->TreeE()->Fill();
834 fRawReader->RewindEvents();
835 fRunLoader->WriteHeader("OVERWRITE");
836 fRunLoader->CdGAFile();
837 fRunLoader->Write(0, TObject::kOverwrite);
838 // AliTracker::SetFieldMap(???);
844 //_____________________________________________________________________________
845 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
847 // get the reconstructor object and the loader for a detector
849 if (fReconstructor[iDet]) return fReconstructor[iDet];
851 // load the reconstructor object
852 TPluginManager* pluginManager = gROOT->GetPluginManager();
853 TString detName = fgkDetectorName[iDet];
854 TString recName = "Ali" + detName + "Reconstructor";
855 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
857 if (detName == "HLT") {
858 if (!gROOT->GetClass("AliLevel3")) {
859 gSystem->Load("libAliL3Src.so");
860 gSystem->Load("libAliL3Misc.so");
861 gSystem->Load("libAliL3Hough.so");
862 gSystem->Load("libAliL3Comp.so");
866 AliReconstructor* reconstructor = NULL;
867 // first check if a plugin is defined for the reconstructor
868 TPluginHandler* pluginHandler =
869 pluginManager->FindHandler("AliReconstructor", detName);
870 // if not, add a plugin for it
871 if (!pluginHandler) {
872 AliDebug(1, Form("defining plugin for %s", recName.Data()));
873 TString libs = gSystem->GetLibraries();
874 if (libs.Contains("lib" + detName + "base.so") ||
875 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
876 pluginManager->AddHandler("AliReconstructor", detName,
877 recName, detName + "rec", recName + "()");
879 pluginManager->AddHandler("AliReconstructor", detName,
880 recName, detName, recName + "()");
882 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
884 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
885 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
888 TObject* obj = fOptions.FindObject(detName.Data());
889 if (obj) reconstructor->SetOption(obj->GetTitle());
890 reconstructor->Init(fRunLoader);
891 fReconstructor[iDet] = reconstructor;
894 // get or create the loader
895 if (detName != "HLT") {
896 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
897 if (!fLoader[iDet]) {
898 AliConfig::Instance()
899 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
901 // first check if a plugin is defined for the loader
902 TPluginHandler* pluginHandler =
903 pluginManager->FindHandler("AliLoader", detName);
904 // if not, add a plugin for it
905 if (!pluginHandler) {
906 TString loaderName = "Ali" + detName + "Loader";
907 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
908 pluginManager->AddHandler("AliLoader", detName,
909 loaderName, detName + "base",
910 loaderName + "(const char*, TFolder*)");
911 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
913 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
915 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
916 fRunLoader->GetEventFolder());
918 if (!fLoader[iDet]) { // use default loader
919 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
921 if (!fLoader[iDet]) {
922 AliWarning(Form("couldn't get loader for %s", detName.Data()));
923 if (fStopOnError) return NULL;
925 fRunLoader->AddLoader(fLoader[iDet]);
926 fRunLoader->CdGAFile();
927 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
928 fRunLoader->Write(0, TObject::kOverwrite);
933 return reconstructor;
936 //_____________________________________________________________________________
937 Bool_t AliReconstruction::CreateVertexer()
939 // create the vertexer
942 AliReconstructor* itsReconstructor = GetReconstructor(0);
943 if (itsReconstructor) {
944 fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
947 AliWarning("couldn't create a vertexer for ITS");
948 if (fStopOnError) return kFALSE;
954 //_____________________________________________________________________________
955 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
957 // create the trackers
959 TString detStr = detectors;
960 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
961 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
962 AliReconstructor* reconstructor = GetReconstructor(iDet);
963 if (!reconstructor) continue;
964 TString detName = fgkDetectorName[iDet];
965 if (detName == "HLT") continue;
967 fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
968 if (!fTracker[iDet] && (iDet < 7)) {
969 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
970 if (fStopOnError) return kFALSE;
977 //_____________________________________________________________________________
978 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
980 // delete trackers and the run loader and close and delete the file
982 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
983 delete fReconstructor[iDet];
984 fReconstructor[iDet] = NULL;
985 fLoader[iDet] = NULL;
986 delete fTracker[iDet];
987 fTracker[iDet] = NULL;
1005 gSystem->Unlink("AliESDs.old.root");
1010 //_____________________________________________________________________________
1011 Bool_t AliReconstruction::ReadESD(AliESD*& esd, const char* recStep) const
1013 // read the ESD event from a file
1015 if (!esd) return kFALSE;
1017 sprintf(fileName, "ESD_%d.%d_%s.root",
1018 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1019 if (gSystem->AccessPathName(fileName)) return kFALSE;
1021 AliDebug(1, Form("reading ESD from file %s", fileName));
1022 TFile* file = TFile::Open(fileName);
1023 if (!file || !file->IsOpen()) {
1024 AliError(Form("opening %s failed", fileName));
1031 esd = (AliESD*) file->Get("ESD");
1037 //_____________________________________________________________________________
1038 void AliReconstruction::WriteESD(AliESD* esd, const char* recStep) const
1040 // write the ESD event to a file
1044 sprintf(fileName, "ESD_%d.%d_%s.root",
1045 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1047 AliDebug(1, Form("writing ESD to file %s", fileName));
1048 TFile* file = TFile::Open(fileName, "recreate");
1049 if (!file || !file->IsOpen()) {
1050 AliError(Form("opening %s failed", fileName));