1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////////////////
20 // class for running the reconstruction //
22 // Clusters and tracks are created for all detectors and all events by //
25 // AliReconstruction rec; //
28 // The Run method returns kTRUE in case of successful execution. //
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 //
33 // rec.SetInput("..."); //
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 //
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 //
45 // rec.SetEventRange(..., ...); //
47 // The index -1 (default) can be used for the last event to indicate no //
48 // upper limit of the event range. //
50 // The name of the galice file can be changed from the default //
51 // "galice.root" by passing it as argument to the AliReconstruction //
52 // constructor or by //
54 // rec.SetGAliceFile("..."); //
56 // The local reconstruction can be switched on or off for individual //
59 // rec.SetRunLocalReconstruction("..."); //
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. //
65 // The reconstruction of the primary vertex position can be switched off by //
67 // rec.SetRunVertexFinder(kFALSE); //
69 // The tracking and the creation of ESD tracks can be switched on for //
70 // selected detectors by //
72 // rec.SetRunTracking("..."); //
74 // Uniform/nonuniform field tracking switches (default: uniform field) //
76 // rec.SetUniformFieldTracking(); ( rec.SetNonuniformFieldTracking(); ) //
78 // The filling of additional ESD information can be steered by //
80 // rec.SetFillESD("..."); //
82 // Again, for both methods the string specifies the list of detectors. //
83 // The default is "ALL". //
85 // The call of the shortcut method //
87 // rec.SetRunReconstruction("..."); //
89 // is equivalent to calling SetRunLocalReconstruction, SetRunTracking and //
90 // SetFillESD with the same detector selecting string as argument. //
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. //
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. //
106 ///////////////////////////////////////////////////////////////////////////////
112 #include <TPluginManager.h>
113 #include <TStopwatch.h>
114 #include <TGeoManager.h>
115 #include <TLorentzVector.h>
117 #include "AliReconstruction.h"
118 #include "AliReconstructor.h"
120 #include "AliRunLoader.h"
122 #include "AliRawReaderFile.h"
123 #include "AliRawReaderDate.h"
124 #include "AliRawReaderRoot.h"
126 #include "AliESDVertex.h"
127 #include "AliTracker.h"
128 #include "AliVertexer.h"
129 #include "AliHeader.h"
130 #include "AliGenEventHeader.h"
132 #include "AliESDpid.h"
133 #include "AliESDtrack.h"
135 #include "AliRunTag.h"
136 //#include "AliLHCTag.h"
137 #include "AliDetectorTag.h"
138 #include "AliEventTag.h"
140 #include "AliTrackPointArray.h"
141 #include "AliCDBManager.h"
143 ClassImp(AliReconstruction)
146 //_____________________________________________________________________________
147 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
149 //_____________________________________________________________________________
150 AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
151 const char* name, const char* title) :
154 fRunLocalReconstruction("ALL"),
155 fUniformField(kTRUE),
156 fRunVertexFinder(kTRUE),
157 fRunHLTTracking(kFALSE),
160 fGAliceFileName(gAliceFilename),
164 fStopOnError(kFALSE),
173 fWriteAlignmentData(kFALSE),
176 // create reconstruction object with default parameters
178 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
179 fReconstructor[iDet] = NULL;
180 fLoader[iDet] = NULL;
181 fTracker[iDet] = NULL;
184 // Import TGeo geometry
185 TString geom(gSystem->DirName(gAliceFilename));
186 geom += "/geometry.root";
187 TGeoManager::Import(geom.Data());
190 //_____________________________________________________________________________
191 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
194 fRunLocalReconstruction(rec.fRunLocalReconstruction),
195 fUniformField(rec.fUniformField),
196 fRunVertexFinder(rec.fRunVertexFinder),
197 fRunHLTTracking(rec.fRunHLTTracking),
198 fRunTracking(rec.fRunTracking),
199 fFillESD(rec.fFillESD),
200 fGAliceFileName(rec.fGAliceFileName),
202 fFirstEvent(rec.fFirstEvent),
203 fLastEvent(rec.fLastEvent),
204 fStopOnError(rec.fStopOnError),
213 fWriteAlignmentData(rec.fWriteAlignmentData),
218 for (Int_t i = 0; i < fOptions.GetEntriesFast(); i++) {
219 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
221 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
222 fReconstructor[iDet] = NULL;
223 fLoader[iDet] = NULL;
224 fTracker[iDet] = NULL;
228 //_____________________________________________________________________________
229 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
231 // assignment operator
233 this->~AliReconstruction();
234 new(this) AliReconstruction(rec);
238 //_____________________________________________________________________________
239 AliReconstruction::~AliReconstruction()
247 //_____________________________________________________________________________
248 void AliReconstruction::InitCDBStorage()
250 // activate a default CDB storage
251 // First check if we have any CDB storage set, because it is used
252 // to retrieve the calibration and alignment constants
254 AliCDBManager* man = AliCDBManager::Instance();
255 if (!man->IsDefaultStorageSet())
257 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
258 AliWarning("Default CDB storage not yet set");
259 AliWarning(Form("Using default storage declared in AliSimulation: %s",fCDBUri.Data()));
260 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
261 SetDefaultStorage(fCDBUri);
266 //_____________________________________________________________________________
267 void AliReconstruction::SetDefaultStorage(const char* uri) {
268 // activate a default CDB storage
270 AliCDBManager::Instance()->SetDefaultStorage(uri);
274 //_____________________________________________________________________________
275 void AliReconstruction::SetSpecificStorage(const char* detName, const char* uri) {
276 // activate a detector-specific CDB storage
278 AliCDBManager::Instance()->SetSpecificStorage(detName, uri);
283 //_____________________________________________________________________________
284 void AliReconstruction::SetGAliceFile(const char* fileName)
286 // set the name of the galice file
288 fGAliceFileName = fileName;
291 //_____________________________________________________________________________
292 void AliReconstruction::SetOption(const char* detector, const char* option)
294 // set options for the reconstruction of a detector
296 TObject* obj = fOptions.FindObject(detector);
297 if (obj) fOptions.Remove(obj);
298 fOptions.Add(new TNamed(detector, option));
302 //_____________________________________________________________________________
303 Bool_t AliReconstruction::Run(const char* input,
304 Int_t firstEvent, Int_t lastEvent)
306 // run the reconstruction
311 if (!input) input = fInput.Data();
312 TString fileName(input);
313 if (fileName.EndsWith("/")) {
314 fRawReader = new AliRawReaderFile(fileName);
315 } else if (fileName.EndsWith(".root")) {
316 fRawReader = new AliRawReaderRoot(fileName);
317 } else if (!fileName.IsNull()) {
318 fRawReader = new AliRawReaderDate(fileName);
319 fRawReader->SelectEvents(7);
322 // get the run loader
323 if (!InitRunLoader()) return kFALSE;
325 // local reconstruction
326 if (!fRunLocalReconstruction.IsNull()) {
327 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
328 if (fStopOnError) {CleanUp(); return kFALSE;}
331 // if (!fRunVertexFinder && fRunTracking.IsNull() &&
332 // fFillESD.IsNull()) return kTRUE;
335 if (fRunVertexFinder && !CreateVertexer()) {
343 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
351 TStopwatch stopwatch;
354 // get the possibly already existing ESD file and tree
355 AliESD* esd = new AliESD; AliESD* hltesd = new AliESD;
356 TFile* fileOld = NULL;
357 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
358 if (!gSystem->AccessPathName("AliESDs.root")){
359 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
360 fileOld = TFile::Open("AliESDs.old.root");
361 if (fileOld && fileOld->IsOpen()) {
362 treeOld = (TTree*) fileOld->Get("esdTree");
363 if (treeOld) treeOld->SetBranchAddress("ESD", &esd);
364 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
365 if (hlttreeOld) hlttreeOld->SetBranchAddress("ESD", &hltesd);
369 // create the ESD output file and tree
370 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
371 if (!file->IsOpen()) {
372 AliError("opening AliESDs.root failed");
373 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
375 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
376 tree->Branch("ESD", "AliESD", &esd);
377 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
378 hlttree->Branch("ESD", "AliESD", &hltesd);
379 delete esd; delete hltesd;
380 esd = NULL; hltesd = NULL;
384 if (fRawReader) fRawReader->RewindEvents();
386 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
387 if (fRawReader) fRawReader->NextEvent();
388 if ((iEvent < firstEvent) || ((lastEvent >= 0) && (iEvent > lastEvent))) {
389 // copy old ESD to the new one
391 treeOld->SetBranchAddress("ESD", &esd);
392 treeOld->GetEntry(iEvent);
396 hlttreeOld->SetBranchAddress("ESD", &hltesd);
397 hlttreeOld->GetEntry(iEvent);
403 AliInfo(Form("processing event %d", iEvent));
404 fRunLoader->GetEvent(iEvent);
407 sprintf(fileName, "ESD_%d.%d_final.root",
408 fRunLoader->GetHeader()->GetRun(),
409 fRunLoader->GetHeader()->GetEventNrInRun());
410 if (!gSystem->AccessPathName(fileName)) continue;
412 // local reconstruction
413 if (!fRunLocalReconstruction.IsNull()) {
414 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
415 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
419 esd = new AliESD; hltesd = new AliESD;
420 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
421 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
422 esd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
423 hltesd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
425 esd->SetMagneticField(AliTracker::GetBz());
426 hltesd->SetMagneticField(AliTracker::GetBz());
432 if (fRunVertexFinder) {
433 if (!ReadESD(esd, "vertex")) {
434 if (!RunVertexFinder(esd)) {
435 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
437 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
442 if (!fRunTracking.IsNull()) {
443 if (fRunHLTTracking) {
444 hltesd->SetVertex(esd->GetVertex());
445 if (!RunHLTTracking(hltesd)) {
446 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
452 if (!fRunTracking.IsNull()) {
453 if (!ReadESD(esd, "tracking")) {
454 if (!RunTracking(esd)) {
455 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
457 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
462 if (!fFillESD.IsNull()) {
463 if (!FillESD(esd, fFillESD)) {
464 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
469 AliESDpid::MakePID(esd);
470 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
477 if (fCheckPointLevel > 0) WriteESD(esd, "final");
479 delete esd; delete hltesd;
480 esd = NULL; hltesd = NULL;
483 AliInfo(Form("Execution time for filling ESD : R:%.2fs C:%.2fs",
484 stopwatch.RealTime(),stopwatch.CpuTime()));
490 // Create tags for the events in the ESD tree (the ESD tree is always present)
491 // In case of empty events the tags will contain dummy values
493 CleanUp(file, fileOld);
499 //_____________________________________________________________________________
500 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
502 // run the local reconstruction
504 TStopwatch stopwatch;
507 TString detStr = detectors;
508 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
509 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
510 AliReconstructor* reconstructor = GetReconstructor(iDet);
511 if (!reconstructor) continue;
512 if (reconstructor->HasLocalReconstruction()) continue;
514 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
515 TStopwatch stopwatchDet;
516 stopwatchDet.Start();
518 fRawReader->RewindEvents();
519 reconstructor->Reconstruct(fRunLoader, fRawReader);
521 reconstructor->Reconstruct(fRunLoader);
523 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
524 fgkDetectorName[iDet],
525 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
528 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
529 AliError(Form("the following detectors were not found: %s",
531 if (fStopOnError) return kFALSE;
534 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
535 stopwatch.RealTime(),stopwatch.CpuTime()));
540 //_____________________________________________________________________________
541 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
543 // run the local reconstruction
545 TStopwatch stopwatch;
548 TString detStr = detectors;
549 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
550 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
551 AliReconstructor* reconstructor = GetReconstructor(iDet);
552 if (!reconstructor) continue;
553 AliLoader* loader = fLoader[iDet];
555 // conversion of digits
556 if (fRawReader && reconstructor->HasDigitConversion()) {
557 AliInfo(Form("converting raw data digits into root objects for %s",
558 fgkDetectorName[iDet]));
559 TStopwatch stopwatchDet;
560 stopwatchDet.Start();
561 loader->LoadDigits("update");
562 loader->CleanDigits();
563 loader->MakeDigitsContainer();
564 TTree* digitsTree = loader->TreeD();
565 reconstructor->ConvertDigits(fRawReader, digitsTree);
566 loader->WriteDigits("OVERWRITE");
567 loader->UnloadDigits();
568 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
569 fgkDetectorName[iDet],
570 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
573 // local reconstruction
574 if (!reconstructor->HasLocalReconstruction()) continue;
575 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
576 TStopwatch stopwatchDet;
577 stopwatchDet.Start();
578 loader->LoadRecPoints("update");
579 loader->CleanRecPoints();
580 loader->MakeRecPointsContainer();
581 TTree* clustersTree = loader->TreeR();
582 if (fRawReader && !reconstructor->HasDigitConversion()) {
583 reconstructor->Reconstruct(fRawReader, clustersTree);
585 loader->LoadDigits("read");
586 TTree* digitsTree = loader->TreeD();
588 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
589 if (fStopOnError) return kFALSE;
591 reconstructor->Reconstruct(digitsTree, clustersTree);
593 loader->UnloadDigits();
595 loader->WriteRecPoints("OVERWRITE");
596 loader->UnloadRecPoints();
597 AliDebug(1,Form("Execution time for %s: R:%.2fs C:%.2fs",
598 fgkDetectorName[iDet],
599 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
602 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
603 AliError(Form("the following detectors were not found: %s",
605 if (fStopOnError) return kFALSE;
608 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
609 stopwatch.RealTime(),stopwatch.CpuTime()));
614 //_____________________________________________________________________________
615 Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
617 // run the barrel tracking
619 TStopwatch stopwatch;
622 AliESDVertex* vertex = NULL;
623 Double_t vtxPos[3] = {0, 0, 0};
624 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
626 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
627 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
628 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
632 AliInfo("running the ITS vertex finder");
633 if (fLoader[0]) fLoader[0]->LoadRecPoints();
634 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
635 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
637 AliWarning("Vertex not found");
638 vertex = new AliESDVertex();
641 vertex->SetTruePos(vtxPos); // store also the vertex from MC
645 AliInfo("getting the primary vertex from MC");
646 vertex = new AliESDVertex(vtxPos, vtxErr);
650 vertex->GetXYZ(vtxPos);
651 vertex->GetSigmaXYZ(vtxErr);
653 AliWarning("no vertex reconstructed");
654 vertex = new AliESDVertex(vtxPos, vtxErr);
656 esd->SetVertex(vertex);
657 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
658 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
662 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
663 stopwatch.RealTime(),stopwatch.CpuTime()));
668 //_____________________________________________________________________________
669 Bool_t AliReconstruction::RunHLTTracking(AliESD*& esd)
671 // run the HLT barrel tracking
673 TStopwatch stopwatch;
677 AliError("Missing runLoader!");
681 AliInfo("running HLT tracking");
683 // Get a pointer to the HLT reconstructor
684 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
685 if (!reconstructor) return kFALSE;
688 for (Int_t iDet = 1; iDet >= 0; iDet--) {
689 TString detName = fgkDetectorName[iDet];
690 AliDebug(1, Form("%s HLT tracking", detName.Data()));
691 reconstructor->SetOption(detName.Data());
692 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
694 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
695 if (fStopOnError) return kFALSE;
699 Double_t vtxErr[3]={0.005,0.005,0.010};
700 const AliESDVertex *vertex = esd->GetVertex();
701 vertex->GetXYZ(vtxPos);
702 tracker->SetVertex(vtxPos,vtxErr);
704 fLoader[iDet]->LoadRecPoints("read");
705 TTree* tree = fLoader[iDet]->TreeR();
707 AliError(Form("Can't get the %s cluster tree", detName.Data()));
710 tracker->LoadClusters(tree);
712 if (tracker->Clusters2Tracks(esd) != 0) {
713 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
717 tracker->UnloadClusters();
722 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
723 stopwatch.RealTime(),stopwatch.CpuTime()));
728 //_____________________________________________________________________________
729 Bool_t AliReconstruction::RunTracking(AliESD*& esd)
731 // run the barrel tracking
733 TStopwatch stopwatch;
736 AliInfo("running tracking");
738 // pass 1: TPC + ITS inwards
739 for (Int_t iDet = 1; iDet >= 0; iDet--) {
740 if (!fTracker[iDet]) continue;
741 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
744 fLoader[iDet]->LoadRecPoints("read");
745 TTree* tree = fLoader[iDet]->TreeR();
747 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
750 fTracker[iDet]->LoadClusters(tree);
753 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
754 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
757 if (fCheckPointLevel > 1) {
758 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
760 // preliminary PID in TPC needed by the ITS tracker
762 GetReconstructor(1)->FillESD(fRunLoader, esd);
763 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
764 AliESDpid::MakePID(esd);
768 // pass 2: ALL backwards
769 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
770 if (!fTracker[iDet]) continue;
771 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
774 if (iDet > 1) { // all except ITS, TPC
776 fLoader[iDet]->LoadRecPoints("read");
777 tree = fLoader[iDet]->TreeR();
779 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
782 fTracker[iDet]->LoadClusters(tree);
786 if (fTracker[iDet]->PropagateBack(esd) != 0) {
787 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
790 if (fCheckPointLevel > 1) {
791 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
795 if (iDet > 2) { // all except ITS, TPC, TRD
796 fTracker[iDet]->UnloadClusters();
797 fLoader[iDet]->UnloadRecPoints();
799 // updated PID in TPC needed by the ITS tracker -MI
801 GetReconstructor(1)->FillESD(fRunLoader, esd);
802 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
803 AliESDpid::MakePID(esd);
807 // write space-points to the ESD in case alignment data output
809 if (fWriteAlignmentData)
810 WriteAlignmentData(esd);
812 // pass 3: TRD + TPC + ITS refit inwards
813 for (Int_t iDet = 2; iDet >= 0; iDet--) {
814 if (!fTracker[iDet]) continue;
815 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
818 if (fTracker[iDet]->RefitInward(esd) != 0) {
819 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
822 if (fCheckPointLevel > 1) {
823 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
827 fTracker[iDet]->UnloadClusters();
828 fLoader[iDet]->UnloadRecPoints();
831 // Propagate track to the vertex - if not done by ITS
833 Int_t ntracks = esd->GetNumberOfTracks();
834 for (Int_t itrack=0; itrack<ntracks; itrack++){
835 const Double_t kRadius = 3; // beam pipe radius
836 const Double_t kMaxStep = 5; // max step
837 const Double_t kMaxD = 123456; // max distance to prim vertex
838 Double_t fieldZ = AliTracker::GetBz(); //
839 AliESDtrack * track = esd->GetTrack(itrack);
840 if (!track) continue;
841 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
842 track->PropagateTo(kRadius, track->GetMass(),kMaxStep,kTRUE);
843 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
846 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
847 stopwatch.RealTime(),stopwatch.CpuTime()));
852 //_____________________________________________________________________________
853 Bool_t AliReconstruction::FillESD(AliESD*& esd, const TString& detectors)
855 // fill the event summary data
857 TStopwatch stopwatch;
859 AliInfo("filling ESD");
861 TString detStr = detectors;
862 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
863 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
864 AliReconstructor* reconstructor = GetReconstructor(iDet);
865 if (!reconstructor) continue;
867 if (!ReadESD(esd, fgkDetectorName[iDet])) {
868 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
869 TTree* clustersTree = NULL;
870 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
871 fLoader[iDet]->LoadRecPoints("read");
872 clustersTree = fLoader[iDet]->TreeR();
874 AliError(Form("Can't get the %s clusters tree",
875 fgkDetectorName[iDet]));
876 if (fStopOnError) return kFALSE;
879 if (fRawReader && !reconstructor->HasDigitConversion()) {
880 reconstructor->FillESD(fRawReader, clustersTree, esd);
882 TTree* digitsTree = NULL;
884 fLoader[iDet]->LoadDigits("read");
885 digitsTree = fLoader[iDet]->TreeD();
887 AliError(Form("Can't get the %s digits tree",
888 fgkDetectorName[iDet]));
889 if (fStopOnError) return kFALSE;
892 reconstructor->FillESD(digitsTree, clustersTree, esd);
893 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
895 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
896 fLoader[iDet]->UnloadRecPoints();
900 reconstructor->FillESD(fRunLoader, fRawReader, esd);
902 reconstructor->FillESD(fRunLoader, esd);
904 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
908 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
909 AliError(Form("the following detectors were not found: %s",
911 if (fStopOnError) return kFALSE;
914 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
915 stopwatch.RealTime(),stopwatch.CpuTime()));
921 //_____________________________________________________________________________
922 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
924 // check whether detName is contained in detectors
925 // if yes, it is removed from detectors
927 // check if all detectors are selected
928 if ((detectors.CompareTo("ALL") == 0) ||
929 detectors.BeginsWith("ALL ") ||
930 detectors.EndsWith(" ALL") ||
931 detectors.Contains(" ALL ")) {
936 // search for the given detector
937 Bool_t result = kFALSE;
938 if ((detectors.CompareTo(detName) == 0) ||
939 detectors.BeginsWith(detName+" ") ||
940 detectors.EndsWith(" "+detName) ||
941 detectors.Contains(" "+detName+" ")) {
942 detectors.ReplaceAll(detName, "");
946 // clean up the detectors string
947 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
948 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
949 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
954 //_____________________________________________________________________________
955 Bool_t AliReconstruction::InitRunLoader()
957 // get or create the run loader
959 if (gAlice) delete gAlice;
962 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
963 // load all base libraries to get the loader classes
964 TString libs = gSystem->GetLibraries();
965 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
966 TString detName = fgkDetectorName[iDet];
967 if (detName == "HLT") continue;
968 if (libs.Contains("lib" + detName + "base.so")) continue;
969 gSystem->Load("lib" + detName + "base.so");
971 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
973 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
977 fRunLoader->CdGAFile();
978 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
979 if (fRunLoader->LoadgAlice() == 0) {
980 gAlice = fRunLoader->GetAliRun();
981 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
984 if (!gAlice && !fRawReader) {
985 AliError(Form("no gAlice object found in file %s",
986 fGAliceFileName.Data()));
991 } else { // galice.root does not exist
993 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
997 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
998 AliConfig::GetDefaultEventFolderName(),
1001 AliError(Form("could not create run loader in file %s",
1002 fGAliceFileName.Data()));
1006 fRunLoader->MakeTree("E");
1008 while (fRawReader->NextEvent()) {
1009 fRunLoader->SetEventNumber(iEvent);
1010 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1012 fRunLoader->MakeTree("H");
1013 fRunLoader->TreeE()->Fill();
1016 fRawReader->RewindEvents();
1017 fRunLoader->WriteHeader("OVERWRITE");
1018 fRunLoader->CdGAFile();
1019 fRunLoader->Write(0, TObject::kOverwrite);
1020 // AliTracker::SetFieldMap(???);
1026 //_____________________________________________________________________________
1027 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1029 // get the reconstructor object and the loader for a detector
1031 if (fReconstructor[iDet]) return fReconstructor[iDet];
1033 // load the reconstructor object
1034 TPluginManager* pluginManager = gROOT->GetPluginManager();
1035 TString detName = fgkDetectorName[iDet];
1036 TString recName = "Ali" + detName + "Reconstructor";
1037 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1039 if (detName == "HLT") {
1040 if (!gROOT->GetClass("AliLevel3")) {
1041 gSystem->Load("libAliL3Src.so");
1042 gSystem->Load("libAliL3Misc.so");
1043 gSystem->Load("libAliL3Hough.so");
1044 gSystem->Load("libAliL3Comp.so");
1048 AliReconstructor* reconstructor = NULL;
1049 // first check if a plugin is defined for the reconstructor
1050 TPluginHandler* pluginHandler =
1051 pluginManager->FindHandler("AliReconstructor", detName);
1052 // if not, add a plugin for it
1053 if (!pluginHandler) {
1054 AliDebug(1, Form("defining plugin for %s", recName.Data()));
1055 TString libs = gSystem->GetLibraries();
1056 if (libs.Contains("lib" + detName + "base.so") ||
1057 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1058 pluginManager->AddHandler("AliReconstructor", detName,
1059 recName, detName + "rec", recName + "()");
1061 pluginManager->AddHandler("AliReconstructor", detName,
1062 recName, detName, recName + "()");
1064 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1066 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1067 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1069 if (reconstructor) {
1070 TObject* obj = fOptions.FindObject(detName.Data());
1071 if (obj) reconstructor->SetOption(obj->GetTitle());
1072 reconstructor->Init(fRunLoader);
1073 fReconstructor[iDet] = reconstructor;
1076 // get or create the loader
1077 if (detName != "HLT") {
1078 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1079 if (!fLoader[iDet]) {
1080 AliConfig::Instance()
1081 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1083 // first check if a plugin is defined for the loader
1084 TPluginHandler* pluginHandler =
1085 pluginManager->FindHandler("AliLoader", detName);
1086 // if not, add a plugin for it
1087 if (!pluginHandler) {
1088 TString loaderName = "Ali" + detName + "Loader";
1089 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1090 pluginManager->AddHandler("AliLoader", detName,
1091 loaderName, detName + "base",
1092 loaderName + "(const char*, TFolder*)");
1093 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1095 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1097 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1098 fRunLoader->GetEventFolder());
1100 if (!fLoader[iDet]) { // use default loader
1101 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1103 if (!fLoader[iDet]) {
1104 AliWarning(Form("couldn't get loader for %s", detName.Data()));
1105 if (fStopOnError) return NULL;
1107 fRunLoader->AddLoader(fLoader[iDet]);
1108 fRunLoader->CdGAFile();
1109 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1110 fRunLoader->Write(0, TObject::kOverwrite);
1115 return reconstructor;
1118 //_____________________________________________________________________________
1119 Bool_t AliReconstruction::CreateVertexer()
1121 // create the vertexer
1124 AliReconstructor* itsReconstructor = GetReconstructor(0);
1125 if (itsReconstructor) {
1126 fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
1129 AliWarning("couldn't create a vertexer for ITS");
1130 if (fStopOnError) return kFALSE;
1136 //_____________________________________________________________________________
1137 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1139 // create the trackers
1141 TString detStr = detectors;
1142 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1143 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1144 AliReconstructor* reconstructor = GetReconstructor(iDet);
1145 if (!reconstructor) continue;
1146 TString detName = fgkDetectorName[iDet];
1147 if (detName == "HLT") {
1148 fRunHLTTracking = kTRUE;
1152 fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
1153 if (!fTracker[iDet] && (iDet < 7)) {
1154 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1155 if (fStopOnError) return kFALSE;
1162 //_____________________________________________________________________________
1163 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1165 // delete trackers and the run loader and close and delete the file
1167 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1168 delete fReconstructor[iDet];
1169 fReconstructor[iDet] = NULL;
1170 fLoader[iDet] = NULL;
1171 delete fTracker[iDet];
1172 fTracker[iDet] = NULL;
1190 gSystem->Unlink("AliESDs.old.root");
1195 //_____________________________________________________________________________
1196 Bool_t AliReconstruction::ReadESD(AliESD*& esd, const char* recStep) const
1198 // read the ESD event from a file
1200 if (!esd) return kFALSE;
1202 sprintf(fileName, "ESD_%d.%d_%s.root",
1203 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1204 if (gSystem->AccessPathName(fileName)) return kFALSE;
1206 AliInfo(Form("reading ESD from file %s", fileName));
1207 AliDebug(1, Form("reading ESD from file %s", fileName));
1208 TFile* file = TFile::Open(fileName);
1209 if (!file || !file->IsOpen()) {
1210 AliError(Form("opening %s failed", fileName));
1217 esd = (AliESD*) file->Get("ESD");
1223 //_____________________________________________________________________________
1224 void AliReconstruction::WriteESD(AliESD* esd, const char* recStep) const
1226 // write the ESD event to a file
1230 sprintf(fileName, "ESD_%d.%d_%s.root",
1231 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1233 AliDebug(1, Form("writing ESD to file %s", fileName));
1234 TFile* file = TFile::Open(fileName, "recreate");
1235 if (!file || !file->IsOpen()) {
1236 AliError(Form("opening %s failed", fileName));
1247 //_____________________________________________________________________________
1248 void AliReconstruction::CreateTag(TFile* file)
1253 Double_t fMUONMASS = 0.105658369;
1256 Double_t fThetaX, fThetaY, fPyz, fChisquare;
1257 Double_t fPxRec,fPyRec, fPzRec, fEnergy;
1259 TLorentzVector fEPvector;
1261 Float_t fZVertexCut = 10.0;
1262 Float_t fRhoVertexCut = 2.0;
1264 Float_t fLowPtCut = 1.0;
1265 Float_t fHighPtCut = 3.0;
1266 Float_t fVeryHighPtCut = 10.0;
1269 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
1271 // Creates the tags for all the events in a given ESD file
1273 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
1274 Int_t nPos, nNeg, nNeutr;
1275 Int_t nK0s, nNeutrons, nPi0s, nGamas;
1276 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
1277 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
1278 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
1279 Float_t maxPt = .0, meanPt = .0, totalP = .0;
1281 AliRunTag *tag = new AliRunTag();
1282 AliEventTag *evTag = new AliEventTag();
1283 TTree ttag("T","A Tree with event tags");
1284 TBranch * btag = ttag.Branch("AliTAG", &tag);
1285 btag->SetCompressionLevel(9);
1287 AliInfo(Form("Creating the tags......."));
1289 if (!file || !file->IsOpen()) {
1290 AliError(Form("opening failed"));
1294 Int_t firstEvent = 0,lastEvent = 0;
1295 TTree *t = (TTree*) file->Get("esdTree");
1296 TBranch * b = t->GetBranch("ESD");
1298 b->SetAddress(&esd);
1300 tag->SetRunId(esd->GetRunNumber());
1302 Int_t iNumberOfEvents = b->GetEntries();
1303 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
1330 b->GetEntry(iEventNumber);
1331 const AliESDVertex * vertexIn = esd->GetVertex();
1333 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
1334 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
1335 UInt_t status = esdTrack->GetStatus();
1337 //select only tracks with ITS refit
1338 if ((status&AliESDtrack::kITSrefit)==0) continue;
1339 //select only tracks with TPC refit
1340 if ((status&AliESDtrack::kTPCrefit)==0) continue;
1342 //select only tracks with the "combined PID"
1343 if ((status&AliESDtrack::kESDpid)==0) continue;
1345 esdTrack->GetPxPyPz(p);
1346 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
1347 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
1350 if(fPt > maxPt) maxPt = fPt;
1352 if(esdTrack->GetSign() > 0) {
1354 if(fPt > fLowPtCut) nCh1GeV++;
1355 if(fPt > fHighPtCut) nCh3GeV++;
1356 if(fPt > fVeryHighPtCut) nCh10GeV++;
1358 if(esdTrack->GetSign() < 0) {
1360 if(fPt > fLowPtCut) nCh1GeV++;
1361 if(fPt > fHighPtCut) nCh3GeV++;
1362 if(fPt > fVeryHighPtCut) nCh10GeV++;
1364 if(esdTrack->GetSign() == 0) nNeutr++;
1368 esdTrack->GetESDpid(prob);
1371 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
1372 if(rcc == 0.0) continue;
1375 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
1378 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
1380 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
1382 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
1384 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
1386 if(fPt > fLowPtCut) nEl1GeV++;
1387 if(fPt > fHighPtCut) nEl3GeV++;
1388 if(fPt > fVeryHighPtCut) nEl10GeV++;
1396 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
1397 // loop over all reconstructed tracks (also first track of combination)
1398 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
1399 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
1400 if (muonTrack == 0x0) continue;
1402 // Coordinates at vertex
1403 fZ = muonTrack->GetZ();
1404 fY = muonTrack->GetBendingCoor();
1405 fX = muonTrack->GetNonBendingCoor();
1407 fThetaX = muonTrack->GetThetaX();
1408 fThetaY = muonTrack->GetThetaY();
1410 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
1411 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
1412 fPxRec = fPzRec * TMath::Tan(fThetaX);
1413 fPyRec = fPzRec * TMath::Tan(fThetaY);
1414 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
1416 //ChiSquare of the track if needed
1417 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
1418 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
1419 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
1421 // total number of muons inside a vertex cut
1422 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
1424 if(fEPvector.Pt() > fLowPtCut) {
1426 if(fEPvector.Pt() > fHighPtCut) {
1428 if (fEPvector.Pt() > fVeryHighPtCut) {
1436 // Fill the event tags
1438 meanPt = meanPt/ntrack;
1440 evTag->SetEventId(iEventNumber+1);
1441 evTag->SetVertexX(vertexIn->GetXv());
1442 evTag->SetVertexY(vertexIn->GetYv());
1443 evTag->SetVertexZ(vertexIn->GetZv());
1445 evTag->SetT0VertexZ(esd->GetT0zVertex());
1447 evTag->SetTrigger(esd->GetTrigger());
1449 evTag->SetZDCNeutronEnergy(esd->GetZDCNEnergy());
1450 evTag->SetZDCProtonEnergy(esd->GetZDCPEnergy());
1451 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
1452 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
1455 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
1456 evTag->SetNumOfPosTracks(nPos);
1457 evTag->SetNumOfNegTracks(nNeg);
1458 evTag->SetNumOfNeutrTracks(nNeutr);
1460 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
1461 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
1462 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
1463 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
1465 evTag->SetNumOfProtons(nProtons);
1466 evTag->SetNumOfKaons(nKaons);
1467 evTag->SetNumOfPions(nPions);
1468 evTag->SetNumOfMuons(nMuons);
1469 evTag->SetNumOfElectrons(nElectrons);
1470 evTag->SetNumOfPhotons(nGamas);
1471 evTag->SetNumOfPi0s(nPi0s);
1472 evTag->SetNumOfNeutrons(nNeutrons);
1473 evTag->SetNumOfKaon0s(nK0s);
1475 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
1476 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
1477 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
1478 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
1479 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
1480 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
1481 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
1482 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
1483 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
1485 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
1486 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
1488 evTag->SetTotalMomentum(totalP);
1489 evTag->SetMeanPt(meanPt);
1490 evTag->SetMaxPt(maxPt);
1492 tag->AddEventTag(*evTag);
1494 lastEvent = iNumberOfEvents;
1500 sprintf(fileName, "Run%d.Event%d_%d.ESD.tag.root",
1501 tag->GetRunId(),firstEvent,lastEvent );
1502 AliInfo(Form("writing tags to file %s", fileName));
1503 AliDebug(1, Form("writing tags to file %s", fileName));
1505 TFile* ftag = TFile::Open(fileName, "recreate");
1514 void AliReconstruction::WriteAlignmentData(AliESD* esd)
1516 // Write space-points which are then used in the alignment procedures
1517 // For the moment only ITS, TRD and TPC
1519 // Load TOF clusters
1521 fLoader[3]->LoadRecPoints("read");
1522 TTree* tree = fLoader[3]->TreeR();
1524 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
1527 fTracker[3]->LoadClusters(tree);
1529 Int_t ntracks = esd->GetNumberOfTracks();
1530 for (Int_t itrack = 0; itrack < ntracks; itrack++)
1532 AliESDtrack *track = esd->GetTrack(itrack);
1535 for (Int_t iDet = 3; iDet >= 0; iDet--)
1536 nsp += track->GetNcls(iDet);
1538 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
1539 track->SetTrackPointArray(sp);
1541 for (Int_t iDet = 3; iDet >= 0; iDet--) {
1542 AliTracker *tracker = fTracker[iDet];
1543 if (!tracker) continue;
1544 Int_t nspdet = track->GetNcls(iDet);
1545 if (nspdet <= 0) continue;
1546 track->GetClusters(iDet,idx);
1550 while (isp < nspdet) {
1551 Bool_t isvalid = tracker->GetTrackPoint(idx[isp2],p); isp2++;
1552 if (!isvalid) continue;
1553 sp->AddPoint(isptrack,&p); isptrack++; isp++;
1559 fTracker[3]->UnloadClusters();
1560 fLoader[3]->UnloadRecPoints();