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"
142 #include "AliCDBEntry.h"
143 #include "AliAlignObj.h"
145 ClassImp(AliReconstruction)
148 //_____________________________________________________________________________
149 const char* AliReconstruction::fgkDetectorName[AliReconstruction::fgkNDetectors] = {"ITS", "TPC", "TRD", "TOF", "PHOS", "RICH", "EMCAL", "MUON", "FMD", "ZDC", "PMD", "START", "VZERO", "CRT", "HLT"};
151 //_____________________________________________________________________________
152 AliReconstruction::AliReconstruction(const char* gAliceFilename, const char* cdbUri,
153 const char* name, const char* title) :
156 fRunLocalReconstruction("ALL"),
157 fUniformField(kTRUE),
158 fRunVertexFinder(kTRUE),
159 fRunHLTTracking(kFALSE),
162 fGAliceFileName(gAliceFilename),
166 fStopOnError(kFALSE),
169 fLoadAlignFromCDB(kTRUE),
170 fLoadAlignData("ALL"),
177 fAlignObjArray(NULL),
178 fWriteAlignmentData(kFALSE),
181 // create reconstruction object with default parameters
183 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
184 fReconstructor[iDet] = NULL;
185 fLoader[iDet] = NULL;
186 fTracker[iDet] = NULL;
191 //_____________________________________________________________________________
192 AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
195 fRunLocalReconstruction(rec.fRunLocalReconstruction),
196 fUniformField(rec.fUniformField),
197 fRunVertexFinder(rec.fRunVertexFinder),
198 fRunHLTTracking(rec.fRunHLTTracking),
199 fRunTracking(rec.fRunTracking),
200 fFillESD(rec.fFillESD),
201 fGAliceFileName(rec.fGAliceFileName),
203 fFirstEvent(rec.fFirstEvent),
204 fLastEvent(rec.fLastEvent),
205 fStopOnError(rec.fStopOnError),
208 fLoadAlignFromCDB(rec.fLoadAlignFromCDB),
209 fLoadAlignData(rec.fLoadAlignData),
216 fAlignObjArray(rec.fAlignObjArray),
217 fWriteAlignmentData(rec.fWriteAlignmentData),
222 for (Int_t i = 0; i < fOptions.GetEntriesFast(); i++) {
223 if (rec.fOptions[i]) fOptions.Add(rec.fOptions[i]->Clone());
225 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
226 fReconstructor[iDet] = NULL;
227 fLoader[iDet] = NULL;
228 fTracker[iDet] = NULL;
232 //_____________________________________________________________________________
233 AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
235 // assignment operator
237 this->~AliReconstruction();
238 new(this) AliReconstruction(rec);
242 //_____________________________________________________________________________
243 AliReconstruction::~AliReconstruction()
251 //_____________________________________________________________________________
252 void AliReconstruction::InitCDBStorage()
254 // activate a default CDB storage
255 // First check if we have any CDB storage set, because it is used
256 // to retrieve the calibration and alignment constants
258 AliCDBManager* man = AliCDBManager::Instance();
259 if (!man->IsDefaultStorageSet())
261 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
262 AliWarning("Default CDB storage not yet set");
263 AliWarning(Form("Using default storage declared in AliSimulation: %s",fCDBUri.Data()));
264 AliWarning("!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
265 SetDefaultStorage(fCDBUri);
267 Int_t cdbRun = AliCDBManager::Instance()->GetRun();
269 AliWarning("AliCDBManager's run number temporarily set to 0!!");
270 AliCDBManager::Instance()->SetRun(0);
276 //_____________________________________________________________________________
277 void AliReconstruction::SetDefaultStorage(const char* uri) {
278 // activate a default CDB storage
280 AliCDBManager::Instance()->SetDefaultStorage(uri);
284 //_____________________________________________________________________________
285 void AliReconstruction::SetSpecificStorage(const char* detName, const char* uri) {
286 // activate a detector-specific CDB storage
288 AliCDBManager::Instance()->SetSpecificStorage(detName, uri);
292 //_____________________________________________________________________________
293 Bool_t AliReconstruction::SetRunNumber()
295 // The method is called in Run() in order
296 // to set a correct run number.
297 // In case of raw data reconstruction the
298 // run number is taken from the raw data header
300 if(AliCDBManager::Instance()->GetRun() < 0) {
302 AliError("No run loader is found !");
305 // read run number from gAlice
306 AliCDBManager::Instance()->SetRun(fRunLoader->GetAliRun()->GetRunNumber());
307 AliInfo(Form("Run number: %d",AliCDBManager::Instance()->GetRun()));
312 //_____________________________________________________________________________
313 Bool_t AliReconstruction::ApplyAlignObjsToGeom(TObjArray* alObjArray)
315 // Read collection of alignment objects (AliAlignObj derived) saved
316 // in the TClonesArray ClArrayName and apply them to the geometry
317 // manager singleton.
320 Int_t nvols = alObjArray->GetEntriesFast();
324 for(Int_t j=0; j<nvols; j++)
326 AliAlignObj* alobj = (AliAlignObj*) alObjArray->UncheckedAt(j);
327 if (alobj->ApplyToGeometry() == kFALSE) flag = kFALSE;
330 if (AliDebugLevelClass() >= 1) {
331 gGeoManager->CheckOverlaps(20);
332 TObjArray* ovexlist = gGeoManager->GetListOfOverlaps();
333 if(ovexlist->GetEntriesFast()){
334 AliError("The application of alignment objects to the geometry caused huge overlaps/extrusions!");
342 //_____________________________________________________________________________
343 Bool_t AliReconstruction::SetAlignObjArraySingleDet(const char* detName)
345 // Fills array of single detector's alignable objects from CDB
347 AliDebug(2, Form("Loading alignment data for detector: %s",detName));
351 AliCDBPath path(detName,"Align","Data");
353 entry=AliCDBManager::Instance()->Get(path.GetPath());
355 AliDebug(2,Form("Couldn't load alignment data for detector %s",detName));
359 TClonesArray *alignArray = (TClonesArray*) entry->GetObject();
360 alignArray->SetOwner(0);
361 AliDebug(2,Form("Found %d alignment objects for %s",
362 alignArray->GetEntries(),detName));
364 AliAlignObj *alignObj=0;
365 TIter iter(alignArray);
367 // loop over align objects in detector
368 while( ( alignObj=(AliAlignObj *) iter.Next() ) ){
369 fAlignObjArray->Add(alignObj);
371 // delete entry --- Don't delete, it is cached!
373 AliDebug(2, Form("fAlignObjArray entries: %d",fAlignObjArray->GetEntries() ));
378 //_____________________________________________________________________________
379 Bool_t AliReconstruction::MisalignGeometry(const TString& detectors)
381 // Read the alignment objects from CDB.
382 // Each detector is supposed to have the
383 // alignment objects in DET/Align/Data CDB path.
384 // All the detector objects are then collected,
385 // sorted by geometry level (starting from ALIC) and
386 // then applied to the TGeo geometry.
387 // Finally an overlaps check is performed.
389 // Load alignment data from CDB and fill fAlignObjArray
390 if(fLoadAlignFromCDB){
391 if(!fAlignObjArray) fAlignObjArray = new TObjArray();
393 //fAlignObjArray->RemoveAll();
394 fAlignObjArray->Clear();
395 fAlignObjArray->SetOwner(0);
397 TString detStr = detectors;
398 TString dataNotLoaded="";
399 TString dataLoaded="";
401 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
402 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
403 if(!SetAlignObjArraySingleDet(fgkDetectorName[iDet])){
404 dataNotLoaded += fgkDetectorName[iDet];
405 dataNotLoaded += " ";
407 dataLoaded += fgkDetectorName[iDet];
410 } // end loop over detectors
412 if ((detStr.CompareTo("ALL") == 0)) detStr = "";
413 dataNotLoaded += detStr;
414 AliInfo(Form("Alignment data loaded for: %s",
416 AliInfo(Form("Didn't/couldn't load alignment data for: %s",
417 dataNotLoaded.Data()));
418 } // fLoadAlignFromCDB flag
420 // Check if the array with alignment objects was
421 // provided by the user. If yes, apply the objects
422 // to the present TGeo geometry
423 if (fAlignObjArray) {
424 if (gGeoManager && gGeoManager->IsClosed()) {
425 if (ApplyAlignObjsToGeom(fAlignObjArray) == kFALSE) {
426 AliError("The misalignment of one or more volumes failed!"
427 "Compare the list of simulated detectors and the list of detector alignment data!");
432 AliError("Can't apply the misalignment! gGeoManager doesn't exist or it is still opened!");
440 //_____________________________________________________________________________
441 void AliReconstruction::SetGAliceFile(const char* fileName)
443 // set the name of the galice file
445 fGAliceFileName = fileName;
448 //_____________________________________________________________________________
449 void AliReconstruction::SetOption(const char* detector, const char* option)
451 // set options for the reconstruction of a detector
453 TObject* obj = fOptions.FindObject(detector);
454 if (obj) fOptions.Remove(obj);
455 fOptions.Add(new TNamed(detector, option));
459 //_____________________________________________________________________________
460 Bool_t AliReconstruction::Run(const char* input,
461 Int_t firstEvent, Int_t lastEvent)
463 // run the reconstruction
468 if (!input) input = fInput.Data();
469 TString fileName(input);
470 if (fileName.EndsWith("/")) {
471 fRawReader = new AliRawReaderFile(fileName);
472 } else if (fileName.EndsWith(".root")) {
473 fRawReader = new AliRawReaderRoot(fileName);
474 } else if (!fileName.IsNull()) {
475 fRawReader = new AliRawReaderDate(fileName);
476 fRawReader->SelectEvents(7);
479 // get the run loader
480 if (!InitRunLoader()) return kFALSE;
482 // Set run number in CDBManager (if it is not already set by the user)
483 if (!SetRunNumber()) if (fStopOnError) return kFALSE;
485 // Import ideal TGeo geometry and apply misalignment
487 TString geom(gSystem->DirName(fGAliceFileName));
488 geom += "/geometry.root";
489 TGeoManager::Import(geom.Data());
490 if (!gGeoManager) if (fStopOnError) return kFALSE;
492 if (!MisalignGeometry(fLoadAlignData)) if (fStopOnError) return kFALSE;
494 // Temporary fix by A.Gheata
495 // Could be removed with the next Root version (>5.11)
497 TIter next(gGeoManager->GetListOfVolumes());
499 while ((vol = (TGeoVolume *)next())) {
500 if (vol->GetVoxels()) {
501 if (vol->GetVoxels()->NeedRebuild()) {
502 vol->GetVoxels()->Voxelize();
509 // local reconstruction
510 if (!fRunLocalReconstruction.IsNull()) {
511 if (!RunLocalReconstruction(fRunLocalReconstruction)) {
512 if (fStopOnError) {CleanUp(); return kFALSE;}
515 // if (!fRunVertexFinder && fRunTracking.IsNull() &&
516 // fFillESD.IsNull()) return kTRUE;
519 if (fRunVertexFinder && !CreateVertexer()) {
527 if (!fRunTracking.IsNull() && !CreateTrackers(fRunTracking)) {
535 TStopwatch stopwatch;
538 // get the possibly already existing ESD file and tree
539 AliESD* esd = new AliESD; AliESD* hltesd = new AliESD;
540 TFile* fileOld = NULL;
541 TTree* treeOld = NULL; TTree *hlttreeOld = NULL;
542 if (!gSystem->AccessPathName("AliESDs.root")){
543 gSystem->CopyFile("AliESDs.root", "AliESDs.old.root", kTRUE);
544 fileOld = TFile::Open("AliESDs.old.root");
545 if (fileOld && fileOld->IsOpen()) {
546 treeOld = (TTree*) fileOld->Get("esdTree");
547 if (treeOld) treeOld->SetBranchAddress("ESD", &esd);
548 hlttreeOld = (TTree*) fileOld->Get("HLTesdTree");
549 if (hlttreeOld) hlttreeOld->SetBranchAddress("ESD", &hltesd);
553 // create the ESD output file and tree
554 TFile* file = TFile::Open("AliESDs.root", "RECREATE");
555 if (!file->IsOpen()) {
556 AliError("opening AliESDs.root failed");
557 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
559 TTree* tree = new TTree("esdTree", "Tree with ESD objects");
560 tree->Branch("ESD", "AliESD", &esd);
561 TTree* hlttree = new TTree("HLTesdTree", "Tree with HLT ESD objects");
562 hlttree->Branch("ESD", "AliESD", &hltesd);
563 delete esd; delete hltesd;
564 esd = NULL; hltesd = NULL;
568 if (fRawReader) fRawReader->RewindEvents();
570 for (Int_t iEvent = 0; iEvent < fRunLoader->GetNumberOfEvents(); iEvent++) {
571 if (fRawReader) fRawReader->NextEvent();
572 if ((iEvent < firstEvent) || ((lastEvent >= 0) && (iEvent > lastEvent))) {
573 // copy old ESD to the new one
575 treeOld->SetBranchAddress("ESD", &esd);
576 treeOld->GetEntry(iEvent);
580 hlttreeOld->SetBranchAddress("ESD", &hltesd);
581 hlttreeOld->GetEntry(iEvent);
587 AliInfo(Form("processing event %d", iEvent));
588 fRunLoader->GetEvent(iEvent);
591 sprintf(fileName, "ESD_%d.%d_final.root",
592 fRunLoader->GetHeader()->GetRun(),
593 fRunLoader->GetHeader()->GetEventNrInRun());
594 if (!gSystem->AccessPathName(fileName)) continue;
596 // local reconstruction
597 if (!fRunLocalReconstruction.IsNull()) {
598 if (!RunLocalEventReconstruction(fRunLocalReconstruction)) {
599 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
603 esd = new AliESD; hltesd = new AliESD;
604 esd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
605 hltesd->SetRunNumber(fRunLoader->GetHeader()->GetRun());
606 esd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
607 hltesd->SetEventNumber(fRunLoader->GetHeader()->GetEventNrInRun());
609 esd->SetMagneticField(AliTracker::GetBz());
610 hltesd->SetMagneticField(AliTracker::GetBz());
616 if (fRunVertexFinder) {
617 if (!ReadESD(esd, "vertex")) {
618 if (!RunVertexFinder(esd)) {
619 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
621 if (fCheckPointLevel > 0) WriteESD(esd, "vertex");
626 if (!fRunTracking.IsNull()) {
627 if (fRunHLTTracking) {
628 hltesd->SetVertex(esd->GetVertex());
629 if (!RunHLTTracking(hltesd)) {
630 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
636 if (!fRunTracking.IsNull()) {
637 if (!ReadESD(esd, "tracking")) {
638 if (!RunTracking(esd)) {
639 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
641 if (fCheckPointLevel > 0) WriteESD(esd, "tracking");
646 if (!fFillESD.IsNull()) {
647 if (!FillESD(esd, fFillESD)) {
648 if (fStopOnError) {CleanUp(file, fileOld); return kFALSE;}
653 AliESDpid::MakePID(esd);
654 if (fCheckPointLevel > 1) WriteESD(esd, "PID");
661 if (fCheckPointLevel > 0) WriteESD(esd, "final");
663 delete esd; delete hltesd;
664 esd = NULL; hltesd = NULL;
667 AliInfo(Form("Execution time for filling ESD : R:%.2fs C:%.2fs",
668 stopwatch.RealTime(),stopwatch.CpuTime()));
674 // Create tags for the events in the ESD tree (the ESD tree is always present)
675 // In case of empty events the tags will contain dummy values
677 CleanUp(file, fileOld);
683 //_____________________________________________________________________________
684 Bool_t AliReconstruction::RunLocalReconstruction(const TString& detectors)
686 // run the local reconstruction
688 TStopwatch stopwatch;
691 TString detStr = detectors;
692 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
693 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
694 AliReconstructor* reconstructor = GetReconstructor(iDet);
695 if (!reconstructor) continue;
696 if (reconstructor->HasLocalReconstruction()) continue;
698 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
699 TStopwatch stopwatchDet;
700 stopwatchDet.Start();
702 fRawReader->RewindEvents();
703 reconstructor->Reconstruct(fRunLoader, fRawReader);
705 reconstructor->Reconstruct(fRunLoader);
707 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
708 fgkDetectorName[iDet],
709 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
712 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
713 AliError(Form("the following detectors were not found: %s",
715 if (fStopOnError) return kFALSE;
718 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
719 stopwatch.RealTime(),stopwatch.CpuTime()));
724 //_____________________________________________________________________________
725 Bool_t AliReconstruction::RunLocalEventReconstruction(const TString& detectors)
727 // run the local reconstruction
729 TStopwatch stopwatch;
732 TString detStr = detectors;
733 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
734 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
735 AliReconstructor* reconstructor = GetReconstructor(iDet);
736 if (!reconstructor) continue;
737 AliLoader* loader = fLoader[iDet];
739 // conversion of digits
740 if (fRawReader && reconstructor->HasDigitConversion()) {
741 AliInfo(Form("converting raw data digits into root objects for %s",
742 fgkDetectorName[iDet]));
743 TStopwatch stopwatchDet;
744 stopwatchDet.Start();
745 loader->LoadDigits("update");
746 loader->CleanDigits();
747 loader->MakeDigitsContainer();
748 TTree* digitsTree = loader->TreeD();
749 reconstructor->ConvertDigits(fRawReader, digitsTree);
750 loader->WriteDigits("OVERWRITE");
751 loader->UnloadDigits();
752 AliInfo(Form("Execution time for %s: R:%.2fs C:%.2fs",
753 fgkDetectorName[iDet],
754 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
757 // local reconstruction
758 if (!reconstructor->HasLocalReconstruction()) continue;
759 AliInfo(Form("running reconstruction for %s", fgkDetectorName[iDet]));
760 TStopwatch stopwatchDet;
761 stopwatchDet.Start();
762 loader->LoadRecPoints("update");
763 loader->CleanRecPoints();
764 loader->MakeRecPointsContainer();
765 TTree* clustersTree = loader->TreeR();
766 if (fRawReader && !reconstructor->HasDigitConversion()) {
767 reconstructor->Reconstruct(fRawReader, clustersTree);
769 loader->LoadDigits("read");
770 TTree* digitsTree = loader->TreeD();
772 AliError(Form("Can't get the %s digits tree", fgkDetectorName[iDet]));
773 if (fStopOnError) return kFALSE;
775 reconstructor->Reconstruct(digitsTree, clustersTree);
777 loader->UnloadDigits();
779 loader->WriteRecPoints("OVERWRITE");
780 loader->UnloadRecPoints();
781 AliDebug(1,Form("Execution time for %s: R:%.2fs C:%.2fs",
782 fgkDetectorName[iDet],
783 stopwatchDet.RealTime(),stopwatchDet.CpuTime()));
786 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
787 AliError(Form("the following detectors were not found: %s",
789 if (fStopOnError) return kFALSE;
792 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
793 stopwatch.RealTime(),stopwatch.CpuTime()));
798 //_____________________________________________________________________________
799 Bool_t AliReconstruction::RunVertexFinder(AliESD*& esd)
801 // run the barrel tracking
803 TStopwatch stopwatch;
806 AliESDVertex* vertex = NULL;
807 Double_t vtxPos[3] = {0, 0, 0};
808 Double_t vtxErr[3] = {0.07, 0.07, 0.1};
810 if (fRunLoader->GetHeader() && fRunLoader->GetHeader()->GenEventHeader()) {
811 fRunLoader->GetHeader()->GenEventHeader()->PrimaryVertex(mcVertex);
812 for (Int_t i = 0; i < 3; i++) vtxPos[i] = mcVertex[i];
816 AliInfo("running the ITS vertex finder");
817 if (fLoader[0]) fLoader[0]->LoadRecPoints();
818 vertex = fVertexer->FindVertexForCurrentEvent(fRunLoader->GetEventNumber());
819 if (fLoader[0]) fLoader[0]->UnloadRecPoints();
821 AliWarning("Vertex not found");
822 vertex = new AliESDVertex();
823 vertex->SetName("default");
826 vertex->SetTruePos(vtxPos); // store also the vertex from MC
827 vertex->SetName("reconstructed");
831 AliInfo("getting the primary vertex from MC");
832 vertex = new AliESDVertex(vtxPos, vtxErr);
836 vertex->GetXYZ(vtxPos);
837 vertex->GetSigmaXYZ(vtxErr);
839 AliWarning("no vertex reconstructed");
840 vertex = new AliESDVertex(vtxPos, vtxErr);
842 esd->SetVertex(vertex);
843 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
844 if (fTracker[iDet]) fTracker[iDet]->SetVertex(vtxPos, vtxErr);
848 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
849 stopwatch.RealTime(),stopwatch.CpuTime()));
854 //_____________________________________________________________________________
855 Bool_t AliReconstruction::RunHLTTracking(AliESD*& esd)
857 // run the HLT barrel tracking
859 TStopwatch stopwatch;
863 AliError("Missing runLoader!");
867 AliInfo("running HLT tracking");
869 // Get a pointer to the HLT reconstructor
870 AliReconstructor *reconstructor = GetReconstructor(fgkNDetectors-1);
871 if (!reconstructor) return kFALSE;
874 for (Int_t iDet = 1; iDet >= 0; iDet--) {
875 TString detName = fgkDetectorName[iDet];
876 AliDebug(1, Form("%s HLT tracking", detName.Data()));
877 reconstructor->SetOption(detName.Data());
878 AliTracker *tracker = reconstructor->CreateTracker(fRunLoader);
880 AliWarning(Form("couldn't create a HLT tracker for %s", detName.Data()));
881 if (fStopOnError) return kFALSE;
885 Double_t vtxErr[3]={0.005,0.005,0.010};
886 const AliESDVertex *vertex = esd->GetVertex();
887 vertex->GetXYZ(vtxPos);
888 tracker->SetVertex(vtxPos,vtxErr);
890 fLoader[iDet]->LoadRecPoints("read");
891 TTree* tree = fLoader[iDet]->TreeR();
893 AliError(Form("Can't get the %s cluster tree", detName.Data()));
896 tracker->LoadClusters(tree);
898 if (tracker->Clusters2Tracks(esd) != 0) {
899 AliError(Form("HLT %s Clusters2Tracks failed", fgkDetectorName[iDet]));
903 tracker->UnloadClusters();
908 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
909 stopwatch.RealTime(),stopwatch.CpuTime()));
914 //_____________________________________________________________________________
915 Bool_t AliReconstruction::RunTracking(AliESD*& esd)
917 // run the barrel tracking
919 TStopwatch stopwatch;
922 AliInfo("running tracking");
924 // pass 1: TPC + ITS inwards
925 for (Int_t iDet = 1; iDet >= 0; iDet--) {
926 if (!fTracker[iDet]) continue;
927 AliDebug(1, Form("%s tracking", fgkDetectorName[iDet]));
930 fLoader[iDet]->LoadRecPoints("read");
931 TTree* tree = fLoader[iDet]->TreeR();
933 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
936 fTracker[iDet]->LoadClusters(tree);
939 if (fTracker[iDet]->Clusters2Tracks(esd) != 0) {
940 AliError(Form("%s Clusters2Tracks failed", fgkDetectorName[iDet]));
943 if (fCheckPointLevel > 1) {
944 WriteESD(esd, Form("%s.tracking", fgkDetectorName[iDet]));
946 // preliminary PID in TPC needed by the ITS tracker
948 GetReconstructor(1)->FillESD(fRunLoader, esd);
949 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
950 AliESDpid::MakePID(esd);
954 // pass 2: ALL backwards
955 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
956 if (!fTracker[iDet]) continue;
957 AliDebug(1, Form("%s back propagation", fgkDetectorName[iDet]));
960 if (iDet > 1) { // all except ITS, TPC
962 fLoader[iDet]->LoadRecPoints("read");
963 tree = fLoader[iDet]->TreeR();
965 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[iDet]));
968 fTracker[iDet]->LoadClusters(tree);
972 if (fTracker[iDet]->PropagateBack(esd) != 0) {
973 AliError(Form("%s backward propagation failed", fgkDetectorName[iDet]));
976 if (fCheckPointLevel > 1) {
977 WriteESD(esd, Form("%s.back", fgkDetectorName[iDet]));
981 if (iDet > 2) { // all except ITS, TPC, TRD
982 fTracker[iDet]->UnloadClusters();
983 fLoader[iDet]->UnloadRecPoints();
985 // updated PID in TPC needed by the ITS tracker -MI
987 GetReconstructor(1)->FillESD(fRunLoader, esd);
988 GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
989 AliESDpid::MakePID(esd);
993 // write space-points to the ESD in case alignment data output
995 if (fWriteAlignmentData)
996 WriteAlignmentData(esd);
998 // pass 3: TRD + TPC + ITS refit inwards
999 for (Int_t iDet = 2; iDet >= 0; iDet--) {
1000 if (!fTracker[iDet]) continue;
1001 AliDebug(1, Form("%s inward refit", fgkDetectorName[iDet]));
1004 if (fTracker[iDet]->RefitInward(esd) != 0) {
1005 AliError(Form("%s inward refit failed", fgkDetectorName[iDet]));
1008 if (fCheckPointLevel > 1) {
1009 WriteESD(esd, Form("%s.refit", fgkDetectorName[iDet]));
1013 fTracker[iDet]->UnloadClusters();
1014 fLoader[iDet]->UnloadRecPoints();
1017 // Propagate track to the vertex - if not done by ITS
1019 Int_t ntracks = esd->GetNumberOfTracks();
1020 for (Int_t itrack=0; itrack<ntracks; itrack++){
1021 const Double_t kRadius = 3; // beam pipe radius
1022 const Double_t kMaxStep = 5; // max step
1023 const Double_t kMaxD = 123456; // max distance to prim vertex
1024 Double_t fieldZ = AliTracker::GetBz(); //
1025 AliESDtrack * track = esd->GetTrack(itrack);
1026 if (!track) continue;
1027 if (track->IsOn(AliESDtrack::kITSrefit)) continue;
1028 track->PropagateTo(kRadius, track->GetMass(),kMaxStep,kTRUE);
1029 track->RelateToVertex(esd->GetVertex(),fieldZ, kMaxD);
1032 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
1033 stopwatch.RealTime(),stopwatch.CpuTime()));
1038 //_____________________________________________________________________________
1039 Bool_t AliReconstruction::FillESD(AliESD*& esd, const TString& detectors)
1041 // fill the event summary data
1043 TStopwatch stopwatch;
1045 AliInfo("filling ESD");
1047 TString detStr = detectors;
1048 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1049 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1050 AliReconstructor* reconstructor = GetReconstructor(iDet);
1051 if (!reconstructor) continue;
1053 if (!ReadESD(esd, fgkDetectorName[iDet])) {
1054 AliDebug(1, Form("filling ESD for %s", fgkDetectorName[iDet]));
1055 TTree* clustersTree = NULL;
1056 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1057 fLoader[iDet]->LoadRecPoints("read");
1058 clustersTree = fLoader[iDet]->TreeR();
1059 if (!clustersTree) {
1060 AliError(Form("Can't get the %s clusters tree",
1061 fgkDetectorName[iDet]));
1062 if (fStopOnError) return kFALSE;
1065 if (fRawReader && !reconstructor->HasDigitConversion()) {
1066 reconstructor->FillESD(fRawReader, clustersTree, esd);
1068 TTree* digitsTree = NULL;
1069 if (fLoader[iDet]) {
1070 fLoader[iDet]->LoadDigits("read");
1071 digitsTree = fLoader[iDet]->TreeD();
1073 AliError(Form("Can't get the %s digits tree",
1074 fgkDetectorName[iDet]));
1075 if (fStopOnError) return kFALSE;
1078 reconstructor->FillESD(digitsTree, clustersTree, esd);
1079 if (fLoader[iDet]) fLoader[iDet]->UnloadDigits();
1081 if (reconstructor->HasLocalReconstruction() && fLoader[iDet]) {
1082 fLoader[iDet]->UnloadRecPoints();
1086 reconstructor->FillESD(fRunLoader, fRawReader, esd);
1088 reconstructor->FillESD(fRunLoader, esd);
1090 if (fCheckPointLevel > 2) WriteESD(esd, fgkDetectorName[iDet]);
1094 if ((detStr.CompareTo("ALL") != 0) && !detStr.IsNull()) {
1095 AliError(Form("the following detectors were not found: %s",
1097 if (fStopOnError) return kFALSE;
1100 AliInfo(Form("Execution time: R:%.2fs C:%.2fs",
1101 stopwatch.RealTime(),stopwatch.CpuTime()));
1107 //_____________________________________________________________________________
1108 Bool_t AliReconstruction::IsSelected(TString detName, TString& detectors) const
1110 // check whether detName is contained in detectors
1111 // if yes, it is removed from detectors
1113 // check if all detectors are selected
1114 if ((detectors.CompareTo("ALL") == 0) ||
1115 detectors.BeginsWith("ALL ") ||
1116 detectors.EndsWith(" ALL") ||
1117 detectors.Contains(" ALL ")) {
1122 // search for the given detector
1123 Bool_t result = kFALSE;
1124 if ((detectors.CompareTo(detName) == 0) ||
1125 detectors.BeginsWith(detName+" ") ||
1126 detectors.EndsWith(" "+detName) ||
1127 detectors.Contains(" "+detName+" ")) {
1128 detectors.ReplaceAll(detName, "");
1132 // clean up the detectors string
1133 while (detectors.Contains(" ")) detectors.ReplaceAll(" ", " ");
1134 while (detectors.BeginsWith(" ")) detectors.Remove(0, 1);
1135 while (detectors.EndsWith(" ")) detectors.Remove(detectors.Length()-1, 1);
1140 //_____________________________________________________________________________
1141 Bool_t AliReconstruction::InitRunLoader()
1143 // get or create the run loader
1145 if (gAlice) delete gAlice;
1148 if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
1149 // load all base libraries to get the loader classes
1150 TString libs = gSystem->GetLibraries();
1151 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1152 TString detName = fgkDetectorName[iDet];
1153 if (detName == "HLT") continue;
1154 if (libs.Contains("lib" + detName + "base.so")) continue;
1155 gSystem->Load("lib" + detName + "base.so");
1157 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
1159 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
1163 fRunLoader->CdGAFile();
1164 if (gFile->GetKey(AliRunLoader::GetGAliceName())) {
1165 if (fRunLoader->LoadgAlice() == 0) {
1166 gAlice = fRunLoader->GetAliRun();
1167 AliTracker::SetFieldMap(gAlice->Field(),fUniformField);
1170 if (!gAlice && !fRawReader) {
1171 AliError(Form("no gAlice object found in file %s",
1172 fGAliceFileName.Data()));
1177 } else { // galice.root does not exist
1179 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
1183 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data(),
1184 AliConfig::GetDefaultEventFolderName(),
1187 AliError(Form("could not create run loader in file %s",
1188 fGAliceFileName.Data()));
1192 fRunLoader->MakeTree("E");
1194 while (fRawReader->NextEvent()) {
1195 fRunLoader->SetEventNumber(iEvent);
1196 fRunLoader->GetHeader()->Reset(fRawReader->GetRunNumber(),
1198 fRunLoader->MakeTree("H");
1199 fRunLoader->TreeE()->Fill();
1202 fRawReader->RewindEvents();
1203 fRunLoader->WriteHeader("OVERWRITE");
1204 fRunLoader->CdGAFile();
1205 fRunLoader->Write(0, TObject::kOverwrite);
1206 // AliTracker::SetFieldMap(???);
1212 //_____________________________________________________________________________
1213 AliReconstructor* AliReconstruction::GetReconstructor(Int_t iDet)
1215 // get the reconstructor object and the loader for a detector
1217 if (fReconstructor[iDet]) return fReconstructor[iDet];
1219 // load the reconstructor object
1220 TPluginManager* pluginManager = gROOT->GetPluginManager();
1221 TString detName = fgkDetectorName[iDet];
1222 TString recName = "Ali" + detName + "Reconstructor";
1223 if (gAlice && !gAlice->GetDetector(detName) && (detName != "HLT")) return NULL;
1225 if (detName == "HLT") {
1226 if (!gROOT->GetClass("AliLevel3")) {
1227 gSystem->Load("libAliL3Src.so");
1228 gSystem->Load("libAliL3Misc.so");
1229 gSystem->Load("libAliL3Hough.so");
1230 gSystem->Load("libAliL3Comp.so");
1234 AliReconstructor* reconstructor = NULL;
1235 // first check if a plugin is defined for the reconstructor
1236 TPluginHandler* pluginHandler =
1237 pluginManager->FindHandler("AliReconstructor", detName);
1238 // if not, add a plugin for it
1239 if (!pluginHandler) {
1240 AliDebug(1, Form("defining plugin for %s", recName.Data()));
1241 TString libs = gSystem->GetLibraries();
1242 if (libs.Contains("lib" + detName + "base.so") ||
1243 (gSystem->Load("lib" + detName + "base.so") >= 0)) {
1244 pluginManager->AddHandler("AliReconstructor", detName,
1245 recName, detName + "rec", recName + "()");
1247 pluginManager->AddHandler("AliReconstructor", detName,
1248 recName, detName, recName + "()");
1250 pluginHandler = pluginManager->FindHandler("AliReconstructor", detName);
1252 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1253 reconstructor = (AliReconstructor*) pluginHandler->ExecPlugin(0);
1255 if (reconstructor) {
1256 TObject* obj = fOptions.FindObject(detName.Data());
1257 if (obj) reconstructor->SetOption(obj->GetTitle());
1258 reconstructor->Init(fRunLoader);
1259 fReconstructor[iDet] = reconstructor;
1262 // get or create the loader
1263 if (detName != "HLT") {
1264 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
1265 if (!fLoader[iDet]) {
1266 AliConfig::Instance()
1267 ->CreateDetectorFolders(fRunLoader->GetEventFolder(),
1269 // first check if a plugin is defined for the loader
1270 TPluginHandler* pluginHandler =
1271 pluginManager->FindHandler("AliLoader", detName);
1272 // if not, add a plugin for it
1273 if (!pluginHandler) {
1274 TString loaderName = "Ali" + detName + "Loader";
1275 AliDebug(1, Form("defining plugin for %s", loaderName.Data()));
1276 pluginManager->AddHandler("AliLoader", detName,
1277 loaderName, detName + "base",
1278 loaderName + "(const char*, TFolder*)");
1279 pluginHandler = pluginManager->FindHandler("AliLoader", detName);
1281 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
1283 (AliLoader*) pluginHandler->ExecPlugin(2, detName.Data(),
1284 fRunLoader->GetEventFolder());
1286 if (!fLoader[iDet]) { // use default loader
1287 fLoader[iDet] = new AliLoader(detName, fRunLoader->GetEventFolder());
1289 if (!fLoader[iDet]) {
1290 AliWarning(Form("couldn't get loader for %s", detName.Data()));
1291 if (fStopOnError) return NULL;
1293 fRunLoader->AddLoader(fLoader[iDet]);
1294 fRunLoader->CdGAFile();
1295 if (gFile && !gFile->IsWritable()) gFile->ReOpen("UPDATE");
1296 fRunLoader->Write(0, TObject::kOverwrite);
1301 return reconstructor;
1304 //_____________________________________________________________________________
1305 Bool_t AliReconstruction::CreateVertexer()
1307 // create the vertexer
1310 AliReconstructor* itsReconstructor = GetReconstructor(0);
1311 if (itsReconstructor) {
1312 fVertexer = itsReconstructor->CreateVertexer(fRunLoader);
1315 AliWarning("couldn't create a vertexer for ITS");
1316 if (fStopOnError) return kFALSE;
1322 //_____________________________________________________________________________
1323 Bool_t AliReconstruction::CreateTrackers(const TString& detectors)
1325 // create the trackers
1327 TString detStr = detectors;
1328 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1329 if (!IsSelected(fgkDetectorName[iDet], detStr)) continue;
1330 AliReconstructor* reconstructor = GetReconstructor(iDet);
1331 if (!reconstructor) continue;
1332 TString detName = fgkDetectorName[iDet];
1333 if (detName == "HLT") {
1334 fRunHLTTracking = kTRUE;
1338 fTracker[iDet] = reconstructor->CreateTracker(fRunLoader);
1339 if (!fTracker[iDet] && (iDet < 7)) {
1340 AliWarning(Form("couldn't create a tracker for %s", detName.Data()));
1341 if (fStopOnError) return kFALSE;
1348 //_____________________________________________________________________________
1349 void AliReconstruction::CleanUp(TFile* file, TFile* fileOld)
1351 // delete trackers and the run loader and close and delete the file
1353 for (Int_t iDet = 0; iDet < fgkNDetectors; iDet++) {
1354 delete fReconstructor[iDet];
1355 fReconstructor[iDet] = NULL;
1356 fLoader[iDet] = NULL;
1357 delete fTracker[iDet];
1358 fTracker[iDet] = NULL;
1376 gSystem->Unlink("AliESDs.old.root");
1381 //_____________________________________________________________________________
1382 Bool_t AliReconstruction::ReadESD(AliESD*& esd, const char* recStep) const
1384 // read the ESD event from a file
1386 if (!esd) return kFALSE;
1388 sprintf(fileName, "ESD_%d.%d_%s.root",
1389 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1390 if (gSystem->AccessPathName(fileName)) return kFALSE;
1392 AliInfo(Form("reading ESD from file %s", fileName));
1393 AliDebug(1, Form("reading ESD from file %s", fileName));
1394 TFile* file = TFile::Open(fileName);
1395 if (!file || !file->IsOpen()) {
1396 AliError(Form("opening %s failed", fileName));
1403 esd = (AliESD*) file->Get("ESD");
1409 //_____________________________________________________________________________
1410 void AliReconstruction::WriteESD(AliESD* esd, const char* recStep) const
1412 // write the ESD event to a file
1416 sprintf(fileName, "ESD_%d.%d_%s.root",
1417 esd->GetRunNumber(), esd->GetEventNumber(), recStep);
1419 AliDebug(1, Form("writing ESD to file %s", fileName));
1420 TFile* file = TFile::Open(fileName, "recreate");
1421 if (!file || !file->IsOpen()) {
1422 AliError(Form("opening %s failed", fileName));
1433 //_____________________________________________________________________________
1434 void AliReconstruction::CreateTag(TFile* file)
1439 Double_t fMUONMASS = 0.105658369;
1442 Double_t fThetaX, fThetaY, fPyz, fChisquare;
1443 Double_t fPxRec,fPyRec, fPzRec, fEnergy;
1445 TLorentzVector fEPvector;
1447 Float_t fZVertexCut = 10.0;
1448 Float_t fRhoVertexCut = 2.0;
1450 Float_t fLowPtCut = 1.0;
1451 Float_t fHighPtCut = 3.0;
1452 Float_t fVeryHighPtCut = 10.0;
1455 Double_t partFrac[5] = {0.01, 0.01, 0.85, 0.10, 0.05};
1457 // Creates the tags for all the events in a given ESD file
1459 Int_t nProtons, nKaons, nPions, nMuons, nElectrons;
1460 Int_t nPos, nNeg, nNeutr;
1461 Int_t nK0s, nNeutrons, nPi0s, nGamas;
1462 Int_t nCh1GeV, nCh3GeV, nCh10GeV;
1463 Int_t nMu1GeV, nMu3GeV, nMu10GeV;
1464 Int_t nEl1GeV, nEl3GeV, nEl10GeV;
1465 Float_t maxPt = .0, meanPt = .0, totalP = .0;
1467 TString fVertexName;
1469 AliRunTag *tag = new AliRunTag();
1470 AliEventTag *evTag = new AliEventTag();
1471 TTree ttag("T","A Tree with event tags");
1472 TBranch * btag = ttag.Branch("AliTAG", &tag);
1473 btag->SetCompressionLevel(9);
1475 AliInfo(Form("Creating the tags......."));
1477 if (!file || !file->IsOpen()) {
1478 AliError(Form("opening failed"));
1482 Int_t firstEvent = 0,lastEvent = 0;
1483 TTree *t = (TTree*) file->Get("esdTree");
1484 TBranch * b = t->GetBranch("ESD");
1486 b->SetAddress(&esd);
1488 tag->SetRunId(esd->GetRunNumber());
1490 Int_t iNumberOfEvents = b->GetEntries();
1491 for (Int_t iEventNumber = 0; iEventNumber < iNumberOfEvents; iEventNumber++) {
1519 b->GetEntry(iEventNumber);
1520 const AliESDVertex * vertexIn = esd->GetVertex();
1521 fVertexName = vertexIn->GetName();
1522 if(fVertexName == "default") fVertexflag = 0;
1523 for (Int_t iTrackNumber = 0; iTrackNumber < esd->GetNumberOfTracks(); iTrackNumber++) {
1524 AliESDtrack * esdTrack = esd->GetTrack(iTrackNumber);
1525 UInt_t status = esdTrack->GetStatus();
1527 //select only tracks with ITS refit
1528 if ((status&AliESDtrack::kITSrefit)==0) continue;
1529 //select only tracks with TPC refit
1530 if ((status&AliESDtrack::kTPCrefit)==0) continue;
1532 //select only tracks with the "combined PID"
1533 if ((status&AliESDtrack::kESDpid)==0) continue;
1535 esdTrack->GetPxPyPz(p);
1536 Double_t momentum = sqrt(pow(p[0],2) + pow(p[1],2) + pow(p[2],2));
1537 Double_t fPt = sqrt(pow(p[0],2) + pow(p[1],2));
1540 if(fPt > maxPt) maxPt = fPt;
1542 if(esdTrack->GetSign() > 0) {
1544 if(fPt > fLowPtCut) nCh1GeV++;
1545 if(fPt > fHighPtCut) nCh3GeV++;
1546 if(fPt > fVeryHighPtCut) nCh10GeV++;
1548 if(esdTrack->GetSign() < 0) {
1550 if(fPt > fLowPtCut) nCh1GeV++;
1551 if(fPt > fHighPtCut) nCh3GeV++;
1552 if(fPt > fVeryHighPtCut) nCh10GeV++;
1554 if(esdTrack->GetSign() == 0) nNeutr++;
1558 esdTrack->GetESDpid(prob);
1561 for(Int_t i = 0; i < AliPID::kSPECIES; i++) rcc += prob[i]*partFrac[i];
1562 if(rcc == 0.0) continue;
1565 for(Int_t i = 0; i < AliPID::kSPECIES; i++) w[i] = prob[i]*partFrac[i]/rcc;
1568 if ((w[4]>w[3])&&(w[4]>w[2])&&(w[4]>w[1])&&(w[4]>w[0])) nProtons++;
1570 if ((w[3]>w[4])&&(w[3]>w[2])&&(w[3]>w[1])&&(w[3]>w[0])) nKaons++;
1572 if ((w[2]>w[4])&&(w[2]>w[3])&&(w[2]>w[1])&&(w[2]>w[0])) nPions++;
1574 if ((w[0]>w[4])&&(w[0]>w[3])&&(w[0]>w[2])&&(w[0]>w[1])) {
1576 if(fPt > fLowPtCut) nEl1GeV++;
1577 if(fPt > fHighPtCut) nEl3GeV++;
1578 if(fPt > fVeryHighPtCut) nEl10GeV++;
1586 Int_t nMuonTracks = esd->GetNumberOfMuonTracks();
1587 // loop over all reconstructed tracks (also first track of combination)
1588 for (Int_t iTrack = 0; iTrack < nMuonTracks; iTrack++) {
1589 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
1590 if (muonTrack == 0x0) continue;
1592 // Coordinates at vertex
1593 fZ = muonTrack->GetZ();
1594 fY = muonTrack->GetBendingCoor();
1595 fX = muonTrack->GetNonBendingCoor();
1597 fThetaX = muonTrack->GetThetaX();
1598 fThetaY = muonTrack->GetThetaY();
1600 fPyz = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
1601 fPzRec = - fPyz / TMath::Sqrt(1.0 + TMath::Tan(fThetaY)*TMath::Tan(fThetaY));
1602 fPxRec = fPzRec * TMath::Tan(fThetaX);
1603 fPyRec = fPzRec * TMath::Tan(fThetaY);
1604 fCharge = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
1606 //ChiSquare of the track if needed
1607 fChisquare = muonTrack->GetChi2()/(2.0 * muonTrack->GetNHit() - 5);
1608 fEnergy = TMath::Sqrt(fMUONMASS * fMUONMASS + fPxRec * fPxRec + fPyRec * fPyRec + fPzRec * fPzRec);
1609 fEPvector.SetPxPyPzE(fPxRec, fPyRec, fPzRec, fEnergy);
1611 // total number of muons inside a vertex cut
1612 if((TMath::Abs(fZ)<fZVertexCut) && (TMath::Sqrt(fY*fY+fX*fX)<fRhoVertexCut)) {
1614 if(fEPvector.Pt() > fLowPtCut) {
1616 if(fEPvector.Pt() > fHighPtCut) {
1618 if (fEPvector.Pt() > fVeryHighPtCut) {
1626 // Fill the event tags
1628 meanPt = meanPt/ntrack;
1630 evTag->SetEventId(iEventNumber+1);
1631 evTag->SetVertexX(vertexIn->GetXv());
1632 evTag->SetVertexY(vertexIn->GetYv());
1633 evTag->SetVertexZ(vertexIn->GetZv());
1634 evTag->SetVertexZError(vertexIn->GetZRes());
1635 evTag->SetVertexFlag(fVertexflag);
1637 evTag->SetT0VertexZ(esd->GetT0zVertex());
1639 evTag->SetTrigger(esd->GetTrigger());
1641 evTag->SetZDCNeutron1Energy(esd->GetZDCN1Energy());
1642 evTag->SetZDCProton1Energy(esd->GetZDCP1Energy());
1643 evTag->SetZDCNeutron2Energy(esd->GetZDCN2Energy());
1644 evTag->SetZDCProton2Energy(esd->GetZDCP2Energy());
1645 evTag->SetZDCEMEnergy(esd->GetZDCEMEnergy());
1646 evTag->SetNumOfParticipants(esd->GetZDCParticipants());
1649 evTag->SetNumOfTracks(esd->GetNumberOfTracks());
1650 evTag->SetNumOfPosTracks(nPos);
1651 evTag->SetNumOfNegTracks(nNeg);
1652 evTag->SetNumOfNeutrTracks(nNeutr);
1654 evTag->SetNumOfV0s(esd->GetNumberOfV0s());
1655 evTag->SetNumOfCascades(esd->GetNumberOfCascades());
1656 evTag->SetNumOfKinks(esd->GetNumberOfKinks());
1657 evTag->SetNumOfPMDTracks(esd->GetNumberOfPmdTracks());
1659 evTag->SetNumOfProtons(nProtons);
1660 evTag->SetNumOfKaons(nKaons);
1661 evTag->SetNumOfPions(nPions);
1662 evTag->SetNumOfMuons(nMuons);
1663 evTag->SetNumOfElectrons(nElectrons);
1664 evTag->SetNumOfPhotons(nGamas);
1665 evTag->SetNumOfPi0s(nPi0s);
1666 evTag->SetNumOfNeutrons(nNeutrons);
1667 evTag->SetNumOfKaon0s(nK0s);
1669 evTag->SetNumOfChargedAbove1GeV(nCh1GeV);
1670 evTag->SetNumOfChargedAbove3GeV(nCh3GeV);
1671 evTag->SetNumOfChargedAbove10GeV(nCh10GeV);
1672 evTag->SetNumOfMuonsAbove1GeV(nMu1GeV);
1673 evTag->SetNumOfMuonsAbove3GeV(nMu3GeV);
1674 evTag->SetNumOfMuonsAbove10GeV(nMu10GeV);
1675 evTag->SetNumOfElectronsAbove1GeV(nEl1GeV);
1676 evTag->SetNumOfElectronsAbove3GeV(nEl3GeV);
1677 evTag->SetNumOfElectronsAbove10GeV(nEl10GeV);
1679 evTag->SetNumOfPHOSClusters(esd->GetNumberOfPHOSClusters());
1680 evTag->SetNumOfEMCALClusters(esd->GetNumberOfEMCALClusters());
1682 evTag->SetTotalMomentum(totalP);
1683 evTag->SetMeanPt(meanPt);
1684 evTag->SetMaxPt(maxPt);
1686 tag->AddEventTag(*evTag);
1688 lastEvent = iNumberOfEvents;
1694 sprintf(fileName, "Run%d.Event%d_%d.ESD.tag.root",
1695 tag->GetRunId(),firstEvent,lastEvent );
1696 AliInfo(Form("writing tags to file %s", fileName));
1697 AliDebug(1, Form("writing tags to file %s", fileName));
1699 TFile* ftag = TFile::Open(fileName, "recreate");
1708 void AliReconstruction::WriteAlignmentData(AliESD* esd)
1710 // Write space-points which are then used in the alignment procedures
1711 // For the moment only ITS, TRD and TPC
1713 // Load TOF clusters
1715 fLoader[3]->LoadRecPoints("read");
1716 TTree* tree = fLoader[3]->TreeR();
1718 AliError(Form("Can't get the %s cluster tree", fgkDetectorName[3]));
1721 fTracker[3]->LoadClusters(tree);
1723 Int_t ntracks = esd->GetNumberOfTracks();
1724 for (Int_t itrack = 0; itrack < ntracks; itrack++)
1726 AliESDtrack *track = esd->GetTrack(itrack);
1729 for (Int_t iDet = 3; iDet >= 0; iDet--)
1730 nsp += track->GetNcls(iDet);
1732 AliTrackPointArray *sp = new AliTrackPointArray(nsp);
1733 track->SetTrackPointArray(sp);
1735 for (Int_t iDet = 3; iDet >= 0; iDet--) {
1736 AliTracker *tracker = fTracker[iDet];
1737 if (!tracker) continue;
1738 Int_t nspdet = track->GetNcls(iDet);
1739 if (nspdet <= 0) continue;
1740 track->GetClusters(iDet,idx);
1744 while (isp < nspdet) {
1745 Bool_t isvalid = tracker->GetTrackPoint(idx[isp2],p); isp2++;
1746 const Int_t kNTPCmax = 159;
1747 if (iDet==1 && isp2>kNTPCmax) break; // to be fixed
1748 if (!isvalid) continue;
1749 sp->AddPoint(isptrack,&p); isptrack++; isp++;
1755 fTracker[3]->UnloadClusters();
1756 fLoader[3]->UnloadRecPoints();