AliPhysics  d014fd0 (d014fd0)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskSOH.cxx
Go to the documentation of this file.
1 // $Id$
2 //
3 // Simulation EMCal task.
4 //
5 // Author: Saehanseul Oh
6 
7 #include "AliAnalysisTaskSOH.h"
8 
9 #include <TH1F.h>
10 #include <TH2F.h>
11 #include <TH3F.h>
12 #include <THnSparse.h>
13 
14 #include "TChain.h"
15 #include "AliAnalysisManager.h"
16 #include "AliAnalysisTask.h"
17 #include "AliEMCALRecoUtils.h"
18 #include "AliEMCALTrack.h"
19 #include "AliESDCaloCluster.h"
20 #include "AliESDEvent.h"
21 #include "AliESDInputHandler.h"
22 #include "AliESDtrack.h"
23 #include "AliESDtrackCuts.h"
24 #include "AliExternalTrackParam.h"
25 #include "AliLog.h"
26 #include "AliMCEvent.h"
27 #include "AliHeader.h"
28 #include "AliGenCocktailEventHeader.h"
29 
30 
31 ClassImp(AliAnalysisTaskSOH)
32 
33 //________________________________________________________________________
34 AliAnalysisTaskSOH::AliAnalysisTaskSOH() :
35  AliAnalysisTaskSE(),
36  fESD(0),
37  fMC(0),
38  fZVtxMax(10),
39  fEsdTrackCuts(0x0),
40  fHybridTrackCuts1(0x0),
41  fHybridTrackCuts2(0x0),
42  fTrackIndices(0x0),
43  fClusterIndices(0x0),
44  fClusterArray(0x0),
45  fMcProcess(kTRUE),
46  fTrackProcess(kTRUE),
47  fSFProcess(kFALSE),
48  fClusterProcess(kFALSE),
49  fOutputList(0x0),
50  fHEventStat(0),
51  fHScaleFactor(0),
52  fHScaleFactor100HC(0),
53  fHEOverPVsPt(0x0),
54  fHEMCalResponsePion(0x0),
55  fHEMCalResponseElec(0x0),
56  fHEMCalResponseProton(0x0),
57  fHEMCalRecdPhidEta(0x0),
58  fHEMCalRecdPhidEtaP(0x0),
59  fHEMCalRecdPhidEtaM(0x0),
60  fHEMCalRecdPhidEta_Truth(0x0),
61  fHEMCalRecdPhidEtaP_Truth(0x0),
62  fHEMCalRecdPhidEtaM_Truth(0x0),
63  fHEMCalRecdPhidEtaposEta(0x0),
64  fHEMCalRecdPhidEtanegEta(0x0),
65  fHPhotonEdiff100HC(0x0),
66  fHPhotonEdiff70HC(0),
67  fHPhotonEdiff30HC(0),
68  fHPhotonEdiff0HC(0x0),
69  fHPhotonEVsClsE(0x0),
70  fHistEsub1Pch(0x0),
71  fHistEsub2Pch(0x0),
72  fHistEsub1PchRat(0x0),
73  fHistEsub2PchRat(0x0),
74  fHClsEoverMcE_All(0x0),
75  fHClsEoverMcE_Photon(0x0),
76  fHClsEoverMcE_Elec(0x0),
77  fHClsEoverMcE_Pion(0x0),
78  fHParGenPion_p(0x0),
79  fHParGenPion_m(0x0),
80  fHParGenPion_rmInj_p(0x0),
81  fHParGenPion_rmInj_m(0x0),
82  fHDetGenFakePion(0x0),
83  fHDetRecFakePion(0x0),
84  fHDetGenSecPion(0x0),
85  fHDetRecSecPion(0x0)
86 {
87  for(Int_t i=0; i<3; i++)
88  {
89  fHDetGenPion_p[i] = 0x0;
90  fHDetRecPion_p[i] = 0x0;
91  fHDetGenPion_m[i] = 0x0;
92  fHDetRecPion_m[i] = 0x0;
93  fHDetGenPion_rmInj_p[i] = 0x0;
94  fHDetRecPion_rmInj_p[i] = 0x0;
95  fHDetGenPion_rmInj_m[i] = 0x0;
96  fHDetRecPion_rmInj_m[i] = 0x0;
97  }
98  DefineInput (0, TChain::Class());
99  DefineOutput(1, TList::Class());
100 
101 }
102 
103 
104 //________________________________________________________________________
105 AliAnalysisTaskSOH::AliAnalysisTaskSOH(const char *name) :
106  AliAnalysisTaskSE(name),
107  fESD(0),
108  fMC(0),
109  fZVtxMax(10),
110  fEsdTrackCuts(0x0),
111  fHybridTrackCuts1(0x0),
112  fHybridTrackCuts2(0x0),
113  fTrackIndices(0x0),
114  fClusterIndices(0x0),
115  fClusterArray(0x0),
116  fMcProcess(kTRUE),
117  fTrackProcess(kTRUE),
118  fSFProcess(kFALSE),
119  fClusterProcess(kFALSE),
120  fOutputList(0x0),
121  fHEventStat(0),
122  fHScaleFactor(0),
123  fHScaleFactor100HC(0),
124  fHEOverPVsPt(0x0),
125  fHEMCalResponsePion(0x0),
126  fHEMCalResponseElec(0x0),
127  fHEMCalResponseProton(0x0),
128  fHEMCalRecdPhidEta(0x0),
129  fHEMCalRecdPhidEtaP(0x0),
130  fHEMCalRecdPhidEtaM(0x0),
131  fHEMCalRecdPhidEta_Truth(0x0),
132  fHEMCalRecdPhidEtaP_Truth(0x0),
133  fHEMCalRecdPhidEtaM_Truth(0x0),
134  fHEMCalRecdPhidEtaposEta(0x0),
135  fHEMCalRecdPhidEtanegEta(0x0),
136  fHPhotonEdiff100HC(0x0),
137  fHPhotonEdiff70HC(0),
138  fHPhotonEdiff30HC(0),
139  fHPhotonEdiff0HC(0x0),
140  fHPhotonEVsClsE(0x0),
141  fHistEsub1Pch(0x0),
142  fHistEsub2Pch(0x0),
143  fHistEsub1PchRat(0x0),
144  fHistEsub2PchRat(0x0),
145  fHClsEoverMcE_All(0x0),
146  fHClsEoverMcE_Photon(0x0),
147  fHClsEoverMcE_Elec(0x0),
148  fHClsEoverMcE_Pion(0x0),
149  fHParGenPion_p(0x0),
150  fHParGenPion_m(0x0),
151  fHParGenPion_rmInj_p(0x0),
152  fHParGenPion_rmInj_m(0x0),
153  fHDetGenFakePion(0x0),
154  fHDetRecFakePion(0x0),
155  fHDetGenSecPion(0x0),
156  fHDetRecSecPion(0x0)
157 {
158  for(Int_t i=0; i<3; i++)
159  {
160  fHDetGenPion_p[i] = 0x0;
161  fHDetRecPion_p[i] = 0x0;
162  fHDetGenPion_m[i] = 0x0;
163  fHDetRecPion_m[i] = 0x0;
164  fHDetGenPion_rmInj_p[i] = 0x0;
165  fHDetRecPion_rmInj_p[i] = 0x0;
166  fHDetGenPion_rmInj_m[i] = 0x0;
167  fHDetRecPion_rmInj_m[i] = 0x0;
168  }
169 
170  // Constructor
171  // Output slot #1 writes into a TH1 container
172  DefineOutput(1, TList::Class());
173 }
174 
175 
176 //________________________________________________________________________
177 AliAnalysisTaskSOH::~AliAnalysisTaskSOH()
178 {
179  // Destructor.
180 
181  if(fEsdTrackCuts) delete fEsdTrackCuts;
182  if(fHybridTrackCuts1) delete fHybridTrackCuts1;
183  if(fHybridTrackCuts2) delete fHybridTrackCuts2;
184  if(fTrackIndices) delete fTrackIndices;
185  if(fClusterIndices) delete fClusterIndices;
186  if(fClusterArray) delete fClusterArray;
187 }
188 
189 
190 //________________________________________________________________________
191 void AliAnalysisTaskSOH::UserCreateOutputObjects()
192 {
193  // Create histograms, called once.
194 
195  OpenFile(1);
196 
197  fOutputList = new TList();
198  fOutputList->SetOwner(1);
199 
200  fHEventStat = new TH1F("fHEventStat","Event statistics for analysis",8,0,8);
201  fHEventStat->GetXaxis()->SetBinLabel(1,"Event");
202  fHEventStat->GetXaxis()->SetBinLabel(2,"cluster");
203  fHEventStat->GetXaxis()->SetBinLabel(3,"good cluster");
204  fHEventStat->GetXaxis()->SetBinLabel(4,"cls/0-truth");
205  fHEventStat->GetXaxis()->SetBinLabel(5,"cls/1-truth");
206  fHEventStat->GetXaxis()->SetBinLabel(6,"cls/2-truth");
207  fHEventStat->GetXaxis()->SetBinLabel(7,"cls/2-goodtruth");
208  fHEventStat->GetXaxis()->SetBinLabel(8,"cls/>3-truth");
209  fOutputList->Add(fHEventStat);
210 
211 
212  if(fSFProcess)
213  {
214  fHScaleFactor = new TH1F("fHScaleFactor", "Scale factor distribution without hadronic correction;Scale factor",100,0,10);
215  fOutputList->Add(fHScaleFactor);
216 
217  fHScaleFactor100HC = new TH1F("fHScaleFactor100HC", "Scale factor distribution with 100% hadronic correction;Scale factor",100,0,10);
218  fOutputList->Add(fHScaleFactor100HC);
219  }
220 
221  if(fClusterProcess)
222  {
223  fHEOverPVsPt = new TH2F("fHEOverPVsPt", "E/P vs track p_{T}; p_{T} (GeV/c); E/P", 200 , 0, 4, 200, 0, 3.2);
224  fOutputList->Add(fHEOverPVsPt);
225 
226  fHEMCalResponsePion = new TH2F("fHEMCalResponsePion", "Pion E/P vs track p_{T}; p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
227  fOutputList->Add(fHEMCalResponsePion);
228 
229  fHEMCalResponseElec = new TH2F("fHEMCalResponseElec", "Electron E/P vs track p_{T}; p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
230  fOutputList->Add(fHEMCalResponseElec);
231 
232  fHEMCalResponseProton = new TH2F("fHEMCalResponseProton", "Proton E/P vs track p_{T}; p_{T} (GeV/c); E/P", 100 , 0, 4, 100, 0, 3.2);
233  fOutputList->Add(fHEMCalResponseProton);
234 
235  fHEMCalRecdPhidEta = new TH2F("fHEMCalRecdPhidEta","EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
236  fOutputList->Add(fHEMCalRecdPhidEta);
237 
238  fHEMCalRecdPhidEtaP = new TH2F("fHEMCalRecdPhidEtaP","EMCAL Charge+ Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
239  fOutputList->Add(fHEMCalRecdPhidEtaP);
240 
241  fHEMCalRecdPhidEtaM = new TH2F("fHEMCalRecdPhidEtaM","EMCAL Charge- Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
242  fOutputList->Add(fHEMCalRecdPhidEtaM);
243 
244  fHEMCalRecdPhidEta_Truth = new TH2F("fHEMCalRecdPhidEta_Truth","EMCAL Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
245  fOutputList->Add(fHEMCalRecdPhidEta_Truth);
246 
247  fHEMCalRecdPhidEtaP_Truth = new TH2F("fHEMCalRecdPhidEtaP_Truth","EMCAL Charge+ Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
248  fOutputList->Add(fHEMCalRecdPhidEtaP_Truth);
249 
250  fHEMCalRecdPhidEtaM_Truth = new TH2F("fHEMCalRecdPhidEtaM_Truth","EMCAL Charge- Cluster-Track(Truth matched) #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
251  fOutputList->Add(fHEMCalRecdPhidEtaM_Truth);
252 
253  fHEMCalRecdPhidEtaposEta = new TH2F("fHEMCalRecdPhidEtaposEta","(+eta track) EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
254  fOutputList->Add(fHEMCalRecdPhidEtaposEta);
255 
256  fHEMCalRecdPhidEtanegEta = new TH2F("fHEMCalRecdPhidEtanegEta","(-eta track) EMCAL Cluster-Track #Delta#phi-#Delta#eta; #Delta#eta; #Delta#phi",1000,-0.1,0.1,1000,-0.5,0.5);
257  fOutputList->Add(fHEMCalRecdPhidEtanegEta);
258 
259  fHPhotonEdiff100HC = new TH2F("fHPhotonEdiff100HC","Photon (E_{Truth}- E_{calc,100% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,100% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
260  fOutputList->Add(fHPhotonEdiff100HC);
261 
262  fHPhotonEdiff70HC = new TH2F("fHPhotonEdiff70HC","Photon (E_{Truth}- E_{calc,70% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,30% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
263  fOutputList->Add(fHPhotonEdiff70HC);
264 
265  fHPhotonEdiff30HC = new TH2F("fHPhotonEdiff30HC","Photon (E_{Truth}- E_{calc,30% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{calc,30% HC})/E_{Truth}",1000,0,10,600,-4.9,1.1);
266  fOutputList->Add(fHPhotonEdiff30HC);
267 
268  fHPhotonEdiff0HC = new TH2F("fHPhotonEdiff0HC","Photon (E_{Truth}- E_{calc,0% HC})/E_{Truth} vs. E_{Truth}; E_{Truth} (GeV); (E_{Truth}- E_{cls})/E_{Truth}",1000,0,10,600,-4.9,1.1);
269  fOutputList->Add(fHPhotonEdiff0HC);
270 
271  fHPhotonEVsClsE = new TH2F("fHPhotonEVsClsE","Cluster E vs. photon E_{Truth}; photon E_{Truth} (GeV); Cluster E (GeV)",500,0,5,500,0,5);
272  fOutputList->Add(fHPhotonEVsClsE);
273 
274  fHistEsub1Pch =new TH2F("fHistEsub1Pch", "(subtracted E in 100% HC) vs. total track P, clusters with 1 matching track; total track P (GeV/c); E_{sub}(GeV)" , 1000, 0., 10, 1000, 0., 10.);
275  fOutputList->Add(fHistEsub1Pch);
276 
277  fHistEsub2Pch =new TH2F("fHistEsub2Pch", "(subtracted E in 100% HC) vs. total track P, clusters with 2 matching tracks; total track P (GeV/c); E_{sub}(GeV)" , 1000, 0., 10, 1000, 0., 10.);
278  fOutputList->Add(fHistEsub2Pch);
279 
280  fHistEsub1PchRat =new TH2F("fHistEsub1PchRat", "(subtracted E in 100% HC)/total track P vs. total track P, clusters with 1 matching track; total track P (GeV/c); E_{sub}/P_{tot}" , 1000, 0., 10, 1100, 0., 1.1);
281  fOutputList->Add(fHistEsub1PchRat);
282 
283  fHistEsub2PchRat =new TH2F("fHistEsub2PchRat", "(subtracted E in 100% HC)/total track P vs. total track P, clusters with 2 matching tracks; total track P (GeV/c); E_{sub}/P_{tot}" , 1000, 0., 10, 1100, 0., 1.1);
284  fOutputList->Add(fHistEsub2PchRat);
285 
286  Int_t bins[4] = {150, 150, 100, 200};
287  Double_t xmin[4] = {1.3, -0.8, 0, 0};
288  Double_t xmax[4] = {3.2, 0.8, 10, 2};
289 
290  fHClsEoverMcE_All = new THnSparseF("fHClsEoverMcE_All", "Cluster E/MC E, clusters with 1 matching particle; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
291  fOutputList->Add(fHClsEoverMcE_All);
292 
293  fHClsEoverMcE_Photon = new THnSparseF("fHClsEoverMcE_Photon", "Cluster E/MC E, clusters with 1 matching photon; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
294  fOutputList->Add(fHClsEoverMcE_Photon);
295 
296  fHClsEoverMcE_Elec = new THnSparseF("fHClsEoverMcE_Elec", "Cluster E/MC E, clusters with 1 matching electron; #phi; #eta; E (GeV); ClsE/McE", 4, bins, xmin, xmax);
297  fOutputList->Add(fHClsEoverMcE_Elec);
298 
299  fHClsEoverMcE_Pion = new THnSparseF("fHClsEoverMcE_Pion", "Cluster E/MC E, clusters with 1 matching pion; #phi; #eta; E (GeV); ClsE/McE",4, bins, xmin, xmax);
300  fOutputList->Add(fHClsEoverMcE_Pion);
301  }
302 
303  fHParGenPion_p = new TH3F("fHParGenPion_p","Particle level truth Phi-Eta-p_{T} distribution of #pi+", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
304  fHParGenPion_p->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
305  fHParGenPion_p->GetYaxis()->SetTitle("#eta");
306  fHParGenPion_p->GetZaxis()->SetTitle("#phi");
307  fOutputList->Add(fHParGenPion_p);
308 
309  fHParGenPion_m = new TH3F("fHParGenPion_m", "Particle level truth Phi-Eta-p_{T} distribution of all #pi-", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
310  fHParGenPion_m->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
311  fHParGenPion_m->GetYaxis()->SetTitle("#eta");
312  fHParGenPion_m->GetZaxis()->SetTitle("#phi");
313  fOutputList->Add(fHParGenPion_m);
314 
315  fHParGenPion_rmInj_p = new TH3F("fHParGenPion_rmInj_p","Particle level truth Phi-Eta-p_{T} distribution of all #pi+ without injected signal", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
316  fHParGenPion_rmInj_p->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
317  fHParGenPion_rmInj_p->GetYaxis()->SetTitle("#eta");
318  fHParGenPion_rmInj_p->GetZaxis()->SetTitle("#phi");
319  fOutputList->Add(fHParGenPion_rmInj_p);
320 
321  fHParGenPion_rmInj_m = new TH3F("fHParGenPion_rmInj_m","Particle level truth Phi-Eta-p_{T} distribution of #pi- without injected signal", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
322  fHParGenPion_rmInj_m->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
323  fHParGenPion_rmInj_m->GetYaxis()->SetTitle("#eta");
324  fHParGenPion_rmInj_m->GetZaxis()->SetTitle("#phi");
325  fOutputList->Add(fHParGenPion_rmInj_m);
326 
327 
328 
329  const char* trackCut[3] = {"cut1","cut2", "cut3"};
330  if(fMcProcess && fTrackProcess)
331  {
332  //Fake
333  fHDetGenFakePion = new TH3F("fHDetGenFakePion", "fake charged pion track Phi-Eta-p_{T} distribution",500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
334  fHDetGenFakePion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
335  fHDetGenFakePion->GetYaxis()->SetTitle("#eta");
336  fHDetGenFakePion->GetZaxis()->SetTitle("#phi");
337  fOutputList->Add(fHDetGenFakePion);
338 
339  fHDetRecFakePion = new TH3F("fHDetRecFakePion", "fake charged pion track Phi-Eta-p_{T} distribution", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
340  fHDetRecFakePion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
341  fHDetRecFakePion->GetYaxis()->SetTitle("#eta");
342  fHDetRecFakePion->GetZaxis()->SetTitle("#phi");
343  fOutputList->Add(fHDetRecFakePion);
344 
345  //Secondary
346  fHDetGenSecPion = new TH3F("fHDetGenSecPion", "secondary charged pion charged pion track Phi-Eta-p_{T} distribution", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
347  fHDetGenSecPion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
348  fHDetGenSecPion->GetYaxis()->SetTitle("#eta");
349  fHDetGenSecPion->GetZaxis()->SetTitle("#phi");
350  fOutputList->Add(fHDetGenSecPion);
351 
352  fHDetRecSecPion = new TH3F("fHDetRecSecPion", "secondary charged pion charged pion track Phi-Eta-p_{T} distribution", 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
353  fHDetRecSecPion->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
354  fHDetRecSecPion->GetYaxis()->SetTitle("#eta");
355  fHDetRecSecPion->GetZaxis()->SetTitle("#phi");
356  fOutputList->Add(fHDetRecSecPion);
357 
358  for(Int_t i=0; i<3; i++)
359  {
360  // pi+
361  fHDetGenPion_p[i] = new TH3F(Form("fHDetGenPion_p_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi+", trackCut[i]), 500, 0, 100, 60, -1.2, 1.2, 128, 0, 6.4);
362  fHDetGenPion_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
363  fHDetGenPion_p[i]->GetYaxis()->SetTitle("#eta");
364  fHDetGenPion_p[i]->GetZaxis()->SetTitle("#phi");
365  fOutputList->Add(fHDetGenPion_p[i]);
366 
367  fHDetRecPion_p[i] = new TH3F(Form("fHDetRecPion_p_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of all #pi+", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
368  fHDetRecPion_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
369  fHDetRecPion_p[i]->GetYaxis()->SetTitle("#eta");
370  fHDetRecPion_p[i]->GetZaxis()->SetTitle("#phi");
371  fOutputList->Add(fHDetRecPion_p[i]);
372 
373  // pi-
374  fHDetGenPion_m[i] = new TH3F(Form("fHDetGenPion_m_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi-", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
375  fHDetGenPion_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
376  fHDetGenPion_m[i]->GetYaxis()->SetTitle("#eta");
377  fHDetGenPion_m[i]->GetZaxis()->SetTitle("#phi");
378  fOutputList->Add(fHDetGenPion_m[i]);
379 
380  fHDetRecPion_m[i] = new TH3F(Form("fHDetRecPion_m_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of #pi-", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
381  fHDetRecPion_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
382  fHDetRecPion_m[i]->GetYaxis()->SetTitle("#eta");
383  fHDetRecPion_m[i]->GetZaxis()->SetTitle("#phi");
384  fOutputList->Add(fHDetRecPion_m[i]);
385 
386  //pi+ without injected signal
387  fHDetGenPion_rmInj_p[i] = new TH3F(Form("fHDetGenPion_rmInj_p_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi+ without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
388  fHDetGenPion_rmInj_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
389  fHDetGenPion_rmInj_p[i]->GetYaxis()->SetTitle("#eta");
390  fHDetGenPion_rmInj_p[i]->GetZaxis()->SetTitle("#phi");
391  fOutputList->Add(fHDetGenPion_rmInj_p[i]);
392 
393  fHDetRecPion_rmInj_p[i] = new TH3F(Form("fHDetRecPion_rmInj_p_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of #pi+ without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
394  fHDetRecPion_rmInj_p[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
395  fHDetRecPion_rmInj_p[i]->GetYaxis()->SetTitle("#eta");
396  fHDetRecPion_rmInj_p[i]->GetZaxis()->SetTitle("#phi");
397  fOutputList->Add(fHDetRecPion_rmInj_p[i]);
398 
399  //pi- charged particle without injected signal
400  fHDetGenPion_rmInj_m[i] = new TH3F(Form("fHDetGenPion_rmInj_m_%s", trackCut[i]), Form("%s: Detector level truth Phi-Eta-p_{T} distribution of #pi- without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
401  fHDetGenPion_rmInj_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
402  fHDetGenPion_rmInj_m[i]->GetYaxis()->SetTitle("#eta");
403  fHDetGenPion_rmInj_m[i]->GetZaxis()->SetTitle("#phi");
404  fOutputList->Add(fHDetGenPion_rmInj_m[i]);
405 
406  fHDetRecPion_rmInj_m[i] = new TH3F(Form("fHDetRecPion_rmInj_m_%s", trackCut[i]), Form("%s: Reconstructed track Phi-Eta-p_{T} distribution of #pi- without injected signal", trackCut[i]), 500, 0, 100, 100, -1.0, 1.0, 120, 0.0,2.*TMath::Pi());
407  fHDetRecPion_rmInj_m[i]->GetXaxis()->SetTitle("p_{T}^{gen} (GeV/c)");
408  fHDetRecPion_rmInj_m[i]->GetYaxis()->SetTitle("#eta");
409  fHDetRecPion_rmInj_m[i]->GetZaxis()->SetTitle("#phi");
410  fOutputList->Add(fHDetRecPion_rmInj_m[i]);
411  }
412  }
413 
414  fTrackIndices = new TArrayI();
415  fClusterIndices = new TArrayI();
416 
417  fClusterArray = new TObjArray();
418  fClusterArray->SetOwner(1);
419 
420  PostData(1, fOutputList);
421 }
422 
423 //________________________________________________________________________
424 void AliAnalysisTaskSOH::UserExec(Option_t *)
425 {
426  // Main loop, called for each event.
427 
428  // Post output data.
429  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
430  if (!fESD) {
431  printf("ERROR: fESD not available\n");
432  return;
433  }
434  if (!EsdVertexOk()) return; // Vetex cut
435 
436  fMC = MCEvent();
437  if (!fMC) {
438  printf("ERROR: fMC not available\n");
439  return;
440  }
441 
442  fHEventStat->Fill(0.5);
443 
444  if(fTrackIndices)
445  fTrackIndices->Reset();
446  if(fClusterIndices)
447  fClusterIndices->Reset();
448  if(fClusterArray)
449  fClusterArray->Delete();
450 
451  if(fTrackProcess)
452  ProcessTrack();
453  if(fClusterProcess)
454  ProcessCluster();
455  if(fMcProcess)
456  ProcessMc();
457  if(fSFProcess)
458  ProcessScaleFactor();
459 
460  PostData(1, fOutputList);
461 }
462 
463 //________________________________________________________________________
464 void AliAnalysisTaskSOH::ProcessTrack()
465 {
466  // Process track.
467 
468  fTrackIndices->Set(fESD->GetNumberOfTracks());
469  AliDebug(3,Form("%s:%d Selecting tracks",(char*)__FILE__,__LINE__));
470 
471  Int_t isMth = 0;
472  Int_t nTracks = 0;
473 
474  Float_t ClsPos[3] = {-999,-999,-999};
475  Double_t emcTrkpos[3] = {-999,-999,-999};
476 
477  for(Int_t itr=0; itr<fESD->GetNumberOfTracks(); itr++)
478  {
479  AliESDtrack *esdtrack = fESD->GetTrack(itr);
480  if(!esdtrack)continue;
481  AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
482  if(!newTrack) continue;
483  if(newTrack->Pt()<0.15 || TMath::Abs(newTrack->Eta())>0.9) {delete newTrack; continue;}
484 
485  if(fClusterProcess)
486  {
487  Double_t clsE = -1;
488  Int_t clsIndex = newTrack->GetEMCALcluster();
489  if(newTrack->GetEMCALcluster()>-1)
490  {
491  AliESDCaloCluster *cluster = fESD->GetCaloCluster(clsIndex);
492  if(IsGoodCluster(cluster))
493  {
494  isMth=1;
495 
496  cluster->GetPosition(ClsPos);
497  TVector3 VClsPos(ClsPos[0], ClsPos[1], ClsPos[2]);
498 
499  AliEMCALTrack EMCTrk(*newTrack);
500  if(!EMCTrk.PropagateToGlobal(ClsPos[0], ClsPos[1], ClsPos[2], 0.0, 0.0)) {continue;}
501  EMCTrk.GetXYZ(emcTrkpos);
502  TVector3 VemcTrkPos(emcTrkpos[0],emcTrkpos[1],emcTrkpos[2]);
503 
504  Double_t dPhi = VClsPos.Phi() - VemcTrkPos.Phi();
505  if (dPhi < -1*TMath::Pi()) dPhi += (2*TMath::Pi());
506  else if (dPhi > TMath::Pi()) dPhi -= (2*TMath::Pi());
507 
508  Double_t dEta = VClsPos.Eta() - VemcTrkPos.Eta();
509 
510  fHEMCalRecdPhidEta->Fill(dEta, dPhi);
511 
512  if((newTrack->GetLabel())>-1 && (newTrack->GetLabel()) < fMC->GetNumberOfTracks())
513  {
514  AliVParticle *vParticle = fMC->GetTrack(newTrack->GetLabel());
515  if(IsGoodMcParticle(vParticle, newTrack->GetLabel()))
516  {
517  fHEMCalRecdPhidEta_Truth->Fill(dEta, dPhi);
518  if(vParticle->Charge() > 0) fHEMCalRecdPhidEtaP_Truth->Fill(dEta, dPhi);
519  if(vParticle->Charge() < 0) fHEMCalRecdPhidEtaM_Truth->Fill(dEta, dPhi);
520  }
521  }
522 
523  if(esdtrack->Charge() > 0) {fHEMCalRecdPhidEtaP->Fill(dEta, dPhi);}
524  if(esdtrack->Charge() < 0) {fHEMCalRecdPhidEtaM->Fill(dEta, dPhi);}
525 
526  if(VemcTrkPos.Eta() > 0) fHEMCalRecdPhidEtaposEta->Fill(dEta, dPhi);
527  if(VemcTrkPos.Eta() < 0) fHEMCalRecdPhidEtanegEta->Fill(dEta, dPhi);
528 
529  clsE = cluster->E();
530  if(newTrack->P()>0) fHEOverPVsPt->Fill(newTrack->Pt(),clsE/newTrack->P());
531  }
532 
533  Int_t ipart = newTrack->GetLabel();
534  if(ipart>-1 && ipart<fMC->GetNumberOfTracks())
535  {
536  AliVParticle* vParticle = fMC->GetTrack(ipart);
537  if(isMth && vParticle)
538  {
539  if(TMath::Abs(vParticle->PdgCode())==211)
540  {
541  fHEMCalResponsePion->Fill(newTrack->Pt(),clsE/newTrack->P());
542  }
543  if(TMath::Abs(vParticle->PdgCode())==11)
544  {
545  fHEMCalResponseElec->Fill(newTrack->Pt(),clsE/newTrack->P());
546  }
547  if(TMath::Abs(vParticle->PdgCode())==2212)
548  {
549  fHEMCalResponseProton->Fill(newTrack->Pt(),clsE/newTrack->P());
550  }
551  }
552  }
553  }
554  }
555  if(newTrack) delete newTrack;
556 
557  // Track Indices
558  fTrackIndices->AddAt(itr,nTracks);
559  nTracks++;
560  }
561 
562  fTrackIndices->Set(nTracks);
563 }
564 
565 //________________________________________________________________________
566 void AliAnalysisTaskSOH::ProcessCluster()
567 {
568  // Process cluster.
569 
570  Int_t nCluster = 0;
571  TLorentzVector gamma;
572  Double_t vertex[3] = {0, 0, 0};
573  fESD->GetVertex()->GetXYZ(vertex);
574  const Int_t nCaloClusters = fESD->GetNumberOfCaloClusters();
575  fClusterIndices->Set(nCaloClusters);
576  Float_t ClsPos[3] = {-999,-999,-999};
577 
578  for(Int_t itr=0; itr<nCaloClusters; itr++)
579  {
580  fHEventStat->Fill(1.5);
581  AliESDCaloCluster *cluster = fESD->GetCaloCluster(itr);
582  if(!IsGoodCluster(cluster)) continue;
583  cluster->GetMomentum(gamma, vertex);
584  if (gamma.Pt() < 0.15) continue;
585  fHEventStat->Fill(2.5);
586 
587  cluster->GetPosition(ClsPos);
588  TVector3 VClsPos(ClsPos[0], ClsPos[1], ClsPos[2]);
589 
590  TArrayI *TrackLabels = cluster->GetTracksMatched();
591 
592  if(TrackLabels->GetSize() == 1)
593  {
594  AliESDtrack *esdtrack = fESD->GetTrack(TrackLabels->operator[](0));
595  AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
596  if(newTrack && TMath::Abs(newTrack->Eta())<0.7)
597  {
598  Double_t Esub = newTrack->P();
599  if (Esub > cluster->E()) Esub = cluster->E();
600  fHistEsub1Pch->Fill(newTrack->P(), Esub);
601  fHistEsub1PchRat->Fill(newTrack->P(), Esub/newTrack->P());
602  }
603  }
604 
605  if(TrackLabels->GetSize() == 2)
606  {
607  AliESDtrack *esdtrack1 = fESD->GetTrack(TrackLabels->operator[](0));
608  AliESDtrack *esdtrack2 = fESD->GetTrack(TrackLabels->operator[](1));
609  AliESDtrack *newTrack1 = GetAcceptTrack(esdtrack1);
610  AliESDtrack *newTrack2 = GetAcceptTrack(esdtrack2);
611  if(newTrack1 && newTrack2 && TMath::Abs(newTrack1->Eta())<0.7 && TMath::Abs(newTrack2->Eta())<0.7)
612  {
613  Double_t Esub = newTrack1->P() + newTrack2->P();
614  if (Esub > cluster->E()) Esub = cluster->E();
615  fHistEsub2Pch->Fill(newTrack1->P() + newTrack2->P(), Esub);
616  fHistEsub2PchRat->Fill(newTrack1->P() + newTrack2->P(), Esub/(newTrack1->P() + newTrack2->P()));
617  }
618  else if(newTrack1 && !(newTrack2) && TMath::Abs(newTrack1->Eta())<0.7)
619  {
620  Double_t Esub = newTrack1->P();
621  if (Esub > cluster->E()) Esub = cluster->E();
622  fHistEsub1Pch->Fill(newTrack1->P(), Esub);
623  fHistEsub1PchRat->Fill(newTrack1->P(), Esub/newTrack1->P());
624  }
625  else if (!(newTrack1) && newTrack2 && TMath::Abs(newTrack2->Eta())<0.7)
626  {
627  Double_t Esub = newTrack2->P();
628  if (Esub > cluster->E()) Esub = cluster->E();
629  fHistEsub1Pch->Fill(newTrack2->P(), Esub);
630  fHistEsub1PchRat->Fill(newTrack2->P(), Esub/newTrack2->P());
631  }
632  else {;}
633  }
634 
635  TArrayI *MCLabels = cluster->GetLabelsArray();
636 
637  if(MCLabels->GetSize() == 0) fHEventStat->Fill(3.5);
638  if(MCLabels->GetSize() == 1)
639  {
640  fHEventStat->Fill(4.5);
641  AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->operator[](0));
642  if(IsGoodMcParticle(vParticle1, MCLabels->operator[](0)))
643  {
644  Double_t Entries[4] = {VClsPos.Phi(), VClsPos.Eta(), vParticle1->E(), cluster->E()/vParticle1->E()};
645  fHClsEoverMcE_All->Fill(Entries);
646  if(vParticle1->PdgCode() == 22)
647  {
648  fHClsEoverMcE_Photon->Fill(Entries);
649  }
650  if(TMath::Abs(vParticle1->PdgCode()) == 11)
651  {
652  fHClsEoverMcE_Elec->Fill(Entries);
653  }
654  if(TMath::Abs(vParticle1->PdgCode()) == 211)
655  {
656  fHClsEoverMcE_Pion->Fill(Entries);
657  }
658  }
659  }
660  if(MCLabels->GetSize() == 2)
661  {
662  fHEventStat->Fill(5.5);
663  AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->operator[](0));
664  AliVParticle* vParticle2 = fMC->GetTrack(MCLabels->operator[](1));
665  if(IsGoodMcParticle(vParticle1, MCLabels->operator[](0)) && IsGoodMcParticle(vParticle2, MCLabels->operator[](1)))
666  {
667  fHEventStat->Fill(6.5);
668  if((vParticle1->PdgCode()==22) && (vParticle2->PdgCode()==22)) {;}
669  else if((vParticle1->PdgCode()!=22) && (vParticle2->PdgCode()!=22)) {;}
670  else
671  {
672  fClusterIndices->AddAt(itr,nCluster);
673  nCluster++;
674  }
675  }
676  }
677  if(MCLabels->GetSize() > 2) fHEventStat->Fill(7.5);
678 
679  AliESDCaloCluster *newCluster = new AliESDCaloCluster(*cluster);
680 
681  Double_t subE = 0;
682  TArrayI arrayTrackMatched(fTrackIndices->GetSize());
683  Int_t nGoodMatch = 0;
684 
685  for(Int_t j=0; j<fTrackIndices->GetSize(); j++)
686  {
687  AliESDtrack *trk = fESD->GetTrack(fTrackIndices->At(j));
688  if(itr==trk->GetEMCALcluster())
689  {
690  arrayTrackMatched[nGoodMatch] = j;
691  nGoodMatch ++;
692  subE += trk->P();
693  }
694  }
695 
696  arrayTrackMatched.Set(nGoodMatch);
697  newCluster->AddTracksMatched(arrayTrackMatched);
698 
699  Double_t clsE = newCluster->E();
700  Double_t newE = clsE-subE;
701  if(newE<0) newE = 0;
702  newCluster->SetDispersion(newE);
703  fClusterArray->Add(newCluster);
704  }
705 
706  fClusterIndices->Set(nCluster);
707 }
708 //________________________________________________________________________
709 void AliAnalysisTaskSOH::ProcessMc()
710 {
711  // Process MC.
712  for(Int_t i=0; i<fTrackIndices->GetSize(); i++)
713  {
714  AliESDtrack *esdtrack = fESD->GetTrack(fTrackIndices->At(i));
715  if(!esdtrack)continue;
716  AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
717  if(!newTrack) continue;
718  if(newTrack->Pt()<0.15 || TMath::Abs(newTrack->Eta())>0.9) {delete newTrack; continue;}
719 
720  Int_t index = esdtrack->GetLabel();
721  if(index < 0)
722  {
723  AliVParticle *vParticle1 = (AliVParticle*)fMC->GetTrack(-1*index);
724  if((TMath::Abs(vParticle1->PdgCode())==211) && IsGoodMcParticle(vParticle1, -1*index))
725  {
726  fHDetRecFakePion->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
727  fHDetGenFakePion->Fill(vParticle1->Pt(), vParticle1->Eta(), vParticle1->Phi());
728  }
729  }
730 
731  AliVParticle* vParticle2 = fMC->GetTrack(TMath::Abs(index));
732  if(!IsGoodMcParticle(vParticle2, TMath::Abs(index)) && (TMath::Abs(vParticle2->PdgCode())==211))
733  {
734  fHDetRecSecPion->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
735  fHDetGenSecPion->Fill(vParticle2->Pt(), vParticle2->Eta(), vParticle2->Phi());
736  }
737 
738  if(newTrack) delete newTrack;
739  }
740 
741  //tracking effciency
742  AliHeader* header = (AliHeader*) fMC->Header();
743  if (!header) AliFatal("fInjectedSignals set but no MC header found");
744 
745  AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
746  if (!cocktailHeader)
747  {
748  header->Dump();
749  AliFatal("fInjectedSignals set but no MC cocktail header found");
750  }
751 
752  AliGenEventHeader* eventHeader = 0;
753  eventHeader = dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->First());
754  if (!eventHeader) AliFatal("First event header not found");
755 
756  for(Int_t ipart=0; ipart<fMC->GetNumberOfTracks(); ipart++)
757  {
758  AliVParticle* vParticle = fMC->GetTrack(ipart);
759  Int_t pdgCode = vParticle->PdgCode();
760  AliMCParticle *McParticle = (AliMCParticle*)fMC->GetTrack(ipart);
761  if(!IsGoodMcParticle(vParticle, ipart)) continue;
762 
763  if(TMath::Abs(vParticle->Eta())<0.9)
764  {
765  if(pdgCode==211)
766  {
767  fHParGenPion_p->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
768  if(McParticle->GetMother() < eventHeader->NProduced()) fHParGenPion_rmInj_p->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
769  }
770 
771  else if(pdgCode==-211)
772  {
773  fHParGenPion_m->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
774  if(McParticle->GetMother() < eventHeader->NProduced()) fHParGenPion_rmInj_m->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
775  }
776  }
777 
778  for(Int_t j=0; j<fTrackIndices->GetSize(); j++)
779  {
780  AliESDtrack *esdtrack = fESD->GetTrack(fTrackIndices->At(j));
781  if(!esdtrack)continue;
782  AliESDtrack *newTrack = GetAcceptTrack(esdtrack);
783  if(!newTrack) continue;
784  if(newTrack->Pt()<0.15 || TMath::Abs(newTrack->Eta())>0.9) {delete newTrack; continue;}
785 
786  Int_t cutType = (Int_t)newTrack->GetTRDQuality();
787 
788  if(newTrack->GetLabel()==ipart && TMath::Abs(vParticle->Eta())<0.9)
789  {
790  if(pdgCode==211)
791  {
792  fHDetGenPion_p[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
793  fHDetRecPion_p[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
794  if(McParticle->GetMother() < eventHeader->NProduced())
795  {
796  fHDetGenPion_rmInj_p[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
797  fHDetRecPion_rmInj_p[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
798  }
799  }
800  else if(pdgCode==-211)
801  {
802  fHDetGenPion_m[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
803  fHDetRecPion_m[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
804  if(McParticle->GetMother() < eventHeader->NProduced())
805  {
806  fHDetGenPion_rmInj_m[cutType]->Fill(vParticle->Pt(), vParticle->Eta(), vParticle->Phi());
807  fHDetRecPion_rmInj_m[cutType]->Fill(newTrack->Pt(), newTrack->Eta(), newTrack->Phi());
808  }
809  }
810 
811 
812  //cluster E vs. truth photon energy
813  if(fClusterProcess)
814  {
815  for(Int_t k=0; k<fClusterIndices->GetSize(); k++)
816  {
817  AliESDCaloCluster *cluster = fESD->GetCaloCluster(fClusterIndices->At(k));
818  Double_t clsE = cluster->E();
819  TArrayI *MCLabels = cluster->GetLabelsArray();
820  AliVParticle* vParticle1 = fMC->GetTrack(MCLabels->operator[](0));
821  AliVParticle* vParticle2 = fMC->GetTrack(MCLabels->operator[](1));
822 
823  if(vParticle1->PdgCode()==22 && vParticle2 == vParticle)
824  {
825  fHPhotonEdiff0HC->Fill(vParticle1->E(), (vParticle1->E() - clsE)/vParticle1->E());
826  fHPhotonEVsClsE->Fill(vParticle1->E(), clsE);
827 
828  if((clsE - 0.3*esdtrack->E())<0) fHPhotonEdiff30HC->Fill(vParticle1->E(), 1);
829  else fHPhotonEdiff30HC->Fill(vParticle1->E(), (vParticle1->E() + 0.3*esdtrack->E() - clsE)/vParticle1->E());
830 
831  if((clsE - 0.7*esdtrack->E())<0) fHPhotonEdiff70HC->Fill(vParticle1->E(), 1);
832  else fHPhotonEdiff70HC->Fill(vParticle1->E(), (vParticle1->E() + 0.7*esdtrack->E() - clsE)/vParticle1->E());
833 
834  if((clsE - esdtrack->E())<0) fHPhotonEdiff100HC->Fill(vParticle1->E(), 1);
835  else fHPhotonEdiff100HC->Fill(vParticle1->E(), (vParticle1->E() + esdtrack->E() - clsE)/vParticle1->E());
836  continue;
837  }
838  if(vParticle2->PdgCode()==22 && vParticle1 == vParticle)
839  {
840  fHPhotonEdiff0HC->Fill(vParticle2->E(), (vParticle2->E() - clsE)/vParticle2->E());
841  fHPhotonEVsClsE->Fill(vParticle2->E(), clsE);
842 
843  if((clsE - 0.3*esdtrack->E())<0) fHPhotonEdiff30HC->Fill(vParticle2->E(), 1);
844  else fHPhotonEdiff30HC->Fill(vParticle2->E(), (vParticle2->E() + 0.3*esdtrack->E() - clsE)/vParticle2->E());
845 
846  if((clsE - 0.7*esdtrack->E())<0) fHPhotonEdiff70HC->Fill(vParticle2->E(), 1);
847  else fHPhotonEdiff70HC->Fill(vParticle2->E(), (vParticle2->E() + 0.7*esdtrack->E() - clsE)/vParticle2->E());
848 
849  if((clsE-esdtrack->E())<0) fHPhotonEdiff100HC->Fill(vParticle2->E(), 1);
850  else fHPhotonEdiff100HC->Fill(vParticle2->E(), (vParticle2->E() + esdtrack->E() - clsE)/vParticle2->E());
851  }
852  }
853  }
854  if(newTrack) delete newTrack;
855  break;
856  }
857  if(newTrack) delete newTrack;
858  }
859  }
860 }
861 
862 
863 //________________________________________________________________________
864 void AliAnalysisTaskSOH::ProcessScaleFactor()
865 {
866  // Scale factor.
867 
868  const Double_t phiMax = 180 * TMath::DegToRad();
869  const Double_t phiMin = 80 * TMath::DegToRad();
870  const Double_t TPCArea= 2*TMath::Pi()*1.8;
871  const Double_t EMCArea = (phiMax-phiMin)*1.4;
872 
873  Double_t PtEMC = 0;
874  Double_t PtTPC = 0;
875 
876  for(Int_t j=0; j<fTrackIndices->GetSize(); j++)
877  {
878  AliESDtrack *trk = fESD->GetTrack(fTrackIndices->At(j));
879  Double_t eta = trk->Eta();
880  Double_t phi = trk->Phi();
881  if(TMath::Abs(eta)<0.9) PtTPC += trk->Pt();
882  if(TMath::Abs(eta)<0.7 && phi > phiMin && phi < phiMax ) PtEMC += trk->Pt();
883  }
884 
885  Double_t EtWithHadCorr = 0;
886  Double_t EtWithoutHadCorr = 0;
887  Double_t vertex[3] = {0, 0, 0};
888  fESD->GetVertex()->GetXYZ(vertex);
889  TLorentzVector gamma;
890 
891  for(Int_t i=0; i<fClusterArray->GetEntriesFast(); i++)
892  {
893  AliESDCaloCluster *cluster = (AliESDCaloCluster*)fClusterArray->At(i);
894  cluster->GetMomentum(gamma, vertex);
895  Double_t sinTheta = TMath::Sqrt(1-TMath::Power(gamma.CosTheta(),2));
896  EtWithoutHadCorr += cluster->E() * sinTheta;
897  EtWithHadCorr += cluster->GetDispersion() * sinTheta;
898  }
899 
900  if(PtTPC>0)
901  {
902  fHScaleFactor->Fill((PtEMC+EtWithoutHadCorr)/EMCArea * TPCArea/PtTPC);
903  fHScaleFactor100HC->Fill((PtEMC+EtWithHadCorr)/EMCArea * TPCArea/PtTPC);
904  }
905 }
906 
907 //________________________________________________________________________
908 AliESDtrack *AliAnalysisTaskSOH::GetAcceptTrack(AliESDtrack *esdtrack)
909 {
910  // Get accepted track.
911  AliESDtrack *newTrack = 0x0;
912  if(fEsdTrackCuts->AcceptTrack(esdtrack))
913  {
914  newTrack = new AliESDtrack(*esdtrack);
915  newTrack->SetTRDQuality(0);
916  }
917  else if(fHybridTrackCuts1->AcceptTrack(esdtrack))
918  {
919  if(esdtrack->GetConstrainedParam())
920  {
921  newTrack = new AliESDtrack(*esdtrack);
922  const AliExternalTrackParam* constrainParam = esdtrack->GetConstrainedParam();
923  newTrack->Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
924  newTrack->SetTRDQuality(1);
925  }
926  else
927  return 0x0;
928  }
929  else if(fHybridTrackCuts2->AcceptTrack(esdtrack))
930  {
931  if(esdtrack->GetConstrainedParam())
932  {
933  newTrack = new AliESDtrack(*esdtrack);
934  const AliExternalTrackParam* constrainParam = esdtrack->GetConstrainedParam();
935  newTrack->Set(constrainParam->GetX(),constrainParam->GetAlpha(),constrainParam->GetParameter(),constrainParam->GetCovariance());
936  newTrack->SetTRDQuality(2);
937  }
938  else
939  return 0x0;
940  }
941  else
942  {
943  return 0x0;
944  }
945 
946  return newTrack;
947 }
948 
949 //________________________________________________________________________
950 Bool_t AliAnalysisTaskSOH::IsGoodMcParticle(AliVParticle* vParticle, Int_t ipart)
951 {
952  // Return true if good MC particle.
953 
954  if(!vParticle) return kFALSE;
955  if(!fMC->IsPhysicalPrimary(ipart)) return kFALSE;
956  if (TMath::Abs(vParticle->Eta())>2) return kFALSE;
957  if(vParticle->Pt()<0.15) return kFALSE;
958  return kTRUE;
959 }
960 
961 //________________________________________________________________________
962 Bool_t AliAnalysisTaskSOH::IsGoodCluster(AliESDCaloCluster *cluster)
963 {
964  // Return true if good cluster.
965 
966  if(!cluster) return kFALSE;
967  if (!cluster->IsEMCAL()) return kFALSE;
968  return kTRUE;
969 }
970 
971 //________________________________________________________________________
972 void AliAnalysisTaskSOH::Terminate(Option_t *)
973 {
974  // Terminate analysis.
975 }
976 
977 //________________________________________________________________________
978 Bool_t AliAnalysisTaskSOH::EsdVertexOk() const
979 {
980  // Modified from AliAnalyseLeadingTrackUE::VertexSelection()
981 
982  const AliESDVertex* vtx = fESD->GetPrimaryVertex();
983  if (!vtx) return kFALSE;
984  Int_t nContributors = vtx->GetNContributors();
985  Double_t zVertex = vtx->GetZ();
986  if( nContributors < 1 || TMath::Abs(zVertex) > fZVtxMax ) return kFALSE;
987  return kTRUE;
988 }
AliAnalysisTaskSOH::fClusterIndices
TArrayI * fClusterIndices
selected track index
Definition: AliAnalysisTaskSOH.h:58
AliAnalysisTaskSOH::EsdVertexOk
Bool_t EsdVertexOk() const
Definition: AliAnalysisTaskSOH.cxx:978
AliAnalysisTaskSOH::fHistEsub1Pch
TH2F * fHistEsub1Pch
cluster E vs. truth photon E
Definition: AliAnalysisTaskSOH.h:87
AliAnalysisTaskSOH::Terminate
void Terminate(Option_t *)
Definition: AliAnalysisTaskSOH.cxx:972
AliAnalysisTaskSOH::fHEMCalResponseElec
TH2F * fHEMCalResponseElec
same as above for pions
Definition: AliAnalysisTaskSOH.h:72
AliAnalysisTaskSOH::fHistEsub1PchRat
TH2F * fHistEsub1PchRat
(subtracted E in 100% HC) vs. total track P, clusters with 2 matching tracks
Definition: AliAnalysisTaskSOH.h:89
AliAnalysisTaskSOH::fHEMCalRecdPhidEtaposEta
TH2F * fHEMCalRecdPhidEtaposEta
same as above with negative truth charge matching
Definition: AliAnalysisTaskSOH.h:80
AliAnalysisTaskSOH::fEsdTrackCuts
AliESDtrackCuts * fEsdTrackCuts
Definition: AliAnalysisTaskSOH.h:54
AliAnalysisTaskSOH::fClusterArray
TObjArray * fClusterArray
cluster with two matched MC track index
Definition: AliAnalysisTaskSOH.h:59
AliAnalysisTaskSOH::fHybridTrackCuts1
AliESDtrackCuts * fHybridTrackCuts1
Definition: AliAnalysisTaskSOH.h:55
AliAnalysisTaskSOH::fHPhotonEdiff70HC
TH2F * fHPhotonEdiff70HC
(truth E - calculated E in 100% HC)/truth E vs. truth E with photon
Definition: AliAnalysisTaskSOH.h:83
AliAnalysisTaskSOH::fHDetGenPion_p
TH3F * fHDetGenPion_p[3]
secondary pion tracks pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:103
AliAnalysisTaskSOH::fHDetRecPion_rmInj_m
TH3F * fHDetRecPion_rmInj_m[3]
minus charged generated detector level particle(pion) without injected signal pt, phi...
Definition: AliAnalysisTaskSOH.h:110
AliAnalysisTaskSOH::fHPhotonEdiff30HC
TH2F * fHPhotonEdiff30HC
(truth E - calculated E in 70% HC)/truth E vs. truth E with photon
Definition: AliAnalysisTaskSOH.h:84
AliAnalysisTaskSOH::fHPhotonEdiff0HC
TH2F * fHPhotonEdiff0HC
(truth E - calculated E in 30% HC)/truth E vs. truth E with photon
Definition: AliAnalysisTaskSOH.h:85
AliAnalysisTaskSOH::fHScaleFactor
TH1F * fHScaleFactor
statistics histo
Definition: AliAnalysisTaskSOH.h:68
AliAnalysisTaskSOH::fHEMCalResponseProton
TH2F * fHEMCalResponseProton
same as above for electrons
Definition: AliAnalysisTaskSOH.h:73
AliAnalysisTaskSOH::IsGoodMcParticle
Bool_t IsGoodMcParticle(AliVParticle *vParticle, Int_t ipart)
Definition: AliAnalysisTaskSOH.cxx:950
AliAnalysisTaskSOH::fHPhotonEVsClsE
TH2F * fHPhotonEVsClsE
(truth E - cluster E)/truth E vs. truth E with photon
Definition: AliAnalysisTaskSOH.h:86
AliAnalysisTaskSOH::fHPhotonEdiff100HC
TH2F * fHPhotonEdiff100HC
same as above for negative eta
Definition: AliAnalysisTaskSOH.h:82
AliAnalysisTaskSOH::fHDetRecPion_m
TH3F * fHDetRecPion_m[3]
minus pion mc detector level pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:106
AliAnalysisTaskSOH::fHDetRecSecPion
TH3F * fHDetRecSecPion
secondary pion tracks pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:102
AliAnalysisTaskSOH::fHDetRecFakePion
TH3F * fHDetRecFakePion
fake pion tracks pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:100
AliAnalysisTaskSOH::fMcProcess
Bool_t fMcProcess
selected cluster array
Definition: AliAnalysisTaskSOH.h:60
AliAnalysisTaskSOH::fZVtxMax
Double_t fZVtxMax
mv event
Definition: AliAnalysisTaskSOH.h:53
ClassImp
ClassImp(AliAnalysisTaskSOH) AliAnalysisTaskSOH
Definition: AliAnalysisTaskSOH.cxx:31
AliAnalysisTaskSOH::fHClsEoverMcE_Pion
THnSparse * fHClsEoverMcE_Pion
above for electron
Definition: AliAnalysisTaskSOH.h:94
AliAnalysisTaskSOH::fHClsEoverMcE_Elec
THnSparse * fHClsEoverMcE_Elec
above for photon
Definition: AliAnalysisTaskSOH.h:93
AliAnalysisTaskSOH::fHEMCalRecdPhidEtanegEta
TH2F * fHEMCalRecdPhidEtanegEta
same as above for positive eta
Definition: AliAnalysisTaskSOH.h:81
AliAnalysisTaskSOH::GetAcceptTrack
AliESDtrack * GetAcceptTrack(AliESDtrack *esdtrack)
Definition: AliAnalysisTaskSOH.cxx:908
AliAnalysisTaskSOH.h
AliAnalysisTaskSOH::fHEMCalResponsePion
TH2F * fHEMCalResponsePion
(cluster energy over reconstructed track p) vs. track pt
Definition: AliAnalysisTaskSOH.h:71
AliAnalysisTaskSOH::fHEMCalRecdPhidEtaM_Truth
TH2F * fHEMCalRecdPhidEtaM_Truth
same as above with positive truth charge matching
Definition: AliAnalysisTaskSOH.h:79
AliAnalysisTaskSOH::fHParGenPion_rmInj_p
TH3F * fHParGenPion_rmInj_p
minus pion mc truth pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:97
AliAnalysisTaskSOH::fHDetGenSecPion
TH3F * fHDetGenSecPion
fake pion tracks pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:101
AliAnalysisTaskSOH::fHEMCalRecdPhidEta_Truth
TH2F * fHEMCalRecdPhidEta_Truth
same as above for negative charge tracks
Definition: AliAnalysisTaskSOH.h:77
AliAnalysisTaskSOH::fHParGenPion_rmInj_m
TH3F * fHParGenPion_rmInj_m
plus charged mc truth(pion) without injected signal pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:98
AliAnalysisTaskSOH::fHEOverPVsPt
TH2F * fHEOverPVsPt
scale factor with 100% HC spectrum
Definition: AliAnalysisTaskSOH.h:70
AliAnalysisTaskSOH::fHEventStat
TH1F * fHEventStat
output list
Definition: AliAnalysisTaskSOH.h:67
AliAnalysisTaskSOH::fHEMCalRecdPhidEtaP
TH2F * fHEMCalRecdPhidEtaP
(EMCal cluster phi - track phi) vs. (EMCal cluster eta - track eta)
Definition: AliAnalysisTaskSOH.h:75
AliAnalysisTaskSOH::fHParGenPion_m
TH3F * fHParGenPion_m
plus pion mc truth pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:96
AliAnalysisTaskSOH::fHDetRecPion_rmInj_p
TH3F * fHDetRecPion_rmInj_p[3]
plus charged generated detector level particle(pion) without injected signal pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:108
AliAnalysisTaskSOH::UserExec
void UserExec(Option_t *option)
Definition: AliAnalysisTaskSOH.cxx:424
AliAnalysisTaskSOH::fHDetGenPion_rmInj_m
TH3F * fHDetGenPion_rmInj_m[3]
plus charged reconstructed detector level pion+ track without injected signal pt, phi...
Definition: AliAnalysisTaskSOH.h:109
AliAnalysisTaskSOH::fHDetGenPion_m
TH3F * fHDetGenPion_m[3]
plus pion reconstructed detector level pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:105
AliAnalysisTaskSOH::fHDetGenPion_rmInj_p
TH3F * fHDetGenPion_rmInj_p[3]
minus pion reconstructed detector level pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:107
AliAnalysisTaskSOH::fHEMCalRecdPhidEtaM
TH2F * fHEMCalRecdPhidEtaM
same as above for positive charge tracks
Definition: AliAnalysisTaskSOH.h:76
AliAnalysisTaskSOH::IsGoodCluster
Bool_t IsGoodCluster(AliESDCaloCluster *cluster)
Definition: AliAnalysisTaskSOH.cxx:962
AliAnalysisTaskSOH::fHDetGenFakePion
TH3F * fHDetGenFakePion
minus charged mc truth(pion) without injected signal pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:99
AliAnalysisTaskSOH::fHScaleFactor100HC
TH1F * fHScaleFactor100HC
scale factor spectrum
Definition: AliAnalysisTaskSOH.h:69
AliAnalysisTaskSOH::fHEMCalRecdPhidEta
TH2F * fHEMCalRecdPhidEta
same as above for protons
Definition: AliAnalysisTaskSOH.h:74
AliAnalysisTaskSOH::fHDetRecPion_p
TH3F * fHDetRecPion_p[3]
plus pion mc detector level pt, phi, eta spectrum
Definition: AliAnalysisTaskSOH.h:104
AliAnalysisTaskSOH::fHEMCalRecdPhidEtaP_Truth
TH2F * fHEMCalRecdPhidEtaP_Truth
same as above with mc truth matching
Definition: AliAnalysisTaskSOH.h:78
AliAnalysisTaskSOH::fHybridTrackCuts2
AliESDtrackCuts * fHybridTrackCuts2
Definition: AliAnalysisTaskSOH.h:56
AliAnalysisTaskSOH::fHClsEoverMcE_All
THnSparse * fHClsEoverMcE_All
(subtracted E in 100% HC)/total track P vs. total track P, clusters with 2 matching tracks ...
Definition: AliAnalysisTaskSOH.h:91
AliAnalysisTaskSOH::fHClsEoverMcE_Photon
THnSparse * fHClsEoverMcE_Photon
cluster E/MC particle E, cluster with only one matching particle
Definition: AliAnalysisTaskSOH.h:92
AliAnalysisTaskSOH::fHParGenPion_p
TH3F * fHParGenPion_p
above for pion
Definition: AliAnalysisTaskSOH.h:95
AliAnalysisTaskSOH::fHistEsub2PchRat
TH2F * fHistEsub2PchRat
(subtracted E in 100% HC)/total track P vs. total track P, clusters with 1 matching track ...
Definition: AliAnalysisTaskSOH.h:90
AliAnalysisTaskSOH::fHistEsub2Pch
TH2F * fHistEsub2Pch
(subtracted E in 100% HC) vs. total track P, clusters with 1 matching track
Definition: AliAnalysisTaskSOH.h:88
AliAnalysisTaskSOH::fMC
AliMCEvent * fMC
esd event
Definition: AliAnalysisTaskSOH.h:52