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 **************************************************************************/
16 //-------------------------------------------------------------------------
18 // ESD reader for jet analysis
19 // Authors: Mercedes Lopez Noriega (mercedes.lopez.noriega@cern.ch)
20 // Magali Estienne <magali.estienne@IReS.in2p3.fr>
21 //-------------------------------------------------------------------------
24 #include <Riostream.h>
26 #include <TLorentzVector.h>
28 #include <TGeoManager.h>
30 #include "AliJetESDReader.h"
31 #include "AliJetESDReaderHeader.h"
32 #include "AliESDEvent.h"
33 #include "AliESDtrack.h"
34 //#include "AliEMCALGeometry.h"
35 #include "AliJetDummyGeo.h"
36 #include "AliJetFillUnitArrayTracks.h"
37 #include "AliJetUnitArray.h"
39 ClassImp(AliJetESDReader)
41 AliJetESDReader::AliJetESDReader():
62 //____________________________________________________________________________
64 AliJetESDReader::~AliJetESDReader()
73 //____________________________________________________________________________
75 void AliJetESDReader::OpenInputFiles()
78 fChain = new TChain("esdTree");
80 // get directory and pattern name from the header
81 const char* dirName=fReaderHeader->GetDirectory();
82 const char* pattern=fReaderHeader->GetPattern();
84 // // Add files matching patters to the chain
86 void *dir = gSystem->OpenDirectory(dirName);
87 const char *name = 0x0;
88 int nesd = ((AliJetESDReaderHeader*) fReaderHeader)->GetNesd();
90 while ((name = gSystem->GetDirEntry(dir))){
91 if (a>=nesd) continue;
93 if (strstr(name,pattern)){
95 sprintf(path,"%s/%s/AliESDs.root",dirName,name);
96 fChain->AddFile(path);
101 gSystem->FreeDirectory(dir);
105 fChain->SetBranchAddress("ESD", &fESD);
107 int nMax = fChain->GetEntries();
109 printf("\n AliJetESDReader: Total number of events in chain= %d \n",nMax);
111 // set number of events in header
112 if (fReaderHeader->GetLastEvent() == -1)
113 fReaderHeader->SetLastEvent(nMax);
115 Int_t nUsr = fReaderHeader->GetLastEvent();
116 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
120 void AliJetESDReader::ConnectTree(TTree* tree, TObject* data) {
122 fChain = (TChain*) tree;
123 fESD = (AliESDEvent*) data;
125 Int_t nMax = fChain->GetEntries();
126 printf("\n AliJetESDReader: Total number of events in chain= %5d \n", nMax);
127 // set number of events in header
128 if (fReaderHeader->GetLastEvent() == -1)
129 fReaderHeader->SetLastEvent(nMax);
131 Int_t nUsr = fReaderHeader->GetLastEvent();
132 fReaderHeader->SetLastEvent(TMath::Min(nMax,nUsr));
136 //____________________________________________________________________________
138 Bool_t AliJetESDReader::FillMomentumArray(Int_t event)
140 // Fill momentum array
147 // clear momentum array
149 fDebug = fReaderHeader->GetDebug();
150 // get event from chain
151 // fChain->GetTree()->GetEntry(event);
157 // get number of tracks in event (for the loop)
158 nt = fESD->GetNumberOfTracks();
159 printf("Fill Momentum Array %5d ", nt);
161 // temporary storage of signal and pt cut flag
162 Int_t* sflag = new Int_t[nt];
163 Int_t* cflag = new Int_t[nt];
165 // get cuts set by user
166 Float_t ptMin = fReaderHeader->GetPtCut();
167 Float_t etaMin = fReaderHeader->GetFiducialEtaMin();
168 Float_t etaMax = fReaderHeader->GetFiducialEtaMax();
171 for (Int_t it = 0; it < nt; it++) {
172 AliESDtrack *track = fESD->GetTrack(it);
173 UInt_t status = track->GetStatus();
176 track->GetPxPyPz(mom);
177 p3.SetXYZ(mom[0],mom[1],mom[2]);
179 if ((status & AliESDtrack::kTPCrefit) == 0) continue; // quality check
180 if ((status & AliESDtrack::kITSrefit) == 0) continue; // quality check
181 if (((AliJetESDReaderHeader*) fReaderHeader)->ReadSignalOnly()
182 && TMath::Abs(track->GetLabel()) > 10000) continue; // quality check
183 if (((AliJetESDReaderHeader*) fReaderHeader)->ReadBkgdOnly()
184 && TMath::Abs(track->GetLabel()) < 10000) continue; // quality check
186 if ( (eta > etaMax) || (eta < etaMin)) continue; // checking eta cut
188 new ((*fMomentumArray)[goodTrack]) TLorentzVector(p3,p3.Mag());
190 if (TMath::Abs(track->GetLabel()) < 10000) sflag[goodTrack]=1;
192 if (pt > ptMin) cflag[goodTrack]=1; // pt cut
195 // set the signal flags
196 fSignalFlag.Set(goodTrack,sflag);
197 fCutFlag.Set(goodTrack,cflag);
201 if (fTpcGrid || fEmcalGrid) {
204 AliJetFillUnitArrayTracks *fillUAFromTracks = new AliJetFillUnitArrayTracks();
205 fillUAFromTracks->SetReaderHeader(fReaderHeader);
206 fillUAFromTracks->SetMomentumArray(fMomentumArray);
207 fillUAFromTracks->SetTPCGrid(fTpcGrid);
208 fillUAFromTracks->SetEMCalGrid(fEmcalGrid);
209 fillUAFromTracks->SetHadCorrection(fHCorrection);
210 fillUAFromTracks->SetHadCorrector(fHadCorr);
211 fNeta = fillUAFromTracks->GetNeta();
212 fNphi = fillUAFromTracks->GetNphi();
213 fillUAFromTracks->SetActive(kFALSE);
214 // TPC only or Digits+TPC or Clusters+TPC
215 if(fOpt%2==!0 && fOpt!=0) {
216 fillUAFromTracks->SetActive(kTRUE);
217 fillUAFromTracks->SetUnitArray(fUnitArray);
218 fillUAFromTracks->ExecuteTask("tpc");
221 delete fillUAFromTracks;
228 void AliJetESDReader::SetEMCALGeometry()
230 // Define EMCAL geometry to be able to read ESDs
231 fGeom = AliJetDummyGeo::GetInstance();
233 fGeom = AliJetDummyGeo::GetInstance("SHISH_77_TRD1_2X2_FINAL_110DEG","EMCAL");
235 // To be setted to run some AliEMCALGeometry functions
236 TGeoManager::Import("geometry.root");
237 fGeom->GetTransformationForSM();
238 printf("\n EMCal Geometry set ! \n");
242 void AliJetESDReader::InitParameters()
244 // Initialise parameters
245 fHCorrection = 0; // For hadron correction
246 fHadCorr = 0; // For hadron correction
247 fNumUnits = fGeom->GetNCells(); // Number of cells in EMCAL
248 if(fDebug>1) printf("\n EMCal parameters initiated ! \n");
251 void AliJetESDReader::InitUnitArray()
253 //Initialises unit arrays
254 Int_t nElements = fTpcGrid->GetNEntries();
255 if(fArrayInitialised) delete [] fUnitArray;
256 if(fTpcGrid->GetGridType()==0)
257 fUnitArray = new AliJetUnitArray[nElements];
258 if(fTpcGrid->GetGridType()==1)
259 fUnitArray = new AliJetUnitArray[fNumUnits+nElements];
260 fArrayInitialised = 1;