7100f1bd91f3d7e4015226cabff9fd5214011398
[u/mrichter/AliRoot.git] / STEER / AliReconstruction.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
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  **************************************************************************/
15
16 /* $Id$ */
17
18 ///////////////////////////////////////////////////////////////////////////////
19 //                                                                           //
20 // class for running the reconstruction                                      //
21 //                                                                           //
22 // Clusters and tracks are created for all detectors and all events by       //
23 // typing:                                                                   //
24 //                                                                           //
25 //   AliReconstruction rec;                                                  //
26 //   rec.Run();                                                              //
27 //                                                                           //
28 // The Run method returns kTRUE in case of successful execution.             //
29 //                                                                           //
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   //
32 //                                                                           //
33 //   rec.SetInput("...");                                                    //
34 //                                                                           //
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                               //
40 //                                                                           //
41 // The name of the galice file can be changed from the default               //
42 // "galice.root" by passing it as argument to the AliReconstruction          //
43 // constructor or by                                                         //
44 //                                                                           //
45 //   rec.SetGAliceFile("...");                                               //
46 //                                                                           //
47 // The local reconstruction can be switched on or off for individual         //
48 // detectors by                                                              //
49 //                                                                           //
50 //   rec.SetRunLocalReconstruction("...");                                   //
51 //                                                                           //
52 // The argument is a (case sensitive) string with the names of the           //
53 // detectors separated by a space. The special string "ALL" selects all      //
54 // available detectors. This is the default.                                 //
55 //                                                                           //
56 // The reconstruction of the primary vertex position can be switched off by  //
57 //                                                                           //
58 //   rec.SetRunVertexFinder(kFALSE);                                         //
59 //                                                                           //
60 // The tracking in ITS, TPC and TRD and the creation of ESD tracks can be    //
61 // switched off by                                                           //
62 //                                                                           //
63 //   rec.SetRunTracking(kFALSE);                                             //
64 //                                                                           //
65 // The filling of additional ESD information can be steered by               //
66 //                                                                           //
67 //   rec.SetFillESD("...");                                                  //
68 //                                                                           //
69 // Again, the string specifies the list of detectors. The default is "ALL".  //
70 //                                                                           //
71 // The reconstruction requires digits or raw data as input. For the creation //
72 // of digits and raw data have a look at the class AliSimulation.            //
73 //                                                                           //
74 // For debug purposes the method SetCheckPointLevel can be used. If the      //
75 // argument is greater than 0, files with ESD events will be written after   //
76 // selected steps of the reconstruction for each event:                      //
77 //   level 1: after tracking and after filling of ESD (final)                //
78 //   level 2: in addition after each tracking step                           //
79 //   level 3: in addition after the filling of ESD for each detector         //
80 // If a final check point file exists for an event, this event will be       //
81 // skipped in the reconstruction. The tracking and the filling of ESD for    //
82 // a detector will be skipped as well, if the corresponding check point      //
83 // file exists. The ESD event will then be loaded from the file instead.     //
84 //                                                                           //
85 ///////////////////////////////////////////////////////////////////////////////
86
87 #include <TArrayF.h>
88 #include <TFile.h>
89 #include <TSystem.h>
90 #include <TROOT.h>
91 #include <TPluginManager.h>
92 #include <TStopwatch.h>
93
94 #include "AliReconstruction.h"
95 #include "AliRunLoader.h"
96 #include "AliRun.h"
97 #include "AliRawReaderFile.h"
98 #include "AliRawReaderDate.h"
99 #include "AliRawReaderRoot.h"
100 #include "AliTracker.h"
101 #include "AliESD.h"
102 #include "AliESDVertex.h"
103 #include "AliVertexer.h"
104 #include "AliHeader.h"
105 #include "AliGenEventHeader.h"
106 #include "AliESDpid.h"
107 #include "AliMagF.h"
108
109 ClassImp(AliReconstruction)
110
111
112 //_____________________________________________________________________________
113 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
114
115 //_____________________________________________________________________________
116 AliReconstruction::AliReconstruction(const char* gAliceFilename,
117                                      const char* name, const char* title) :
118   TNamed(name, title),
119
120   fRunLocalReconstruction("ALL"),
121   fRunVertexFinder(kTRUE),
122   fRunTracking(kTRUE),
123   fFillESD("ALL"),
124   fGAliceFileName(gAliceFilename),
125   fInput(""),
126   fStopOnError(kFALSE),
127   fCheckPointLevel(0),
128
129   fRunLoader(NULL),
130   fRawReader(NULL),
131   fITSLoader(NULL),
132   fITSVertexer(NULL),
133   fITSTracker(NULL),
134   fTPCLoader(NULL),
135   fTPCTracker(NULL),
136   fTRDLoader(NULL),
137   fTRDTracker(NULL),
138   fTOFLoader(NULL),
139   fTOFTracker(NULL),
140
141   fReconstructors(),
142   fOptions()
143 {
144 // create reconstruction object with default parameters
145
146 }
147
148 //_____________________________________________________________________________
149 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
150   TNamed(rec),
151
152   fRunLocalReconstruction(rec.fRunLocalReconstruction),
153   fRunVertexFinder(rec.fRunVertexFinder),
154   fRunTracking(rec.fRunTracking),
155   fFillESD(rec.fFillESD),
156   fGAliceFileName(rec.fGAliceFileName),
157   fInput(rec.fInput),
158   fStopOnError(rec.fStopOnError),
159   fCheckPointLevel(0),
160
161   fRunLoader(NULL),
162   fRawReader(NULL),
163   fITSLoader(NULL),
164   fITSVertexer(NULL),
165   fITSTracker(NULL),
166   fTPCLoader(NULL),
167   fTPCTracker(NULL),
168   fTRDLoader(NULL),
169   fTRDTracker(NULL),
170   fTOFLoader(NULL),
171   fTOFTracker(NULL),
172
173   fReconstructors(),
174   fOptions()
175 {
176 // copy constructor
177
178   for (Int_t i = 0; i < fOptions.GetEntriesFast(); i++) {
179     if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
180   }
181 }
182
183 //_____________________________________________________________________________
184 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
185 {
186 // assignment operator
187
188   this->~AliReconstruction();
189   new(this) AliReconstruction(rec);
190   return *this;
191 }
192
193 //_____________________________________________________________________________
194 AliReconstruction::~AliReconstruction()
195 {
196 // clean up
197
198   CleanUp();
199   fOptions.Delete();
200 }
201
202
203 //_____________________________________________________________________________
204 void AliReconstruction::SetGAliceFile(const char* fileName)
205 {
206 // set the name of the galice file
207
208   fGAliceFileName = fileName;
209 }
210
211 //_____________________________________________________________________________
212 void AliReconstruction::SetOption(const char* detector, const char* option)
213 {
214 // set options for the reconstruction of a detector
215
216   TObject* obj = fOptions.FindObject(detector);
217   if (obj) fOptions.Remove(obj);
218   fOptions.Add(new TNamed(detector, option));
219 }
220
221
222 //_____________________________________________________________________________
223 Bool_t AliReconstruction::Run(const char* input)
224 {
225 // run the reconstruction
226
227   // set the input
228   if (!input) input = fInput.Data();
229   TString fileName(input);
230   if (fileName.EndsWith("/")) {
231     fRawReader = new AliRawReaderFile(fileName);
232   } else if (fileName.EndsWith(".root")) {
233     fRawReader = new AliRawReaderRoot(fileName);
234   } else if (!fileName.IsNull()) {
235     fRawReader = new AliRawReaderDate(fileName);
236     fRawReader->SelectEvents(7);
237   }
238
239   // open the run loader
240   fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
241   if (!fRunLoader) {
242     Error("Run", "no run loader found in file %s", 
243           fGAliceFileName.Data());
244     CleanUp();
245     return kFALSE;
246   }
247   fRunLoader->LoadgAlice();
248   AliRun* aliRun = fRunLoader->GetAliRun();
249   if (!aliRun) {
250     Error("Run", "no gAlice object found in file %s", 
251           fGAliceFileName.Data());
252     CleanUp();
253     return kFALSE;
254   }
255   gAlice = aliRun;
256   AliTracker::SetFieldMap(gAlice->Field());
257
258   // load the reconstructor objects
259   TPluginManager* pluginManager = gROOT->GetPluginManager();
260   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
261     TString detName = fgkDetectorName[iDet];
262     TString recName = "Ali" + detName + "Reconstructor";
263     if (!gAlice->GetDetector(detName) && detName != "HLT") continue;
264
265     if(detName == "HLT") {
266       if (!gROOT->GetClass("AliLevel3")) {
267         gSystem->Load("libAliL3Src.so");
268         gSystem->Load("libAliL3Misc.so");
269         gSystem->Load("libAliL3Hough.so");
270         gSystem->Load("libAliL3Comp.so");
271       }
272     }
273
274     AliReconstructor* reconstructor = NULL;
275     // first check if a plugin is defined for the reconstructor
276     TPluginHandler* pluginHandler = 
277       pluginManager->FindHandler("AliReconstructor", detName);
278     // if not, but the reconstructor class is implemented, add a plugin for it
279     if (!pluginHandler && gROOT->GetClass(recName.Data())) {
280       Info("Run", "defining plugin for %s", recName.Data());
281       pluginManager->AddHandler("AliReconstructor", detName, 
282                                 recName, detName, recName + "()");
283       pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
284     }
285     if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
286       reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
287     }
288     // if there is no reconstructor class for the detector use the dummy one
289     if (!reconstructor && gAlice->GetDetector(detName)) {
290       Info("Run", "using dummy reconstructor for %s", detName.Data());
291       reconstructor = new AliDummyReconstructor(gAlice->GetDetector(detName));
292     }
293     if (reconstructor) {
294       TObject* obj = fOptions.FindObject(detName.Data());
295       if (obj) reconstructor->SetOption(obj->GetTitle());
296       fReconstructors.Add(reconstructor);
297     }
298   }
299
300   // local reconstruction
301   if (!fRunLocalReconstruction.IsNull()) {
302     if (!RunLocalReconstruction(fRunLocalReconstruction)) {
303       if (fStopOnError) {CleanUp(); return kFALSE;}
304     }
305   }
306   if (!fRunVertexFinder && !fRunTracking && fFillESD.IsNull()) return kTRUE;
307
308   // get vertexer
309   if (fRunVertexFinder && !CreateVertexer()) {
310     if (fStopOnError) {
311       CleanUp(); 
312       return kFALSE;
313     }
314   }
315
316   // get loaders and trackers
317   if (fRunTracking && !CreateTrackers()) {
318     if (fStopOnError) {
319       CleanUp(); 
320       return kFALSE;
321     }      
322   }
323
324   // create the ESD output file and tree
325   TFile* file = TFile::Open("AliESDs.root", "RECREATE");
326   if (!file->IsOpen()) {
327     Error("Run", "opening AliESDs.root failed");
328     if (fStopOnError) {CleanUp(file); return kFALSE;}    
329   }
330   AliESD* esd = new AliESD;
331   TTree* tree = new TTree("esdTree", "Tree with ESD objects");
332   tree->Branch("ESD", "AliESD", &esd);
333   delete esd;
334   gROOT->cd();
335
336   // loop over events
337   if (fRawReader) fRawReader->RewindEvents();
338   for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
339     Info("Run", "processing event %d", iEvent);
340     fRunLoader->GetEvent(iEvent);
341     if (fRawReader) fRawReader->NextEvent();
342
343     char fileName[256];
344     sprintf(fileName, "ESD_%d.%d_final.root", 
345             aliRun->GetRunNumber(), aliRun->GetEvNumber());
346     if (!gSystem->AccessPathName(fileName)) continue;
347
348     esd = new AliESD;
349     esd->SetRunNumber(aliRun->GetRunNumber());
350     esd->SetEventNumber(aliRun->GetEvNumber());
351     esd->SetMagneticField(aliRun->Field()->SolenoidField());
352
353     // vertex finder
354     if (fRunVertexFinder) {
355       if (!ReadESD(esd, "vertex")) {
356         if (!RunVertexFinder(esd)) {
357           if (fStopOnError) {CleanUp(file); return kFALSE;}
358         }
359         if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
360       }
361     }
362
363     // barrel tracking
364     if (fRunTracking) {
365       if (!ReadESD(esd, "tracking")) {
366         if (!RunTracking(esd)) {
367           if (fStopOnError) {CleanUp(file); return kFALSE;}
368         }
369         if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
370       }
371     }
372
373     // fill ESD
374     if (!fFillESD.IsNull()) {
375       if (!FillESD(esd, fFillESD)) {
376         if (fStopOnError) {CleanUp(file); return kFALSE;}
377       }
378     }
379
380     // combined PID
381     AliESDpid::MakePID(esd);
382     if (fCheckPointLevel > 1) WriteESD(esd, "PID");
383
384     // write ESD
385     tree->Fill();
386
387     if (fCheckPointLevel > 0) WriteESD(esd, "final");
388     delete esd;
389   }
390
391   file->cd();
392   tree->Write();
393   CleanUp(file);
394
395   return kTRUE;
396 }
397
398
399 //_____________________________________________________________________________
400 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
401 {
402 // run the local reconstruction
403
404   TStopwatch stopwatch;
405   stopwatch.Start();
406
407   TString detStr = detectors;
408   for (Int_t iDet = 0; iDet < fReconstructors.GetEntriesFast(); iDet++) {
409     AliReconstructor* reconstructor = 
410       (AliReconstructor*) fReconstructors[iDet];
411     TString detName = reconstructor->GetDetectorName();
412     if (IsSelected(detName, detStr)) {
413       Info("RunLocalReconstruction", "running reconstruction for %s", 
414            detName.Data());
415       TStopwatch stopwatchDet;
416       stopwatchDet.Start();
417       if (fRawReader) {
418         fRawReader->RewindEvents();
419         reconstructor->Reconstruct(fRunLoader, fRawReader);
420       } else {
421         reconstructor->Reconstruct(fRunLoader);
422       }
423       Info("RunLocalReconstruction", "execution time for %s:", detName.Data());
424       stopwatchDet.Print();
425     }
426   }
427
428   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
429     Error("RunLocalReconstruction", 
430           "the following detectors were not found: %s", detStr.Data());
431     if (fStopOnError) return kFALSE;
432   }
433
434   Info("RunLocalReconstruction", "execution time:");
435   stopwatch.Print();
436
437   return kTRUE;
438 }
439
440 //_____________________________________________________________________________
441 Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
442 {
443 // run the barrel tracking
444
445   TStopwatch stopwatch;
446   stopwatch.Start();
447
448   AliESDVertex* vertex = NULL;
449   Double_t vtxPos[3] = {0, 0, 0};
450   Double_t vtxErr[3] = {0.07, 0.07, 0.1};
451   TArrayF mcVertex(3); 
452   if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
453     fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
454     for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
455   }
456
457   if (fITSVertexer) {
458     Info("RunVertexFinder", "running the ITS vertex finder");
459     fITSVertexer->SetDebug(1);
460     vertex = fITSVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
461     if(!vertex){
462       Warning("RunVertexFinder","Vertex not found \n");
463       vertex = new AliESDVertex();
464     }
465     else {
466       vertex->SetTruePos(vtxPos);  // store also the vertex from MC
467     }
468
469   } else {
470     Info("RunVertexFinder", "getting the primary vertex from MC");
471     vertex = new AliESDVertex(vtxPos, vtxErr);
472   }
473
474   if (vertex) {
475     vertex->GetXYZ(vtxPos);
476     vertex->GetSigmaXYZ(vtxErr);
477   } else {
478     Warning("RunVertexFinder", "no vertex reconstructed");
479     vertex = new AliESDVertex(vtxPos, vtxErr);
480   }
481   esd->SetVertex(vertex);
482   if (fITSTracker) fITSTracker->SetVertex(vtxPos, vtxErr);
483   if (fTPCTracker) fTPCTracker->SetVertex(vtxPos, vtxErr);
484   if (fTRDTracker) fTRDTracker->SetVertex(vtxPos, vtxErr);
485   delete vertex;
486
487   Info("RunVertexFinder", "execution time:");
488   stopwatch.Print();
489
490   return kTRUE;
491 }
492
493 //_____________________________________________________________________________
494 Bool_t AliReconstruction::RunTracking(AliESD*& esd)
495 {
496 // run the barrel tracking
497
498   TStopwatch stopwatch;
499   stopwatch.Start();
500
501   if (!fTPCTracker) {
502     Error("RunTracking", "no TPC tracker");
503     return kFALSE;
504   }
505
506   // TPC tracking
507   Info("RunTracking", "TPC tracking");
508   fTPCLoader->LoadRecPoints("read");
509   TTree* tpcTree = fTPCLoader->TreeR();
510   if (!tpcTree) {
511     Error("RunTracking", "Can't get the TPC cluster tree");
512     return kFALSE;
513   }
514   fTPCTracker->LoadClusters(tpcTree);
515   if (fTPCTracker->Clusters2Tracks(esd) != 0) {
516     Error("RunTracking", "TPC Clusters2Tracks failed");
517     return kFALSE;
518   }
519   if (fCheckPointLevel > 1) WriteESD(esd, "TPC.tracking");
520
521   if (!fITSTracker) {
522     Warning("RunTracking", "no ITS tracker");
523   } else {
524
525     GetReconstructor("TPC")->FillESD(fRunLoader, esd); // preliminary
526     AliESDpid::MakePID(esd);                  // PID for the ITS tracker
527
528     // ITS tracking
529     Info("RunTracking", "ITS tracking");
530     fITSLoader->LoadRecPoints("read");
531     TTree* itsTree = fITSLoader->TreeR();
532     if (!itsTree) {
533       Error("RunTracking", "Can't get the ITS cluster tree");
534       return kFALSE;
535     }
536     fITSTracker->LoadClusters(itsTree);
537     if (fITSTracker->Clusters2Tracks(esd) != 0) {
538       Error("RunTracking", "ITS Clusters2Tracks failed");
539       return kFALSE;
540     }
541     if (fCheckPointLevel > 1) WriteESD(esd, "ITS.tracking");
542
543     if (!fTRDTracker) {
544       Warning("RunTracking", "no TRD tracker");
545     } else {
546       // ITS back propagation
547       Info("RunTracking", "ITS back propagation");
548       if (fITSTracker->PropagateBack(esd) != 0) {
549         Error("RunTracking", "ITS backward propagation failed");
550         return kFALSE;
551       }
552       if (fCheckPointLevel > 1) WriteESD(esd, "ITS.back");
553
554       // TPC back propagation
555       Info("RunTracking", "TPC back propagation");
556       if (fTPCTracker->PropagateBack(esd) != 0) {
557         Error("RunTracking", "TPC backward propagation failed");
558         return kFALSE;
559       }
560       if (fCheckPointLevel > 1) WriteESD(esd, "TPC.back");
561
562       // TRD back propagation
563       Info("RunTracking", "TRD back propagation");
564       fTRDLoader->LoadRecPoints("read");
565       TTree* trdTree = fTRDLoader->TreeR();
566       if (!trdTree) {
567         Error("RunTracking", "Can't get the TRD cluster tree");
568         return kFALSE;
569       }
570       fTRDTracker->LoadClusters(trdTree);
571       if (fTRDTracker->PropagateBack(esd) != 0) {
572         Error("RunTracking", "TRD backward propagation failed");
573         return kFALSE;
574       }
575       if (fCheckPointLevel > 1) WriteESD(esd, "TRD.back");
576
577       if (!fTOFTracker) {
578         Warning("RunTracking", "no TOF tracker");
579       } else {
580         // TOF back propagation
581         Info("RunTracking", "TOF back propagation");
582         fTOFLoader->LoadDigits("read");
583         TTree* tofTree = fTOFLoader->TreeD();
584         if (!tofTree) {
585           Error("RunTracking", "Can't get the TOF digits tree");
586           return kFALSE;
587         }
588         fTOFTracker->LoadClusters(tofTree);
589         if (fTOFTracker->PropagateBack(esd) != 0) {
590           Error("RunTracking", "TOF backward propagation failed");
591           return kFALSE;
592         }
593         if (fCheckPointLevel > 1) WriteESD(esd, "TOF.back");
594         fTOFTracker->UnloadClusters();
595         fTOFLoader->UnloadDigits();
596       }
597
598       // TRD inward refit
599       Info("RunTracking", "TRD inward refit");
600       if (fTRDTracker->RefitInward(esd) != 0) {
601         Error("RunTracking", "TRD inward refit failed");
602         return kFALSE;
603       }
604       if (fCheckPointLevel > 1) WriteESD(esd, "TRD.refit");
605       fTRDTracker->UnloadClusters();
606       fTRDLoader->UnloadRecPoints();
607     
608       // TPC inward refit
609       Info("RunTracking", "TPC inward refit");
610       if (fTPCTracker->RefitInward(esd) != 0) {
611         Error("RunTracking", "TPC inward refit failed");
612         return kFALSE;
613       }
614       if (fCheckPointLevel > 1) WriteESD(esd, "TPC.refit");
615     
616       // ITS inward refit
617       Info("RunTracking", "ITS inward refit");
618       if (fITSTracker->RefitInward(esd) != 0) {
619         Error("RunTracking", "ITS inward refit failed");
620         return kFALSE;
621       }
622       if (fCheckPointLevel > 1) WriteESD(esd, "ITS.refit");
623
624     }  // if TRD tracker
625     fITSTracker->UnloadClusters();
626     fITSLoader->UnloadRecPoints();
627
628   }  // if ITS tracker
629   fTPCTracker->UnloadClusters();
630   fTPCLoader->UnloadRecPoints();
631
632   Info("RunTracking", "execution time:");
633   stopwatch.Print();
634
635   return kTRUE;
636 }
637
638 //_____________________________________________________________________________
639 Bool_t AliReconstruction::FillESD(AliESD*& esd, const TString& detectors)
640 {
641 // fill the event summary data
642
643   TStopwatch stopwatch;
644   stopwatch.Start();
645
646   TString detStr = detectors;
647   for (Int_t iDet = 0; iDet < fReconstructors.GetEntriesFast(); iDet++) {
648     AliReconstructor* reconstructor = 
649       (AliReconstructor*) fReconstructors[iDet];
650     TString detName = reconstructor->GetDetectorName();
651     if (IsSelected(detName, detStr)) {
652       if (!ReadESD(esd, detName.Data())) {
653         Info("FillESD", "filling ESD for %s", detName.Data());
654         if (fRawReader) {
655           reconstructor->FillESD(fRunLoader, fRawReader, esd);
656         } else {
657           reconstructor->FillESD(fRunLoader, esd);
658         }
659         if (fCheckPointLevel > 2) WriteESD(esd, detName.Data());
660       }
661     }
662   }
663
664   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
665     Error("FillESD", "the following detectors were not found: %s", 
666           detStr.Data());
667     if (fStopOnError) return kFALSE;
668   }
669
670   Info("FillESD", "execution time:");
671   stopwatch.Print();
672
673   return kTRUE;
674 }
675
676
677 //_____________________________________________________________________________
678 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
679 {
680 // check whether detName is contained in detectors
681 // if yes, it is removed from detectors
682
683   // check if all detectors are selected
684   if ((detectors.CompareTo("ALL") == 0) ||
685       detectors.BeginsWith("ALL ") ||
686       detectors.EndsWith(" ALL") ||
687       detectors.Contains(" ALL ")) {
688     detectors = "ALL";
689     return kTRUE;
690   }
691
692   // search for the given detector
693   Bool_t result = kFALSE;
694   if ((detectors.CompareTo(detName) == 0) ||
695       detectors.BeginsWith(detName+" ") ||
696       detectors.EndsWith(" "+detName) ||
697       detectors.Contains(" "+detName+" ")) {
698     detectors.ReplaceAll(detName, "");
699     result = kTRUE;
700   }
701
702   // clean up the detectors string
703   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
704   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
705   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
706
707   return result;
708 }
709
710 //_____________________________________________________________________________
711 AliReconstructor* AliReconstruction::GetReconstructor(const char* detName) const
712 {
713 // get the reconstructor object for a detector
714
715   for (Int_t iDet = 0; iDet < fReconstructors.GetEntriesFast(); iDet++) {
716     AliReconstructor* reconstructor = 
717       (AliReconstructor*) fReconstructors[iDet];
718     if (strcmp(reconstructor->GetDetectorName(), detName) == 0) {
719       return reconstructor;
720     }
721   }
722   return NULL;
723 }
724
725 //_____________________________________________________________________________
726 Bool_t AliReconstruction::CreateVertexer()
727 {
728 // create the vertexer
729
730   fITSVertexer = NULL;
731   AliReconstructor* itsReconstructor = GetReconstructor("ITS");
732   if (itsReconstructor) {
733     fITSVertexer = itsReconstructor->CreateVertexer(fRunLoader);
734   }
735   if (!fITSVertexer) {
736     Warning("CreateVertexer", "couldn't create a vertexer for ITS");
737     if (fStopOnError) return kFALSE;
738   }
739
740   return kTRUE;
741 }
742
743 //_____________________________________________________________________________
744 Bool_t AliReconstruction::CreateTrackers()
745 {
746 // get the loaders and create the trackers
747
748   fITSTracker = NULL;
749   fITSLoader = fRunLoader->GetLoader("ITSLoader");
750   if (!fITSLoader) {
751     Warning("CreateTrackers", "no ITS loader found");
752     if (fStopOnError) return kFALSE;
753   } else {
754     AliReconstructor* itsReconstructor = GetReconstructor("ITS");
755     if (itsReconstructor) {
756       fITSTracker = itsReconstructor->CreateTracker(fRunLoader);
757     }
758     if (!fITSTracker) {
759       Warning("CreateTrackers", "couldn't create a tracker for ITS");
760       if (fStopOnError) return kFALSE;
761     }
762   }
763     
764   fTPCTracker = NULL;
765   fTPCLoader = fRunLoader->GetLoader("TPCLoader");
766   if (!fTPCLoader) {
767     Error("CreateTrackers", "no TPC loader found");
768     if (fStopOnError) return kFALSE;
769   } else {
770     AliReconstructor* tpcReconstructor = GetReconstructor("TPC");
771     if (tpcReconstructor) {
772       fTPCTracker = tpcReconstructor->CreateTracker(fRunLoader);
773     }
774     if (!fTPCTracker) {
775       Error("CreateTrackers", "couldn't create a tracker for TPC");
776       if (fStopOnError) return kFALSE;
777     }
778   }
779     
780   fTRDTracker = NULL;
781   fTRDLoader = fRunLoader->GetLoader("TRDLoader");
782   if (!fTRDLoader) {
783     Warning("CreateTrackers", "no TRD loader found");
784     if (fStopOnError) return kFALSE;
785   } else {
786     AliReconstructor* trdReconstructor = GetReconstructor("TRD");
787     if (trdReconstructor) {
788       fTRDTracker = trdReconstructor->CreateTracker(fRunLoader);
789     }
790     if (!fTRDTracker) {
791       Warning("CreateTrackers", "couldn't create a tracker for TRD");
792       if (fStopOnError) return kFALSE;
793     }
794   }
795     
796   fTOFTracker = NULL;
797   fTOFLoader = fRunLoader->GetLoader("TOFLoader");
798   if (!fTOFLoader) {
799     Warning("CreateTrackers", "no TOF loader found");
800     if (fStopOnError) return kFALSE;
801   } else {
802     AliReconstructor* tofReconstructor = GetReconstructor("TOF");
803     if (tofReconstructor) {
804       fTOFTracker = tofReconstructor->CreateTracker(fRunLoader);
805     }
806     if (!fTOFTracker) {
807       Warning("CreateTrackers", "couldn't create a tracker for TOF");
808       if (fStopOnError) return kFALSE;
809     }
810   }
811
812   return kTRUE;
813 }
814
815 //_____________________________________________________________________________
816 void AliReconstruction::CleanUp(TFile* file)
817 {
818 // delete trackers and the run loader and close and delete the file
819
820   fReconstructors.Delete();
821
822   delete fITSVertexer;
823   fITSVertexer = NULL;
824   delete fITSTracker;
825   fITSTracker = NULL;
826   delete fTPCTracker;
827   fTPCTracker = NULL;
828   delete fTRDTracker;
829   fTRDTracker = NULL;
830   delete fTOFTracker;
831   fTOFTracker = NULL;
832
833   delete fRunLoader;
834   fRunLoader = NULL;
835   delete fRawReader;
836   fRawReader = NULL;
837
838   if (file) {
839     file->Close();
840     delete file;
841   }
842 }
843
844
845 //_____________________________________________________________________________
846 Bool_t AliReconstruction::ReadESD(AliESD*& esd, const char* recStep) const
847 {
848 // read the ESD event from a file
849
850   if (!esd) return kFALSE;
851   char fileName[256];
852   sprintf(fileName, "ESD_%d.%d_%s.root", 
853           esd->GetRunNumber(), esd->GetEventNumber(), recStep);
854   if (gSystem->AccessPathName(fileName)) return kFALSE;
855
856   Info("ReadESD", "reading ESD from file %s", fileName);
857   TFile* file = TFile::Open(fileName);
858   if (!file || !file->IsOpen()) {
859     Error("ReadESD", "opening %s failed", fileName);
860     delete file;
861     return kFALSE;
862   }
863
864   gROOT->cd();
865   delete esd;
866   esd = (AliESD*) file->Get("ESD");
867   file->Close();
868   delete file;
869   return kTRUE;
870 }
871
872 //_____________________________________________________________________________
873 void AliReconstruction::WriteESD(AliESD* esd, const char* recStep) const
874 {
875 // write the ESD event to a file
876
877   if (!esd) return;
878   char fileName[256];
879   sprintf(fileName, "ESD_%d.%d_%s.root", 
880           esd->GetRunNumber(), esd->GetEventNumber(), recStep);
881
882   Info("WriteESD", "writing ESD to file %s", fileName);
883   TFile* file = TFile::Open(fileName, "recreate");
884   if (!file || !file->IsOpen()) {
885     Error("WriteESD", "opening %s failed", fileName);
886   } else {
887     esd->Write("ESD");
888     file->Close();
889   }
890   delete file;
891 }