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