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