]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrDep/AliAnaParticleIsolation.cxx
Fixed access to MC-only analysis
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaParticleIsolation.cxx
CommitLineData
a3aebfff 1/**************************************************************************
477d6cee 2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15/* $Id: AliAnaParticleIsolation.cxx 28688 2008-09-11 15:04:07Z gconesab $ */
16
17//_________________________________________________________________________
18// Class for analysis of particle isolation
19// Input is selected particles put in AOD branch (AliAODPWG4ParticleCorrelation)
20//
21// Class created from old AliPHOSGammaJet
22// (see AliRoot versions previous Release 4-09)
23//
24// -- Author: Gustavo Conesa (LNF-INFN)
25//////////////////////////////////////////////////////////////////////////////
26
27
28// --- ROOT system ---
29#include <TClonesArray.h>
30#include <TList.h>
0c1383b5 31#include <TObjString.h>
477d6cee 32#include <TH2F.h>
a3aebfff 33//#include <Riostream.h>
9415d854 34#include <TClass.h>
477d6cee 35
36// --- Analysis system ---
37#include "AliAnaParticleIsolation.h"
38#include "AliCaloTrackReader.h"
39#include "AliIsolationCut.h"
40#include "AliNeutralMesonSelection.h"
41#include "AliAODPWG4ParticleCorrelation.h"
42#include "AliMCAnalysisUtils.h"
43#include "AliAODTrack.h"
44#include "AliAODCaloCluster.h"
45
46ClassImp(AliAnaParticleIsolation)
47
48//____________________________________________________________________________
49 AliAnaParticleIsolation::AliAnaParticleIsolation() :
233e0df8 50 AliAnaPartCorrBaseClass(), fCalorimeter(""),
477d6cee 51 fReMakeIC(0), fMakeSeveralIC(0), fMakeInvMass(0),
52 fhPtIso(0),fhPhiIso(0),fhEtaIso(0), fhConeSumPt(0), fhPtInCone(0),
53 //Several IC
54 fNCones(0),fNPtThresFrac(0), fConeSizes(), fPtThresholds(), fPtFractions(),
55 //MC
56 fhPtIsoPrompt(0),fhPhiIsoPrompt(0),fhEtaIsoPrompt(0),
57 fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(),
58 fhPtIsoFragmentation(0),fhPhiIsoFragmentation(0),fhEtaIsoFragmentation(0),
59 fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(),
60 fhPtIsoPi0Decay(0),fhPhiIsoPi0Decay(0),fhEtaIsoPi0Decay(0),
61 fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(),
62 fhPtIsoOtherDecay(0),fhPhiIsoOtherDecay(0),fhEtaIsoOtherDecay(0),
63 fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(),
64 fhPtIsoConversion(0),fhPhiIsoConversion(0),fhEtaIsoConversion(0),
65 fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(),
66 fhPtIsoUnknown(0),fhPhiIsoUnknown(0),fhEtaIsoUnknown(0),
67 fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown(),
68 //Histograms settings
69 fHistoNPtSumBins(0), fHistoPtSumMax(0.), fHistoPtSumMin(0.),
70 fHistoNPtInConeBins(0), fHistoPtInConeMax(0.), fHistoPtInConeMin(0.)
71{
72 //default ctor
73
74 //Initialize parameters
75 InitParameters();
233e0df8 76
477d6cee 77 for(Int_t i = 0; i < 5 ; i++){
78 fConeSizes[i] = 0 ;
79 fhPtSumIsolated[i] = 0 ;
80
81 fhPtSumIsolatedPrompt[i] = 0 ;
82 fhPtSumIsolatedFragmentation[i] = 0 ;
83 fhPtSumIsolatedPi0Decay[i] = 0 ;
84 fhPtSumIsolatedOtherDecay[i] = 0 ;
85 fhPtSumIsolatedConversion[i] = 0 ;
86 fhPtSumIsolatedUnknown[i] = 0 ;
87
88 for(Int_t j = 0; j < 5 ; j++){
89 fhPtThresIsolated[i][j] = 0 ;
90 fhPtFracIsolated[i][j] = 0 ;
91
92 fhPtThresIsolatedPrompt[i][j] = 0 ;
93 fhPtThresIsolatedFragmentation[i][j] = 0 ;
94 fhPtThresIsolatedPi0Decay[i][j] = 0 ;
95 fhPtThresIsolatedOtherDecay[i][j] = 0 ;
96 fhPtThresIsolatedConversion[i][j] = 0 ;
97 fhPtThresIsolatedUnknown[i][j] = 0 ;
98
99 fhPtFracIsolatedPrompt[i][j] = 0 ;
100 fhPtFracIsolatedFragmentation[i][j] = 0 ;
101 fhPtFracIsolatedPi0Decay[i][j] = 0 ;
102 fhPtFracIsolatedOtherDecay[i][j] = 0 ;
103 fhPtFracIsolatedConversion[i][j] = 0 ;
104 fhPtFracIsolatedUnknown[i][j] = 0 ;
105
106 }
107 }
108
109 for(Int_t i = 0; i < 5 ; i++){
110 fPtFractions[i]= 0 ;
111 fPtThresholds[i]= 0 ;
112 }
113
114
115}
78219bac 116/*
477d6cee 117//____________________________________________________________________________
118AliAnaParticleIsolation::AliAnaParticleIsolation(const AliAnaParticleIsolation & g) :
233e0df8 119 AliAnaPartCorrBaseClass(g), fCalorimeter(g.fCalorimeter),
477d6cee 120 fReMakeIC(g.fReMakeIC), fMakeSeveralIC(g.fMakeSeveralIC), fMakeInvMass(g.fMakeInvMass),
121 fhPtIso(g.fhPtIso),fhPhiIso(g.fhPhiIso),fhEtaIso(g.fhEtaIso),
122 fhConeSumPt(g.fhConeSumPt), fhPtInCone(g.fhPtInCone),
123 //Several IC
124 fNCones(g.fNCones),fNPtThresFrac(g.fNPtThresFrac), fConeSizes(), fPtThresholds(), fPtFractions(),
125 fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(),
126 //MC
127 fhPtIsoPrompt(g.fhPtIsoPrompt),fhPhiIsoPrompt(g.fhPhiIsoPrompt),fhEtaIsoPrompt(g.fhEtaIsoPrompt),
128 fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(),
129 fhPtIsoFragmentation(g.fhPtIsoFragmentation),fhPhiIsoFragmentation(g.fhPhiIsoFragmentation),fhEtaIsoFragmentation(g.fhEtaIsoFragmentation),
130 fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(),
131 fhPtIsoPi0Decay(g.fhPtIsoPi0Decay),fhPhiIsoPi0Decay(g.fhPhiIsoPi0Decay),fhEtaIsoPi0Decay(g.fhEtaIsoPi0Decay),
132 fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(),
133 fhPtIsoOtherDecay(g.fhPtIsoOtherDecay),fhPhiIsoOtherDecay(g.fhPhiIsoOtherDecay),fhEtaIsoOtherDecay(g.fhEtaIsoOtherDecay),
134 fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(),
135 fhPtIsoConversion(g. fhPtIsoConversion),fhPhiIsoConversion(g.fhPhiIsoConversion),fhEtaIsoConversion(g.fhEtaIsoConversion),
136 fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(),
137 fhPtIsoUnknown(g.fhPtIsoUnknown),fhPhiIsoUnknown(g.fhPhiIsoUnknown),fhEtaIsoUnknown(g.fhEtaIsoUnknown),
138 fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown(),
139 //Histograms
140 fHistoNPtSumBins(g.fHistoNPtSumBins), fHistoPtSumMax(g.fHistoPtSumMax), fHistoPtSumMin(g.fHistoPtSumMax),
141 fHistoNPtInConeBins(g.fHistoNPtInConeBins), fHistoPtInConeMax(g.fHistoPtInConeMax), fHistoPtInConeMin(g.fHistoPtInConeMin)
142{
143 // cpy ctor
233e0df8 144
477d6cee 145 //Several IC
146 for(Int_t i = 0; i < fNCones ; i++){
147 fConeSizes[i] = g.fConeSizes[i];
148 fhPtSumIsolated[i] = g.fhPtSumIsolated[i];
149
150 fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i];
151 fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i];
152 fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i];
153 fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i];
154 fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i];
155 fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i];
156
157 for(Int_t j = 0; j < fNPtThresFrac ; j++){
158 fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j];
159 fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j];
160
161 fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j];
162 fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j];
163 fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j];
164 fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j];
165 fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j];
166 fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j];
167
168 fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j];
169 fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j];
170 fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j];
171 fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j];
172 fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j];
173 fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j];
174
175 }
176 }
177
178 for(Int_t i = 0; i < fNPtThresFrac ; i++){
179 fPtFractions[i]= g.fPtFractions[i];
180 fPtThresholds[i]= g.fPtThresholds[i];
181 }
182
183}
184
185//_________________________________________________________________________
186AliAnaParticleIsolation & AliAnaParticleIsolation::operator = (const AliAnaParticleIsolation & g)
187{
188 // assignment operator
189
190 if(&g == this) return *this;
191
233e0df8 192 fReMakeIC = g.fReMakeIC ;
477d6cee 193 fMakeSeveralIC = g.fMakeSeveralIC ;
233e0df8 194 fMakeInvMass = g.fMakeInvMass ;
195 fCalorimeter = g.fCalorimeter ;
196
197 fhConeSumPt = g.fhConeSumPt ;
198 fhPtInCone = g.fhPtInCone;
477d6cee 199
233e0df8 200 fhPtIso = g.fhPtIso ;
477d6cee 201 fhPhiIso = g.fhPhiIso ;
202 fhEtaIso = g.fhEtaIso ;
203
233e0df8 204 fhPtIsoPrompt = g.fhPtIsoPrompt;
477d6cee 205 fhPhiIsoPrompt = g.fhPhiIsoPrompt;
206 fhEtaIsoPrompt = g.fhEtaIsoPrompt;
233e0df8 207 fhPtIsoFragmentation = g.fhPtIsoFragmentation;
477d6cee 208 fhPhiIsoFragmentation = g.fhPhiIsoFragmentation;
209 fhEtaIsoFragmentation = g.fhEtaIsoFragmentation;
233e0df8 210 fhPtIsoPi0Decay = g.fhPtIsoPi0Decay;
477d6cee 211 fhPhiIsoPi0Decay = g.fhPhiIsoPi0Decay;
212 fhEtaIsoPi0Decay = g.fhEtaIsoPi0Decay;
233e0df8 213 fhPtIsoOtherDecay = g.fhPtIsoOtherDecay;
477d6cee 214 fhPhiIsoOtherDecay = g.fhPhiIsoOtherDecay;
215 fhEtaIsoOtherDecay = g.fhEtaIsoOtherDecay;
233e0df8 216 fhPtIsoConversion = g. fhPtIsoConversion;
477d6cee 217 fhPhiIsoConversion = g.fhPhiIsoConversion;
218 fhEtaIsoConversion = g.fhEtaIsoConversion;
233e0df8 219 fhPtIsoUnknown = g.fhPtIsoUnknown;
477d6cee 220 fhPhiIsoUnknown = g.fhPhiIsoUnknown;
221 fhEtaIsoUnknown = g.fhEtaIsoUnknown;
222
223 //Several IC
224 fNCones = g.fNCones ;
225 fNPtThresFrac = g.fNPtThresFrac ;
226
227 for(Int_t i = 0; i < fNCones ; i++){
228 fConeSizes[i] = g.fConeSizes[i];
229 fhPtSumIsolated[i] = g.fhPtSumIsolated[i] ;
230
231 fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i];
232 fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i];
233 fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i];
234 fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i];
235 fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i];
236 fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i];
237
238 for(Int_t j = 0; j < fNPtThresFrac ; j++){
239 fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j] ;
240 fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j] ;
241
242 fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j];
243 fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j];
244 fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j];
245 fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j];
246 fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j];
247 fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j];
248
249 fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j];
250 fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j];
251 fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j];
252 fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j];
253 fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j];
254 fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j];
255
256 }
257 }
258
259 for(Int_t i = 0; i < fNPtThresFrac ; i++){
260 fPtThresholds[i]= g.fPtThresholds[i];
261 fPtFractions[i]= g.fPtFractions[i];
262 }
263
264
265 fHistoNPtSumBins = g.fHistoNPtSumBins;
266 fHistoPtSumMax = g.fHistoPtSumMax;
267 fHistoPtSumMin = g.fHistoPtSumMax;
268 fHistoNPtInConeBins = g.fHistoNPtInConeBins;
269 fHistoPtInConeMax = g.fHistoPtInConeMax;
270 fHistoPtInConeMin = g.fHistoPtInConeMin;
271
272 return *this;
273
274}
78219bac 275*/
477d6cee 276
277//____________________________________________________________________________
278AliAnaParticleIsolation::~AliAnaParticleIsolation()
279{
280 //dtor
281 //do not delete histograms
282
7787a778 283 //delete [] fConeSizes ;
284 //delete [] fPtThresholds ;
285 //delete [] fPtFractions ;
233e0df8 286
477d6cee 287}
288
289//_________________________________________________________________________
decca433 290Bool_t AliAnaParticleIsolation::CheckInvMass(const Int_t iaod, const AliAODPWG4Particle * part1)
477d6cee 291{
292 // Search if there is a companion decay photon to the candidate
293 // and discard it in such case
294 // Use it only if isolation candidates are photons
295 // Make sure that no selection on photon pt is done in the input aod photon list.
296 TLorentzVector mom1 = *(part1->Momentum());
297 TLorentzVector mom2 ;
298 for(Int_t jaod = 0; jaod < GetInputAODBranch()->GetEntriesFast(); jaod++){
299 if(iaod == jaod) continue ;
300 AliAODPWG4ParticleCorrelation * part2 = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(jaod));
301 mom2 = *(part2->Momentum());
302 //Select good pair (good phi, pt cuts, aperture and invariant mass)
303 if(GetNeutralMesonSelection()->SelectPair(mom1, mom2)){
1e1befe8 304 if(GetDebug() > 1)printf("AliAnaParticleIsolation::CheckInvMass() - Selected gamma pair: pt %f, phi %f, eta%f\n",(mom1+mom2).Pt(), (mom1+mom2).Phi()*180./3.1416, (mom1+mom2).Eta());
477d6cee 305 return kTRUE ;
306 }
307 }//loop
308
309 return kFALSE;
310}
311
0c1383b5 312//________________________________________________________________________
313TObjString * AliAnaParticleIsolation::GetAnalysisCuts()
314{
315 //Save parameters used for analysis
316 TString parList ; //this will be list of parameters used for this analysis.
317 char onePar[255] ;
318
319 sprintf(onePar,"--- AliAnaParticleIsolation ---\n") ;
320 parList+=onePar ;
321 sprintf(onePar,"Calorimeter: %s\n",fCalorimeter.Data()) ;
322 parList+=onePar ;
323 sprintf(onePar,"fReMakeIC =%d (Flag for reisolation during histogram filling) \n",fReMakeIC) ;
324 parList+=onePar ;
325 sprintf(onePar,"fMakeSeveralIC=%d (Flag for isolation with several cuts at the same time ) \n",fMakeSeveralIC) ;
326 parList+=onePar ;
327 sprintf(onePar,"fMakeInvMass=%d (Flag for rejection of candidates with a pi0 inv mass pair) \n",fMakeInvMass) ;
328 parList+=onePar ;
329
330 if(fMakeSeveralIC){
331 sprintf(onePar,"fNCones =%d (Number of cone sizes) \n",fNCones) ;
332 parList+=onePar ;
333 sprintf(onePar,"fNPtThresFrac=%d (Flag for isolation with several cuts at the same time ) \n",fNPtThresFrac) ;
334 parList+=onePar ;
335
336 for(Int_t icone = 0; icone < fNCones ; icone++){
337 sprintf(onePar,"fConeSizes[%d]=%1.2f (isolation cone size) \n",icone, fConeSizes[icone]) ;
338 parList+=onePar ;
339 }
340 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++){
341 sprintf(onePar,"fPtThresholds[%d]=%1.2f (isolation pt threshold) \n",ipt, fPtThresholds[ipt]) ;
342 parList+=onePar ;
343 }
344 for(Int_t ipt = 0; ipt < fNPtThresFrac ; ipt++){
345 sprintf(onePar,"fPtFractions[%d]=%1.2f (isolation pt fraction threshold) \n",ipt, fPtFractions[ipt]) ;
346 parList+=onePar ;
347 }
348 }
349
350 //Get parameters set in base class.
351 parList += GetBaseParametersList() ;
352
353 //Get parameters set in IC class.
354 if(!fMakeSeveralIC)parList += GetIsolationCut()->GetICParametersList() ;
355
356 return new TObjString(parList) ;
357
358}
359
477d6cee 360//________________________________________________________________________
361TList * AliAnaParticleIsolation::GetCreateOutputObjects()
362{
363 // Create histograms to be saved in output file and
364 // store them in outputContainer
365 TList * outputContainer = new TList() ;
366 outputContainer->SetName("IsolatedParticleHistos") ;
367
5a2dbc3c 368 Int_t nptbins = GetHistoPtBins();
369 Int_t nphibins = GetHistoPhiBins();
370 Int_t netabins = GetHistoEtaBins();
477d6cee 371 Float_t ptmax = GetHistoPtMax();
372 Float_t phimax = GetHistoPhiMax();
373 Float_t etamax = GetHistoEtaMax();
374 Float_t ptmin = GetHistoPtMin();
375 Float_t phimin = GetHistoPhiMin();
376 Float_t etamin = GetHistoEtaMin();
377
378 Int_t nptsumbins = fHistoNPtSumBins;
379 Float_t ptsummax = fHistoPtSumMax;
380 Float_t ptsummin = fHistoPtSumMin;
381 Int_t nptinconebins = fHistoNPtInConeBins;
382 Float_t ptinconemax = fHistoPtInConeMax;
383 Float_t ptinconemin = fHistoPtInConeMin;
384
385 if(!fMakeSeveralIC){
386
387 fhConeSumPt = new TH2F
a3aebfff 388 ("hConePtSum","#Sigma p_{T} in isolation cone ",nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
477d6cee 389 fhConeSumPt->SetYTitle("#Sigma p_{T}");
a3aebfff 390 fhConeSumPt->SetXTitle("p_{T} (GeV/c)");
477d6cee 391 outputContainer->Add(fhConeSumPt) ;
392
393 fhPtInCone = new TH2F
a3aebfff 394 ("hPtInCone","p_{T} in isolation cone ",nptbins,ptmin,ptmax,nptinconebins,ptinconemin,ptinconemax);
477d6cee 395 fhPtInCone->SetYTitle("p_{T in cone} (GeV/c)");
a3aebfff 396 fhPtInCone->SetXTitle("p_{T} (GeV/c)");
477d6cee 397 outputContainer->Add(fhPtInCone) ;
398
a3aebfff 399 fhPtIso = new TH1F("hPt","Isolated Number of particles",nptbins,ptmin,ptmax);
477d6cee 400 fhPtIso->SetYTitle("N");
a3aebfff 401 fhPtIso->SetXTitle("p_{T}(GeV/c)");
477d6cee 402 outputContainer->Add(fhPtIso) ;
403
404 fhPhiIso = new TH2F
a3aebfff 405 ("hPhi","Isolated Number of particles",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 406 fhPhiIso->SetYTitle("#phi");
a3aebfff 407 fhPhiIso->SetXTitle("p_{T} (GeV/c)");
477d6cee 408 outputContainer->Add(fhPhiIso) ;
409
410 fhEtaIso = new TH2F
a3aebfff 411 ("hEta","Isolated Number of particles",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 412 fhEtaIso->SetYTitle("#eta");
a3aebfff 413 fhEtaIso->SetXTitle("p_{T} (GeV/c)");
477d6cee 414 outputContainer->Add(fhEtaIso) ;
415
416 if(IsDataMC()){
417
a3aebfff 418 fhPtIsoPrompt = new TH1F("hPtMCPrompt","Isolated Number of #gamma prompt",nptbins,ptmin,ptmax);
477d6cee 419 fhPtIsoPrompt->SetYTitle("N");
420 fhPtIsoPrompt->SetXTitle("p_{T #gamma}(GeV/c)");
421 outputContainer->Add(fhPtIsoPrompt) ;
422
423 fhPhiIsoPrompt = new TH2F
a3aebfff 424 ("hPhiMCPrompt","Isolated Number of #gamma prompt",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 425 fhPhiIsoPrompt->SetYTitle("#phi");
426 fhPhiIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
427 outputContainer->Add(fhPhiIsoPrompt) ;
428
429 fhEtaIsoPrompt = new TH2F
a3aebfff 430 ("hEtaMCPrompt","Isolated Number of #gamma prompt ",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 431 fhEtaIsoPrompt->SetYTitle("#eta");
432 fhEtaIsoPrompt->SetXTitle("p_{T #gamma} (GeV/c)");
433 outputContainer->Add(fhEtaIsoPrompt) ;
434
a3aebfff 435 fhPtIsoFragmentation = new TH1F("hPtMCFragmentation","Isolated Number of #gamma",nptbins,ptmin,ptmax);
477d6cee 436 fhPtIsoFragmentation->SetYTitle("N");
437 fhPtIsoFragmentation->SetXTitle("p_{T #gamma}(GeV/c)");
438 outputContainer->Add(fhPtIsoFragmentation) ;
439
440 fhPhiIsoFragmentation = new TH2F
a3aebfff 441 ("hPhiMCFragmentation","Isolated Number of #gamma fragmentation",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 442 fhPhiIsoFragmentation->SetYTitle("#phi");
443 fhPhiIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
444 outputContainer->Add(fhPhiIsoFragmentation) ;
445
446 fhEtaIsoFragmentation = new TH2F
a3aebfff 447 ("hEtaMCFragmentation","Isolated Number of #gamma fragmentation",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 448 fhEtaIsoFragmentation->SetYTitle("#eta");
449 fhEtaIsoFragmentation->SetXTitle("p_{T #gamma} (GeV/c)");
450 outputContainer->Add(fhEtaIsoFragmentation) ;
451
a3aebfff 452 fhPtIsoPi0Decay = new TH1F("hPtMCPi0Decay","Isolated Number of #gamma from #pi^{0} decay",nptbins,ptmin,ptmax);
477d6cee 453 fhPtIsoPi0Decay->SetYTitle("N");
454 fhPtIsoPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)");
455 outputContainer->Add(fhPtIsoPi0Decay) ;
456
457 fhPhiIsoPi0Decay = new TH2F
a3aebfff 458 ("hPhiMCPi0Decay","Isolated Number of #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 459 fhPhiIsoPi0Decay->SetYTitle("#phi");
460 fhPhiIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
461 outputContainer->Add(fhPhiIsoPi0Decay) ;
462
463 fhEtaIsoPi0Decay = new TH2F
a3aebfff 464 ("hEtaMCPi0Decay","Isolated Number of #gamma from #pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 465 fhEtaIsoPi0Decay->SetYTitle("#eta");
466 fhEtaIsoPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)");
467 outputContainer->Add(fhEtaIsoPi0Decay) ;
468
a3aebfff 469 fhPtIsoOtherDecay = new TH1F("hPtMCOtherDecay","Isolated Number of #gamma from non #pi^{0} decay",nptbins,ptmin,ptmax);
477d6cee 470 fhPtIsoOtherDecay->SetYTitle("N");
471 fhPtIsoOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)");
472 outputContainer->Add(fhPtIsoOtherDecay) ;
473
474 fhPhiIsoOtherDecay = new TH2F
a3aebfff 475 ("hPhiMCOtherDecay","Isolated Number of #gamma from non #pi^{0} decay",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 476 fhPhiIsoOtherDecay->SetYTitle("#phi");
477 fhPhiIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
478 outputContainer->Add(fhPhiIsoOtherDecay) ;
479
480 fhEtaIsoOtherDecay = new TH2F
a3aebfff 481 ("hEtaMCOtherDecay","Isolated Number of #gamma non #pi^{0} decay",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 482 fhEtaIsoOtherDecay->SetYTitle("#eta");
483 fhEtaIsoOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)");
484 outputContainer->Add(fhEtaIsoOtherDecay) ;
485
a3aebfff 486 fhPtIsoConversion = new TH1F("hPtMCConversion","Isolated Number of #gamma converted",nptbins,ptmin,ptmax);
477d6cee 487 fhPtIsoConversion->SetYTitle("N");
488 fhPtIsoConversion->SetXTitle("p_{T #gamma}(GeV/c)");
489 outputContainer->Add(fhPtIsoConversion) ;
490
491 fhPhiIsoConversion = new TH2F
a3aebfff 492 ("hPhiMCConversion","Isolated Number of #gamma converted",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 493 fhPhiIsoConversion->SetYTitle("#phi");
494 fhPhiIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
495 outputContainer->Add(fhPhiIsoConversion) ;
496
497 fhEtaIsoConversion = new TH2F
a3aebfff 498 ("hEtaMCConversion","Isolated Number of #gamma converted",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 499 fhEtaIsoConversion->SetYTitle("#eta");
500 fhEtaIsoConversion->SetXTitle("p_{T #gamma} (GeV/c)");
501 outputContainer->Add(fhEtaIsoConversion) ;
502
a3aebfff 503 fhPtIsoUnknown = new TH1F("hPtMCUnknown","Isolated Number of non #gamma particles",nptbins,ptmin,ptmax);
477d6cee 504 fhPtIsoUnknown->SetYTitle("N");
a3aebfff 505 fhPtIsoUnknown->SetXTitle("p_{T}(GeV/c)");
477d6cee 506 outputContainer->Add(fhPtIsoUnknown) ;
507
508 fhPhiIsoUnknown = new TH2F
a3aebfff 509 ("hPhiMCUnknown","Isolated Number of non #gamma particles",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
477d6cee 510 fhPhiIsoUnknown->SetYTitle("#phi");
a3aebfff 511 fhPhiIsoUnknown->SetXTitle("p_{T} (GeV/c)");
477d6cee 512 outputContainer->Add(fhPhiIsoUnknown) ;
513
514 fhEtaIsoUnknown = new TH2F
a3aebfff 515 ("hEtaMCUnknown","Isolated Number of non #gamma particles",nptbins,ptmin,ptmax,netabins,etamin,etamax);
477d6cee 516 fhEtaIsoUnknown->SetYTitle("#eta");
a3aebfff 517 fhEtaIsoUnknown->SetXTitle("p_{T} (GeV/c)");
477d6cee 518 outputContainer->Add(fhEtaIsoUnknown) ;
519 }//Histos with MC
520
521 }
522
523 if(fMakeSeveralIC){
524 char name[128];
525 char title[128];
526 for(Int_t icone = 0; icone<fNCones; icone++){
a3aebfff 527 sprintf(name,"hPtSum_Cone_%d",icone);
477d6cee 528 sprintf(title,"Candidate cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
529 fhPtSumIsolated[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
530 fhPtSumIsolated[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
531 fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)");
532 outputContainer->Add(fhPtSumIsolated[icone]) ;
533
534 if(IsDataMC()){
a3aebfff 535 sprintf(name,"hPtSumPrompt_Cone_%d",icone);
477d6cee 536 sprintf(title,"Candidate Prompt cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
537 fhPtSumIsolatedPrompt[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
538 fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
539 fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)");
540 outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ;
541
a3aebfff 542 sprintf(name,"hPtSumFragmentation_Cone_%d",icone);
477d6cee 543 sprintf(title,"Candidate Fragmentation cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
544 fhPtSumIsolatedFragmentation[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
545 fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
546 fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)");
547 outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ;
548
a3aebfff 549 sprintf(name,"hPtSumPi0Decay_Cone_%d",icone);
477d6cee 550 sprintf(title,"Candidate Pi0Decay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
551 fhPtSumIsolatedPi0Decay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
552 fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
553 fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)");
554 outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ;
555
a3aebfff 556 sprintf(name,"hPtSumOtherDecay_Cone_%d",icone);
477d6cee 557 sprintf(title,"Candidate OtherDecay cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
558 fhPtSumIsolatedOtherDecay[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
559 fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
560 fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)");
561 outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ;
562
a3aebfff 563 sprintf(name,"hPtSumConversion_Cone_%d",icone);
477d6cee 564 sprintf(title,"Candidate Conversion cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
565 fhPtSumIsolatedConversion[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
566 fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
567 fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)");
568 outputContainer->Add(fhPtSumIsolatedConversion[icone]) ;
569
a3aebfff 570 sprintf(name,"hPtSumUnknown_Cone_%d",icone);
477d6cee 571 sprintf(title,"Candidate Unknown cone sum p_{T} for cone size %d vs candidate p_{T}",icone);
572 fhPtSumIsolatedUnknown[icone] = new TH2F(name, title,nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
573 fhPtSumIsolatedUnknown[icone]->SetYTitle("#Sigma p_{T} (GeV/c)");
574 fhPtSumIsolatedUnknown[icone]->SetXTitle("p_{T} (GeV/c)");
575 outputContainer->Add(fhPtSumIsolatedUnknown[icone]) ;
576
577 }//Histos with MC
578
579 for(Int_t ipt = 0; ipt<fNPtThresFrac;ipt++){
a3aebfff 580 sprintf(name,"hPtThres_Cone_%d_Pt%d",icone,ipt);
477d6cee 581 sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
582 fhPtThresIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
583 fhPtThresIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
584 outputContainer->Add(fhPtThresIsolated[icone][ipt]) ;
585
a3aebfff 586 sprintf(name,"hPtFrac_Cone_%d_Pt%d",icone,ipt);
477d6cee 587 sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
588 fhPtFracIsolated[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
589 fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
590 outputContainer->Add(fhPtFracIsolated[icone][ipt]) ;
591
592 if(IsDataMC()){
a3aebfff 593 sprintf(name,"hPtThresMCPrompt_Cone_%d_Pt%d",icone,ipt);
477d6cee 594 sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
595 fhPtThresIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
596 fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
597 outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ;
598
a3aebfff 599 sprintf(name,"hPtFracMCPrompt_Cone_%d_Pt%d",icone,ipt);
477d6cee 600 sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
601 fhPtFracIsolatedPrompt[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
602 fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
603 outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ;
604
a3aebfff 605 sprintf(name,"hPtThresMCFragmentation_Cone_%d_Pt%d",icone,ipt);
477d6cee 606 sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
607 fhPtThresIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
608 fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
609 outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ;
610
a3aebfff 611 sprintf(name,"hPtFracMCFragmentation_Cone_%d_Pt%d",icone,ipt);
477d6cee 612 sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
613 fhPtFracIsolatedFragmentation[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
614 fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
615 outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ;
616
a3aebfff 617 sprintf(name,"hPtThresMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
477d6cee 618 sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
619 fhPtThresIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
620 fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
621 outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ;
622
a3aebfff 623 sprintf(name,"hPtFracMCPi0Decay_Cone_%d_Pt%d",icone,ipt);
477d6cee 624 sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
625 fhPtFracIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
626 fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
627 outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ;
628
a3aebfff 629 sprintf(name,"hPtThresMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
477d6cee 630 sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
631 fhPtThresIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
632 fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
633 outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ;
634
a3aebfff 635 sprintf(name,"hPtFracMCOtherDecay_Cone_%d_Pt%d",icone,ipt);
477d6cee 636 sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
637 fhPtFracIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
638 fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
639 outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ;
640
a3aebfff 641 sprintf(name,"hPtThresMCConversion_Cone_%d_Pt%d",icone,ipt);
477d6cee 642 sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
643 fhPtThresIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
644 fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
645 outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ;
646
a3aebfff 647 sprintf(name,"hPtFracMCConversion_Cone_%d_Pt%d",icone,ipt);
477d6cee 648 sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
649 fhPtFracIsolatedConversion[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
650 fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
651 outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ;
652
a3aebfff 653 sprintf(name,"hPtThresMCUnknown_Cone_%d_Pt%d",icone,ipt);
477d6cee 654 sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
655 fhPtThresIsolatedUnknown[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
656 fhPtThresIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
657 outputContainer->Add(fhPtThresIsolatedUnknown[icone][ipt]) ;
658
a3aebfff 659 sprintf(name,"hPtFracMCUnknown_Cone_%d_Pt%d",icone,ipt);
477d6cee 660 sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt);
661 fhPtFracIsolatedUnknown[icone][ipt] = new TH1F(name, title,nptbins,ptmin,ptmax);
662 fhPtFracIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)");
663 outputContainer->Add(fhPtFracIsolatedUnknown[icone][ipt]) ;
664
665 }//Histos with MC
666
667 }//icone loop
668 }//ipt loop
669 }
670
671 //Keep neutral meson selection histograms if requiered
672 //Setting done in AliNeutralMesonSelection
673 if(fMakeInvMass && GetNeutralMesonSelection()){
674 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
675 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
676 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
a14fd526 677 delete nmsHistos;
477d6cee 678 }
679
477d6cee 680 return outputContainer ;
681
682}
683
684//__________________________________________________________________
685void AliAnaParticleIsolation::MakeAnalysisFillAOD()
686{
687 //Do analysis and fill aods
688 //Search for the isolated photon in fCalorimeter with pt > GetMinPt()
689
690 if(!GetInputAODBranch()){
a3aebfff 691 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 692 abort();
693 }
a3aebfff 694
9415d854 695 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
696 printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
697 abort();
698 }
699
477d6cee 700 Int_t n = 0, nfrac = 0;
701 Bool_t isolated = kFALSE ;
702 Float_t coneptsum = 0 ;
4a745797 703 TObjArray * pl = 0x0; ;
477d6cee 704
705 //Select the calorimeter for candidate isolation with neutral particles
706 if(fCalorimeter == "PHOS")
707 pl = GetAODPHOS();
708 else if (fCalorimeter == "EMCAL")
709 pl = GetAODEMCAL();
710
477d6cee 711 //Loop on AOD branch, filled previously in AliAnaPhoton
712 TLorentzVector mom ;
1e1befe8 713 Int_t naod = GetInputAODBranch()->GetEntriesFast();
714 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - Input aod branch entries %d\n", naod);
715 for(Int_t iaod = 0; iaod < naod; iaod++){
a3aebfff 716 AliAODPWG4ParticleCorrelation * aodinput = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
477d6cee 717
718 //If too small or too large pt, skip
a3aebfff 719 if(aodinput->Pt() < GetMinPt() || aodinput->Pt() > GetMaxPt() ) continue ;
477d6cee 720
721 //Check invariant mass, if pi0, skip.
722 Bool_t decay = kFALSE ;
a3aebfff 723 if(fMakeInvMass) decay = CheckInvMass(iaod,aodinput);
477d6cee 724 if(decay) continue ;
a3aebfff 725
726 //After cuts, study isolation
477d6cee 727 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
233e0df8 728 GetIsolationCut()->MakeIsolationCut(GetAODCTS(),pl,GetReader(), kTRUE, aodinput, GetAODObjArrayName(), n,nfrac,coneptsum, isolated);
a3aebfff 729 aodinput->SetIsolated(isolated);
1e1befe8 730 if(GetDebug() > 1 && isolated) printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() : Particle %d IS ISOLATED \n",iaod);
731
477d6cee 732 }//loop
733
1e1befe8 734 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillAOD() - End fill AODs \n");
477d6cee 735
736}
737
738//__________________________________________________________________
739void AliAnaParticleIsolation::MakeAnalysisFillHistograms()
740{
741 //Do analysis and fill histograms
742 Int_t n = 0, nfrac = 0;
743 Bool_t isolated = kFALSE ;
744 Float_t coneptsum = 0 ;
477d6cee 745 //Loop on stored AOD
746 Int_t naod = GetInputAODBranch()->GetEntriesFast();
1e1befe8 747 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Histo aod branch entries %d\n", naod);
233e0df8 748
749 //Get vertex for photon momentum calculation
750 Double_t vertex[]={0,0,0} ; //vertex ;
751 Double_t vertex2[]={0,0,0} ; //vertex ;
edd59991 752 if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {
233e0df8 753 GetReader()->GetVertex(vertex);
754 if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex2);
755 }
756
477d6cee 757 for(Int_t iaod = 0; iaod < naod ; iaod++){
758 AliAODPWG4ParticleCorrelation* aod = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
a3aebfff 759
760 Bool_t isolation = aod->IsIsolated();
477d6cee 761 Float_t ptcluster = aod->Pt();
762 Float_t phicluster = aod->Phi();
763 Float_t etacluster = aod->Eta();
a3aebfff 764 //Recover reference arrays with clusters and tracks
591cc579 765 TObjArray * refclusters = aod->GetObjArray(GetAODObjArrayName()+"Clusters");
b6991fc4 766 TObjArray * reftracks = aod->GetObjArray(GetAODObjArrayName()+"Tracks");
a3aebfff 767
b6991fc4 768 //If too small or too large pt, skip
769 if(aod->Pt() < GetMinPt() || aod->Pt() > GetMaxPt() ) continue ;
770
477d6cee 771 if(fMakeSeveralIC) {
772 //Analysis of multiple IC at same time
773 MakeSeveralICAnalysis(aod);
774 continue;
775 }
776 else if(fReMakeIC){
777 //In case a more strict IC is needed in the produced AOD
778 n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0;
233e0df8 779 GetIsolationCut()->MakeIsolationCut(reftracks, refclusters, GetReader(), kFALSE, aod, "", n,nfrac,coneptsum, isolated);
1e1befe8 780 fhConeSumPt->Fill(ptcluster,coneptsum);
781 if(GetDebug() > 0) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Energy Sum in Isolation Cone %2.2f\n", coneptsum);
782 }
477d6cee 783
784 //Fill pt distribution of particles in cone
785 //Tracks
786 coneptsum=0;
a3aebfff 787 if(reftracks){
788 for(Int_t itrack=0; itrack < reftracks->GetEntriesFast(); itrack++){
789 AliAODTrack* track = (AliAODTrack *) reftracks->At(itrack);
790 fhPtInCone->Fill(ptcluster,TMath::Sqrt(track->Px()*track->Px()+track->Py()*track->Py()));
791 coneptsum+=track->Pt();
792 }
793 }
794
477d6cee 795 //CaloClusters
a3aebfff 796 if(refclusters){
797 TLorentzVector mom ;
798 for(Int_t icalo=0; icalo < refclusters->GetEntriesFast(); icalo++){
799 AliAODCaloCluster* calo = (AliAODCaloCluster *) refclusters->At(icalo);
233e0df8 800 Int_t input = 0;
801 if (fCalorimeter == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= icalo) input = 1 ;
802 else if(fCalorimeter == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= icalo) input = 1;
803
804 //Get Momentum vector,
805 if (input == 0) calo->GetMomentum(mom,vertex) ;//Assume that come from vertex in straight line
806 else if(input == 1) calo->GetMomentum(mom,vertex2);//Assume that come from vertex in straight line
807
a3aebfff 808 fhPtInCone->Fill(ptcluster, mom.Pt());
809 coneptsum+=mom.Pt();
810 }
477d6cee 811 }
a3aebfff 812
1e1befe8 813 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d Energy Sum in Isolation Cone %2.2f\n", iaod, coneptsum);
814
815 if(!fReMakeIC) fhConeSumPt->Fill(ptcluster,coneptsum);
477d6cee 816
817 if(isolation){
1e1befe8 818
819 if(GetDebug() > 1) printf("AliAnaParticleIsolation::MakeAnalysisFillHistograms() - Particle %d ISOLATED, fill histograms\n", iaod);
820
477d6cee 821 fhPtIso ->Fill(ptcluster);
822 fhPhiIso ->Fill(ptcluster,phicluster);
823 fhEtaIso ->Fill(ptcluster,etacluster);
824
825 if(IsDataMC()){
826 Int_t tag =aod->GetTag();
827
591cc579 828 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)){
477d6cee 829 fhPtIsoPrompt ->Fill(ptcluster);
830 fhPhiIsoPrompt ->Fill(ptcluster,phicluster);
831 fhEtaIsoPrompt ->Fill(ptcluster,etacluster);
832 }
591cc579 833 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
477d6cee 834 {
835 fhPtIsoFragmentation ->Fill(ptcluster);
836 fhPhiIsoFragmentation ->Fill(ptcluster,phicluster);
837 fhEtaIsoFragmentation ->Fill(ptcluster,etacluster);
838 }
591cc579 839 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
477d6cee 840 {
841 fhPtIsoPi0Decay ->Fill(ptcluster);
842 fhPhiIsoPi0Decay ->Fill(ptcluster,phicluster);
843 fhEtaIsoPi0Decay ->Fill(ptcluster,etacluster);
844 }
591cc579 845 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) || GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay))
477d6cee 846 {
847 fhPtIsoOtherDecay ->Fill(ptcluster);
848 fhPhiIsoOtherDecay ->Fill(ptcluster,phicluster);
849 fhEtaIsoOtherDecay ->Fill(ptcluster,etacluster);
850 }
591cc579 851 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion))
477d6cee 852 {
853 fhPtIsoConversion ->Fill(ptcluster);
854 fhPhiIsoConversion ->Fill(ptcluster,phicluster);
855 fhEtaIsoConversion ->Fill(ptcluster,etacluster);
856 }
857 else
858 {
859 fhPtIsoUnknown ->Fill(ptcluster);
860 fhPhiIsoUnknown ->Fill(ptcluster,phicluster);
861 fhEtaIsoUnknown ->Fill(ptcluster,etacluster);
862 }
863 }//Histograms with MC
864
865 }//Isolated histograms
866
867 }// aod loop
868
869}
870
871//____________________________________________________________________________
872void AliAnaParticleIsolation::InitParameters()
873{
874
875 //Initialize the parameters of the analysis.
a3aebfff 876 SetInputAODName("PWG4Particle");
591cc579 877 SetAODObjArrayName("IsolationCone");
a3aebfff 878 AddToHistogramsName("AnaIsolation_");
879
477d6cee 880 fCalorimeter = "PHOS" ;
881 fReMakeIC = kFALSE ;
882 fMakeSeveralIC = kFALSE ;
883 fMakeInvMass = kFALSE ;
884
885 //----------- Several IC-----------------
886 fNCones = 5 ;
887 fNPtThresFrac = 5 ;
888 fConeSizes[0] = 0.1; fConeSizes[1] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4; fConeSizes[4] = 0.5;
889 fPtThresholds[0] = 1.; fPtThresholds[1] = 2.; fPtThresholds[2] = 3.; fPtThresholds[3] = 4.; fPtThresholds[4] = 5.;
890 fPtFractions[0] = 0.05; fPtFractions[1] = 0.075; fPtFractions[2] = 0.1; fPtFractions[3] = 1.25; fPtFractions[4] = 1.5;
891
892//------------- Histograms settings -------
893 fHistoNPtSumBins = 100 ;
894 fHistoPtSumMax = 50 ;
895 fHistoPtSumMin = 0. ;
896
897 fHistoNPtInConeBins = 100 ;
898 fHistoPtInConeMax = 50 ;
899 fHistoPtInConeMin = 0. ;
900
901}
902
903//__________________________________________________________________
904void AliAnaParticleIsolation::MakeSeveralICAnalysis(AliAODPWG4ParticleCorrelation* ph)
905{
906 //Isolation Cut Analysis for both methods and different pt cuts and cones
907 Float_t ptC = ph->Pt();
908 Int_t tag = ph->GetTag();
909
910 //Keep original setting used when filling AODs, reset at end of analysis
911 Float_t ptthresorg = GetIsolationCut()->GetPtThreshold();
912 Float_t ptfracorg = GetIsolationCut()->GetPtFraction();
913 Float_t rorg = GetIsolationCut()->GetConeSize();
914
915 Float_t coneptsum = 0 ;
916 Int_t n[10][10];//[fNCones][fNPtThresFrac];
917 Int_t nfrac[10][10];//[fNCones][fNPtThresFrac];
918 Bool_t isolated = kFALSE;
919
920 //Loop on cone sizes
921 for(Int_t icone = 0; icone<fNCones; icone++){
922 GetIsolationCut()->SetConeSize(fConeSizes[icone]);
923 coneptsum = 0 ;
a3aebfff 924
477d6cee 925 //Loop on ptthresholds
926 for(Int_t ipt = 0; ipt<fNPtThresFrac ;ipt++){
927 n[icone][ipt]=0;
928 nfrac[icone][ipt]=0;
929 GetIsolationCut()->SetPtThreshold(fPtThresholds[ipt]);
591cc579 930 GetIsolationCut()->MakeIsolationCut(ph->GetObjArray(GetAODObjArrayName()+"Tracks"),
931 ph->GetObjArray(GetAODObjArrayName()+"Clusters"),
233e0df8 932 GetReader(), kFALSE, ph, "",n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated);
a3aebfff 933
477d6cee 934 //Normal ptThreshold cut
935 if(n[icone][ipt] == 0) {
936 fhPtThresIsolated[icone][ipt]->Fill(ptC);
937 if(IsDataMC()){
591cc579 938 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtThresIsolatedPrompt[icone][ipt]->Fill(ptC) ;
939 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtThresIsolatedConversion[icone][ipt]->Fill(ptC) ;
940 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
941 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtThresIsolatedPi0Decay[icone][ipt]->Fill(ptC) ;
942 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) || GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtThresIsolatedOtherDecay[icone][ipt]->Fill(ptC) ;
477d6cee 943 else fhPtThresIsolatedUnknown[icone][ipt]->Fill(ptC) ;
944 }
945 }
946
947 //Pt threshold on pt cand/ pt in cone fraction
948 if(nfrac[icone][ipt] == 0) {
949 fhPtFracIsolated[icone][ipt]->Fill(ptC);
950 if(IsDataMC()){
591cc579 951 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtFracIsolatedPrompt[icone][ipt]->Fill(ptC) ;
952 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtFracIsolatedConversion[icone][ipt]->Fill(ptC) ;
953 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC) ;
954 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtFracIsolatedPi0Decay[icone][ipt]->Fill(ptC) ;
955 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) || GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtFracIsolatedOtherDecay[icone][ipt]->Fill(ptC) ;
477d6cee 956 else fhPtFracIsolatedUnknown[icone][ipt]->Fill(ptC) ;
957 }
958 }
959 }//pt thresh loop
960
961 //Sum in cone histograms
962 fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ;
963 if(IsDataMC()){
591cc579 964 if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt)) fhPtSumIsolatedPrompt[icone]->Fill(ptC,coneptsum) ;
965 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) fhPtSumIsolatedConversion[icone]->Fill(ptC,coneptsum) ;
966 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation)) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ;
967 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)) fhPtSumIsolatedPi0Decay[icone]->Fill(ptC,coneptsum) ;
968 else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) || GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay)) fhPtSumIsolatedOtherDecay[icone]->Fill(ptC,coneptsum) ;
477d6cee 969 else fhPtSumIsolatedUnknown[icone]->Fill(ptC,coneptsum) ;
970 }
971
972 }//cone size loop
973
974 //Reset original parameters for AOD analysis
975 GetIsolationCut()->SetPtThreshold(ptthresorg);
976 GetIsolationCut()->SetPtFraction(ptfracorg);
977 GetIsolationCut()->SetConeSize(rorg);
978
979}
980
981//__________________________________________________________________
982void AliAnaParticleIsolation::Print(const Option_t * opt) const
983{
984
985 //Print some relevant parameters set for the analysis
986 if(! opt)
987 return;
988
989 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
a3aebfff 990 AliAnaPartCorrBaseClass::Print(" ");
477d6cee 991
992 printf("ReMake Isolation = %d \n", fReMakeIC) ;
993 printf("Make Several Isolation = %d \n", fMakeSeveralIC) ;
994 printf("Calorimeter for isolation = %s \n", fCalorimeter.Data()) ;
995
996 if(fMakeSeveralIC){
997 printf("N Cone Sizes = %d\n", fNCones) ;
998 printf("Cone Sizes = \n") ;
999 for(Int_t i = 0; i < fNCones; i++)
1000 printf(" %1.2f;", fConeSizes[i]) ;
1001 printf(" \n") ;
1002
1003 printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ;
1004 printf(" pT thresholds = \n") ;
1005 for(Int_t i = 0; i < fNPtThresFrac; i++)
1006 printf(" %2.2f;", fPtThresholds[i]) ;
1007
1008 printf(" \n") ;
1009
1010 printf(" pT fractions = \n") ;
1011 for(Int_t i = 0; i < fNPtThresFrac; i++)
1012 printf(" %2.2f;", fPtFractions[i]) ;
1013
1014 }
1015
1016 printf("Histograms: %3.1f < pT sum < %3.1f, Nbin = %d\n", fHistoPtSumMin, fHistoPtSumMax, fHistoNPtSumBins);
1017 printf("Histograms: %3.1f < pT in cone < %3.1f, Nbin = %d\n", fHistoPtInConeMin, fHistoPtInConeMax, fHistoNPtInConeBins);
1018
1019 printf(" \n") ;
1020
1021}
1022