Corrected memory leaks
[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 ///////////////////////////////////////////////////////////////////////////////
58
59
60 #include "AliReconstruction.h"
61 #include "AliRunLoader.h"
62 #include "AliRun.h"
63 #include "AliModule.h"
64 #include "AliDetector.h"
65 #include "AliTracker.h"
66 #include "AliESD.h"
67 #include "AliHeader.h"
68 #include "AliGenEventHeader.h"
69 #include "AliESDpid.h"
70 #include "AliMagF.h"
71 #include <TArrayF.h>
72
73
74 ClassImp(AliReconstruction)
75
76
77 //_____________________________________________________________________________
78 AliReconstruction::AliReconstruction(const char* gAliceFilename,
79                                      const char* name, const char* title) :
80   TNamed(name, title),
81
82   fRunReconstruction("ALL"),
83   fRunTracking(kTRUE),
84   fFillESD("ALL"),
85   fGAliceFileName(gAliceFilename),
86   fStopOnError(kFALSE),
87
88   fRunLoader(NULL),
89   fITSLoader(NULL),
90   fITSTracker(NULL),
91   fTPCLoader(NULL),
92   fTPCTracker(NULL),
93   fTRDLoader(NULL),
94   fTRDTracker(NULL),
95   fTOFLoader(NULL),
96   fTOFTracker(NULL)
97 {
98 // create reconstruction object with default parameters
99
100 }
101
102 //_____________________________________________________________________________
103 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
104   TNamed(rec),
105
106   fRunReconstruction(rec.fRunReconstruction),
107   fRunTracking(rec.fRunTracking),
108   fFillESD(rec.fFillESD),
109   fGAliceFileName(rec.fGAliceFileName),
110   fStopOnError(rec.fStopOnError),
111
112   fRunLoader(NULL),
113   fITSLoader(NULL),
114   fITSTracker(NULL),
115   fTPCLoader(NULL),
116   fTPCTracker(NULL),
117   fTRDLoader(NULL),
118   fTRDTracker(NULL),
119   fTOFLoader(NULL),
120   fTOFTracker(NULL)
121 {
122 // copy constructor
123
124 }
125
126 //_____________________________________________________________________________
127 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
128 {
129 // assignment operator
130
131   this->~AliReconstruction();
132   new(this) AliReconstruction(rec);
133   return *this;
134 }
135
136 //_____________________________________________________________________________
137 AliReconstruction::~AliReconstruction()
138 {
139 // clean up
140
141   CleanUp();
142 }
143
144
145 //_____________________________________________________________________________
146 void AliReconstruction::SetGAliceFile(const char* fileName)
147 {
148 // set the name of the galice file
149
150   fGAliceFileName = fileName;
151 }
152
153
154 //_____________________________________________________________________________
155 Bool_t AliReconstruction::Run()
156 {
157 // run the reconstruction
158
159   // open the run loader
160   fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
161   if (!fRunLoader) {
162     Error("Run", "no run loader found in file %s", 
163           fGAliceFileName.Data());
164     CleanUp();
165     return kFALSE;
166   }
167   fRunLoader->LoadgAlice();
168   AliRun* aliRun = fRunLoader->GetAliRun();
169   if (!aliRun) {
170     Error("Run", "no gAlice object found in file %s", 
171           fGAliceFileName.Data());
172     CleanUp();
173     return kFALSE;
174   }
175   gAlice = aliRun;
176
177   // local reconstruction
178   if (!fRunReconstruction.IsNull()) {
179     if (!RunReconstruction(fRunReconstruction)) {
180       if (fStopOnError) {CleanUp(); return kFALSE;}
181     }
182   }
183   if (!fRunTracking && fFillESD.IsNull()) return kTRUE;
184
185   if (fRunTracking) {
186     // get loaders and trackers
187     fITSLoader = fRunLoader->GetLoader("ITSLoader");
188     if (!fITSLoader) {
189       Error("Run", "no ITS loader found");
190       if (fStopOnError) {CleanUp(); return kFALSE;}
191     }
192     fITSTracker = NULL;
193     if (aliRun->GetDetector("ITS")) {
194       fITSTracker = aliRun->GetDetector("ITS")->CreateTracker();
195     }
196     if (!fITSTracker) {
197       Error("Run", "couldn't create a tracker for ITS");
198       if (fStopOnError) {CleanUp(); return kFALSE;}
199     }
200     
201     fTPCLoader = fRunLoader->GetLoader("TPCLoader");
202     if (!fTPCLoader) {
203       Error("Run", "no TPC loader found");
204       if (fStopOnError) {CleanUp(); return kFALSE;}
205     }
206     fTPCTracker = NULL;
207     if (aliRun->GetDetector("TPC")) {
208       fTPCTracker = aliRun->GetDetector("TPC")->CreateTracker();
209     }
210     if (!fTPCTracker) {
211       Error("Run", "couldn't create a tracker for TPC");
212       if (fStopOnError) {CleanUp(); return kFALSE;}
213     }
214     
215     fTRDLoader = fRunLoader->GetLoader("TRDLoader");
216     if (!fTRDLoader) {
217       Error("Run", "no TRD loader found");
218       if (fStopOnError) {CleanUp(); return kFALSE;}
219     }
220     fTRDTracker = NULL;
221     if (aliRun->GetDetector("TRD")) {
222       fTRDTracker = aliRun->GetDetector("TRD")->CreateTracker();
223     }
224     if (!fTRDTracker) {
225       Error("Run", "couldn't create a tracker for TRD");
226       if (fStopOnError) {CleanUp(); return kFALSE;}
227     }
228     
229     fTOFLoader = fRunLoader->GetLoader("TOFLoader");
230     if (!fTOFLoader) {
231       Error("Run", "no TOF loader found");
232       if (fStopOnError) {CleanUp(); return kFALSE;}
233     }
234     fTOFTracker = NULL;
235     if (aliRun->GetDetector("TOF")) {
236       fTOFTracker = aliRun->GetDetector("TOF")->CreateTracker();
237     }
238     if (!fTOFTracker) {
239       Error("Run", "couldn't create a tracker for TOF");
240       if (fStopOnError) {CleanUp(); return kFALSE;}
241     }
242   }
243   // create the ESD output file
244   TFile* file = TFile::Open("AliESDs.root", "RECREATE");
245   if (!file->IsOpen()) {
246     Error("Run", "opening AliESDs.root failed");
247     if (fStopOnError) {CleanUp(file); return kFALSE;}    
248   }
249
250   // loop over events
251   for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
252     Info("Run", "processing event %d", iEvent);
253     AliESD* esd = new AliESD;
254     fRunLoader->GetEvent(iEvent);
255     esd->SetRunNumber(aliRun->GetRunNumber());
256     esd->SetEventNumber(aliRun->GetEvNumber());
257     esd->SetMagneticField(aliRun->Field()->SolenoidField());
258
259     // barrel tracking
260     if (fRunTracking) {
261       if (!RunTracking(esd)) {
262         if (fStopOnError) {CleanUp(file); return kFALSE;}
263       }
264     }
265
266     // fill ESD
267     if (!fFillESD.IsNull()) {
268       if (!FillESD(esd, fFillESD)) {
269         if (fStopOnError) {CleanUp(file); return kFALSE;}
270       }
271     }
272
273     // combined PID
274     AliESDpid::MakePID(esd);
275
276     // write ESD
277     char name[100]; 
278     sprintf(name, "ESD%d", iEvent);
279     file->cd();
280     if (!esd->Write(name)) {
281       Error("Run", "writing ESD failed");
282       if (fStopOnError) {CleanUp(file); return kFALSE;}
283     }
284     delete esd;
285   }
286
287   CleanUp(file);
288
289   return kTRUE;
290 }
291
292
293 //_____________________________________________________________________________
294 Bool_t AliReconstruction::RunReconstruction(const TString& detectors)
295 {
296 // run the reconstruction
297
298   TStopwatch stopwatch;
299   stopwatch.Start();
300
301   TString detStr = detectors;
302   TObjArray* detArray = fRunLoader->GetAliRun()->Detectors();
303   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
304     AliModule* det = (AliModule*) detArray->At(iDet);
305     if (!det || !det->IsActive()) continue;
306     if (IsSelected(det->GetName(), detStr)) {
307       Info("RunReconstruction", "running reconstruction for %s", 
308            det->GetName());
309       TStopwatch stopwatchDet;
310       stopwatchDet.Start();
311       det->Reconstruct();
312       Info("RunReconstruction", "execution time for %s:", det->GetName());
313       stopwatchDet.Print();
314     }
315   }
316
317   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
318     Error("RunReconstruction", "the following detectors were not found: %s", 
319           detStr.Data());
320     if (fStopOnError) return kFALSE;
321   }
322
323   Info("RunReconstruction", "execution time:");
324   stopwatch.Print();
325
326   return kTRUE;
327 }
328
329 //_____________________________________________________________________________
330 Bool_t AliReconstruction::RunTracking(AliESD* esd)
331 {
332 // run the barrel tracking
333
334   TStopwatch stopwatch;
335   stopwatch.Start();
336
337   // get the primary vertex (from MC for the moment)
338   TArrayF vertex(3);     
339   fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vertex);
340   Double_t vtxPos[3] = {vertex[0], vertex[1], vertex[2]};
341   Double_t vtxCov[6] = {
342     0.005,
343     0.000, 0.005,
344     0.000, 0.000, 0.010
345   };
346   Double_t vtxErr[3] = {vtxCov[0], vtxCov[2], vtxCov[5]}; // diag. elements
347   esd->SetVertex(vtxPos, vtxCov);
348   fITSTracker->SetVertex(vtxPos, vtxErr);
349   fTPCTracker->SetVertex(vtxPos, vtxErr);
350   fTRDTracker->SetVertex(vtxPos, vtxErr);
351
352   // TPC tracking
353   Info("RunTracking", "TPC tracking");
354   fTPCLoader->LoadRecPoints("read");
355   TTree* tpcTree = fTPCLoader->TreeR();
356   if (!tpcTree) {
357     Error("RunTracking", "Can't get the TPC cluster tree");
358     return kFALSE;
359   }     
360   fTPCTracker->LoadClusters(tpcTree);
361   if (fTPCTracker->Clusters2Tracks(esd) != 0) {
362     Error("RunTracking", "TPC Clusters2Tracks failed");
363     return kFALSE;
364   }
365
366   fRunLoader->GetAliRun()->GetDetector("TPC")->FillESD(esd); // preliminary PID
367   AliESDpid::MakePID(esd);                  // for the ITS tracker
368
369   // ITS tracking
370   Info("RunTracking", "ITS tracking");
371   fITSLoader->LoadRecPoints("read");
372   TTree* itsTree = fITSLoader->TreeR();
373   if (!itsTree) {
374     Error("RunTracking", "Can't get the ITS cluster tree");
375     return kFALSE;
376   }     
377   fITSTracker->LoadClusters(itsTree);
378   if (fITSTracker->Clusters2Tracks(esd) != 0) {
379     Error("RunTracking", "ITS Clusters2Tracks failed");
380     return kFALSE;
381   }
382
383   // ITS back propagation
384   Info("RunTracking", "ITS back propagation");
385   if (fITSTracker->PropagateBack(esd) != 0) {
386     Error("RunTracking", "ITS backward propagation failed");
387     return kFALSE;
388   }
389
390   // TPC back propagation
391   Info("RunTracking", "TPC back propagation");
392   if (fTPCTracker->PropagateBack(esd) != 0) {
393     Error("RunTracking", "TPC backward propagation failed");
394     return kFALSE;
395   }
396
397   // TRD back propagation
398   Info("RunTracking", "TRD back propagation");
399   fTRDLoader->LoadRecPoints("read");
400   TTree* trdTree = fTRDLoader->TreeR();
401   if (!trdTree) {
402     Error("RunTracking", "Can't get the TRD cluster tree");
403     return kFALSE;
404   }     
405   fTRDTracker->LoadClusters(trdTree);
406   if (fTRDTracker->PropagateBack(esd) != 0) {
407     Error("RunTracking", "TRD backward propagation failed");
408     return kFALSE;
409   }
410
411   // TOF back propagation
412   Info("RunTracking", "TOF back propagation");
413   fTOFLoader->LoadDigits("read");
414   TTree* tofTree = fTOFLoader->TreeD();
415   if (!tofTree) {
416     Error("RunTracking", "Can't get the TOF digits tree");
417     return kFALSE;
418   }     
419   fTOFTracker->LoadClusters(tofTree);
420   if (fTOFTracker->PropagateBack(esd) != 0) {
421     Error("RunTracking", "TOF backward propagation failed");
422     return kFALSE;
423   }
424   fTOFTracker->UnloadClusters();
425   fTOFLoader->UnloadDigits();
426
427   // TRD inward refit
428   Info("RunTracking", "TRD inward refit");
429   if (fTRDTracker->RefitInward(esd) != 0) {
430     Error("RunTracking", "TRD inward refit failed");
431     return kFALSE;
432   }
433   fTRDTracker->UnloadClusters();
434   fTRDLoader->UnloadRecPoints();
435     
436   // TPC inward refit
437   Info("RunTracking", "TPC inward refit");
438   if (fTPCTracker->RefitInward(esd) != 0) {
439     Error("RunTracking", "TPC inward refit failed");
440     return kFALSE;
441   }
442   fTPCTracker->UnloadClusters();
443   fTPCLoader->UnloadRecPoints();
444     
445   // ITS inward refit
446   Info("RunTracking", "ITS inward refit");
447   if (fITSTracker->RefitInward(esd) != 0) {
448     Error("RunTracking", "ITS inward refit failed");
449     return kFALSE;
450   }
451   fITSTracker->UnloadClusters();
452   fITSLoader->UnloadRecPoints();
453
454   Info("RunTracking", "execution time:");
455   stopwatch.Print();
456
457   return kTRUE;
458 }
459
460 //_____________________________________________________________________________
461 Bool_t AliReconstruction::FillESD(AliESD* esd, const TString& detectors)
462 {
463 // fill the event summary data
464
465   TStopwatch stopwatch;
466   stopwatch.Start();
467
468   TString detStr = detectors;
469   TObjArray* detArray = fRunLoader->GetAliRun()->Detectors();
470   for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
471     AliModule* det = (AliModule*) detArray->At(iDet);
472     if (!det || !det->IsActive()) continue;
473     if (IsSelected(det->GetName(), detStr)) {
474       Info("FillESD", "filling ESD for %s", 
475            det->GetName());
476       det->FillESD(esd);
477     }
478   }
479
480   if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
481     Error("FillESD", "the following detectors were not found: %s", 
482           detStr.Data());
483     if (fStopOnError) return kFALSE;
484   }
485
486   Info("FillESD", "execution time:");
487   stopwatch.Print();
488
489   return kTRUE;
490 }
491
492
493 //_____________________________________________________________________________
494 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
495 {
496 // check whether detName is contained in detectors
497 // if yes, it is removed from detectors
498
499   // check if all detectors are selected
500   if ((detectors.CompareTo("ALL") == 0) ||
501       detectors.BeginsWith("ALL ") ||
502       detectors.EndsWith(" ALL") ||
503       detectors.Contains(" ALL ")) {
504     detectors = "ALL";
505     return kTRUE;
506   }
507
508   // search for the given detector
509   Bool_t result = kFALSE;
510   if ((detectors.CompareTo(detName) == 0) ||
511       detectors.BeginsWith(detName+" ") ||
512       detectors.EndsWith(" "+detName) ||
513       detectors.Contains(" "+detName+" ")) {
514     detectors.ReplaceAll(detName, "");
515     result = kTRUE;
516   }
517
518   // clean up the detectors string
519   while (detectors.Contains("  ")) detectors.ReplaceAll("  ", " ");
520   while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
521   while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
522
523   return result;
524 }
525
526 //_____________________________________________________________________________
527 void AliReconstruction::CleanUp(TFile* file)
528 {
529 // delete trackers and the run loader and close and delete the file
530
531   delete fITSTracker;
532   fITSTracker = NULL;
533   delete fTPCTracker;
534   fTPCTracker = NULL;
535   delete fTRDTracker;
536   fTRDTracker = NULL;
537   delete fTOFTracker;
538   fTOFTracker = NULL;
539
540   delete fRunLoader;
541   fRunLoader = NULL;
542
543   if (file) {
544     file->Close();
545     delete file;
546   }
547 }