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