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