]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGLF/totEt/AliAnalysisEt.cxx
Commit for Simone
[u/mrichter/AliRoot.git] / PWGLF / totEt / AliAnalysisEt.cxx
1 //_________________________________________________________________________
2 //  Utility Class for transverse energy studies
3 //  Base class for ESD & MC analysis
4 //  - reconstruction and MonteCarlo output
5 // implementation file
6 //
7 //*-- Authors: Oystein Djuvsland (Bergen), David Silvermyr (ORNL)
8 //_________________________________________________________________________
9
10 #include "AliAnalysisEt.h"
11 #include "TMath.h"
12 #include "TList.h"
13 #include "TH1F.h"
14 #include "TH2F.h"
15 #include "TH1I.h" 
16 #include "TTree.h"
17 #include <iostream>
18 #include "AliAnalysisEtCuts.h"
19 #include "AliESDtrackCuts.h"
20 #include "AliESDCaloCluster.h"
21 #include "AliVEvent.h"
22 #include "Rtypes.h"
23 #include "TString.h"
24 #include "AliCentrality.h"
25 #include "AliAnalysisEtSelector.h"
26 #include "AliAnalysisEtTrackMatchCorrections.h"
27 #include "AliAnalysisEtRecEffCorrection.h"
28 #include "TFile.h"
29 #include "TVector3.h"
30 #include "AliPIDResponse.h"
31 #include "AliTPCPIDResponse.h" 
32 #include "AliInputEventHandler.h"
33 #include "AliAnalysisManager.h"
34
35 using namespace std;
36 ClassImp(AliAnalysisEt);
37
38 /* Auxiliary Histogram variables */
39 Int_t AliAnalysisEt::fgnumOfEtaBins = 16;
40 Float_t AliAnalysisEt::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};
41
42 Int_t AliAnalysisEt::fgNumOfPtBins = 111;
43 Float_t AliAnalysisEt::fgPtAxis[117]=
44    {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,
45         0.2, .22, .24, .26, .28, 0.30, 0.32, .34, .36, .38, 0.40, .42, .44, .46, .48,
46         0.5, .52, .54, .56, .58, 0.60, 0.62, .64, .66, .68, 0.70, .72, .74, .76, .78,
47         .80, .82, .84, .86, .88, 0.90, 0.92, .94, .96, .98, 1.00,1.05, 1.1,1.15, 1.2,
48         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,
49         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,
50         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,
51         20.0,25.0,30.0,35.0,40.0, 45.0, 50.0}; 
52
53 Int_t AliAnalysisEt::fgNumOfEBins = 78;
54 Float_t AliAnalysisEt::fgEAxis[79]={0., 0.05, 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
55                                                                         1.,1.5,2.,2.5,3.,3.5,4.,4.5,5.,5.5,6.,6.5,7.,7.5,8.,8.5,9.,9.5,10.,11.,
56                                                                         12.,13.,14.,15.,16.,17.,18.,19.,20.,21.,22.,23.,24.,25.,26.,27.,28.,29.,30.,31.,
57                                                                         32.,33.,34.,35.,36.,37.,38.,39.,40.,41.,42.,43.,44.,45.,46.,47.,48.,49.,50.};
58
59 Int_t AliAnalysisEt::fgNumOfRBins = 47;
60 Float_t AliAnalysisEt::fgRAxis[48]={-2.,-1.,0.,0.0005,0.001,0.0015,0.002,0.0025,0.003,0.0035,0.004,0.0045,0.005,0.0055,0.006,0.0065,0.007,0.0075,0.008,0.0085,0.009,0.0095,0.01,
61                                                                         0.011,0.012,0.013,0.014,0.015,0.016,0.017,0.018,0.019,0.020,0.022,0.024,0.026,0.028,0.03,0.032,0.034,0.036,0.038,0.04,0.042,0.044,0.046,0.048,0.05};
62
63
64 AliAnalysisEt::AliAnalysisEt() : AliAnalysisEtCommon()
65                                ,fTmCorrections(0)
66                                ,fReCorrections(0)
67                                ,fEventSummaryTree(0)
68                                ,fAcceptedTree(0)
69                                ,fDepositTree(0)
70                                ,fTotEt(0)
71                                ,fTotEtAcc(0)
72                                ,fTotNeutralEt(0)
73                                ,fTotNeutralEtAcc(0)
74                                ,fTotChargedEt(0)
75                                ,fTotChargedEtAcc(0)
76                                ,fMultiplicity(0)
77                                ,fChargedMultiplicity(0)
78                                ,fNeutralMultiplicity(0)
79                                ,fProtonEt(0)
80                                ,fAntiProtonEt(0)
81                                ,fNeutronEt(0)
82                                ,fAntiNeutronEt(0)
83                                ,fPi0Et(0)
84                                ,fPiPlusEt(0)
85                                ,fPiMinusEt(0)
86                                ,fKPlusEt(0)
87                                ,fKMinusEt(0)
88                                ,fK0sEt(0)
89                                ,fK0lEt(0)
90                                ,fMuMinusEt(0)
91                                ,fMuPlusEt(0)
92                                ,fEMinusEt(0)
93                                ,fEPlusEt(0)
94                                ,fGammaEt(0)
95                                ,fProtonRemovedEt(0)
96                                ,fAntiProtonRemovedEt(0)
97                                ,fNeutronRemovedEt(0)
98                                ,fAntiNeutronRemovedEt(0)
99                                ,fPi0RemovedEt(0)
100                                ,fPiPlusRemovedEt(0)
101                                ,fPiMinusRemovedEt(0)
102                                ,fKPlusRemovedEt(0)
103                                ,fKMinusRemovedEt(0)
104                                ,fK0sRemovedEt(0)
105                                ,fK0lRemovedEt(0)
106                                ,fMuMinusRemovedEt(0)
107                                ,fMuPlusRemovedEt(0)
108                                ,fEMinusRemovedEt(0)
109                                ,fEPlusRemovedEt(0)
110                                ,fGammaRemovedEt(0)
111                                ,fProtonMult(0)
112                                ,fAntiProtonMult(0)
113                                ,fNeutronMult(0)
114                                ,fAntiNeutronMult(0)
115                                ,fPi0Mult(0)
116                                ,fPiPlusMult(0)
117                                ,fPiMinusMult(0)
118                                ,fKPlusMult(0)
119                                ,fKMinusMult(0)
120                                ,fK0sMult(0)
121                                ,fK0lMult(0)
122                                ,fMuMinusMult(0)
123                                ,fMuPlusMult(0)
124                                ,fEMinusMult(0)
125                                ,fEPlusMult(0)
126                                ,fGammaMult(0)
127                                ,fProtonRemovedMult(0)
128                                ,fAntiProtonRemovedMult(0)
129                                ,fNeutronRemovedMult(0)
130                                ,fAntiNeutronRemovedMult(0)
131                                ,fPi0RemovedMult(0)
132                                ,fPiPlusRemovedMult(0)
133                                ,fPiMinusRemovedMult(0)
134                                ,fKPlusRemovedMult(0)
135                                ,fKMinusRemovedMult(0)
136                                ,fK0sRemovedMult(0)
137                                ,fK0lRemovedMult(0)
138                                ,fMuMinusRemovedMult(0)
139                                ,fMuPlusRemovedMult(0)
140                                ,fEMinusRemovedMult(0)
141                                ,fEPlusRemovedMult(0)
142                                ,fGammaRemovedMult(0)
143                                ,fEnergyDeposited(0)
144                                ,fMomentumTPC(0)
145                                ,fCharge(0)
146                                ,fParticlePid(0)
147                                ,fPidProb(0)
148                                ,fTrackPassedCut(kFALSE)
149                                ,fCentClass(0)
150                                ,fDetectorRadius(0)
151                                ,fSingleCellEnergyCut(0)
152                                ,fChargedEnergyRemoved(0)
153                                ,fNeutralEnergyRemoved(0)
154                                ,fGammaEnergyAdded(0)
155                                ,fHistEt(0)
156                                ,fHistNeutralMult(0)
157                                ,fHistPhivsPtPos(0)
158                                ,fHistPhivsPtNeg(0)
159                                ,fCentrality(0)
160                                ,fMakeSparse(kFALSE)
161                                ,fCutFlow(0)
162                                ,fSelector(0)
163                                ,fPIDResponse(0)
164                
165 {}
166
167 AliAnalysisEt::~AliAnalysisEt()
168 {//Destructor
169   delete fTmCorrections;
170   delete fReCorrections;
171   if(fDepositTree){
172     fDepositTree->Clear();
173     delete fDepositTree; // optional TTree
174   }
175   if(fEventSummaryTree)
176   {
177     fEventSummaryTree->Clear();
178     delete fEventSummaryTree;
179   }
180   if(fAcceptedTree)
181   {
182     fAcceptedTree->Clear();
183     delete fAcceptedTree;
184   }
185   delete fHistEt; //Et spectrum
186   delete fHistNeutralMult; //Neutral multiplicity
187   delete fHistPhivsPtPos; //phi vs pT plot for positive tracks
188   delete fHistPhivsPtNeg; //phi vs pT Moplot for negative tracks
189   //delete fCentrality;//this code does not actually own AliCentrality so we don't have to worry about deleting it...  we just borrow it...
190   delete fCutFlow;
191   delete fSelector;
192   delete fPIDResponse;
193 }
194
195 void AliAnalysisEt::FillOutputList(TList *list)
196 { // histograms to be added to output
197     list->Add(fHistEt);
198     list->Add(fHistNeutralMult);
199
200     list->Add(fHistPhivsPtPos);
201     list->Add(fHistPhivsPtNeg);
202
203     if (fCuts) {
204       if (fCuts->GetHistMakeTree()) {
205         //list->Add(fTree);
206         list->Add(fEventSummaryTree);
207       }
208       if (fCuts->GetHistMakeTreeDeposit()) {
209         list->Add(fDepositTree);
210       }
211     }
212     
213     list->Add(fCutFlow);
214
215     AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
216     if (!man) {
217       AliFatal("Analysis manager needed");
218       return;
219     }
220
221     AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
222     if (!inputHandler) {
223       AliFatal("Input handler needed");
224       return;
225     }
226
227     //pid response object
228     fPIDResponse=inputHandler->GetPIDResponse();
229     if (!fPIDResponse) AliError("PIDResponse object was not created");
230
231 }
232
233 void AliAnalysisEt::Init()
234 {// clear variables, set up cuts and PDG info
235   AliAnalysisEtCommon::Init();
236   if(ReadCorrections("calocorrections.root") != 0)
237   {
238     // Shouldn't do this, why oh why are exceptions not allowed?
239     exit(-1);
240   }
241   ResetEventValues();
242
243 }
244
245 void AliAnalysisEt::CreateHistograms()
246 { // create histograms..
247   // histogram binning for E_T, p_T and Multiplicity: defaults for p+p
248   Int_t nbinsEt = 10000;
249   Double_t minEt = 0.0;
250   Double_t maxEt = 1000;
251   Int_t nbinsPt = 200;
252   Double_t minPt = 0;
253   Double_t maxPt = 20;
254   Int_t nbinsMult = 200;
255   Double_t minMult = -0.5; // offset -0.5 to have integer bins centered around 0
256   Double_t maxMult = nbinsMult + minMult; // 1 bin per integer value
257
258   // see if we should change histogram limits etc, and possibly create a tree
259   if (fCuts) {
260     if (fCuts->GetHistMakeTree()) {
261       CreateTrees();
262     }
263
264     nbinsMult = fCuts->GetHistNbinsMult();
265     minMult = fCuts->GetHistMinMult();
266     maxMult = fCuts->GetHistMaxMult();
267
268     nbinsEt = fCuts->GetHistNbinsTotEt();
269     minEt = fCuts->GetHistMinTotEt();
270     maxEt = fCuts->GetHistMaxTotEt();
271
272     nbinsPt = fCuts->GetHistNbinsParticlePt();
273     minPt = fCuts->GetHistMinParticlePt();
274     maxPt = fCuts->GetHistMaxParticlePt();
275   }
276
277     TString histname = "fHistEt" + fHistogramNameSuffix;
278     fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", nbinsEt, minEt, maxEt);
279     fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
280     fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
281
282     histname = "fHistNeutralMult" + fHistogramNameSuffix;
283     fHistNeutralMult = new TH1F(histname.Data(), "Neutral Multiplicity", nbinsMult, minMult, maxMult);
284     fHistNeutralMult->GetXaxis()->SetTitle("N");
285     fHistNeutralMult->GetYaxis()->SetTitle("Multiplicity");
286
287     histname = "fHistPhivsPtPos" + fHistogramNameSuffix;
288     fHistPhivsPtPos = new TH2F(histname.Data(), "Phi vs pT of positively charged tracks hitting the calorimeter",       200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
289
290     histname = "fHistPhivsPtNeg" + fHistogramNameSuffix;
291     fHistPhivsPtNeg = new TH2F(histname.Data(), "Phi vs pT of negatively charged tracks hitting the calorimeter",       200, 0, 2*TMath::Pi(), nbinsPt, minPt, maxPt);
292
293     histname = "fCutFlow" + fHistogramNameSuffix;
294     fCutFlow = new TH1I(histname, histname, 20, -0.5, 19.5);
295 }
296
297 TH2F* AliAnalysisEt::CreateEtaEHisto2D(TString name, TString title, TString ztitle)
298 {     //creates a 2-d histogram in eta and E and adds it to the list of histograms to be saved
299         TString histoname   = name + fHistogramNameSuffix;
300         
301         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgnumOfEtaBins, fgEtaAxis);
302         histo->SetYTitle("#eta");
303         histo->SetXTitle("E (GeV)");
304         histo->SetZTitle(ztitle.Data());
305         histo->Sumw2();
306         
307         return histo; 
308 }
309
310
311 TH2F* AliAnalysisEt::CreateEtaPtHisto2D(TString name, TString title, TString ztitle)
312 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
313         TString histoname   = name + fHistogramNameSuffix;
314         
315         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfPtBins, fgPtAxis, fgnumOfEtaBins, fgEtaAxis);
316         histo->SetYTitle("#eta");
317         histo->SetXTitle("p_{T}");
318         histo->SetZTitle(ztitle.Data());
319         histo->Sumw2();
320         
321         return histo; 
322 }
323
324 TH2F* AliAnalysisEt::CreateEtaEtHisto2D(TString name, TString title, TString ztitle)
325 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
326         TString histoname   = name + fHistogramNameSuffix;
327         
328         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgnumOfEtaBins, fgEtaAxis);
329         histo->SetYTitle("#eta");
330         histo->SetXTitle("E_{T}");
331         histo->SetZTitle(ztitle.Data());
332         histo->Sumw2();
333         
334         return histo; 
335 }
336
337 TH2F* AliAnalysisEt::CreateResEHisto2D(TString name, TString title, TString ztitle)
338 {     //creates a 2-d histogram in eta and E and adds it to the list of histograms to be saved
339         TString histoname   = name + fHistogramNameSuffix;
340         
341         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfEBins, fgEAxis, fgNumOfRBins, fgRAxis);
342         histo->SetYTitle("R");
343         histo->SetXTitle("E (GeV)");
344         histo->SetZTitle(ztitle.Data());
345         histo->Sumw2();
346         
347         return histo; 
348 }
349
350
351 TH2F* AliAnalysisEt::CreateResPtHisto2D(TString name, TString title, TString ztitle)
352 {     //creates a 2-d histogram in eta and phi and adds it to the list of histograms to be saved
353         TString histoname   = name + fHistogramNameSuffix;
354         
355         TH2F *histo = new TH2F(histoname.Data(),title.Data(),fgNumOfPtBins, fgPtAxis, fgNumOfRBins, fgRAxis);
356         histo->SetYTitle("R");
357         histo->SetXTitle("p_{T}");
358         histo->SetZTitle(ztitle.Data());
359         histo->Sumw2();
360         
361         return histo; 
362 }
363
364 THnSparseF* AliAnalysisEt::CreateClusterHistoSparse(TString name, TString title)
365 {     //creates a 2D sparse histogram
366         TString histoname   = name + fHistogramNameSuffix;
367         
368         Int_t nBins[4] = {200,200,240,20};
369         Double_t min[4] = {0.,-1.,70.,0.5};
370         Double_t max[4] = {50.,1.,190,20.5};
371         
372         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),4,nBins,min, max);
373     histo->GetAxis(0)->SetTitle("E");
374     histo->GetAxis(1)->SetTitle("#eta_{cluster}");
375     histo->GetAxis(2)->SetTitle("#phi_{cluster}");
376     histo->GetAxis(3)->SetTitle("n_{cells}");
377         histo->Sumw2();
378         
379         return histo; 
380 }
381
382 THnSparseF* AliAnalysisEt::CreateNeutralPartHistoSparse(TString name, TString title)
383 {     //creates a sparse neutral particle histogram
384         TString histoname   = name + fHistogramNameSuffix;
385         
386         Int_t nBins[6] = {20,200,200,200,240,20};
387         Double_t min[6] = {-1.,0.,0.,-1.,70.,0.5};
388         Double_t max[6] = {1.,50.,50.,1.,190,20.5};
389         
390         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),6,nBins,min, max);
391     histo->GetAxis(0)->SetTitle("#eta");
392     histo->GetAxis(1)->SetTitle("p_{T}");
393     histo->GetAxis(2)->SetTitle("E");
394     histo->GetAxis(3)->SetTitle("#eta_{cluster}");
395     histo->GetAxis(4)->SetTitle("#phi_{cluster}");
396     histo->GetAxis(5)->SetTitle("n_{cells}");
397         histo->Sumw2();
398         
399         return histo; 
400 }
401
402 THnSparseF* AliAnalysisEt::CreateChargedPartHistoSparse(TString name, TString title)
403 {     //creates a sparse charged particle histogram
404         TString histoname   = name + fHistogramNameSuffix;
405         
406         Int_t nBins[7] = {20,200,200,200,240,20,100};
407         Double_t min[7] = {-1.,0.,0.,-1.,70.,0.5,0.};
408         Double_t max[7] = {1.,50.,50.,1.,190,20.5,0.1};
409         
410         THnSparseF *histo = new THnSparseF(histoname.Data(),title.Data(),7,nBins,min, max);
411     histo->GetAxis(0)->SetTitle("#eta");
412     histo->GetAxis(1)->SetTitle("p_{T}");
413     histo->GetAxis(2)->SetTitle("E");
414     histo->GetAxis(3)->SetTitle("#eta_{cluster}");
415     histo->GetAxis(4)->SetTitle("#phi_{cluster}");
416     histo->GetAxis(5)->SetTitle("n_{cells}");
417     histo->GetAxis(6)->SetTitle("R_{match}");
418         histo->Sumw2();
419         
420         return histo; 
421 }
422
423
424 void AliAnalysisEt::CreateTrees()
425 { // create tree..
426   TString treename = "fEventSummaryTree" + fHistogramNameSuffix;
427   if(fCuts->GetHistMakeTree())
428   {
429   
430     fEventSummaryTree = new TTree(treename, treename);
431     fEventSummaryTree->Branch("fTotEt",&fTotEt,"fTotEt/D");
432     fEventSummaryTree->Branch("fTotEtAcc",&fTotEtAcc,"fTotEtAcc/D");
433     fEventSummaryTree->Branch("fTotNeutralEt",&fTotNeutralEt,"fTotNeutralEt/D");
434     fEventSummaryTree->Branch("fTotNeutralEtAcc",&fTotNeutralEtAcc,"fTotNeutralEtAcc/D");
435     fEventSummaryTree->Branch("fTotChargedEt",&fTotChargedEt,"fTotChargedEt/D");
436     fEventSummaryTree->Branch("fTotChargedEtAcc",&fTotChargedEtAcc,"fTotChargedEtAcc/D");
437     fEventSummaryTree->Branch("fMultiplicity",&fMultiplicity,"fMultiplicity/I");
438     fEventSummaryTree->Branch("fChargedMultiplicity",&fChargedMultiplicity,"fChargedMultiplicity/I");
439     fEventSummaryTree->Branch("fNeutralMultiplicity",&fNeutralMultiplicity,"fNeutralMultiplicity/I");
440     fEventSummaryTree->Branch("fCentClass",&fCentClass,"fCentClass/I");
441     fEventSummaryTree->Branch("fChargedEnergyRemoved", &fChargedEnergyRemoved, "fChargedEnergyRemoved/D");
442     fEventSummaryTree->Branch("fNeutralEnergyRemoved", &fNeutralEnergyRemoved, "fNeutralEnergyRemoved/D");
443     fEventSummaryTree->Branch("fGammaEnergyAdded", &fGammaEnergyAdded, "fGammaEnergyAdded/D");
444     
445
446     fEventSummaryTree->Branch("fProtonEt",&fProtonEt,"fProtonEt/D");
447     fEventSummaryTree->Branch("fAntiProtonEt",&fAntiProtonEt,"fAntiProtonEt/D");
448
449     fEventSummaryTree->Branch("fNeutronEt",&fNeutronEt,"fNeutronEt/D");
450     fEventSummaryTree->Branch("fAntiNeutronEt",&fAntiNeutronEt,"fAntiNeutronEt/D");
451
452     fEventSummaryTree->Branch("fPi0Et",&fPi0Et,"fPi0Et/D");
453     fEventSummaryTree->Branch("fPiPlusEt",&fPiPlusEt,"fPiPlusEt/D");
454     fEventSummaryTree->Branch("fPiMinusEt",&fPiMinusEt,"fPiMinusEt/D");
455
456     fEventSummaryTree->Branch("fKPlusEt",&fKPlusEt,"fKPlusEt/D");
457     fEventSummaryTree->Branch("fKMinusEt",&fKMinusEt,"fKMinusEt/D");
458     fEventSummaryTree->Branch("fK0sEt",&fK0sEt,"fK0sEt/D");
459     fEventSummaryTree->Branch("fK0lEt",&fK0lEt,"fK0lEt/D");
460
461     fEventSummaryTree->Branch("fMuMinusEt",&fMuMinusEt,"fMuMinusEt/D");
462     fEventSummaryTree->Branch("fMuPlusEt",&fMuPlusEt,"fMuPlusEt/D");
463     
464     fEventSummaryTree->Branch("fEMinusEt",&fEMinusEt,"fEMinusEt/D");
465     fEventSummaryTree->Branch("fEPlusEt",&fEPlusEt,"fEPlusEt/D");
466     
467     fEventSummaryTree->Branch("fGammaEt", &fGammaEt, "fGammaEt/D");
468     
469     fEventSummaryTree->Branch("fProtonRemovedEt",&fProtonRemovedEt,"fProtonRemovedEt/D");
470     fEventSummaryTree->Branch("fAntiProtonRemovedEt",&fAntiProtonRemovedEt,"fAntiProtonRemovedEt/D");
471
472     fEventSummaryTree->Branch("fNeutronRemovedEt",&fNeutronRemovedEt,"fNeutronRemovedEt/D");
473     fEventSummaryTree->Branch("fAntiNeutronRemovedEt",&fAntiNeutronRemovedEt,"fAntiNeutronRemovedEt/D");
474
475     fEventSummaryTree->Branch("fPi0RemovedEt",&fPi0RemovedEt,"fPi0RemovedEt/D");
476     fEventSummaryTree->Branch("fPiPlusRemovedEt",&fPiPlusRemovedEt,"fPiPlusRemovedEt/D");
477     fEventSummaryTree->Branch("fPiMinusRemovedEt",&fPiMinusRemovedEt,"fPiMinusRemovedEt/D");
478
479     fEventSummaryTree->Branch("fKPlusRemovedEt",&fKPlusRemovedEt,"fKPlusRemovedEt/D");
480     fEventSummaryTree->Branch("fKMinusRemovedEt",&fKMinusEt,"fKMinusRemovedEt/D");
481     fEventSummaryTree->Branch("fK0sRemovedEt",&fK0sEt,"fK0sRemovedEt/D");
482     fEventSummaryTree->Branch("fK0lRemovedEt",&fK0lRemovedEt,"fK0lRemovedEt/D");
483
484     fEventSummaryTree->Branch("fMuMinusRemovedEt",&fMuMinusRemovedEt,"fMuMinusRemovedEt/D");
485     fEventSummaryTree->Branch("fMuPlusRemovedEt",&fMuPlusRemovedEt,"fMuPlusRemovedEt/D");
486     
487     fEventSummaryTree->Branch("fEMinusRemovedEt",&fEMinusRemovedEt,"fEMinusRemovedEt/D");
488     fEventSummaryTree->Branch("fEPlusRemovedEt",&fEPlusRemovedEt,"fEPlusRemovedEt/D");
489     
490     fEventSummaryTree->Branch("fGammaRemovedEt", &fGammaRemovedEt, "fGammaEtRemoved/D");
491
492     fEventSummaryTree->Branch("fProtonMult",&fProtonMult,"fProtonMult/D");
493     fEventSummaryTree->Branch("fAntiProtonMult",&fAntiProtonMult,"fAntiProtonMult/D");
494
495     fEventSummaryTree->Branch("fNeutronMult",&fNeutronMult,"fNeutronMult/D");
496     fEventSummaryTree->Branch("fAntiNeutronMult",&fAntiNeutronMult,"fAntiNeutronMult/D");
497
498     fEventSummaryTree->Branch("fPi0Mult",&fPi0Mult,"fPi0Mult/D");
499     fEventSummaryTree->Branch("fPiPlusMult",&fPiPlusMult,"fPiPlusMult/D");
500     fEventSummaryTree->Branch("fPiMinusMult",&fPiMinusMult,"fPiMinusMult/D");
501
502     fEventSummaryTree->Branch("fKPlusMult",&fKPlusMult,"fKPlusMult/D");
503     fEventSummaryTree->Branch("fKMinusMult",&fKMinusMult,"fKMinusMult/D");
504     fEventSummaryTree->Branch("fK0sMult",&fK0sMult,"fK0sMult/D");
505     fEventSummaryTree->Branch("fK0lMult",&fK0lMult,"fK0lMult/D");
506
507     fEventSummaryTree->Branch("fMuMinusMult",&fMuMinusMult,"fMuMinusMult/D");
508     fEventSummaryTree->Branch("fMuPlusMult",&fMuPlusMult,"fMuPlusMult/D");
509     
510     fEventSummaryTree->Branch("fEMinusMult",&fEMinusMult,"fEMinusMult/D");
511     fEventSummaryTree->Branch("fEPlusMult",&fEPlusMult,"fEPlusMult/D");
512     
513     fEventSummaryTree->Branch("fGammaMult", &fGammaMult, "fGammaMult/D");
514     
515     fEventSummaryTree->Branch("fProtonRemovedMult",&fProtonRemovedMult,"fProtonRemovedMult/D");
516     fEventSummaryTree->Branch("fAntiProtonRemovedMult",&fAntiProtonRemovedMult,"fAntiProtonRemovedMult/D");
517
518     fEventSummaryTree->Branch("fNeutronRemovedMult",&fNeutronRemovedMult,"fNeutronRemovedMult/D");
519     fEventSummaryTree->Branch("fAntiNeutronRemovedMult",&fAntiNeutronRemovedMult,"fAntiNeutronRemovedMult/D");
520
521     fEventSummaryTree->Branch("fPi0RemovedMult",&fPi0RemovedMult,"fPi0RemovedMult/D");
522     fEventSummaryTree->Branch("fPiPlusRemovedMult",&fPiPlusRemovedMult,"fPiPlusRemovedMult/D");
523     fEventSummaryTree->Branch("fPiMinusRemovedMult",&fPiMinusRemovedMult,"fPiMinusRemovedMult/D");
524
525     fEventSummaryTree->Branch("fKPlusRemovedMult",&fKPlusRemovedMult,"fKPlusRemovedMult/D");
526     fEventSummaryTree->Branch("fKMinusRemovedMult",&fKMinusMult,"fKMinusRemovedMult/D");
527     fEventSummaryTree->Branch("fK0sRemovedMult",&fK0sMult,"fK0sRemovedMult/D");
528     fEventSummaryTree->Branch("fK0lRemovedMult",&fK0lRemovedMult,"fK0lRemovedMult/D");
529
530     fEventSummaryTree->Branch("fMuMinusRemovedMult",&fMuMinusRemovedMult,"fMuMinusRemovedMult/D");
531     fEventSummaryTree->Branch("fMuPlusRemovedMult",&fMuPlusRemovedMult,"fMuPlusRemovedMult/D");
532     
533     fEventSummaryTree->Branch("fEMinusRemovedMult",&fEMinusRemovedMult,"fEMinusRemovedMult/D");
534     fEventSummaryTree->Branch("fEPlusRemovedMult",&fEPlusRemovedMult,"fEPlusRemovedMult/D");
535     
536     fEventSummaryTree->Branch("fGammaRemovedMult", &fGammaRemovedMult, "fGammaMultRemoved/D");
537
538     
539     
540   }
541   
542   if(fCuts->GetHistMakeTreeDeposit())
543   {
544     treename = "fTreeDeposit" + fHistogramNameSuffix;
545     fDepositTree = new TTree(treename, treename);
546   
547     fDepositTree->Branch("fEnergyDeposited", &fEnergyDeposited, "fEnergyDeposited/F");
548     fDepositTree->Branch("fMomentumTPC", &fMomentumTPC, "fMomentumTPC/F");
549     fDepositTree->Branch("fCharge", &fCharge, "fCharge/S");
550     fDepositTree->Branch("fParticlePid", &fParticlePid, "fParticlePid/S");
551     fDepositTree->Branch("fPidProb", &fPidProb, "fPidProb/F");
552     fDepositTree->Branch("fTrackPassedCut", &fTrackPassedCut, "fTrackPassedCut/B");
553   }
554
555   return;
556 }
557 void AliAnalysisEt::FillHistograms()
558 { // fill histograms..
559     fHistEt->Fill(fTotEt);
560
561     fHistNeutralMult->Fill(fNeutralMultiplicity);
562
563     if (fCuts) {
564       if (fCuts->GetHistMakeTree()) {
565         fEventSummaryTree->Fill();
566       }
567     }
568     if(fCuts)
569     {
570       if(fCuts->GetHistMakeTreeDeposit())
571       {
572         fDepositTree->Fill();
573       }
574     }
575       
576 }
577
578 Int_t AliAnalysisEt::AnalyseEvent(AliVEvent *event)
579 { //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.
580   AliAnalysisEtCommon::AnalyseEvent(event);
581   //fSelector->SetEvent(event);
582   ResetEventValues();
583   return 0;
584 }
585
586 void AliAnalysisEt::ResetEventValues()
587 { // clear
588   AliAnalysisEtCommon::ResetEventValues();
589   fTotEt = 0;
590   fTotEtAcc = 0;
591   fTotNeutralEt = 0;
592   fTotChargedEt = 0;
593   fTotNeutralEtAcc = 0;
594   fTotChargedEtAcc  = 0;
595   fMultiplicity = 0;
596   fChargedMultiplicity = 0;
597   fNeutralMultiplicity = 0;
598   
599   fProtonEt = 0;
600   fAntiProtonEt = 0;
601   
602   fNeutronEt = 0;
603   fAntiNeutronEt = 0;
604   
605   fPi0Et = 0;
606   fPiPlusEt = 0;
607   fPiMinusEt = 0;
608   
609   fKPlusEt = 0;
610   fKMinusEt = 0;
611   fK0sEt = 0;
612   fK0lEt = 0;
613   
614   fMuMinusEt = 0;
615   fMuPlusEt = 0;
616   
617   fEMinusEt = 0;
618   fEPlusEt = 0;
619   
620   fGammaEt = 0;
621   
622   fProtonRemovedEt = 0;
623   fAntiProtonRemovedEt = 0;
624   
625   fNeutronRemovedEt = 0;
626   fAntiNeutronRemovedEt = 0;
627   
628   fPi0RemovedEt = 0;
629   fPiPlusRemovedEt = 0;
630   fPiMinusRemovedEt = 0;
631   
632   fKPlusRemovedEt = 0;
633   fKMinusRemovedEt = 0;
634   fK0sRemovedEt = 0;
635   fK0lRemovedEt = 0;
636   
637   fMuMinusRemovedEt = 0;
638   fMuPlusRemovedEt = 0;
639   
640   fEMinusRemovedEt = 0;
641   fEPlusRemovedEt = 0;
642   
643   fGammaRemovedEt = 0;
644   
645   
646   fProtonMult = 0;
647   fAntiProtonMult = 0;
648   
649   fNeutronMult = 0;
650   fAntiNeutronMult = 0;
651   
652   fPi0Mult = 0;
653   fPiPlusMult = 0;
654   fPiMinusMult = 0;
655   
656   fKPlusMult = 0;
657   fKMinusMult = 0;
658   fK0sMult = 0;
659   fK0lMult = 0;
660   
661   fMuMinusMult = 0;
662   fMuPlusMult = 0;
663   
664   fEMinusMult = 0;
665   fEPlusMult = 0;
666   
667   fGammaMult = 0;
668   
669   fProtonRemovedMult = 0;
670   fAntiProtonRemovedMult = 0;
671   
672   fNeutronRemovedMult = 0;
673   fAntiNeutronRemovedMult = 0;
674   
675   fPi0RemovedMult = 0;
676   fPiPlusRemovedMult = 0;
677   fPiMinusRemovedMult = 0;
678   
679   fKPlusRemovedMult = 0;
680   fKMinusRemovedMult = 0;
681   fK0sRemovedMult = 0;
682   fK0lRemovedMult = 0;
683   
684   fMuMinusRemovedMult = 0;
685   fMuPlusRemovedMult = 0;
686   
687   fEMinusRemovedMult = 0;
688   fEPlusRemovedMult = 0;
689   
690   fGammaRemovedMult = 0;
691   
692   return;
693 }
694
695 Int_t AliAnalysisEt::ReadCorrections(TString filename)
696 {
697   TFile *f = TFile::Open(filename, "READ");
698   if(f)
699   {
700     TString det = "Phos";
701     if(fHistogramNameSuffix.Contains("Emcal"))
702     {
703       det = "Emcal";
704     }
705     cout<<"Histo name suffix "<<fHistogramNameSuffix<<endl;
706     TString name = "TmCorrections" + det;
707     std::cout << name << std::endl;
708     fTmCorrections = dynamic_cast<AliAnalysisEtTrackMatchCorrections*>(f->Get(name));
709     if(!fTmCorrections)
710     {
711       cout<<"No corrections with name "<<name<<endl;
712       Printf("Could not load TM corrections");
713       return -1;
714     }
715     name = "ReCorrections" + det;
716     fReCorrections = dynamic_cast<AliAnalysisEtRecEffCorrection*>(f->Get(name));
717     if(!fReCorrections)
718     {
719       Printf("Could not load rec eff corrections");
720       return -1;
721     }
722     return 0;
723   }
724  return -1; 
725 }
726
727 Double_t AliAnalysisEt::CorrectForReconstructionEfficiency(const AliESDCaloCluster& cluster, Int_t cent)
728 {
729   Float_t pos[3];
730   cluster.GetPosition(pos);
731   TVector3 cp(pos);
732   Double_t corrEnergy = fReCorrections->CorrectedEnergy(cluster.E(), cent);
733   
734   //std::cout << "Original energy: " << cluster.E() << ", corrected energy: " << corrEnergy << std::endl;
735   return TMath::Sin(cp.Theta())*corrEnergy;
736 }
737