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