AliRoot Core  edcc906 (edcc906)
MUONRefit.C
Go to the documentation of this file.
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 
23 
24 #if !defined(__CINT__) || defined(__MAKECINT__)
25 #include <TStopwatch.h>
26 #include <TFile.h>
27 #include <TTree.h>
28 #include <TString.h>
29 #include <Riostream.h>
30 #include <TGeoManager.h>
31 #include <TRandom.h>
32 #include <TROOT.h>
33 
34 // STEER includes
35 #include "AliESDEvent.h"
36 #include "AliESDMuonTrack.h"
37 #include "AliCDBManager.h"
38 #include "AliGeomManager.h"
39 
40 // MUON includes
41 #include "AliMUONCDB.h"
42 #include "AliMUONRecoParam.h"
43 #include "AliMUONESDInterface.h"
44 #include "AliMUONRefitter.h"
45 #include "AliMUONVDigit.h"
46 #include "AliMUONTrack.h"
47 #include "AliMUONVTrackStore.h"
48 #include "AliMUONTrackParam.h"
49 #endif
50 
51 const Int_t printLevel = 1;
52 const Bool_t reconstructFromDigits = kTRUE; // kFALSE = reconstruct from clusters
53 
54 TTree* GetESDTree(TFile *esdFile);
55 
56 //-----------------------------------------------------------------------
57 void MUONRefit(Int_t nevents = -1, const char* esdFileNameIn = "AliESDs.root", const char* esdFileNameOut = "AliESDs_New.root",
58  const char* geoFilename = "geometry.root", const char* ocdbPath = "local://$ALICE_ROOT/OCDB")
59 {
65 
66  // open the ESD file and tree
67  TFile* esdFile = TFile::Open(esdFileNameIn);
68  TTree* esdTree = GetESDTree(esdFile);
69 
70  // connect ESD event to the ESD tree
71  AliESDEvent* esd = new AliESDEvent();
72  esd->ReadFromTree(esdTree);
73 
74  // create the ESD output file and tree
75  TFile* newESDFile = TFile::Open(esdFileNameOut, "RECREATE");
76  newESDFile->SetCompressionLevel(2);
77 
78  // connect ESD event to the ESD tree (recreate track branch for backward compatibility)
79  esdTree->SetBranchStatus("*MuonTracks*",0);
80  TTree* newESDTree = esdTree->CloneTree(0);
81  esdTree->SetBranchStatus("*MuonTracks*",1);
82  esd->WriteToTree(newESDTree);
83 
84  // get run number
85  if (esdTree->GetEvent(0) <= 0) {
86  Error("MUONRefit", "no ESD object found for event 0");
87  return;
88  }
89  Int_t runNumber = esd->GetRunNumber();
90 
91  // Import TGeo geometry
92  if (!gGeoManager) {
93  AliGeomManager::LoadGeometry(geoFilename);
94  if (!gGeoManager) {
95  Error("MUONRefit", "getting geometry from file %s failed", "generated/galice.root");
96  exit(-1);
97  }
98  }
99 
100  // load necessary data from OCDB
102  AliCDBManager::Instance()->SetRun(runNumber);
103  if (!AliMUONCDB::LoadField()) return;
104  if (!AliMUONCDB::LoadMapping(kTRUE)) return;
105 
106  // reconstruction parameters for the refitting
108  Info("MUONRefit", "\n Reconstruction parameters for refitting:");
109  recoParam->Print("FULL");
110 
111  AliMUONESDInterface esdInterface;
112  AliMUONRefitter refitter(recoParam);
113  refitter.Connect(&esdInterface);
114  gRandom->SetSeed(1);
115 
116  // timer start...
117  TStopwatch timer;
118 
119  // Loop over ESD events
120  if (nevents > 0) nevents = TMath::Min(nevents,(Int_t)esdTree->GetEntries());
121  else nevents = (Int_t)esdTree->GetEntries();
122  for (Int_t iEvent = 0; iEvent < nevents; iEvent++) {
123  if (printLevel>0) cout<<endl<<" ****************event #"<<iEvent+1<<"****************"<<endl;
124 
125  //----------------------------------------------//
126  // -------------- process event --------------- //
127  //----------------------------------------------//
128  // get the ESD of current event
129  esdTree->GetEvent(iEvent);
130  if (!esd) {
131  Error("MUONRefit", "no ESD object found for event %d", iEvent);
132  return;
133  }
134  Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks();
135  if (nTracks < 1) continue;
136 
137  // load the current event
138  esdInterface.LoadEvent(*esd, kFALSE);
139 
140  // remove digits and clusters from ESD
141  esd->FindListObject("MuonClusters")->Clear("C");
142  esd->FindListObject("MuonPads")->Clear("C");
143 
144  AliMUONVTrackStore* newTrackStore = 0x0;
145  if (reconstructFromDigits) {
146 
147  // loop over digits to modify their charge
148  AliMUONVDigit *digit;
149  TIter next(esdInterface.CreateDigitIterator());
150  while ((digit = static_cast<AliMUONVDigit*>(next()))) {
151  digit->SetCharge(digit->ADC());
152  digit->Calibrated(kFALSE);
153  }
154 
155  // refit the tracks from digits
156  refitter.SetFirstClusterIndex(0);
157  newTrackStore = refitter.ReconstructFromDigits();
158 
159  } else {
160 
161  // refit the tracks from clusters
162  newTrackStore = refitter.ReconstructFromClusters();
163 
164  }
165 
166  //----------------------------------------------//
167  // ------ fill new ESD and print results ------ //
168  //----------------------------------------------//
169  // loop over the list of ESD tracks
170  TClonesArray *esdTracks = (TClonesArray*) esd->FindListObject("MuonTracks");
171  for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
172 
173  // get the ESD track
174  AliESDMuonTrack* esdTrack = (AliESDMuonTrack*) esdTracks->UncheckedAt(iTrack);
175 
176  // skip ghost tracks (leave them unchanged in the new ESD file)
177  if (!esdTrack->ContainTrackerData()) continue;
178 
179  // get the corresponding MUON track
180  AliMUONTrack* track = esdInterface.FindTrack(esdTrack->GetUniqueID());
181 
182  // Find the corresponding re-fitted MUON track
183  AliMUONTrack* newTrack = (AliMUONTrack*) newTrackStore->FindObject(esdTrack->GetUniqueID());
184 
185  // replace the content of the current ESD track or remove it if the refitting has failed
186  if (newTrack && (!recoParam->ImproveTracks() || newTrack->IsImproved())) {
187 
188  // fill the track info
189  Double_t vertex[3] = {esdTrack->GetNonBendingCoor(), esdTrack->GetBendingCoor(), esdTrack->GetZ()};
190  AliMUONESDInterface::MUONToESD(*newTrack, *esdTrack, vertex);
191 
192  // add the clusters (and the digits if previously there)
193  for (Int_t i = 0; i < newTrack->GetNClusters(); i++) {
194  AliMUONTrackParam *param = static_cast<AliMUONTrackParam*>(newTrack->GetTrackParamAtCluster()->UncheckedAt(i));
195  AliMUONESDInterface::MUONToESD(*(param->GetClusterPtr()), *esd, esdInterface.GetDigits());
196  }
197 
198  } else esdTracks->Remove(esdTrack);
199 
200  // print initial and re-fitted track parameters at first cluster if any
201  if (printLevel>0) {
202  cout<<" ----------------track #"<<iTrack+1<<"----------------"<<endl;
203  cout<<"before refit:"<<endl;
204  AliMUONTrackParam *param = (AliMUONTrackParam*) track->GetTrackParamAtCluster()->First();
205  param->Print("FULL");
206  if (printLevel>1) param->GetCovariances().Print();
207  if (!newTrack) continue;
208  cout<<"after refit:"<<endl;
209  param = (AliMUONTrackParam*) newTrack->GetTrackParamAtCluster()->First();
210  param->Print("FULL");
211  if (printLevel>1) param->GetCovariances().Print();
212  cout<<" ----------------------------------------"<<endl;
213  }
214 
215  }
216 
217  // free memory
218  delete newTrackStore;
219 
220  // fill new ESD tree with new tracks
221  esdTracks->Compress();
222  newESDTree->Fill();
223 
224  if (printLevel>0) cout<<" ****************************************"<<endl;
225  }
226 
227  // ...timer stop
228  timer.Stop();
229 
230  // write output ESD tree
231  newESDFile->cd();
232  newESDTree->Write();
233  delete newESDTree;
234  newESDFile->Close();
235 
236  // free memory
237  esdFile->Close();
238  delete esd;
239  delete recoParam;
240 
241  cout<<endl<<"time to refit: R:"<<timer.RealTime()<<" C:"<<timer.CpuTime()<<endl<<endl;
242 }
243 
244 //-----------------------------------------------------------------------
245 TTree* GetESDTree(TFile *esdFile)
246 {
249 
250  if (!esdFile || !esdFile->IsOpen()) {
251  Error("GetESDTree", "opening ESD file failed");
252  exit(-1);
253  }
254 
255  TTree* tree = (TTree*) esdFile->Get("esdTree");
256  if (!tree) {
257  Error("GetESDTree", "no ESD tree found");
258  exit(-1);
259  }
260 
261  return tree;
262 
263 }
264 
Double_t GetBendingCoor(void) const
Base class of a track container.
const Bool_t reconstructFromDigits
Definition: MUONRefit.C:52
TFile * Open(const char *filename, Long64_t &nevents)
static void MUONToESD(const AliMUONTrack &track, AliESDEvent &esd, const Double_t vertex[3], const AliMUONVDigitStore *digits=0x0, const AliMUONLocalTrigger *locTrg=0x0)
Int_t GetNClusters() const
return the number of clusters attached to the track
Definition: AliMUONTrack.h:46
void ImproveTracks(Bool_t flag)
switch on/off the track improvement and keep the default cut in sigma to apply on cluster (local chi2...
Class to describe the MUON tracks in the Event Summary Data class.
TStopwatch timer
Definition: kNNSpeed.C:15
void LoadEvent(AliESDEvent &esdEvent, Bool_t refit=kTRUE)
virtual Int_t ADC() const =0
Raw ADC value of this digit.
AliMUONVTrackStore * ReconstructFromClusters()
Track parameters in ALICE dimuon spectrometer.
Class with MUON reconstruction parameters.
AliTPCfastTrack * track
virtual void Calibrated(Bool_t value)=0
Set the calibrated status (see note 1 in AliMUONVDigit.cxx)
Bool_t LoadMapping(Bool_t segmentationOnly)
Definition: AliMUONCDB.cxx:502
virtual void SetCharge(Float_t q)=0
Set the charge of this digit.
TTree * tree
Bool_t ContainTrackerData() const
class to refit the ESD clusters/tracks
Bool_t LoadField()
Definition: AliMUONCDB.cxx:478
Double_t GetZ(void) const
Int_t GetNumberOfMuonTracks() const
Definition: AliESDEvent.h:543
Int_t GetRunNumber() const
Definition: AliESDEvent.h:141
TIterator * CreateDigitIterator() const
Converter between MUON track/cluster/digit and ESDMuon track/cluster/pad.
void WriteToTree(TTree *tree) const
TGeoManager * gGeoManager
AliMUONVTrackStore * ReconstructFromDigits()
static AliMUONRecoParam * GetLowFluxParam()
TObject * FindListObject(const char *name) const
void SetRun(Int_t run)
virtual void Print(Option_t *opt="") const
Double_t GetNonBendingCoor(void) const
void SetDefaultStorage(const char *dbString)
void SetFirstClusterIndex(Int_t index)
void ReadFromTree(TTree *tree, Option_t *opt="")
AliMUONTrack * FindTrack(UInt_t trackId) const
ABC of a MUON digit.
Definition: AliMUONVDigit.h:18
AliMUONVCluster * GetClusterPtr() const
get pointeur to associated cluster
static Int_t runNumber
Definition: pdc06_config.C:126
AliMUONVDigitStore * GetDigits() const
Return internal track store.
Reconstructed track in ALICE dimuon spectrometer.
Definition: AliMUONTrack.h:24
Bool_t IsImproved() const
return kTRUE if the track has been improved
Definition: AliMUONTrack.h:74
static AliCDBManager * Instance(TMap *entryCache=NULL, Int_t run=-1)
TTree * GetESDTree(TFile *esdFile)
Definition: MUONRefit.C:245
void MUONRefit(Int_t nevents=-1, const char *esdFileNameIn="AliESDs.root", const char *esdFileNameOut="AliESDs_New.root", const char *geoFilename="geometry.root", const char *ocdbPath="local://$ALICE_ROOT/OCDB")
Definition: MUONRefit.C:57
void Connect(const AliMUONESDInterface *esdInterface)
connect to the ESD interface containing MUON data to refit
const Int_t printLevel
Definition: MUONRefit.C:51
const TMatrixD & GetCovariances() const
static void LoadGeometry(const char *geomFileName=NULL)
virtual TObject * FindObject(const char *name) const
Find an object by name.
TObjArray * GetTrackParamAtCluster() const
virtual void Print(Option_t *option="") const