Changes suggested by Effective C++ (F.Carminati)
[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 // By default all events are reconstructed. The reconstruction can be        //
42 // limited to a range of events by giving the index of the first and the     //
43 // last event as an argument to the Run method or by calling                 //
44 //                                                                           //
45 //   rec.SetEventRange(..., ...);                                            //
46 //                                                                           //
47 // The index -1 (default) can be used for the last event to indicate no      //
48 // upper limit of the event range.                                           //
49 //                                                                           //
50 // The name of the galice file can be changed from the default               //
51 // "galice.root" by passing it as argument to the AliReconstruction          //
52 // constructor or by                                                         //
53 //                                                                           //
54 //   rec.SetGAliceFile("...");                                               //
55 //                                                                           //
56 // The local reconstruction can be switched on or off for individual         //
57 // detectors by                                                              //
58 //                                                                           //
59 //   rec.SetRunLocalReconstruction("...");                                   //
60 //                                                                           //
61 // The argument is a (case sensitive) string with the names of the           //
62 // detectors separated by a space. The special string "ALL" selects all      //
63 // available detectors. This is the default.                                 //
64 //                                                                           //
65 // The reconstruction of the primary vertex position can be switched off by  //
66 //                                                                           //
67 //   rec.SetRunVertexFinder(kFALSE);                                         //
68 //                                                                           //
69 // The tracking and the creation of ESD tracks can be switched on for        //
70 // selected detectors by                                                     //
71 //                                                                           //
72 //   rec.SetRunTracking("...");                                              //
73 //                                                                           //
74 // The filling of additional ESD information can be steered by               //
75 //                                                                           //
76 //   rec.SetFillESD("...");                                                  //
77 //                                                                           //
78 // Again, for both methods the string specifies the list of detectors.       //
79 // The default is "ALL".                                                     //
80 //                                                                           //
81 // The call of the shortcut method                                           //
82 //                                                                           //
83 //   rec.SetRunReconstruction("...");                                        //
84 //                                                                           //
85 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and    //
86 // SetFillESD with the same detector selecting string as argument.           //
87 //                                                                           //
88 // The reconstruction requires digits or raw data as input. For the creation //
89 // of digits and raw data have a look at the class AliSimulation.            //
90 //                                                                           //
91 // For debug purposes the method SetCheckPointLevel can be used. If the      //
92 // argument is greater than 0, files with ESD events will be written after   //
93 // selected steps of the reconstruction for each event:                      //
94 //   level 1: after tracking and after filling of ESD (final)                //
95 //   level 2: in addition after each tracking step                           //
96 //   level 3: in addition after the filling of ESD for each detector         //
97 // If a final check point file exists for an event, this event will be       //
98 // skipped in the reconstruction. The tracking and the filling of ESD for    //
99 // a detector will be skipped as well, if the corresponding check point      //
100 // file exists. The ESD event will then be loaded from the file instead.     //
101 //                                                                           //
102 ///////////////////////////////////////////////////////////////////////////////
103
104 #include <TArrayF.h>
105 #include <TFile.h>
106 #include <TSystem.h>
107 #include <TROOT.h>
108 #include <TPluginManager.h>
109 #include <TStopwatch.h>
110
111 #include "AliReconstruction.h"
112 #include "AliReconstructor.h"
113 #include "AliLog.h"
114 #include "AliRunLoader.h"
115 #include "AliRun.h"
116 #include "AliRawReaderFile.h"
117 #include "AliRawReaderDate.h"
118 #include "AliRawReaderRoot.h"
119 #include "AliTracker.h"
120 #include "AliESD.h"
121 #include "AliESDVertex.h"
122 #include "AliVertexer.h"
123 #include "AliHeader.h"
124 #include "AliGenEventHeader.h"
125 #include "AliPID.h"
126 #include "AliESDpid.h"
127 #include "AliMagF.h"
128
129 ClassImp(AliReconstruction)
130
131
132 //_____________________________________________________________________________
133 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
134
135 //_____________________________________________________________________________
136 AliReconstruction::AliReconstruction(const char* gAliceFilename,
137                                      const char* name, const char* title) :
138   TNamed(name, title),
139
140   fRunLocalReconstruction("ALL"),
141   fRunVertexFinder(kTRUE),
142   fRunTracking("ALL"),
143   fFillESD("ALL"),
144   fGAliceFileName(gAliceFilename),
145   fInput(""),
146   fFirstEvent(0),
147   fLastEvent(-1),
148   fStopOnError(kFALSE),
149   fCheckPointLevel(0),
150   fOptions(),
151
152   fRunLoader(NULL),
153   fRawReader(NULL),
154
155   fVertexer(NULL)
156 {
157 // create reconstruction object with default parameters
158   
159   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
160     fReconstructor[iDet] = NULL;
161     fLoader[iDet] = NULL;
162     fTracker[iDet] = NULL;
163   }
164   //  AliPID::Init();
165 }
166
167 //_____________________________________________________________________________
168 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
169   TNamed(rec),
170
171   fRunLocalReconstruction(rec.fRunLocalReconstruction),
172   fRunVertexFinder(rec.fRunVertexFinder),
173   fRunTracking(rec.fRunTracking),
174   fFillESD(rec.fFillESD),
175   fGAliceFileName(rec.fGAliceFileName),
176   fInput(rec.fInput),
177   fFirstEvent(rec.fFirstEvent),
178   fLastEvent(rec.fLastEvent),
179   fStopOnError(rec.fStopOnError),
180   fCheckPointLevel(0),
181   fOptions(),
182
183   fRunLoader(NULL),
184   fRawReader(NULL),
185
186   fVertexer(NULL)
187 {
188 // copy constructor
189
190   for (Int_t i = 0; i < fOptions.GetEntriesFast(); i++) {
191     if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
192   }
193   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
194     fReconstructor[iDet] = NULL;
195     fLoader[iDet] = NULL;
196     fTracker[iDet] = NULL;
197   }
198 }
199
200 //_____________________________________________________________________________
201 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
202 {
203 // assignment operator
204
205   this->~AliReconstruction();
206   new(this) AliReconstruction(rec);
207   return *this;
208 }
209
210 //_____________________________________________________________________________
211 AliReconstruction::~AliReconstruction()
212 {
213 // clean up
214
215   CleanUp();
216   fOptions.Delete();
217 }
218
219
220 //_____________________________________________________________________________
221 void AliReconstruction::SetGAliceFile(const char* fileName)
222 {
223 // set the name of the galice file
224
225   fGAliceFileName = fileName;
226 }
227
228 //_____________________________________________________________________________
229 void AliReconstruction::SetOption(const char* detector, const char* option)
230 {
231 // set options for the reconstruction of a detector
232
233   TObject* obj = fOptions.FindObject(detector);
234   if (obj) fOptions.Remove(obj);
235   fOptions.Add(new TNamed(detector, option));
236 }
237
238
239 //_____________________________________________________________________________
240 Bool_t AliReconstruction::Run(const char* input,
241                               Int_t firstEvent, Int_t lastEvent)
242 {
243 // run the reconstruction
244
245   // set the input
246   if (!input) input = fInput.Data();
247   TString fileName(input);
248   if (fileName.EndsWith("/")) {
249     fRawReader = new AliRawReaderFile(fileName);
250   } else if (fileName.EndsWith(".root")) {
251     fRawReader = new AliRawReaderRoot(fileName);
252   } else if (!fileName.IsNull()) {
253     fRawReader = new AliRawReaderDate(fileName);
254     fRawReader->SelectEvents(7);
255   }
256
257   // get the run loader
258   if (!InitRunLoader()) return kFALSE;
259
260   // local reconstruction
261   if (!fRunLocalReconstruction.IsNull()) {
262     if (!RunLocalReconstruction(fRunLocalReconstruction)) {
263       if (fStopOnError) {CleanUp(); return kFALSE;}
264     }
265   }
266 //  if (!fRunVertexFinder && fRunTracking.IsNull() && 
267 //      fFillESD.IsNull()) return kTRUE;
268
269   // get vertexer
270   if (fRunVertexFinder && !CreateVertexer()) {
271     if (fStopOnError) {
272       CleanUp(); 
273       return kFALSE;
274     }
275   }
276
277   // get trackers
278   if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
279     if (fStopOnError) {
280       CleanUp(); 
281       return kFALSE;
282     }      
283   }
284
285   // get the possibly already existing ESD file and tree
286   AliESD* esd = new AliESD;
287   TFile* fileOld = NULL;
288   TTree* treeOld = NULL;
289   if (!gSystem->AccessPathName("AliESDs.root")){
290     gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
291     fileOld = TFile::Open("AliESDs.old.root");
292     if (fileOld && fileOld->IsOpen()) {
293       treeOld = (TTree*) fileOld->Get("esdTree");
294       if (treeOld) treeOld->SetBranchAddress("ESD", &esd);
295     }
296   }
297
298   // create the ESD output file and tree
299   TFile* file = TFile::Open("AliESDs.root", "RECREATE");
300   if (!file->IsOpen()) {
301     AliError("opening AliESDs.root failed");
302     if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}    
303   }
304   TTree* tree = new TTree("esdTree", "Tree with ESD objects");
305   tree->Branch("ESD", "AliESD", &esd);
306   delete esd;
307   esd = NULL;
308   gROOT->cd();
309
310   // loop over events
311   if (fRawReader) fRawReader->RewindEvents();
312   
313   for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
314     if (fRawReader) fRawReader->NextEvent();
315     if ((iEvent < firstEvent) || ((lastEvent >= 0) && (iEvent > lastEvent))) {
316       // copy old ESD to the new one
317       if (treeOld) {
318         treeOld->SetBranchAddress("ESD", &esd);
319         treeOld->GetEntry(iEvent);
320       }
321       tree->Fill();
322       continue;
323     }
324
325     AliInfo(Form("processing event %d", iEvent));
326     fRunLoader->GetEvent(iEvent);
327
328     char fileName[256];
329     sprintf(fileName, "ESD_%d.%d_final.root", 
330             fRunLoader->GetHeader()->GetRun(), 
331             fRunLoader->GetHeader()->GetEventNrInRun());
332     if (!gSystem->AccessPathName(fileName)) continue;
333
334     // local reconstruction
335     if (!fRunLocalReconstruction.IsNull()) {
336       if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
337         if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
338       }
339     }
340
341     esd = new AliESD;
342     esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
343     esd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
344     if (gAlice) {
345       esd->SetMagneticField(gAlice->Field()->SolenoidField());
346     } else {
347       // ???
348     }
349
350     // vertex finder
351     if (fRunVertexFinder) {
352       if (!ReadESD(esd, "vertex")) {
353         if (!RunVertexFinder(esd)) {
354           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
355         }
356         if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
357       }
358     }
359
360     // barrel tracking
361     if (!fRunTracking.IsNull()) {
362       if (!ReadESD(esd, "tracking")) {
363         if (!RunTracking(esd)) {
364           if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
365         }
366         if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
367       }
368     }
369
370     // fill ESD
371     if (!fFillESD.IsNull()) {
372       if (!FillESD(esd, fFillESD)) {
373         if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
374       }
375     }
376
377     // combined PID
378     AliESDpid::MakePID(esd);
379     if (fCheckPointLevel > 1) WriteESD(esd, "PID");
380
381     // write ESD
382     tree->Fill();
383
384     if (fCheckPointLevel > 0) WriteESD(esd, "final");
385     delete esd;
386     esd = NULL;
387   }
388
389   file->cd();
390   tree->Write();
391   CleanUp(file, fileOld);
392
393   return kTRUE;
394 }
395
396
397 //_____________________________________________________________________________
398 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
399 {
400 // run the local reconstruction
401
402   TStopwatch stopwatch;
403   stopwatch.Start();
404
405   TString detStr = detectors;
406   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
407     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
408     AliReconstructor* reconstructor = GetReconstructor(iDet);
409     if (!reconstructor) continue;
410     if (reconstructor->HasLocalReconstruction()) continue;
411
412     AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
413     TStopwatch stopwatchDet;
414     stopwatchDet.Start();
415     if (fRawReader) {
416       fRawReader->RewindEvents();
417       reconstructor->Reconstruct(fRunLoader, fRawReader);
418     } else {
419       reconstructor->Reconstruct(fRunLoader);
420     }
421     AliInfo(Form("execution time for %s:", fgkDetectorName[iDet]));
422     ToAliInfo(stopwatchDet.Print());
423   }
424
425   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
426     AliError(Form("the following detectors were not found: %s",
427                   detStr.Data()));
428     if (fStopOnError) return kFALSE;
429   }
430
431   AliInfo("execution time:");
432   ToAliInfo(stopwatch.Print());
433
434   return kTRUE;
435 }
436
437 //_____________________________________________________________________________
438 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
439 {
440 // run the local reconstruction
441
442   TStopwatch stopwatch;
443   stopwatch.Start();
444
445   TString detStr = detectors;
446   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
447     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
448     AliReconstructor* reconstructor = GetReconstructor(iDet);
449     if (!reconstructor) continue;
450     AliLoader* loader = fLoader[iDet];
451
452     // conversion of digits
453     if (fRawReader && reconstructor->HasDigitConversion()) {
454       AliInfo(Form("converting raw data digits into root objects for %s", 
455                    fgkDetectorName[iDet]));
456       TStopwatch stopwatchDet;
457       stopwatchDet.Start();
458       loader->LoadDigits("update");
459       loader->CleanDigits();
460       loader->MakeDigitsContainer();
461       TTree* digitsTree = loader->TreeD();
462       reconstructor->ConvertDigits(fRawReader, digitsTree);
463       loader->WriteDigits("OVERWRITE");
464       loader->UnloadDigits();
465       AliDebug(1, Form("execution time for %s:", fgkDetectorName[iDet]));
466       ToAliDebug(1, stopwatchDet.Print());
467     }
468
469     // local reconstruction
470     if (!reconstructor->HasLocalReconstruction()) continue;
471     AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
472     TStopwatch stopwatchDet;
473     stopwatchDet.Start();
474     loader->LoadRecPoints("update");
475     loader->CleanRecPoints();
476     loader->MakeRecPointsContainer();
477     TTree* clustersTree = loader->TreeR();
478     if (fRawReader && !reconstructor->HasDigitConversion()) {
479       reconstructor->Reconstruct(fRawReader, clustersTree);
480     } else {
481       loader->LoadDigits("read");
482       TTree* digitsTree = loader->TreeD();
483       if (!digitsTree) {
484         AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
485         if (fStopOnError) return kFALSE;
486       } else {
487         reconstructor->Reconstruct(digitsTree, clustersTree);
488       }
489       loader->UnloadDigits();
490     }
491     loader->WriteRecPoints("OVERWRITE");
492     loader->UnloadRecPoints();
493     AliDebug(1, Form("execution time for %s:", fgkDetectorName[iDet]));
494     ToAliDebug(1, stopwatchDet.Print());
495   }
496
497   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
498     AliError(Form("the following detectors were not found: %s",
499                   detStr.Data()));
500     if (fStopOnError) return kFALSE;
501   }
502
503   AliInfo("execution time:");
504   ToAliInfo(stopwatch.Print());
505
506   return kTRUE;
507 }
508
509 //_____________________________________________________________________________
510 Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
511 {
512 // run the barrel tracking
513
514   TStopwatch stopwatch;
515   stopwatch.Start();
516
517   AliESDVertex* vertex = NULL;
518   Double_t vtxPos[3] = {0, 0, 0};
519   Double_t vtxErr[3] = {0.07, 0.07, 0.1};
520   TArrayF mcVertex(3); 
521   if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
522     fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
523     for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
524   }
525
526   if (fVertexer) {
527     AliInfo("running the ITS vertex finder");
528     if (fLoader[0]) fLoader[0]->LoadRecPoints();
529     vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
530     if (fLoader[0]) fLoader[0]->UnloadRecPoints();
531     if(!vertex){
532       AliWarning("Vertex not found");
533       vertex = new AliESDVertex();
534     }
535     else {
536       vertex->SetTruePos(vtxPos);  // store also the vertex from MC
537     }
538
539   } else {
540     AliInfo("getting the primary vertex from MC");
541     vertex = new AliESDVertex(vtxPos, vtxErr);
542   }
543
544   if (vertex) {
545     vertex->GetXYZ(vtxPos);
546     vertex->GetSigmaXYZ(vtxErr);
547   } else {
548     AliWarning("no vertex reconstructed");
549     vertex = new AliESDVertex(vtxPos, vtxErr);
550   }
551   esd->SetVertex(vertex);
552   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
553     if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
554   }  
555   delete vertex;
556
557   AliInfo("execution time:");
558   ToAliInfo(stopwatch.Print());
559
560   return kTRUE;
561 }
562
563 //_____________________________________________________________________________
564 Bool_t AliReconstruction::RunTracking(AliESD*& esd)
565 {
566 // run the barrel tracking
567
568   TStopwatch stopwatch;
569   stopwatch.Start();
570
571   AliInfo("running tracking");
572
573   // pass 1: TPC + ITS inwards
574   for (Int_t iDet = 1; iDet >= 0; iDet--) {
575     if (!fTracker[iDet]) continue;
576     AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
577
578     // load clusters
579     fLoader[iDet]->LoadRecPoints("read");
580     TTree* tree = fLoader[iDet]->TreeR();
581     if (!tree) {
582       AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
583       return kFALSE;
584     }
585     fTracker[iDet]->LoadClusters(tree);
586
587     // run tracking
588     if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
589       AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
590       return kFALSE;
591     }
592     if (fCheckPointLevel > 1) {
593       WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
594     }
595     // preliminary PID in TPC needed by the ITS tracker
596     if (iDet == 1) {
597       GetReconstructor(1)->FillESD(fRunLoader, esd);
598       GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
599       AliESDpid::MakePID(esd);
600     }
601   }
602
603   // pass 2: ALL backwards
604   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
605     if (!fTracker[iDet]) continue;
606     AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
607
608     // load clusters
609     if (iDet > 1) {     // all except ITS, TPC
610       TTree* tree = NULL;
611       if (iDet == 3) {   // TOF
612         fLoader[iDet]->LoadDigits("read");
613         tree = fLoader[iDet]->TreeD();
614       } else {
615         fLoader[iDet]->LoadRecPoints("read");
616         tree = fLoader[iDet]->TreeR();
617       }
618       if (!tree) {
619         AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
620         return kFALSE;
621       }
622       fTracker[iDet]->LoadClusters(tree);
623     }
624
625     // run tracking
626     if (fTracker[iDet]->PropagateBack(esd) != 0) {
627       AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
628       return kFALSE;
629     }
630     if (fCheckPointLevel > 1) {
631       WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
632     }
633
634     // unload clusters
635     if (iDet > 2) {     // all except ITS, TPC, TRD
636       fTracker[iDet]->UnloadClusters();
637       if (iDet == 3) {   // TOF
638         fLoader[iDet]->UnloadDigits();
639       } else {
640         fLoader[iDet]->UnloadRecPoints();
641       }
642     }
643   }
644
645   // pass 3: TRD + TPC + ITS refit inwards
646   for (Int_t iDet = 2; iDet >= 0; iDet--) {
647     if (!fTracker[iDet]) continue;
648     AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
649
650     // run tracking
651     if (fTracker[iDet]->RefitInward(esd) != 0) {
652       AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
653       return kFALSE;
654     }
655     if (fCheckPointLevel > 1) {
656       WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
657     }
658
659     // unload clusters
660     fTracker[iDet]->UnloadClusters();
661     fLoader[iDet]->UnloadRecPoints();
662   }
663
664   AliInfo("execution time:");
665   ToAliInfo(stopwatch.Print());
666
667   return kTRUE;
668 }
669
670 //_____________________________________________________________________________
671 Bool_t AliReconstruction::FillESD(AliESD*& esd, const TString& detectors)
672 {
673 // fill the event summary data
674
675   TStopwatch stopwatch;
676   stopwatch.Start();
677   AliInfo("filling ESD");
678
679   TString detStr = detectors;
680   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
681     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
682     AliReconstructor* reconstructor = GetReconstructor(iDet);
683     if (!reconstructor) continue;
684
685     if (!ReadESD(esd, fgkDetectorName[iDet])) {
686       AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
687       TTree* clustersTree = NULL;
688       if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
689         fLoader[iDet]->LoadRecPoints("read");
690         clustersTree = fLoader[iDet]->TreeR();
691         if (!clustersTree) {
692           AliError(Form("Can't get the %s clusters tree", 
693                         fgkDetectorName[iDet]));
694           if (fStopOnError) return kFALSE;
695         }
696       }
697       if (fRawReader && !reconstructor->HasDigitConversion()) {
698         reconstructor->FillESD(fRawReader, clustersTree, esd);
699       } else {
700         TTree* digitsTree = NULL;
701         if (fLoader[iDet]) {
702           fLoader[iDet]->LoadDigits("read");
703           digitsTree = fLoader[iDet]->TreeD();
704           if (!digitsTree) {
705             AliError(Form("Can't get the %s digits tree", 
706                           fgkDetectorName[iDet]));
707             if (fStopOnError) return kFALSE;
708           }
709         }
710         reconstructor->FillESD(digitsTree, clustersTree, esd);
711         if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
712       }
713       if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
714         fLoader[iDet]->UnloadRecPoints();
715       }
716
717       if (fRawReader) {
718         reconstructor->FillESD(fRunLoader, fRawReader, esd);
719       } else {
720         reconstructor->FillESD(fRunLoader, esd);
721       }
722       if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
723     }
724   }
725
726   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
727     AliError(Form("the following detectors were not found: %s", 
728                   detStr.Data()));
729     if (fStopOnError) return kFALSE;
730   }
731
732   AliInfo("execution time:");
733   ToAliInfo(stopwatch.Print());
734
735   return kTRUE;
736 }
737
738
739 //_____________________________________________________________________________
740 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
741 {
742 // check whether detName is contained in detectors
743 // if yes, it is removed from detectors
744
745   // check if all detectors are selected
746   if ((detectors.CompareTo("ALL") == 0) ||
747       detectors.BeginsWith("ALL ") ||
748       detectors.EndsWith(" ALL") ||
749       detectors.Contains(" ALL ")) {
750     detectors = "ALL";
751     return kTRUE;
752   }
753
754   // search for the given detector
755   Bool_t result = kFALSE;
756   if ((detectors.CompareTo(detName) == 0) ||
757       detectors.BeginsWith(detName+" ") ||
758       detectors.EndsWith(" "+detName) ||
759       detectors.Contains(" "+detName+" ")) {
760     detectors.ReplaceAll(detName, "");
761     result = kTRUE;
762   }
763
764   // clean up the detectors string
765   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
766   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
767   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
768
769   return result;
770 }
771
772 //_____________________________________________________________________________
773 Bool_t AliReconstruction::InitRunLoader()
774 {
775 // get or create the run loader
776
777   if (gAlice) delete gAlice;
778   gAlice = NULL;
779
780   if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
781     // load all base libraries to get the loader classes
782     TString libs = gSystem->GetLibraries();
783     for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
784       TString detName = fgkDetectorName[iDet];
785       if (detName == "HLT") continue;
786       if (libs.Contains("lib" + detName + "base.so")) continue;
787       gSystem->Load("lib" + detName + "base.so");
788     }
789     fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
790     if (!fRunLoader) {
791       AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
792       CleanUp();
793       return kFALSE;
794     }
795     fRunLoader->CdGAFile();
796     if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
797       if (fRunLoader->LoadgAlice() == 0) {
798         gAlice = fRunLoader->GetAliRun();
799         AliTracker::SetFieldMap(gAlice->Field());
800       }
801     }
802     if (!gAlice && !fRawReader) {
803       AliError(Form("no gAlice object found in file %s",
804                     fGAliceFileName.Data()));
805       CleanUp();
806       return kFALSE;
807     }
808
809   } else {               // galice.root does not exist
810     if (!fRawReader) {
811       AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
812       CleanUp();
813       return kFALSE;
814     }
815     fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
816                                     AliConfig::GetDefaultEventFolderName(),
817                                     "recreate");
818     if (!fRunLoader) {
819       AliError(Form("could not create run loader in file %s", 
820                     fGAliceFileName.Data()));
821       CleanUp();
822       return kFALSE;
823     }
824     fRunLoader->MakeTree("E");
825     Int_t iEvent = 0;
826     while (fRawReader->NextEvent()) {
827       fRunLoader->SetEventNumber(iEvent);
828       fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(), 
829                                      iEvent, iEvent);
830       fRunLoader->MakeTree("H");
831       fRunLoader->TreeE()->Fill();
832       iEvent++;
833     }
834     fRawReader->RewindEvents();
835     fRunLoader->WriteHeader("OVERWRITE");
836     fRunLoader->CdGAFile();
837     fRunLoader->Write(0, TObject::kOverwrite);
838 //    AliTracker::SetFieldMap(???);
839   }
840
841   return kTRUE;
842 }
843
844 //_____________________________________________________________________________
845 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
846 {
847 // get the reconstructor object and the loader for a detector
848
849   if (fReconstructor[iDet]) return fReconstructor[iDet];
850
851   // load the reconstructor object
852   TPluginManager* pluginManager = gROOT->GetPluginManager();
853   TString detName = fgkDetectorName[iDet];
854   TString recName = "Ali" + detName + "Reconstructor";
855   if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
856
857   if (detName == "HLT") {
858     if (!gROOT->GetClass("AliLevel3")) {
859       gSystem->Load("libAliL3Src.so");
860       gSystem->Load("libAliL3Misc.so");
861       gSystem->Load("libAliL3Hough.so");
862       gSystem->Load("libAliL3Comp.so");
863     }
864   }
865
866   AliReconstructor* reconstructor = NULL;
867   // first check if a plugin is defined for the reconstructor
868   TPluginHandler* pluginHandler = 
869     pluginManager->FindHandler("AliReconstructor", detName);
870   // if not, add a plugin for it
871   if (!pluginHandler) {
872     AliDebug(1, Form("defining plugin for %s", recName.Data()));
873     TString libs = gSystem->GetLibraries();
874     if (libs.Contains("lib" + detName + "base.so") ||
875         (gSystem->Load("lib" + detName + "base.so") >= 0)) {
876       pluginManager->AddHandler("AliReconstructor", detName, 
877                                 recName, detName + "rec", recName + "()");
878     } else {
879       pluginManager->AddHandler("AliReconstructor", detName, 
880                                 recName, detName, recName + "()");
881     }
882     pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
883   }
884   if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
885     reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
886   }
887   if (reconstructor) {
888     TObject* obj = fOptions.FindObject(detName.Data());
889     if (obj) reconstructor->SetOption(obj->GetTitle());
890     reconstructor->Init(fRunLoader);
891     fReconstructor[iDet] = reconstructor;
892   }
893
894   // get or create the loader
895   if (detName != "HLT") {
896     fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
897     if (!fLoader[iDet]) {
898       AliConfig::Instance()
899         ->CreateDetectorFolders(fRunLoader->GetEventFolder(), 
900                                 detName, detName);
901       // first check if a plugin is defined for the loader
902       TPluginHandler* pluginHandler = 
903         pluginManager->FindHandler("AliLoader", detName);
904       // if not, add a plugin for it
905       if (!pluginHandler) {
906         TString loaderName = "Ali" + detName + "Loader";
907         AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
908         pluginManager->AddHandler("AliLoader", detName, 
909                                   loaderName, detName + "base", 
910                                   loaderName + "(const char*, TFolder*)");
911         pluginHandler = pluginManager->FindHandler("AliLoader", detName);
912       }
913       if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
914         fLoader[iDet] = 
915           (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(), 
916                                                  fRunLoader->GetEventFolder());
917       }
918       if (!fLoader[iDet]) {   // use default loader
919         fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
920       }
921       if (!fLoader[iDet]) {
922         AliWarning(Form("couldn't get loader for %s", detName.Data()));
923         if (fStopOnError) return NULL;
924       } else {
925         fRunLoader->AddLoader(fLoader[iDet]);
926         fRunLoader->CdGAFile();
927         if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
928         fRunLoader->Write(0, TObject::kOverwrite);
929       }
930     }
931   }
932       
933   return reconstructor;
934 }
935
936 //_____________________________________________________________________________
937 Bool_t AliReconstruction::CreateVertexer()
938 {
939 // create the vertexer
940
941   fVertexer = NULL;
942   AliReconstructor* itsReconstructor = GetReconstructor(0);
943   if (itsReconstructor) {
944     fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
945   }
946   if (!fVertexer) {
947     AliWarning("couldn't create a vertexer for ITS");
948     if (fStopOnError) return kFALSE;
949   }
950
951   return kTRUE;
952 }
953
954 //_____________________________________________________________________________
955 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
956 {
957 // create the trackers
958
959   TString detStr = detectors;
960   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
961     if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
962     AliReconstructor* reconstructor = GetReconstructor(iDet);
963     if (!reconstructor) continue;
964     TString detName = fgkDetectorName[iDet];
965     if (detName == "HLT") continue;
966
967     fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
968     if (!fTracker[iDet] && (iDet < 7)) {
969       AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
970       if (fStopOnError) return kFALSE;
971     }
972   }
973
974   return kTRUE;
975 }
976
977 //_____________________________________________________________________________
978 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
979 {
980 // delete trackers and the run loader and close and delete the file
981
982   for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
983     delete fReconstructor[iDet];
984     fReconstructor[iDet] = NULL;
985     fLoader[iDet] = NULL;
986     delete fTracker[iDet];
987     fTracker[iDet] = NULL;
988   }
989   delete fVertexer;
990   fVertexer = NULL;
991
992   delete fRunLoader;
993   fRunLoader = NULL;
994   delete fRawReader;
995   fRawReader = NULL;
996
997   if (file) {
998     file->Close();
999     delete file;
1000   }
1001
1002   if (fileOld) {
1003     fileOld->Close();
1004     delete fileOld;
1005     gSystem->Unlink("AliESDs.old.root");
1006   }
1007 }
1008
1009
1010 //_____________________________________________________________________________
1011 Bool_t AliReconstruction::ReadESD(AliESD*& esd, const char* recStep) const
1012 {
1013 // read the ESD event from a file
1014
1015   if (!esd) return kFALSE;
1016   char fileName[256];
1017   sprintf(fileName, "ESD_%d.%d_%s.root", 
1018           esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1019   if (gSystem->AccessPathName(fileName)) return kFALSE;
1020
1021   AliDebug(1, Form("reading ESD from file %s", fileName));
1022   TFile* file = TFile::Open(fileName);
1023   if (!file || !file->IsOpen()) {
1024     AliError(Form("opening %s failed", fileName));
1025     delete file;
1026     return kFALSE;
1027   }
1028
1029   gROOT->cd();
1030   delete esd;
1031   esd = (AliESD*) file->Get("ESD");
1032   file->Close();
1033   delete file;
1034   return kTRUE;
1035 }
1036
1037 //_____________________________________________________________________________
1038 void AliReconstruction::WriteESD(AliESD* esd, const char* recStep) const
1039 {
1040 // write the ESD event to a file
1041
1042   if (!esd) return;
1043   char fileName[256];
1044   sprintf(fileName, "ESD_%d.%d_%s.root", 
1045           esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1046
1047   AliDebug(1, Form("writing ESD to file %s", fileName));
1048   TFile* file = TFile::Open(fileName, "recreate");
1049   if (!file || !file->IsOpen()) {
1050     AliError(Form("opening %s failed", fileName));
1051   } else {
1052     esd->Write("ESD");
1053     file->Close();
1054   }
1055   delete file;
1056 }