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