]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/EMCALTasks/AliAnalysisTaskEMCALIsoPhoton.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
changing the output list setter argument to be kTRUE
[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskEMCALIsoPhoton.cxx
CommitLineData
bd092f0f 1// $Id$
30159e6f 2
bd092f0f 3#include "AliAnalysisTaskEMCALIsoPhoton.h"
30159e6f 4
bd092f0f 5#include <TCanvas.h>
6#include <TChain.h>
7#include <TFile.h>
8#include <TH1F.h>
9#include <TH2F.h>
10#include <TLorentzVector.h>
11#include <TTree.h>
12
13#include "AliAnalysisManager.h"
14#include "AliAnalysisTask.h"
15#include "AliEMCALGeometry.h"
16#include "AliEMCALRecoUtils.h"
17#include "AliESDCaloCells.h"
18#include "AliESDCaloCluster.h"
30159e6f 19#include "AliESDEvent.h"
20#include "AliESDHeader.h"
30159e6f 21#include "AliESDInputHandler.h"
bd092f0f 22#include "AliESDUtils.h"
30159e6f 23#include "AliESDpid.h"
30159e6f 24#include "AliESDtrack.h"
25#include "AliESDtrackCuts.h"
26#include "AliESDv0.h"
bd092f0f 27#include "AliKFParticle.h"
28#include "AliMCEvent.h"
29#include "AliMCEventHandler.h"
30#include "AliStack.h"
30159e6f 31#include "AliV0vertexer.h"
30159e6f 32#include "AliVCluster.h"
33
30159e6f 34ClassImp(AliAnalysisTaskEMCALIsoPhoton)
35
36//________________________________________________________________________
bd092f0f 37AliAnalysisTaskEMCALIsoPhoton::AliAnalysisTaskEMCALIsoPhoton() :
38 AliAnalysisTaskSE(),
39 fCaloClusters(0),
40 fSelPrimTracks(0),
41 fEMCalCells(0),
42 fPrTrCuts(0),
43 fGeom(0x0),
9148fa89 44 fGeoName("EMCAL_COMPLETEV1"),
45 fPeriod("LHC11c"),
46 fTrigBit("kEMC7"),
47 fIsTrain(0),
48 fExoticCut(0.97),
49 fIsoConeR(0.4),
bd092f0f 50 fESD(0),
51 fOutputList(0),
52 fEvtSel(0),
53 fPVtxZ(0),
54 fCellAbsIdVsAmpl(0),
55 fNClusHighClusE(0),
56 fM02Et(0),
57 fM02EtTM(0),
58 fM02EtCeIso1(0),
59 fM02EtCeIso2(0),
60 fM02EtCeIso5(0),
61 fM02EtTrIso1(0),
62 fM02EtTrIso2(0),
63 fM02EtTrIso5(0),
64 fM02EtAllIso1(0),
65 fM02EtAllIso2(0),
66 fM02EtAllIso5(0),
67 fCeIsoVsEtPho(0),
68 fTrIsoVsEtPho(0),
69 fAllIsoVsEtPho(0),
c6ee6f73 70 fCeIsoVsEtPi0(0),
71 fTrIsoVsEtPi0(0),
72 fAllIsoVsEtPi0(0),
bd092f0f 73 fM02EtCeIso1TM(0),
74 fM02EtCeIso2TM(0),
75 fM02EtCeIso5TM(0),
76 fM02EtTrIso1TM(0),
77 fM02EtTrIso2TM(0),
78 fM02EtTrIso5TM(0),
79 fM02EtAllIso1TM(0),
80 fM02EtAllIso2TM(0),
81 fM02EtAllIso5TM(0),
82 fCeIsoVsEtPhoTM(0),
83 fTrIsoVsEtPhoTM(0),
c6ee6f73 84 fAllIsoVsEtPhoTM(0),
85 fCeIsoVsEtPi0TM(0),
86 fTrIsoVsEtPi0TM(0),
87 fAllIsoVsEtPi0TM(0)
bd092f0f 88{
89 // Default constructor.
90}
91
92//________________________________________________________________________
93AliAnalysisTaskEMCALIsoPhoton::AliAnalysisTaskEMCALIsoPhoton(const char *name) :
94 AliAnalysisTaskSE(name),
30159e6f 95 fCaloClusters(0),
96 fSelPrimTracks(0),
97 fEMCalCells(0),
98 fPrTrCuts(0),
99 fGeom(0x0),
100 fGeoName("EMCAL_COMPLETEV1"),
101 fPeriod("LHC11c"),
751194e8 102 fTrigBit("kEMC7"),
30159e6f 103 fIsTrain(0),
104 fExoticCut(0.97),
105 fIsoConeR(0.4),
30159e6f 106 fESD(0),
30159e6f 107 fOutputList(0),
30159e6f 108 fEvtSel(0),
bd092f0f 109 fPVtxZ(0),
110 fCellAbsIdVsAmpl(0),
111 fNClusHighClusE(0),
112 fM02Et(0),
113 fM02EtTM(0),
114 fM02EtCeIso1(0),
115 fM02EtCeIso2(0),
116 fM02EtCeIso5(0),
117 fM02EtTrIso1(0),
118 fM02EtTrIso2(0),
119 fM02EtTrIso5(0),
120 fM02EtAllIso1(0),
121 fM02EtAllIso2(0),
122 fM02EtAllIso5(0),
123 fCeIsoVsEtPho(0),
124 fTrIsoVsEtPho(0),
125 fAllIsoVsEtPho(0),
c6ee6f73 126 fCeIsoVsEtPi0(0),
127 fTrIsoVsEtPi0(0),
128 fAllIsoVsEtPi0(0),
bd092f0f 129 fM02EtCeIso1TM(0),
130 fM02EtCeIso2TM(0),
131 fM02EtCeIso5TM(0),
132 fM02EtTrIso1TM(0),
133 fM02EtTrIso2TM(0),
134 fM02EtTrIso5TM(0),
135 fM02EtAllIso1TM(0),
136 fM02EtAllIso2TM(0),
137 fM02EtAllIso5TM(0),
138 fCeIsoVsEtPhoTM(0),
139 fTrIsoVsEtPhoTM(0),
c6ee6f73 140 fAllIsoVsEtPhoTM(0),
141 fCeIsoVsEtPi0TM(0),
142 fTrIsoVsEtPi0TM(0),
143 fAllIsoVsEtPi0TM(0)
30159e6f 144{
145 // Constructor
146
147 // Define input and output slots here
148 // Input slot #0 works with a TChain
149 DefineInput(0, TChain::Class());
150 // Output slot #0 id reserved by the base class for AOD
151 // Output slot #1 writes into a TH1 container
152 DefineOutput(1, TList::Class());
153}
bd092f0f 154
30159e6f 155//________________________________________________________________________
156void AliAnalysisTaskEMCALIsoPhoton::UserCreateOutputObjects()
157{
bd092f0f 158 // Create histograms, called once.
30159e6f 159
160 fCaloClusters = new TRefArray();
161 fSelPrimTracks = new TObjArray();
162
163
164 fOutputList = new TList();
78c932b4 165 fOutputList->SetOwner(kTRUE);// Container cleans up all histos (avoids leaks in merging)
30159e6f 166
167 fGeom = AliEMCALGeometry::GetInstance(fGeoName);
168
169 fEvtSel = new TH1F("hEvtSel","Event selection counter (0=all trg, 1=pvz cut) ;evt cut ;dN/dcut}",2,0,2);
170 fOutputList->Add(fEvtSel);
171
172
173 fPVtxZ = new TH1F("hPVtxZ","primary vertex Z before cut;prim-vz(cm) ;",200,-100,100);
174 fOutputList->Add(fPVtxZ);
175
176 fCellAbsIdVsAmpl = new TH2F("hCellAbsIdVsAmpl","cell abs id vs cell amplitude (energy);E (GeV);ID",200,0,100,24*48*10,-0.5,24*48*10-0.5);
177 fOutputList->Add(fPVtxZ);
178
179 fNClusHighClusE = new TH2F("hNClusHighClusE","total number of clusters vs. highest clus energy in the event;E (GeV);NClus",200,0,100,301,-0.5,300.5);
180 fOutputList->Add(fNClusHighClusE);
181
182 fM02Et = new TH2F("fM02Et","M02 vs Et for all clusters;E_{T} (GeV);M02",1000,0,100,400,0,4);
183 fM02Et->Sumw2();
184 fOutputList->Add(fM02Et);
185
186 fM02EtTM = new TH2F("fM02EtTM","M02 vs Et for all track-matched clusters;E_{T} (GeV);M02",1000,0,100,400,0,4);
187 fM02EtTM->Sumw2();
188 fOutputList->Add(fM02EtTM);
189
190 fM02EtCeIso1 = new TH2F("fM02EtCeIso1","M02 vs Et for all clusters (ISO_{EMC}<1GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
191 fM02EtCeIso1->Sumw2();
192 fOutputList->Add(fM02EtCeIso1);
193
194 fM02EtCeIso2 = new TH2F("fM02EtCeIso2","M02 vs Et for all clusters (ISO_{EMC}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
195 fM02EtCeIso2->Sumw2();
196 fOutputList->Add(fM02EtCeIso2);
197
198 fM02EtCeIso5 = new TH2F("fM02EtCeIso5","M02 vs Et for all clusters (ISO_{EMC}<5GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
199 fM02EtCeIso5->Sumw2();
200 fOutputList->Add(fM02EtCeIso5);
201
202 fM02EtTrIso1 = new TH2F("fM02EtTrIso1","M02 vs Et for all clusters (ISO_{TRK}<1GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
203 fM02EtTrIso1->Sumw2();
204 fOutputList->Add(fM02EtTrIso1);
205
206 fM02EtTrIso2 = new TH2F("fM02EtTrIso2","M02 vs Et for all clusters (ISO_{TRK}<1GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
207 fM02EtTrIso2->Sumw2();
208 fOutputList->Add(fM02EtTrIso2);
209
210 fM02EtTrIso5 = new TH2F("fM02EtTrIso5","M02 vs Et for all clusters (ISO_{TRK}<5GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
211 fM02EtTrIso5->Sumw2();
212 fOutputList->Add(fM02EtTrIso5);
213
214 fM02EtAllIso1 = new TH2F("fM02EtAllIso1","M02 vs Et for all clusters (ISO_{EMC+TRK}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
215 fM02EtAllIso1->Sumw2();
216 fOutputList->Add(fM02EtAllIso1);
217
218 fM02EtAllIso2 = new TH2F("fM02EtAllIso2","M02 vs Et for all clusters (ISO_{EMC+TRK}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
219 fM02EtAllIso2->Sumw2();
220 fOutputList->Add(fM02EtAllIso2);
221
222 fM02EtAllIso5 = new TH2F("fM02EtAllIso5","M02 vs Et for all clusters (ISO_{EMC+TRK}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
223 fM02EtAllIso5->Sumw2();
224 fOutputList->Add(fM02EtAllIso5);
225
226 fCeIsoVsEtPho = new TH2F("fCeIsoVsEtPho","ISO_{EMC} vs. E_{T}^{clus} (0.1<#lambda_{0}^{2}<0.3);E_{T} (GeV);ISO_{EMC}",1000,0,100,1000,-10,190);
227 fCeIsoVsEtPho->Sumw2();
228 fOutputList->Add(fCeIsoVsEtPho);
229
230 fTrIsoVsEtPho = new TH2F("fTrIsoVsEtPho","ISO_{TRK} vs. E_{T}^{clus} (0.1<#lambda_{0}^{2}<0.3);E_{T} (GeV);ISO_{TRK}",1000,0,100,1000,-10,190);
231 fTrIsoVsEtPho->Sumw2();
232 fOutputList->Add(fTrIsoVsEtPho);
233
234 fAllIsoVsEtPho = new TH2F("fAllIsoVsEtPho","ISO_{EMC+TRK} vs. E_{T}^{clus} (0.1<#lambda_{0}^{2}<0.3);E_{T} (GeV);ISO_{EMC+TRK}",1000,0,100,1000,-10,190);
235 fAllIsoVsEtPho->Sumw2();
236 fOutputList->Add(fAllIsoVsEtPho);
237
c6ee6f73 238 fCeIsoVsEtPi0 = new TH2F("fCeIsoVsEtPi0","ISO_{EMC} vs. E_{T}^{clus} (#pi^{0} selection for BG);E_{T} (GeV);ISO_{EMC}",1000,0,100,1000,-10,190);
239 fCeIsoVsEtPi0->Sumw2();
240 fOutputList->Add(fCeIsoVsEtPi0);
241
242 fTrIsoVsEtPi0 = new TH2F("fTrIsoVsEtPi0","ISO_{TRK} vs. E_{T}^{clus} (#pi^{0} selection for BG);E_{T} (GeV);ISO_{TRK}",1000,0,100,1000,-10,190);
243 fTrIsoVsEtPi0->Sumw2();
244 fOutputList->Add(fTrIsoVsEtPi0);
245
246 fAllIsoVsEtPi0 = new TH2F("fAllIsoVsEtPi0","ISO_{EMC+TRK} vs. E_{T}^{clus} (#pi^{0} selection for BG);E_{T} (GeV);ISO_{EMC+TRK}",1000,0,100,1000,-10,190);
247 fAllIsoVsEtPi0->Sumw2();
248 fOutputList->Add(fAllIsoVsEtPi0);
249
30159e6f 250 fM02EtCeIso1TM = new TH2F("fM02EtCeIso1TM","M02 vs Et for all clusters (ISO_{EMC}<1GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
251 fM02EtCeIso1TM->Sumw2();
252 fOutputList->Add(fM02EtCeIso1TM);
253
254 fM02EtCeIso2TM = new TH2F("fM02EtCeIso2TM","M02 vs Et for all clusters (ISO_{EMC}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
255 fM02EtCeIso2TM->Sumw2();
256 fOutputList->Add(fM02EtCeIso2TM);
257
258 fM02EtCeIso5TM = new TH2F("fM02EtCeIso5TM","M02 vs Et for all clusters (ISO_{EMC}<5GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
259 fM02EtCeIso5TM->Sumw2();
260 fOutputList->Add(fM02EtCeIso5TM);
261
262 fM02EtTrIso1TM = new TH2F("fM02EtTrIso1TM","M02 vs Et for all clusters (ISO_{TRK}<1GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
263 fM02EtTrIso1TM->Sumw2();
264 fOutputList->Add(fM02EtTrIso1TM);
265
266 fM02EtTrIso2TM = new TH2F("fM02EtTrIso2TM","M02 vs Et for all clusters (ISO_{TRK}<1GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
267 fM02EtTrIso2TM->Sumw2();
268 fOutputList->Add(fM02EtTrIso2TM);
269
270 fM02EtTrIso5TM = new TH2F("fM02EtTrIso5TM","M02 vs Et for all clusters (ISO_{TRK}<5GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
271 fM02EtTrIso5TM->Sumw2();
272 fOutputList->Add(fM02EtTrIso5TM);
273
274 fM02EtAllIso1TM = new TH2F("fM02EtAllIso1TM","M02 vs Et for all clusters (ISO_{EMC+TRK}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
275 fM02EtAllIso1TM->Sumw2();
276 fOutputList->Add(fM02EtAllIso1TM);
277
278 fM02EtAllIso2TM = new TH2F("fM02EtAllIso2TM","M02 vs Et for all clusters (ISO_{EMC+TRK}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
279 fM02EtAllIso2TM->Sumw2();
280 fOutputList->Add(fM02EtAllIso2TM);
281
282 fM02EtAllIso5TM = new TH2F("fM02EtAllIso5TM","M02 vs Et for all clusters (ISO_{EMC+TRK}<2GeV);E_{T} (GeV);M02",1000,0,100,400,0,4);
283 fM02EtAllIso5TM->Sumw2();
284 fOutputList->Add(fM02EtAllIso5TM);
285
286 fCeIsoVsEtPhoTM = new TH2F("fCeIsoVsEtPhoTM","ISO_{EMC} vs. E_{T}^{clus} (0.1<#lambda_{0}^{2}<0.3);E_{T} (GeV);ISO_{EMC}",1000,0,100,1000,-10,190);
287 fCeIsoVsEtPhoTM->Sumw2();
288 fOutputList->Add(fCeIsoVsEtPhoTM);
289
290 fTrIsoVsEtPhoTM = new TH2F("fTrIsoVsEtPhoTM","ISO_{TRK} vs. E_{T}^{clus} (0.1<#lambda_{0}^{2}<0.3);E_{T} (GeV);ISO_{TRK}",1000,0,100,1000,-10,190);
291 fTrIsoVsEtPhoTM->Sumw2();
292 fOutputList->Add(fTrIsoVsEtPhoTM);
293
294 fAllIsoVsEtPhoTM = new TH2F("fAllIsoVsEtPhoTM","ISO_{EMC+TRK} vs. E_{T}^{clus} (0.1<#lambda_{0}^{2}<0.3);E_{T} (GeV);ISO_{EMC+TRK}",1000,0,100,1000,-10,190);
295 fAllIsoVsEtPhoTM->Sumw2();
296 fOutputList->Add(fAllIsoVsEtPhoTM);
297
c6ee6f73 298 fCeIsoVsEtPi0TM = new TH2F("fCeIsoVsEtPi0TM","ISO_{EMC} vs. E_{T}^{clus} (#pi^{0} selection for BG);E_{T} (GeV);ISO_{EMC}",1000,0,100,1000,-10,190);
299 fCeIsoVsEtPi0TM->Sumw2();
300 fOutputList->Add(fCeIsoVsEtPi0TM);
301
302 fTrIsoVsEtPi0TM = new TH2F("fTrIsoVsEtPi0TM","ISO_{TRK} vs. E_{T}^{clus} (#pi^{0} selection for BG);E_{T} (GeV);ISO_{TRK}",1000,0,100,1000,-10,190);
303 fTrIsoVsEtPi0TM->Sumw2();
304 fOutputList->Add(fTrIsoVsEtPi0TM);
305
306 fAllIsoVsEtPi0TM = new TH2F("fAllIsoVsEtPi0TM","ISO_{EMC+TRK} vs. E_{T}^{clus} (#pi^{0} selection for BG);E_{T} (GeV);ISO_{EMC+TRK}",1000,0,100,1000,-10,190);
307 fAllIsoVsEtPi0TM->Sumw2();
308 fOutputList->Add(fAllIsoVsEtPi0TM);
309
30159e6f 310 PostData(1, fOutputList);
311}
312
313//________________________________________________________________________
314void AliAnalysisTaskEMCALIsoPhoton::UserExec(Option_t *)
315{
bd092f0f 316 // Main loop, called for each event.
317
318 // event trigger selection
30159e6f 319 Bool_t isSelected = 0;
751194e8 320 if(fPeriod.Contains("11a")){
321 if(fTrigBit.Contains("kEMC"))
322 isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kEMC1);
323 else
324 isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
325 }
326 else{
327 if(fTrigBit.Contains("kEMC"))
328 isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kEMC7);
329 else
330 isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kINT7);
331 }
30159e6f 332 if(!isSelected )
333 return;
30159e6f 334
30159e6f 335 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
336 if (!fESD) {
337 printf("ERROR: fESD not available\n");
338 return;
339 }
340
341 fEvtSel->Fill(0);
342
343 AliESDVertex *pv = (AliESDVertex*)fESD->GetPrimaryVertex();
344 if(!pv)
345 return;
346 fPVtxZ->Fill(pv->GetZ());
347 if(TMath::Abs(pv->GetZ())>15)
348 return;
349
350 fEvtSel->Fill(1);
bd092f0f 351
30159e6f 352 // Track loop to fill a pT spectrum
353 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
354 AliESDtrack* track = (AliESDtrack*)fESD->GetTrack(iTracks);
355 if (!track)
356 continue;
357 if (fPrTrCuts && fPrTrCuts->IsSelected(track)){
358 fSelPrimTracks->Add(track);
359 //printf("pt,eta,phi:%1.1f,%1.1f,%1.1f \n",track->Pt(),track->Eta(), track->Phi());
360 }
bd092f0f 361 }
362
30159e6f 363 if(!fIsTrain){
364 for(Int_t mod=0; mod < (fGeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
365 if(fGeoName=="EMCAL_FIRSTYEARV1" && mod>3)
366 break;
367 fGeom->SetMisalMatrix(fESD->GetEMCALMatrix(mod), mod);
368 }
369 }
370 fESD->GetEMCALClusters(fCaloClusters);
371 fEMCalCells = fESD->GetEMCALCells();
372
373
374 FillClusHists();
375
376 fCaloClusters->Clear();
377 fSelPrimTracks->Clear();
378 PostData(1, fOutputList);
379}
bd092f0f 380
30159e6f 381//________________________________________________________________________
382void AliAnalysisTaskEMCALIsoPhoton::FillClusHists()
383{
bd092f0f 384 // Fill cluster histograms.
385
30159e6f 386 if(!fCaloClusters)
387 return;
388 const Int_t nclus = fCaloClusters->GetEntries();
389 if(nclus==0)
390 return;
e681d9ce 391 Double_t maxE = 0;
30159e6f 392 for(Int_t ic=0;ic<nclus;ic++){
393 maxE=0;
394 AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(fCaloClusters->At(ic));
395 if(!c)
396 continue;
397 if(!c->IsEMCAL())
398 continue;
399 Short_t id;
400 Double_t Emax = GetMaxCellEnergy( c, id);
401 Double_t Ecross = GetCrossEnergy( c, id);
402 if((1-Ecross/Emax)>fExoticCut)
403 continue;
404 Float_t clsPos[3] = {0,0,0};
405 c->GetPosition(clsPos);
406 TVector3 clsVec(clsPos);
3dc2825e 407 Double_t Et = c->E()*TMath::Sin(clsVec.Theta());
30159e6f 408 Float_t ceiso, cephiband, cecore;
409 Float_t triso, trphiband, trcore;
410 Float_t alliso, allphiband, allcore;
411 Float_t phibandArea = (1.4 - 2*fIsoConeR)*2*fIsoConeR;
412 Float_t netConeArea = TMath::Pi()*(fIsoConeR*fIsoConeR - 0.04*0.04);
413 GetCeIso(clsVec, ceiso, cephiband, cecore);
414 GetTrIso(clsVec, triso, trphiband, trcore);
415 alliso = ceiso + triso;
416 allphiband = cephiband + trphiband;
417 allcore = cecore + trcore;
418 Float_t ceisoue = cephiband/phibandArea*netConeArea;
419 Float_t trisoue = trphiband/phibandArea*netConeArea;
420 Float_t allisoue = allphiband/phibandArea*netConeArea;
30159e6f 421 Double_t dR = TMath::Sqrt(pow(c->GetTrackDx(),2)+pow(c->GetTrackDz(),2));
c6ee6f73 422 Double_t M02u;
423 if(Et<12)
424 M02u = 0.02486*Et*Et - 0.7289*Et + 6.266;
425 else
426 M02u = 14.32/Et - 0.09863;
427 if(M02u<0.65)
428 M02u = 0.65;
429 Double_t M02l = 12.88/Et - 0.3194;
430 if(M02l<0.4)
431 M02l = 0.4;
432 if(dR>0.028){
433 fM02Et->Fill(Et, c->GetM02());
434 if(ceiso-cecore<1)
435 fM02EtCeIso1->Fill(Et, c->GetM02());
436 if(ceiso-cecore<2)
437 fM02EtCeIso2->Fill(Et, c->GetM02());
438 if(ceiso-cecore<5)
439 fM02EtCeIso5->Fill(Et, c->GetM02());
440 if(triso-trcore<1)
441 fM02EtTrIso1->Fill(Et, c->GetM02());
442 if(triso-trcore<2)
443 fM02EtTrIso2->Fill(Et, c->GetM02());
444 if(triso-trcore<5)
445 fM02EtTrIso5->Fill(Et, c->GetM02());
446 if(alliso-allcore<1)
447 fM02EtAllIso1->Fill(Et, c->GetM02());
448 if(alliso-allcore<2)
449 fM02EtAllIso2->Fill(Et, c->GetM02());
450 if(alliso-allcore<5)
451 fM02EtAllIso5->Fill(Et, c->GetM02());
452 if(c->GetM02()>0.1 && c->GetM02()<0.3){
453 fCeIsoVsEtPho->Fill(Et, ceiso - cecore - ceisoue);
454 fTrIsoVsEtPho->Fill(Et, triso - trcore - trisoue);
455 fAllIsoVsEtPho->Fill(Et, alliso - allcore - allisoue);
456 }
457 if(c->GetM02()>M02l && c->GetM02()<M02u){
458 fCeIsoVsEtPi0->Fill(Et, ceiso - cecore - ceisoue);
459 fTrIsoVsEtPi0->Fill(Et, triso - trcore - trisoue);
460 fAllIsoVsEtPi0->Fill(Et, alliso - allcore - allisoue);
461 }
462 }
463 else
464 {
30159e6f 465 fM02EtTM->Fill(Et, c->GetM02());
70484e97 466 if(ceiso-cecore<1)
30159e6f 467 fM02EtCeIso1TM->Fill(Et, c->GetM02());
70484e97 468 if(ceiso-cecore<2)
30159e6f 469 fM02EtCeIso2TM->Fill(Et, c->GetM02());
70484e97 470 if(ceiso-cecore<5)
30159e6f 471 fM02EtCeIso5TM->Fill(Et, c->GetM02());
70484e97 472 if(triso-trcore<1)
30159e6f 473 fM02EtTrIso1TM->Fill(Et, c->GetM02());
70484e97 474 if(triso-trcore<2)
30159e6f 475 fM02EtTrIso2TM->Fill(Et, c->GetM02());
70484e97 476 if(triso-trcore<5)
30159e6f 477 fM02EtTrIso5TM->Fill(Et, c->GetM02());
70484e97 478 if(alliso-allcore<1)
30159e6f 479 fM02EtAllIso1TM->Fill(Et, c->GetM02());
70484e97 480 if(alliso-allcore<2)
30159e6f 481 fM02EtAllIso2TM->Fill(Et, c->GetM02());
70484e97 482 if(alliso-allcore<5)
30159e6f 483 fM02EtAllIso5TM->Fill(Et, c->GetM02());
484 if(c->GetM02()>0.1 && c->GetM02()<0.3){
c6ee6f73 485 fCeIsoVsEtPhoTM->Fill(Et, ceiso - cecore - ceisoue);
486 fTrIsoVsEtPhoTM->Fill(Et, triso - trcore - trisoue);
487 fAllIsoVsEtPhoTM->Fill(Et, alliso - allcore - allisoue);
488 }
489 if(c->GetM02()>M02l && c->GetM02()<M02u){
490 fCeIsoVsEtPi0TM->Fill(Et, ceiso - cecore - ceisoue);
491 fTrIsoVsEtPi0TM->Fill(Et, triso - trcore - trisoue);
492 fAllIsoVsEtPi0TM->Fill(Et, alliso - allcore - allisoue);
30159e6f 493 }
494 }
495 if(c->E()>maxE)
496 maxE = c->E();
497 }
498 fNClusHighClusE->Fill(maxE,nclus);
499}
bd092f0f 500
30159e6f 501//________________________________________________________________________
502void AliAnalysisTaskEMCALIsoPhoton::GetCeIso(TVector3 vec, Float_t &iso, Float_t &phiband, Float_t &core)
503{
bd092f0f 504 // Get cell isolation.
505
30159e6f 506 if(!fEMCalCells)
507 return;
508 const Int_t ncells = fEMCalCells->GetNumberOfCells();
509 Float_t totiso=0;
510 Float_t totband=0;
511 Float_t totcore=0;
512 Float_t etacl = vec.Eta();
513 Float_t phicl = vec.Phi();
514 Float_t thetacl = vec.Theta();
515 if(phicl<0)
516 phicl+=TMath::TwoPi();
517 Int_t absid = -1;
518 Float_t eta=-1, phi=-1;
519 for(int icell=0;icell<ncells;icell++){
520 absid = TMath::Abs(fEMCalCells->GetCellNumber(icell));
521 if(!fGeom)
522 return;
523 fGeom->EtaPhiFromIndex(absid,eta,phi);
524 Float_t dphi = TMath::Abs(phi-phicl);
525 Float_t deta = TMath::Abs(eta-etacl);
526 Float_t R = TMath::Sqrt(deta*deta + dphi*dphi);
527 Float_t etcell = fEMCalCells->GetCellAmplitude(absid)*TMath::Sin(thetacl);
528 if(R<fIsoConeR){
529 totiso += etcell;
530 if(R<0.04)
531 totcore += etcell;
532 }
533 else{
534 if(dphi>fIsoConeR)
535 continue;
536 if(deta<fIsoConeR)
537 continue;
538 totband += fEMCalCells->GetCellAmplitude(absid)*TMath::Sin(thetacl);
539 }
540 }
541 iso = totiso;
542 phiband = totband;
543 core = totcore;
544}
bd092f0f 545
30159e6f 546//________________________________________________________________________
547void AliAnalysisTaskEMCALIsoPhoton::GetTrIso(TVector3 vec, Float_t &iso, Float_t &phiband, Float_t &core)
548{
bd092f0f 549 // Get track isolation.
550
30159e6f 551 if(!fSelPrimTracks)
552 return;
553 const Int_t ntracks = fSelPrimTracks->GetEntries();
554 Double_t totiso=0;
555 Double_t totband=0;
556 Double_t totcore=0;
557 Double_t etacl = vec.Eta();
558 Double_t phicl = vec.Phi();
559 if(phicl<0)
560 phicl+=TMath::TwoPi();
561 for(int itrack=0;itrack<ntracks;itrack++){
562 AliESDtrack *track = static_cast<AliESDtrack*> (fSelPrimTracks->At(itrack));
563 if(!track)
564 continue;
565 Double_t dphi = TMath::Abs(track->Phi()-phicl);
566 Double_t deta = TMath::Abs(track->Eta()-etacl);
567 Double_t R = TMath::Sqrt(deta*deta + dphi*dphi);
568 Double_t pt = track->Pt();
569 if(R<fIsoConeR){
570 totiso += track->Pt();
571 if(R<0.04)
572 totcore += pt;
573 }
574 else{
575 if(dphi>fIsoConeR)
576 continue;
577 if(deta<fIsoConeR)
578 continue;
579 totband += track->Pt();
580 }
581 }
582 iso = totiso;
583 phiband = totband;
584 core = totcore;
585}
bd092f0f 586
30159e6f 587//________________________________________________________________________
588Double_t AliAnalysisTaskEMCALIsoPhoton::GetCrossEnergy(const AliVCluster *cluster, Short_t &idmax)
589{
590 // Calculate the energy of cross cells around the leading cell.
591
592 AliVCaloCells *cells = 0;
593 cells = fESD->GetEMCALCells();
594 if (!cells)
595 return 0;
596
30159e6f 597 if (!fGeom)
598 return 0;
599
600 Int_t iSupMod = -1;
601 Int_t iTower = -1;
602 Int_t iIphi = -1;
603 Int_t iIeta = -1;
604 Int_t iphi = -1;
605 Int_t ieta = -1;
606 Int_t iphis = -1;
607 Int_t ietas = -1;
608
609 Double_t crossEnergy = 0;
610
611 fGeom->GetCellIndex(idmax,iSupMod,iTower,iIphi,iIeta);
612 fGeom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphis,ietas);
613
614 Int_t ncells = cluster->GetNCells();
615 for (Int_t i=0; i<ncells; i++) {
616 Int_t cellAbsId = cluster->GetCellAbsId(i);
617 fGeom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
618 fGeom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
619 Int_t aphidiff = TMath::Abs(iphi-iphis);
620 if (aphidiff>1)
621 continue;
622 Int_t aetadiff = TMath::Abs(ieta-ietas);
623 if (aetadiff>1)
624 continue;
625 if ( (aphidiff==1 && aetadiff==0) ||
626 (aphidiff==0 && aetadiff==1) ) {
627 crossEnergy += cells->GetCellAmplitude(cellAbsId);
628 }
629 }
630
631 return crossEnergy;
632}
633
30159e6f 634//________________________________________________________________________
635Double_t AliAnalysisTaskEMCALIsoPhoton ::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const
636{
637 // Get maximum energy of attached cell.
638
639 id = -1;
640
641 AliVCaloCells *cells = 0;
642 cells = fESD->GetEMCALCells();
643 if (!cells)
644 return 0;
645
646 Double_t maxe = 0;
647 Int_t ncells = cluster->GetNCells();
648 for (Int_t i=0; i<ncells; i++) {
649 Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
650 if (e>maxe) {
651 maxe = e;
652 id = cluster->GetCellAbsId(i);
653 }
654 }
655 return maxe;
656}
bd092f0f 657
30159e6f 658//________________________________________________________________________
659void AliAnalysisTaskEMCALIsoPhoton::Terminate(Option_t *)
660{
bd092f0f 661 // Called once at the end of the query.
30159e6f 662}