]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/totEt/AliAnalysisHadEt.cxx
ALIROOT-5420 Missing include
[u/mrichter/AliRoot.git] / PWGLF / totEt / AliAnalysisHadEt.cxx
CommitLineData
cf6522d1 1//Created by Christine Nattrass, Rebecca Scott, Irakli Martashvili
641e1e0c 2//University of Tennessee at Knoxville
cf6522d1 3//
4// This class is designed for the analysis of the hadronic component of
5// transverse energy. It is used by AliAnalysisTaskHadEt.
6// This gets information about the hadronic component of the transverse energy
7// from tracks reconstructed in an event
8// it has daughters, AliAnalysisHadEtMonteCarlo and
9// AliAnalysisHadEtReconstructed which loop over either Monte Carlo data or
10// real data to get Et
3ce6b879 11
641e1e0c 12#include "AliAnalysisHadEt.h"
13#include "TMath.h"
14#include "TList.h"
15#include "TH1F.h"
16#include "TH2F.h"
17#include <iostream>
18#include "AliAnalysisEtCuts.h"
020abff0 19#include "AliMCEvent.h"
641e1e0c 20#include "AliVEvent.h"
020abff0 21#include "AliStack.h"
22#include "AliESDtrackCuts.h"
23#include "TDatabasePDG.h"
24#include "TParticle.h"
25#include "Rtypes.h"
4b40b2b1 26#include "AliPDG.h"
9a365626 27#include "AliPIDResponse.h"
28#include "AliTPCPIDResponse.h"
29#include "AliInputEventHandler.h"
30#include "AliAnalysisManager.h"
31#include "AliLog.h"
641e1e0c 32
33using namespace std;
34
35ClassImp(AliAnalysisHadEt);
36
37
4b40b2b1 38Int_t AliAnalysisHadEt::fgnumOfEtaBins = 16;
39Float_t AliAnalysisHadEt::fgEtaAxis[17]={-0.78, -0.7, -0.58, -0.46, -0.34, -0.22, -0.12, -0.06, -0.0, 0.06, 0.12, 0.22, 0.34, 0.46, 0.58, 0.7, 0.78};
020abff0 40Int_t AliAnalysisHadEt::fgNumOfPtBins = 111;
41Float_t AliAnalysisHadEt::fgPtAxis[117]=
641e1e0c 42 {0.0,0.01,0.02,0.03,0.04, 0.05, 0.06,0.07,0.08,0.09, 0.10,0.11, .12,0.13, .14,0.15, .16,0.17, .18,0.19,
43 0.2, .22, .24, .26, .28, 0.30, 0.32, .34, .36, .38, 0.40, .42, .44, .46, .48,
44 0.5, .52, .54, .56, .58, 0.60, 0.62, .64, .66, .68, 0.70, .72, .74, .76, .78,
45 .80, .82, .84, .86, .88, 0.90, 0.92, .94, .96, .98, 1.00,1.05, 1.1,1.15, 1.2,
46 1.25, 1.3,1.35,1.40,1.45, 1.50, 1.55, 1.6,1.65, 1.7, 1.75, 1.8,1.85, 1.9,1.95,
47 2.0, 2.2, 2.4, 2.6, 2.8, 3.00, 3.20, 3.4, 3.6, 3.8, 4.00, 4.2, 4.4, 4.6, 4.8,
48 5.0, 5.5, 6.0, 6.5, 7.0, 7.50, 8.00, 8.5, 9.0, 9.5, 10.0,12.0,14.0,16.0,18.0,
66e4ff36 49 20.0,25.0,30.0,35.0,40.0, 45.0, 50.0};
50Float_t AliAnalysisHadEt::fgResAxis[81] = {-0.150,-0.140,-0.130,-0.120,-0.110,-0.100,-0.090,-0.080,-0.070,-0.060,
51 -0.050,-0.045,-0.040,-0.035,-0.030,-0.025,-0.024,-0.023,-0.022,-0.021,
52 -0.020,-0.019,-0.018,-0.017,-0.016,-0.015,-0.014,-0.013,-0.012,-0.011,
53 -0.010,-0.009,-0.008,-0.007,-0.006,-0.005,-0.004,-0.003,-0.002,-0.001,
54 -0.000,0.001,0.002,0.003,0.004,0.005,0.006,0.007,0.008,0.009,
55 0.010,0.011,0.012,0.013,0.014,0.015,0.016,0.017,0.018,0.019,
56 0.020,0.021,0.022,0.023,0.024,0.025,0.030,0.035,0.040,0.045,
57 0.050,0.060,0.070,0.080,0.090,0.100,0.110,0.120,0.130,0.140,
58 0.150,};
59// Float_t AliAnalysisHadEt::fgResAxis[31] = {-0.15,-0.14,-0.13,-0.12,-0.11,-0.1,-0.09,-0.08,-0.07,-0.06,-0.05,-0.04,-0.03,-0.02,-0.01,0.0,0.01,0.02,0.03,0.04,0.05,0.06,0.07,0.08,0.09,0.1,0.11,0.12,0.13,0.14,0.15};
60Int_t AliAnalysisHadEt::fgNumOfResBins = 80;
641e1e0c 61
d6214a64 62
7d2d1773 63AliAnalysisHadEt::AliAnalysisHadEt() : AliAnalysisEtCommon()
9a365626 64 ,fPIDResponse(0)
7d2d1773 65 ,fSumEt(0)
66 ,fSumEtAcc(0)
67 ,fTotEt(0)
68 ,fTotEtAcc(0)
69 ,fTotNeutralEt(0)
70 ,fTotNeutralEtAcc(0)
71 ,fTotChargedEt(0)
72 ,fTotChargedEtAcc(0)
73 ,fMultiplicity(0)
74 ,fChargedMultiplicity(0)
75 ,fNeutralMultiplicity(0)
76 ,fhistoList(0)
66e4ff36 77 ,fGoodEvent(0)
020abff0 78{//default constructor
641e1e0c 79
80}
81
82AliAnalysisHadEt::~AliAnalysisHadEt()
020abff0 83{//destructor
951efd81 84 delete fEsdtrackCutsITSTPC;
85 delete fEsdtrackCutsITS;
86 delete fEsdtrackCutsTPC;
9a365626 87 delete fPIDResponse;
641e1e0c 88}
89
90Int_t AliAnalysisHadEt::AnalyseEvent(AliVEvent *event)
020abff0 91{ //this line is basically here to eliminate a compiler warning that event is not used. Making it a virtual function did not work with the plugin.
7d2d1773 92 AliAnalysisEtCommon::AnalyseEvent(event);
0065531a 93 ResetEventValues();
641e1e0c 94 return 0;
95}
96
641e1e0c 97
98void AliAnalysisHadEt::Init()
0065531a 99{// clear variables, set up cuts and PDG info
7d2d1773 100 AliAnalysisEtCommon::Init();
9a365626 101
641e1e0c 102}
103
104void AliAnalysisHadEt::ResetEventValues()
020abff0 105{//Resets event values of et to zero
7d2d1773 106 AliAnalysisEtCommon::ResetEventValues();
0065531a 107 fTotEt = 0;
108 fTotEtAcc = 0;
109 fTotNeutralEt = 0;
110 fTotNeutralEtAcc = 0;
111 fTotChargedEt = 0;
112 fTotChargedEtAcc = 0;
113 fMultiplicity = 0;
114 fChargedMultiplicity = 0;
115 fNeutralMultiplicity = 0;
116
4998becf 117}
641e1e0c 118void AliAnalysisHadEt::CreateEtaPtHisto2D(TString name, TString title)
020abff0 119{ //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
641e1e0c 120 TString *histoname = new TString();
121 TString *histotitle = new TString();
122
123 histoname->Append(name);
124 histotitle->Append(title);
641e1e0c 125
020abff0 126 TH2F *histo = new TH2F(histoname->Data(),histotitle->Data(),fgNumOfPtBins, fgPtAxis, fgnumOfEtaBins, fgEtaAxis);
641e1e0c 127 histo->SetYTitle("#eta");
128 histo->SetXTitle("p_{T}");
129 histo->SetZTitle("E_{T}");
130 histo->Sumw2();
020abff0 131 fhistoList->Add(histo);
641e1e0c 132 delete histoname;
133 delete histotitle;
134
66e4ff36 135}
136void AliAnalysisHadEt::CreateResolutionPtHisto2D(TString name, TString title, TString xtitle, TString ytitle)
137{ //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
138 TString *histoname = new TString();
139 TString *histotitle = new TString();
140
141 histoname->Append(name);
142 histotitle->Append(title);
143
144 //TH2F *histo = new TH2F(histoname->Data(),histotitle->Data(),fgNumOfPtBins, fgPtAxis, fgNumOfResBins, fgResAxis);
145 TH2F *histo = new TH2F(histoname->Data(),histotitle->Data(),fgNumOfPtBins, fgPtAxis, fgNumOfResBins, fgResAxis);
146 histo->SetYTitle(ytitle);
147 histo->SetXTitle(xtitle);
148 histo->SetZTitle("Number of entries");
149 histo->Sumw2();
150 fhistoList->Add(histo);
151 delete histoname;
152 delete histotitle;
153
154}
155void AliAnalysisHadEt::CreatePtHisto1D(TString name, TString title, TString xtitle, TString ytitle)
156{ //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
157 TString *histoname = new TString();
158 TString *histotitle = new TString();
159
160 histoname->Append(name);
161 histotitle->Append(title);
162
163 TH1F *histo = new TH1F(histoname->Data(),histotitle->Data(),fgNumOfPtBins, fgPtAxis);
164 histo->SetYTitle(ytitle);
165 histo->SetXTitle(xtitle);
166 histo->Sumw2();
167 fhistoList->Add(histo);
168 delete histoname;
169 delete histotitle;
170
641e1e0c 171}
172
173void AliAnalysisHadEt::CreateHisto1D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Float_t xlow,Float_t xhigh)
020abff0 174{ //creates a 1d histogram of the given dimensions and adds it to the list of histograms to be saved
641e1e0c 175 TString *histoname = new TString();
176 TString *histotitle = new TString();
641e1e0c 177 histoname->Append(name);
178 histotitle->Append(title);
641e1e0c 179 TH1F *histo = new TH1F(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh);
180 histo->SetYTitle(ytitle);
181 histo->SetXTitle(xtitle);
182 histo->Sumw2();
020abff0 183 fhistoList->Add(histo);
641e1e0c 184 delete histoname;
185 delete histotitle;
641e1e0c 186}
187void AliAnalysisHadEt::CreateIntHisto1D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Int_t xlow,Int_t xhigh)
020abff0 188{ //creates a 1d integer histogram and adds it to the list of histograms to be saved
641e1e0c 189 TString *histoname = new TString();
190 TString *histotitle = new TString();
641e1e0c 191 histoname->Append(name);
192 histotitle->Append(title);
641e1e0c 193 TH1I *histo = new TH1I(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh);
194 histo->SetYTitle(ytitle);
195 histo->SetXTitle(xtitle);
196 histo->Sumw2();
020abff0 197 fhistoList->Add(histo);
641e1e0c 198 delete histoname;
199 delete histotitle;
200
201}
202void AliAnalysisHadEt::CreateHisto2D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Float_t xlow,Float_t xhigh,Int_t ybins,Float_t ylow,Float_t yhigh)
020abff0 203{ //creates a 2d histogram and adds it to the list of histograms to be saved
641e1e0c 204 TString *histoname = new TString();
205 TString *histotitle = new TString();
641e1e0c 206 histoname->Append(name);
207 histotitle->Append(title);
641e1e0c 208 TH2F *histo = new TH2F(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh,ybins,ylow,yhigh);
209 histo->SetYTitle(ytitle);
210 histo->SetXTitle(xtitle);
211 histo->Sumw2();
020abff0 212 fhistoList->Add(histo);
641e1e0c 213 delete histoname;
214 delete histotitle;
215
216}
217void AliAnalysisHadEt::CreateIntHisto2D(TString name, TString title, TString xtitle, TString ytitle,Int_t xbins, Int_t xlow,Int_t xhigh,Int_t ybins,Int_t ylow,Int_t yhigh)
020abff0 218{ //creates a 2-d integer histogram and adds it to the list of histograms to be saved
641e1e0c 219 TString *histoname = new TString();
220 TString *histotitle = new TString();
641e1e0c 221 histoname->Append(name);
222 histotitle->Append(title);
641e1e0c 223 TH2I *histo = new TH2I(histoname->Data(),histotitle->Data(),xbins,xlow,xhigh,ybins,ylow,yhigh);
224 histo->SetYTitle(ytitle);
225 histo->SetXTitle(xtitle);
226 histo->Sumw2();
020abff0 227 fhistoList->Add(histo);
641e1e0c 228 delete histoname;
229 delete histotitle;
230
231}
232
233void AliAnalysisHadEt::CreateEtaHisto1D(TString name, TString title)
020abff0 234{ //creates 1d histogram in eta and adds it to the list of histograms to be saved
641e1e0c 235 TString *histoname = new TString();
236 TString *histotitle = new TString();
641e1e0c 237 histoname->Append(name);
238 histotitle->Append(title);
020abff0 239 TH1F *histo = new TH1F(histoname->Data(),histotitle->Data(),fgnumOfEtaBins, fgEtaAxis);
641e1e0c 240 histo->SetYTitle("E_{T}");
241 histo->SetXTitle("#eta");
242 histo->Sumw2();
020abff0 243 fhistoList->Add(histo);
641e1e0c 244 delete histoname;
245 delete histotitle;
641e1e0c 246}
247void AliAnalysisHadEt::FillHisto1D(TString histname, Float_t x, Float_t weight)
020abff0 248{//fills a 1d histogram with the name histoname with the value x and the weight "weight"
641e1e0c 249 TH1F *histo;
250 TString *name = new TString();
641e1e0c 251 name->Append(histname);
020abff0 252 histo = (TH1F *)fhistoList->FindObject(name->Data());
641e1e0c 253 if(histo){
254 histo->Fill((Double_t)x, weight);
255 }
256 else{cerr<<"CorrelationMaker::FillHisto1D: no histogram "<< name->Data()<<endl;}
257 delete name;
258}
259void AliAnalysisHadEt::FillHisto2D(TString histname, Float_t x, Float_t y, Float_t weight)
020abff0 260{//fills a 2d histogram with the name histoname with the value x and the weight "weight"
641e1e0c 261 TH2F *histo;
262 TString *name = new TString();
641e1e0c 263 name->Append(histname);
020abff0 264 histo = (TH2F *)fhistoList->FindObject(name->Data());
641e1e0c 265 if(histo){
266 histo->Fill((Double_t)x,(Double_t)y, weight);
267 }
268 else{cerr<<"CorrelationMaker::FillHisto2D: no histogram "<< name->Data()<<endl;}
269 delete name;
270}
271
272
020abff0 273Float_t AliAnalysisHadEt::Et(TParticle *part, float mass){//function to calculate et in the same way as it would be calculated in a calorimeter
274 if(mass+1000<0.01){//if no mass given return default. The default argument is -1000
641e1e0c 275 if(TMath::Abs(part->GetPDG(0)->PdgCode())==2212 || TMath::Abs(part->GetPDG(0)->PdgCode())==2112){
276 if(part->GetPDG(0)->PdgCode()==-2212 || part->GetPDG(0)->PdgCode()==-2112){//antiproton or antineutron
277 //for antinucleons we specifically want to return the kinetic energy plus twice the rest mass
278 return (part->Energy()+part->GetMass())*TMath::Sin(part->Theta());
279 }
0e866ddc 280 if(part->GetPDG(0)->PdgCode()==2212 || part->GetPDG(0)->PdgCode()==2112){//proton or neutron
641e1e0c 281 //for nucleons we specifically want to return the kinetic energy only
282 return (part->Energy()-part->GetMass())*TMath::Sin(part->Theta());
283 }
284 }
285 else{//otherwise go to the default
286 return part->Energy()*TMath::Sin(part->Theta());
287 }
288 }
289 else{//otherwise use the mass that was given
290 return (TMath::Sqrt(TMath::Power(part->P(),2.0)+TMath::Power(mass,2.0)))*TMath::Sin(part->Theta());
291 }
292 return 0.0;
293}
464aa50c 294Float_t AliAnalysisHadEt::Et(Float_t p, Float_t theta, Int_t pid, Short_t charge) const {//function to calculate et in the same way as it would be calculated in a calorimeter
7d2d1773 295 if(pid==AliAnalysisEtCommon::fgPiPlusCode || pid==AliAnalysisEtCommon::fgPiMinusCode){//Nothing special for pions
4b40b2b1 296 return TMath::Sqrt(p*p + fgPionMass*fgPionMass) * TMath::Sin(theta);
3ce6b879 297 }
7d2d1773 298 if(pid==AliAnalysisEtCommon::fgKPlusCode || pid==AliAnalysisEtCommon::fgKMinusCode){//Nothing special for kaons
299 return TMath::Sqrt(p*p + AliAnalysisEtCommon::fgKaonMass*AliAnalysisEtCommon::fgKaonMass) * TMath::Sin(theta);
3ce6b879 300 }
7d2d1773 301 if(pid==AliAnalysisEtCommon::fgEPlusCode || pid==AliAnalysisEtCommon::fgEMinusCode){//Nothing special for electrons
302 return TMath::Sqrt(p*p + AliAnalysisEtCommon::fgElectronMass*AliAnalysisEtCommon::fgElectronMass) * TMath::Sin(theta);
3ce6b879 303 }
7d2d1773 304 if(pid==AliAnalysisEtCommon::fgProtonCode || pid==AliAnalysisEtCommon::fgAntiProtonCode){//But for protons we must be careful...
3ce6b879 305 if(charge<0.0){//antiprotns: kinetic energy plus twice the rest mass
7d2d1773 306 return (TMath::Sqrt(p*p + AliAnalysisEtCommon::fgProtonMass*AliAnalysisEtCommon::fgProtonMass) + AliAnalysisEtCommon::fgProtonMass) * TMath::Sin(theta);
3ce6b879 307 }
308 if(charge>0.0){//antiprotns: kinetic energy only
7d2d1773 309 return (TMath::Sqrt(p*p + AliAnalysisEtCommon::fgProtonMass*AliAnalysisEtCommon::fgProtonMass) - AliAnalysisEtCommon::fgProtonMass) * TMath::Sin(theta);
3ce6b879 310 }
311 }
66e4ff36 312 //cerr<<"Uh-oh! Et not set properly!"<<endl;
3ce6b879 313 return 0.0;
314}
7d2d1773 315
d263b18f 316Float_t AliAnalysisHadEt::TrueP(float pTrec) const {
66e4ff36 317 if(pTrec>1.0) return pTrec;
318 return pTrec/(1-599.334*pTrec+7285.15*pTrec*pTrec)+pTrec;
319}