Additional protection (Y.Schutz)
[u/mrichter/AliRoot.git] / STEER / AliReconstruction.cxx
CommitLineData
596a855f 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 //
e583c30d 30// "galice.root" by passing it as argument to the AliReconstruction //
31// constructor or by //
596a855f 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"
a866ac60 70#include "AliMagF.h"
596a855f 71#include <TArrayF.h>
72
73
74ClassImp(AliReconstruction)
75
76
77//_____________________________________________________________________________
e583c30d 78AliReconstruction::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)
596a855f 97{
98// create reconstruction object with default parameters
99
596a855f 100}
101
102//_____________________________________________________________________________
103AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
e583c30d 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)
596a855f 121{
122// copy constructor
123
596a855f 124}
125
126//_____________________________________________________________________________
127AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
128{
129// assignment operator
130
131 this->~AliReconstruction();
132 new(this) AliReconstruction(rec);
133 return *this;
134}
135
136//_____________________________________________________________________________
137AliReconstruction::~AliReconstruction()
138{
139// clean up
140
e583c30d 141 CleanUp();
596a855f 142}
143
144
145//_____________________________________________________________________________
146void AliReconstruction::SetGAliceFile(const char* fileName)
147{
148// set the name of the galice file
149
150 fGAliceFileName = fileName;
151}
152
153
154//_____________________________________________________________________________
155Bool_t AliReconstruction::Run()
156{
157// run the reconstruction
158
159 // open the run loader
596a855f 160 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
161 if (!fRunLoader) {
162 Error("Run", "no run loader found in file %s",
163 fGAliceFileName.Data());
e583c30d 164 CleanUp();
596a855f 165 return kFALSE;
166 }
167 fRunLoader->LoadgAlice();
e583c30d 168 AliRun* aliRun = fRunLoader->GetAliRun();
169 if (!aliRun) {
596a855f 170 Error("Run", "no gAlice object found in file %s",
171 fGAliceFileName.Data());
e583c30d 172 CleanUp();
596a855f 173 return kFALSE;
174 }
e583c30d 175 gAlice = aliRun;
596a855f 176
177 // local reconstruction
178 if (!fRunReconstruction.IsNull()) {
179 if (!RunReconstruction(fRunReconstruction)) {
e583c30d 180 if (fStopOnError) {CleanUp(); return kFALSE;}
596a855f 181 }
182 }
183 if (!fRunTracking && fFillESD.IsNull()) return kTRUE;
184
96ae9b65 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 }
596a855f 242 }
596a855f 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");
e583c30d 247 if (fStopOnError) {CleanUp(file); return kFALSE;}
596a855f 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);
e583c30d 255 esd->SetRunNumber(aliRun->GetRunNumber());
256 esd->SetEventNumber(aliRun->GetEvNumber());
257 esd->SetMagneticField(aliRun->Field()->SolenoidField());
596a855f 258
259 // barrel tracking
260 if (fRunTracking) {
261 if (!RunTracking(esd)) {
e583c30d 262 if (fStopOnError) {CleanUp(file); return kFALSE;}
596a855f 263 }
264 }
265
266 // fill ESD
267 if (!fFillESD.IsNull()) {
268 if (!FillESD(esd, fFillESD)) {
e583c30d 269 if (fStopOnError) {CleanUp(file); return kFALSE;}
596a855f 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");
e583c30d 282 if (fStopOnError) {CleanUp(file); return kFALSE;}
596a855f 283 }
284 }
285
e583c30d 286 CleanUp(file);
596a855f 287
288 return kTRUE;
289}
290
291
292//_____________________________________________________________________________
293Bool_t AliReconstruction::RunReconstruction(const TString& detectors)
294{
295// run the reconstruction
296
030b532d 297 TStopwatch stopwatch;
298 stopwatch.Start();
299
596a855f 300 TString detStr = detectors;
e583c30d 301 TObjArray* detArray = fRunLoader->GetAliRun()->Detectors();
596a855f 302 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
303 AliModule* det = (AliModule*) detArray->At(iDet);
304 if (!det || !det->IsActive()) continue;
305 if (IsSelected(det->GetName(), detStr)) {
306 Info("RunReconstruction", "running reconstruction for %s",
307 det->GetName());
030b532d 308 TStopwatch stopwatchDet;
309 stopwatchDet.Start();
596a855f 310 det->Reconstruct();
030b532d 311 Info("RunReconstruction", "execution time for %s:", det->GetName());
312 stopwatchDet.Print();
596a855f 313 }
314 }
315
316 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
317 Error("RunReconstruction", "the following detectors were not found: %s",
318 detStr.Data());
319 if (fStopOnError) return kFALSE;
320 }
321
030b532d 322 Info("RunReconstruction", "execution time:");
323 stopwatch.Print();
324
596a855f 325 return kTRUE;
326}
327
328//_____________________________________________________________________________
329Bool_t AliReconstruction::RunTracking(AliESD* esd)
330{
331// run the barrel tracking
332
030b532d 333 TStopwatch stopwatch;
334 stopwatch.Start();
335
596a855f 336 // get the primary vertex (from MC for the moment)
337 TArrayF vertex(3);
338 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vertex);
339 Double_t vtxPos[3] = {vertex[0], vertex[1], vertex[2]};
340 Double_t vtxCov[6] = {
341 0.005,
342 0.000, 0.005,
343 0.000, 0.000, 0.010
344 };
345 Double_t vtxErr[3] = {vtxCov[0], vtxCov[2], vtxCov[5]}; // diag. elements
346 esd->SetVertex(vtxPos, vtxCov);
347 fITSTracker->SetVertex(vtxPos, vtxErr);
348 fTPCTracker->SetVertex(vtxPos, vtxErr);
349 fTRDTracker->SetVertex(vtxPos, vtxErr);
350
85555ca8 351 // TPC tracking
352 Info("RunTracking", "TPC tracking");
596a855f 353 fTPCLoader->LoadRecPoints("read");
354 TTree* tpcTree = fTPCLoader->TreeR();
355 if (!tpcTree) {
356 Error("RunTracking", "Can't get the TPC cluster tree");
357 return kFALSE;
358 }
359 fTPCTracker->LoadClusters(tpcTree);
360 if (fTPCTracker->Clusters2Tracks(esd) != 0) {
85555ca8 361 Error("RunTracking", "TPC Clusters2Tracks failed");
596a855f 362 return kFALSE;
363 }
364
e583c30d 365 fRunLoader->GetAliRun()->GetDetector("TPC")->FillESD(esd); // preliminary PID
85555ca8 366 AliESDpid::MakePID(esd); // for the ITS tracker
367
368 // ITS tracking
369 Info("RunTracking", "ITS tracking");
596a855f 370 fITSLoader->LoadRecPoints("read");
371 TTree* itsTree = fITSLoader->TreeR();
372 if (!itsTree) {
373 Error("RunTracking", "Can't get the ITS cluster tree");
374 return kFALSE;
375 }
376 fITSTracker->LoadClusters(itsTree);
377 if (fITSTracker->Clusters2Tracks(esd) != 0) {
85555ca8 378 Error("RunTracking", "ITS Clusters2Tracks failed");
596a855f 379 return kFALSE;
380 }
381
382 // ITS back propagation
383 Info("RunTracking", "ITS back propagation");
384 if (fITSTracker->PropagateBack(esd) != 0) {
385 Error("RunTracking", "ITS backward propagation failed");
386 return kFALSE;
387 }
388
389 // TPC back propagation
390 Info("RunTracking", "TPC back propagation");
391 if (fTPCTracker->PropagateBack(esd) != 0) {
392 Error("RunTracking", "TPC backward propagation failed");
393 return kFALSE;
394 }
395
396 // TRD back propagation
397 Info("RunTracking", "TRD back propagation");
398 fTRDLoader->LoadRecPoints("read");
399 TTree* trdTree = fTRDLoader->TreeR();
400 if (!trdTree) {
401 Error("RunTracking", "Can't get the TRD cluster tree");
402 return kFALSE;
403 }
404 fTRDTracker->LoadClusters(trdTree);
405 if (fTRDTracker->PropagateBack(esd) != 0) {
406 Error("RunTracking", "TRD backward propagation failed");
407 return kFALSE;
408 }
409
410 // TOF back propagation
411 Info("RunTracking", "TOF back propagation");
412 fTOFLoader->LoadDigits("read");
413 TTree* tofTree = fTOFLoader->TreeD();
414 if (!tofTree) {
415 Error("RunTracking", "Can't get the TOF digits tree");
416 return kFALSE;
417 }
418 fTOFTracker->LoadClusters(tofTree);
419 if (fTOFTracker->PropagateBack(esd) != 0) {
420 Error("RunTracking", "TOF backward propagation failed");
421 return kFALSE;
422 }
423 fTOFTracker->UnloadClusters();
424 fTOFLoader->UnloadDigits();
425
426 // TRD inward refit
427 Info("RunTracking", "TRD inward refit");
428 if (fTRDTracker->RefitInward(esd) != 0) {
429 Error("RunTracking", "TRD inward refit failed");
430 return kFALSE;
431 }
432 fTRDTracker->UnloadClusters();
433 fTRDLoader->UnloadRecPoints();
434
435 // TPC inward refit
436 Info("RunTracking", "TPC inward refit");
437 if (fTPCTracker->RefitInward(esd) != 0) {
438 Error("RunTracking", "TPC inward refit failed");
439 return kFALSE;
440 }
441 fTPCTracker->UnloadClusters();
442 fTPCLoader->UnloadRecPoints();
443
444 // ITS inward refit
445 Info("RunTracking", "ITS inward refit");
446 if (fITSTracker->RefitInward(esd) != 0) {
447 Error("RunTracking", "ITS inward refit failed");
448 return kFALSE;
449 }
450 fITSTracker->UnloadClusters();
451 fITSLoader->UnloadRecPoints();
452
030b532d 453 Info("RunTracking", "execution time:");
454 stopwatch.Print();
455
596a855f 456 return kTRUE;
457}
458
459//_____________________________________________________________________________
460Bool_t AliReconstruction::FillESD(AliESD* esd, const TString& detectors)
461{
462// fill the event summary data
463
030b532d 464 TStopwatch stopwatch;
465 stopwatch.Start();
466
596a855f 467 TString detStr = detectors;
e583c30d 468 TObjArray* detArray = fRunLoader->GetAliRun()->Detectors();
596a855f 469 for (Int_t iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
470 AliModule* det = (AliModule*) detArray->At(iDet);
471 if (!det || !det->IsActive()) continue;
472 if (IsSelected(det->GetName(), detStr)) {
473 Info("FillESD", "filling ESD for %s",
474 det->GetName());
475 det->FillESD(esd);
476 }
477 }
478
479 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
480 Error("FillESD", "the following detectors were not found: %s",
481 detStr.Data());
482 if (fStopOnError) return kFALSE;
483 }
484
030b532d 485 Info("FillESD", "execution time:");
486 stopwatch.Print();
487
596a855f 488 return kTRUE;
489}
490
491
492//_____________________________________________________________________________
493Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
494{
495// check whether detName is contained in detectors
496// if yes, it is removed from detectors
497
498 // check if all detectors are selected
499 if ((detectors.CompareTo("ALL") == 0) ||
500 detectors.BeginsWith("ALL ") ||
501 detectors.EndsWith(" ALL") ||
502 detectors.Contains(" ALL ")) {
503 detectors = "ALL";
504 return kTRUE;
505 }
506
507 // search for the given detector
508 Bool_t result = kFALSE;
509 if ((detectors.CompareTo(detName) == 0) ||
510 detectors.BeginsWith(detName+" ") ||
511 detectors.EndsWith(" "+detName) ||
512 detectors.Contains(" "+detName+" ")) {
513 detectors.ReplaceAll(detName, "");
514 result = kTRUE;
515 }
516
517 // clean up the detectors string
518 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
519 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
520 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
521
522 return result;
523}
e583c30d 524
525//_____________________________________________________________________________
526void AliReconstruction::CleanUp(TFile* file)
527{
528// delete trackers and the run loader and close and delete the file
529
530 delete fITSTracker;
531 fITSTracker = NULL;
532 delete fTPCTracker;
533 fTPCTracker = NULL;
534 delete fTRDTracker;
535 fTRDTracker = NULL;
536 delete fTOFTracker;
537 fTOFTracker = NULL;
538
539 delete fRunLoader;
540 fRunLoader = NULL;
541
542 if (file) {
543 file->Close();
544 delete file;
545 }
546}