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