added low and high flux parameter options
[u/mrichter/AliRoot.git] / EMCAL / beamtest07 / TestESDCaloCluster.C
CommitLineData
35eca021 1#if !defined(__CINT__) || defined(__MAKECINT__)
2
3//Root include files
4#include <Riostream.h>
5#include <TFile.h>
6#include <TChain.h>
7#include <TParticle.h>
8#include <TNtuple.h>
9#include <TCanvas.h>
10#include <TObjArray.h>
11#include <TSystem.h>
12#include <TString.h>
13#include <TH1F.h>
14#include <TVector.h>
15#include <TParticle.h>
16
17//AliRoot include files
18#include "AliESD.h"
19#include "AliESDEvent.h"
20#include "AliEMCALLoader.h"
21#include "AliESDCaloCluster.h"
22#include "AliEMCALRecPoint.h"
23#include "AliPID.h"
24#include "AliLog.h"
25
26//macros
27
28#endif
29
30TChain * AliReadESDfromdisk(const UInt_t eventsToRead,
31 const TString dirName,
32 const TString esdTreeName,
33 const char * pattern)
34{
35 // Reads ESDs from Disk
36 TChain * rv = 0;
37
38 // create a TChain of all the files
39 TChain * cESDTree = new TChain(esdTreeName) ;
40
41 // read from the directory file until the require number of events are collected
42 void * from = gSystem->OpenDirectory(dirName) ;
43 if (!from)
44 rv = 0 ;
45 else{ // reading file names from directory
46 const char * subdir ;
47 // search all subdirectories witch matching pattern
48 while( (subdir = gSystem->GetDirEntry(from)) &&
49 (cESDTree->GetEntries() < eventsToRead)) {
50 if ( strstr(subdir, pattern) != 0 ) {
51 char file[200] ;
52 sprintf(file, "%s%s/AliESDs.root", dirName.Data(), subdir);
53 cESDTree->Add(file) ;
54 }
55 } // while file
56
57 rv = cESDTree ;
58
59 } // reading file names from directory
60 return rv ;
61}
62
63//======================================================================
64TChain * AliReadESD(const UInt_t eventsToRead,
65 const TString dirName,
66 const TString esdTreeName,
67 const char * pattern)
68{
69 // Read AliESDs files and return a Chain of events
70
71 if ( dirName == "" )
72 return 0 ;
73 if ( esdTreeName == "" )
74 return AliReadESDfromdisk(eventsToRead, dirName,"","") ;//Last 2 arguments are not necessary but pdsf compiler complains "","'
75 else if ( strcmp(pattern, "") == 0 )
76 return AliReadESDfromdisk(eventsToRead, dirName, esdTreeName,"") ;//Last argument is not necessary but pdsf compiler complains "","'
77 else
78 return AliReadESDfromdisk(eventsToRead, dirName, esdTreeName, pattern) ;
79}
80
81//=====================================================================
82// Do:
83// .L TestESDCaloCluster.C++
84// TestESDCaloCluster(number of events to process)
85//=====================================================================
86void TestESDCaloCluster(const UInt_t eventsToProcess = 5,
87 TString dirName = "./",
88 const TString esdTreeName = "esdTree",
89 const char * pattern = ".")
90{
91
92 //Create chain of esd trees
93 //AliReadESD(eventsToProcess, directoryName,esdTreeName,patternOfDirectory) ;
94 //By default the root files are in the same directory
95 TChain * t = AliReadESD(eventsToProcess, dirName,esdTreeName,pattern) ;
96
97
98 // ESD
99 AliESDEvent * esd = new AliESDEvent();
100 esd->ReadFromTree(t);
101
102 //Define few variables to be used in macro
103 TString alirunName = "" ;
104
105 //Define example histograms
106 TH1F * hEnergy = new TH1F("hEnergy","Energy Distribution",100,0.,100.);
107 TH1F * hEta = new TH1F("hEta","Eta Distribution",100,-0.7,0.7);
108 TH1F * hPhi = new TH1F("hPhi","Phi Distribution",100,0,2*TMath::Pi());
109
110 Int_t beg, end ;
111 UInt_t event ;
112 Float_t pos[3] ;
113
114 for (event = 0; event < eventsToProcess; event++) {//event loop
115 //AliInfo( Form("Event %d \n",event) );
116 Int_t nbytes = t->GetEntry(event); // store event in esd
117 //cout<<"nbytes "<<nbytes<<endl;
118 if ( nbytes == 0 ) //If nothing in ESD
119 break ;
120
121 // Check that name of file is correct
122 if (alirunName != t->GetFile()->GetName()) {
123 alirunName = t->GetFile()->GetName() ;
124 alirunName.ReplaceAll("galice.root", "AliESDs.root") ;
125 }
126
127 //get reconstructed vertex position
128
129 //Double_t vertex_position[3] ;
130 // esd->GetVertex()->GetXYZ(vertex_position) ;
131
132 cout<<"Event >>>>>>>>>>> "<<event<<endl;
133
134
135 //select EMCAL tracks only
136 end = esd->GetNumberOfCaloClusters() ;
137 beg = 0;
138 Int_t nphP ;
139 //cout<<"begin "<<beg<<" end "<<end<<endl;
140 Int_t iclus = 0;
141 for (nphP = beg; nphP < end; nphP++) {//////////////EMCAL cluster loop
142 AliESDCaloCluster * clus = esd->GetCaloCluster(nphP) ; // retrieve cluster from esd
143
144 //Get the cluster parameters
145 Float_t energy = clus->E() ;
146 cout << "Cluster " << nphP << " energy = " << energy << endl;
147
148 // Int_t mult = clus->GetNumberOfDigits() ;
149 // Float_t disp = clus->GetClusterDisp() ;
150 // UShort_t * amp = clus->GetDigitAmplitude() ;
151 // UShort_t * time = clus->GetDigitTime() ;
152 // UShort_t * index = clus->GetDigitIndex() ;
153 // Int_t iprim = clus->GetPrimaryIndex();
154
155 /*
156 clus->GetPosition(pos) ;
157 TVector3 vpos(pos[0],pos[1],pos[2]) ;
158 //Print values on screen
159 cout<<"ESD cluster "<<iclus <<"; Energy "<<energy
160 <<"; Phi "<<vpos.Phi()*180/TMath::Pi()<<"; Eta "<<vpos.Eta()<<endl;
161 //cout<<"Dispersion "<<disp<<"; multiplicity "<<mult<<endl;
162 //Get the parent main values
163 */
164 //Fill histograms
165 hEnergy->Fill(energy) ;
166 // hEta->Fill(vpos.Eta()) ;
167 //hPhi->Fill(vpos.Phi()) ;
168
169 }// track loop
170 }//////////////////////////////////////////event loop
171 hEnergy->Draw();//Draw histogram on screen
172 //Write histograms in Root file
173 TFile outf("histo.root","recreate") ;
174 hEnergy->Write() ;
175 hPhi->Write() ;
176 hEta->Write() ;
177 outf.Close() ;
178}
179
180