]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrDep/AliAnaParticleHadronCorrelation.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Update of the Azimuthal isotropic expansion analysis by C.Andrei
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaParticleHadronCorrelation.cxx
CommitLineData
1c5acb87 1/**************************************************************************
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 is 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: $ */
16
17//_________________________________________________________________________
18// Class for the analysis of particle - hadron correlations
19// Particle (for example direct gamma) must be found in a previous analysis
20//-- Author: Gustavo Conesa (LNF-INFN)
dd4871cd 21
22// Modified by Yaxian Mao:
23// 1. add the UE subtraction for corrlation study
24// 2. change the correlation variable
1c5acb87 25//////////////////////////////////////////////////////////////////////////////
26
27
28// --- ROOT system ---
29#include "TH2F.h"
477d6cee 30#include "TClonesArray.h"
9415d854 31#include "TClass.h"
dd4871cd 32#include "TMath.h"
1c5acb87 33
34//---- ANALYSIS system ----
1c5acb87 35#include "AliNeutralMesonSelection.h"
36#include "AliAnaParticleHadronCorrelation.h"
37#include "AliCaloTrackReader.h"
38#include "AliCaloPID.h"
39#include "AliAODPWG4ParticleCorrelation.h"
ff45398a 40#include "AliFiducialCut.h"
477d6cee 41#include "AliAODTrack.h"
42#include "AliAODCaloCluster.h"
9415d854 43#include "AliMCAnalysisUtils.h"
44#include "TParticle.h"
45#include "AliStack.h"
591cc579 46#include "AliAODMCParticle.h"
1c5acb87 47
48ClassImp(AliAnaParticleHadronCorrelation)
49
50
51//____________________________________________________________________________
dde5a268 52 AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation():
1c5acb87 53 AliAnaPartCorrBaseClass(),
dde5a268 54 fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fSelectIsolated(0),
55 fMakeSeveralUE(0), fUeDeltaPhiMaxCut(0.), fUeDeltaPhiMinCut(0.),
477d6cee 56 fhPhiCharged(0), fhPhiNeutral(0), fhEtaCharged(0), fhEtaNeutral(0),
1c5acb87 57 fhDeltaPhiCharged(0), fhDeltaPhiNeutral(0),
58 fhDeltaEtaCharged(0), fhDeltaEtaNeutral(0),
123fc3bd 59 fhDeltaPhiChargedPt(0), fhDeltaPhiNeutralPt(0),
60 fhDeltaPhiUeChargedPt(0), fhDeltaPhiUeNeutralPt(0),
dde5a268 61 fhPtImbalanceNeutral(0), fhPtImbalanceCharged(0),
123fc3bd 62 fhPtImbalanceUeCharged(0),fhPtImbalanceUeNeutral(0),
dd4871cd 63 fhPtHbpCharged(0), fhPtHbpUeCharged(0),
64 fhPtHbpNeutral(0), fhPtHbpUeNeutral(0),
123fc3bd 65 fhDeltaPhiUeLeftCharged(0),fhDeltaPhiUeRightCharged(0),
66 fhDeltaPhiUeLeftNeutral(0),fhDeltaPhiUeRightNeutral(0),
67 fhPtImbalanceUeLeftCharged(0),fhPtImbalanceUeRightCharged(0),
dd4871cd 68 fhPtImbalanceUeLeftNeutral(0),fhPtImbalanceUeRightNeutral(0),
69 fhPtHbpUeLeftCharged(0),fhPtHbpUeRightCharged(0),
70 fhPtHbpUeLeftNeutral(0),fhPtHbpUeRightNeutral(0)
1c5acb87 71{
72 //Default Ctor
477d6cee 73
1c5acb87 74 //Initialize parameters
75 InitParameters();
76}
77
78//____________________________________________________________________________
dde5a268 79AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation(const AliAnaParticleHadronCorrelation & g):
1c5acb87 80 AliAnaPartCorrBaseClass(g),
81 fDeltaPhiMaxCut(g.fDeltaPhiMaxCut), fDeltaPhiMinCut(g.fDeltaPhiMinCut),
6639984f 82 fSelectIsolated(g.fSelectIsolated),
dde5a268 83 fMakeSeveralUE(g.fMakeSeveralUE), fUeDeltaPhiMaxCut(g.fUeDeltaPhiMaxCut),
84 fUeDeltaPhiMinCut(g.fUeDeltaPhiMinCut),
1c5acb87 85 fhPhiCharged(g.fhPhiCharged), fhPhiNeutral(g.fhPhiNeutral),
86 fhEtaCharged(g.fhEtaCharged), fhEtaNeutral(g.fhEtaNeutral),
87 fhDeltaPhiCharged(g.fhDeltaPhiCharged),
88 fhDeltaPhiNeutral(g.fhDeltaPhiNeutral),
89 fhDeltaEtaCharged(g.fhDeltaEtaCharged),
90 fhDeltaEtaNeutral(g.fhDeltaEtaNeutral),
dde5a268 91 fhDeltaPhiChargedPt(g.fhDeltaPhiChargedPt),
1c5acb87 92 fhDeltaPhiNeutralPt(g.fhDeltaPhiNeutralPt),
dde5a268 93 fhDeltaPhiUeChargedPt(g.fhDeltaPhiUeChargedPt),
123fc3bd 94 fhDeltaPhiUeNeutralPt(g.fhDeltaPhiUeNeutralPt),
1c5acb87 95 fhPtImbalanceNeutral(g.fhPtImbalanceNeutral),
dde5a268 96 fhPtImbalanceCharged(g.fhPtImbalanceCharged),
97 fhPtImbalanceUeCharged(g.fhPtImbalanceUeCharged),
123fc3bd 98 fhPtImbalanceUeNeutral(g.fhPtImbalanceUeNeutral),
dd4871cd 99 fhPtHbpCharged(g.fhPtHbpCharged),
100 fhPtHbpUeCharged(g.fhPtHbpUeCharged),
101 fhPtHbpNeutral(g.fhPtHbpNeutral),
102 fhPtHbpUeNeutral(g.fhPtHbpUeNeutral),
123fc3bd 103 fhDeltaPhiUeLeftCharged(g.fhDeltaPhiUeLeftCharged),
104 fhDeltaPhiUeRightCharged(g.fhDeltaPhiUeRightCharged),
105 fhDeltaPhiUeLeftNeutral(g.fhDeltaPhiUeLeftNeutral),
106 fhDeltaPhiUeRightNeutral(g.fhDeltaPhiUeRightNeutral),
dde5a268 107 fhPtImbalanceUeLeftCharged(g.fhPtImbalanceUeLeftCharged),
123fc3bd 108 fhPtImbalanceUeRightCharged(g.fhPtImbalanceUeRightCharged),
109 fhPtImbalanceUeLeftNeutral(g.fhPtImbalanceUeLeftNeutral),
dd4871cd 110 fhPtImbalanceUeRightNeutral(g.fhPtImbalanceUeRightNeutral),
111 fhPtHbpUeLeftCharged(g.fhPtHbpUeLeftCharged),
112 fhPtHbpUeRightCharged(g.fhPtHbpUeRightCharged),
113 fhPtHbpUeLeftNeutral(g.fhPtHbpUeLeftNeutral),
114 fhPtHbpUeRightNeutral(g.fhPtHbpUeRightNeutral)
1c5acb87 115{
116 // cpy ctor
477d6cee 117
1c5acb87 118}
119
120//_________________________________________________________________________
121AliAnaParticleHadronCorrelation & AliAnaParticleHadronCorrelation::operator = (const AliAnaParticleHadronCorrelation & source)
122{
123 // assignment operator
477d6cee 124
1c5acb87 125 if(this == &source)return *this;
126 ((AliAnaPartCorrBaseClass *)this)->operator=(source);
127
128 fDeltaPhiMaxCut = source.fDeltaPhiMaxCut ;
129 fDeltaPhiMinCut = source.fDeltaPhiMinCut ;
6639984f 130 fSelectIsolated = source.fSelectIsolated ;
dde5a268 131 fMakeSeveralUE = source.fMakeSeveralUE ;
132 fUeDeltaPhiMaxCut = source.fUeDeltaPhiMaxCut ;
133 fUeDeltaPhiMinCut = source.fUeDeltaPhiMinCut ;
1c5acb87 134 fhPhiCharged = source.fhPhiCharged ; fhPhiNeutral = source.fhPhiNeutral ;
135 fhEtaCharged = source.fhEtaCharged ; fhEtaNeutral = source.fhEtaNeutral ;
136 fhDeltaPhiCharged = source.fhDeltaPhiCharged ;
137 fhDeltaPhiNeutral = source.fhDeltaPhiNeutral ;
1c5acb87 138 fhDeltaEtaCharged = source.fhDeltaEtaCharged ;
139 fhDeltaEtaNeutral = source.fhDeltaEtaNeutral ;
140 fhDeltaPhiChargedPt = source.fhDeltaPhiChargedPt ;
dde5a268 141 fhDeltaPhiNeutralPt = source.fhDeltaPhiNeutralPt ;
142 fhDeltaPhiUeChargedPt = source.fhDeltaPhiUeChargedPt ;
123fc3bd 143 fhDeltaPhiUeNeutralPt = source.fhDeltaPhiUeNeutralPt ;
1c5acb87 144 fhPtImbalanceNeutral = source.fhPtImbalanceNeutral ;
145 fhPtImbalanceCharged = source.fhPtImbalanceCharged ;
dd4871cd 146 fhPtImbalanceUeCharged = source.fhPtImbalanceUeCharged ;
123fc3bd 147 fhPtImbalanceUeNeutral = source.fhPtImbalanceUeNeutral ;
dd4871cd 148 fhPtHbpCharged = source.fhPtHbpCharged ;
149 fhPtHbpUeCharged = source.fhPtHbpUeCharged;
150 fhPtHbpNeutral = source.fhPtHbpNeutral ;
151 fhPtHbpUeNeutral = source.fhPtHbpUeNeutral;
123fc3bd 152 fhDeltaPhiUeLeftCharged = source.fhDeltaPhiUeLeftCharged ;
153 fhDeltaPhiUeRightCharged = source.fhDeltaPhiUeRightCharged ;
154 fhDeltaPhiUeLeftNeutral = source.fhDeltaPhiUeLeftNeutral ;
155 fhDeltaPhiUeRightNeutral = source.fhDeltaPhiUeRightNeutral ;
dde5a268 156 fhPtImbalanceUeLeftCharged = source.fhPtImbalanceUeLeftCharged ;
157 fhPtImbalanceUeRightCharged = source.fhPtImbalanceUeRightCharged ;
123fc3bd 158 fhPtImbalanceUeLeftNeutral = source.fhPtImbalanceUeLeftNeutral ;
159 fhPtImbalanceUeRightNeutral = source.fhPtImbalanceUeRightNeutral ;
dd4871cd 160 fhPtHbpUeLeftCharged = source.fhPtHbpUeLeftCharged ;
161 fhPtHbpUeRightCharged = source.fhPtHbpUeRightCharged ;
162 fhPtHbpUeLeftNeutral = source.fhPtHbpUeLeftNeutral ;
163 fhPtHbpUeRightNeutral = source.fhPtHbpUeRightNeutral ;
164
1c5acb87 165 return *this;
166
167}
168
169//________________________________________________________________________
170TList * AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
171{
477d6cee 172
173 // Create histograms to be saved in output file and
174 // store them in fOutputContainer
175 TList * outputContainer = new TList() ;
176 outputContainer->SetName("CorrelationHistos") ;
177
178 Int_t nptbins = GetHistoNPtBins();
179 Int_t nphibins = GetHistoNPhiBins();
180 Int_t netabins = GetHistoNEtaBins();
181 Float_t ptmax = GetHistoPtMax();
182 Float_t phimax = GetHistoPhiMax();
183 Float_t etamax = GetHistoEtaMax();
184 Float_t ptmin = GetHistoPtMin();
185 Float_t phimin = GetHistoPhiMin();
186 Float_t etamin = GetHistoEtaMin();
187
188 //Correlation with charged hadrons
189 if(GetReader()->IsCTSSwitchedOn()) {
190 fhPhiCharged = new TH2F
dd4871cd 191 ("PhiCharged","#phi_{h^{#pm}} vs p_{T #pm}",
477d6cee 192 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
193 fhPhiCharged->SetYTitle("#phi_{h^{#pm}} (rad)");
dd4871cd 194 fhPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
477d6cee 195
196 fhEtaCharged = new TH2F
dd4871cd 197 ("EtaCharged","#eta_{h^{#pm}} vs p_{T #pm}",
477d6cee 198 nptbins,ptmin,ptmax,netabins,etamin,etamax);
199 fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
dd4871cd 200 fhEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
477d6cee 201
202 fhDeltaPhiCharged = new TH2F
203 ("DeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
dde5a268 204 nptbins,ptmin,ptmax,700,-2.,5.);
477d6cee 205 fhDeltaPhiCharged->SetYTitle("#Delta #phi");
206 fhDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
207
208 fhDeltaPhiChargedPt = new TH2F
9415d854 209 ("DeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
dde5a268 210 nptbins,ptmin,ptmax,700,-2.,5.);
477d6cee 211 fhDeltaPhiChargedPt->SetYTitle("#Delta #phi");
212 fhDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
dde5a268 213
214 fhDeltaPhiUeChargedPt = new TH2F
215 ("DeltaPhiUeChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}}",
216 nptbins,ptmin,ptmax,700,-2.,5.);
217 fhDeltaPhiUeChargedPt->SetYTitle("#Delta #phi");
218 fhDeltaPhiUeChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
477d6cee 219
220 fhDeltaEtaCharged = new TH2F
221 ("DeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger}",
222 nptbins,ptmin,ptmax,200,-2,2);
223 fhDeltaEtaCharged->SetYTitle("#Delta #eta");
224 fhDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
225
226 fhPtImbalanceCharged =
227 new TH2F("CorrelationCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T trigger}",
dde5a268 228 nptbins,ptmin,ptmax,1000,0.,2.);
477d6cee 229 fhPtImbalanceCharged->SetYTitle("z_{trigger h^{#pm}}");
230 fhPtImbalanceCharged->SetXTitle("p_{T trigger}");
231
dde5a268 232 fhPtImbalanceUeCharged =
233 new TH2F("CorrelationUeCharged","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger}",
234 nptbins,ptmin,ptmax,1000,0.,2.);
235 fhPtImbalanceUeCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
236 fhPtImbalanceUeCharged->SetXTitle("p_{T trigger}");
237
dd4871cd 238 fhPtHbpCharged =
239 new TH2F("HbpCharged","#xi = ln(1/x_{E}) with charged hadrons",
240 nptbins,ptmin,ptmax,1000,0.,10.);
241 fhPtHbpCharged->SetYTitle("ln(1/x_{E})");
242 fhPtHbpCharged->SetXTitle("p_{T trigger}");
243
244 fhPtHbpUeCharged =
245 new TH2F("HbpUeCharged","#xi = ln(1/x_{E}) with charged hadrons",
246 nptbins,ptmin,ptmax,1000,0.,10.);
247 fhPtHbpUeCharged->SetYTitle("ln(1/x_{E})");
248 fhPtHbpUeCharged->SetXTitle("p_{T trigger}");
249
477d6cee 250 outputContainer->Add(fhPhiCharged) ;
251 outputContainer->Add(fhEtaCharged) ;
252 outputContainer->Add(fhDeltaPhiCharged) ;
253 outputContainer->Add(fhDeltaEtaCharged) ;
477d6cee 254 outputContainer->Add(fhDeltaPhiChargedPt) ;
dde5a268 255 outputContainer->Add(fhDeltaPhiUeChargedPt) ;
256 outputContainer->Add(fhPtImbalanceCharged) ;
257 outputContainer->Add(fhPtImbalanceUeCharged) ;
dd4871cd 258 outputContainer->Add(fhPtHbpCharged) ;
259 outputContainer->Add(fhPtHbpUeCharged) ;
dde5a268 260 if(fMakeSeveralUE){
123fc3bd 261 fhDeltaPhiUeLeftCharged = new TH2F
dde5a268 262 ("DeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE left side range of trigger particles",
263 nptbins,ptmin,ptmax,700,-2.,5.);
123fc3bd 264 fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #phi");
265 fhDeltaPhiUeLeftCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
266 outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
dde5a268 267
123fc3bd 268 fhDeltaPhiUeRightCharged = new TH2F
dde5a268 269 ("DeltaPhiUeRightChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE right side range of trigger particles",
270 nptbins,ptmin,ptmax,700,-2.,5.);
123fc3bd 271 fhDeltaPhiUeRightCharged->SetYTitle("#Delta #phi");
272 fhDeltaPhiUeRightCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
273 outputContainer->Add(fhDeltaPhiUeRightCharged) ;
dde5a268 274
275 fhPtImbalanceUeLeftCharged =
276 new TH2F("CorrelationUeChargedLeft","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE left side of trigger",
277 nptbins,ptmin,ptmax,1000,0.,2.);
278 fhPtImbalanceUeLeftCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
279 fhPtImbalanceUeLeftCharged->SetXTitle("p_{T trigger}");
280 outputContainer->Add(fhPtImbalanceUeLeftCharged) ;
281
282 fhPtImbalanceUeRightCharged =
283 new TH2F("CorrelationUeChargedRight","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE right side of trigger",
284 nptbins,ptmin,ptmax,1000,0.,2.);
285 fhPtImbalanceUeRightCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
286 fhPtImbalanceUeRightCharged->SetXTitle("p_{T trigger}");
287 outputContainer->Add(fhPtImbalanceUeRightCharged) ;
dd4871cd 288
289 fhPtHbpUeLeftCharged =
290 new TH2F("HbpUeChargedLeft","#xi = ln(1/x_{E}) with charged UE left side of trigger",
291 nptbins,ptmin,ptmax,1000,0.,10.);
292 fhPtHbpUeLeftCharged->SetYTitle("ln(1/x_{E})");
293 fhPtHbpUeLeftCharged->SetXTitle("p_{T trigger}");
294 outputContainer->Add(fhPtHbpUeLeftCharged) ;
295
296 fhPtHbpUeRightCharged =
297 new TH2F("HbpUeChargedRight","#xi = ln(1/x_{E}) with charged UE right side of trigger",
298 nptbins,ptmin,ptmax,1000,0.,10.);
299 fhPtHbpUeRightCharged->SetYTitle("ln(1/x_{E})");
300 fhPtHbpUeRightCharged->SetXTitle("p_{T trigger}");
301 outputContainer->Add(fhPtHbpUeRightCharged) ;
302
dde5a268 303 }
477d6cee 304 } //Correlation with charged hadrons
305
306 //Correlation with neutral hadrons
307 if(GetReader()->IsEMCALSwitchedOn() || GetReader()->IsPHOSSwitchedOn()){
308
309 fhPhiNeutral = new TH2F
dd4871cd 310 ("PhiNeutral","#phi_{#pi^{0}} vs p_{T #pi^{0}}",
477d6cee 311 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
312 fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
dd4871cd 313 fhPhiNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
477d6cee 314
315 fhEtaNeutral = new TH2F
dd4871cd 316 ("EtaNeutral","#eta_{#pi^{0}} vs p_{T #pi^{0}}",
477d6cee 317 nptbins,ptmin,ptmax,netabins,etamin,etamax);
318 fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
dd4871cd 319 fhEtaNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
477d6cee 320
321 fhDeltaPhiNeutral = new TH2F
322 ("DeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T trigger}",
323 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
324 fhDeltaPhiNeutral->SetYTitle("#Delta #phi");
325 fhDeltaPhiNeutral->SetXTitle("p_{T trigger} (GeV/c)");
326
327 fhDeltaPhiNeutralPt = new TH2F
328 ("DeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
dde5a268 329 nptbins,ptmin,ptmax,700,-2.,5.);
477d6cee 330 fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi");
123fc3bd 331 fhDeltaPhiNeutralPt->SetXTitle("p_{T h^{0} (GeV/c)");
332
333 fhDeltaPhiUeNeutralPt = new TH2F
334 ("DeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
335 nptbins,ptmin,ptmax,700,-2.,5.);
336 fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi");
337 fhDeltaPhiUeNeutralPt->SetXTitle("p_{T h^{0} (GeV/c)");
477d6cee 338
339 fhDeltaEtaNeutral = new TH2F
340 ("DeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs p_{T trigger}",
341 nptbins,ptmin,ptmax,200,-2,2);
342 fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
343 fhDeltaEtaNeutral->SetXTitle("p_{T trigger} (GeV/c)");
344
345 fhPtImbalanceNeutral =
346 new TH2F("CorrelationNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
dde5a268 347 nptbins,ptmin,ptmax,1000,0.,2.);
477d6cee 348 fhPtImbalanceNeutral->SetYTitle("z_{trigger #pi^{0}}");
349 fhPtImbalanceNeutral->SetXTitle("p_{T trigger}");
123fc3bd 350
351 fhPtImbalanceUeNeutral =
352 new TH2F("CorrelationUeNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
353 nptbins,ptmin,ptmax,1000,0.,2.);
354 fhPtImbalanceUeNeutral->SetYTitle("z_{trigger #pi^{0}}");
355 fhPtImbalanceUeNeutral->SetXTitle("p_{T trigger}");
dd4871cd 356
357 fhPtHbpNeutral =
358 new TH2F("HbpNeutral","#xi = ln(1/x_{E}) with neutral particles",
359 nptbins,ptmin,ptmax,1000,0.,10.);
360 fhPtHbpNeutral->SetYTitle("ln(1/x_{E})");
361 fhPtHbpNeutral->SetXTitle("p_{T trigger}");
362
363 fhPtHbpUeNeutral =
364 new TH2F("HbpUeNeutral","#xi = ln(1/x_{E}) with neutral particles",
365 nptbins,ptmin,ptmax,1000,0.,10.);
366 fhPtHbpUeNeutral->SetYTitle("ln(1/x_{E})");
367 fhPtHbpUeNeutral->SetXTitle("p_{T trigger}");
368
369
123fc3bd 370
477d6cee 371 outputContainer->Add(fhPhiNeutral) ;
372 outputContainer->Add(fhEtaNeutral) ;
373 outputContainer->Add(fhDeltaPhiNeutral) ;
123fc3bd 374 outputContainer->Add(fhDeltaPhiNeutralPt) ;
375 outputContainer->Add(fhDeltaPhiUeNeutralPt) ;
477d6cee 376 outputContainer->Add(fhDeltaEtaNeutral) ;
377 outputContainer->Add(fhPtImbalanceNeutral) ;
dd4871cd 378 outputContainer->Add(fhPtImbalanceUeNeutral) ;
379 outputContainer->Add(fhPtHbpNeutral) ;
380 outputContainer->Add(fhPtHbpUeNeutral) ;
123fc3bd 381
382 if(fMakeSeveralUE){
383 fhDeltaPhiUeLeftNeutral = new TH2F
dd4871cd 384 ("DeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T h^{0}} with neutral UE left side range of trigger particles",
123fc3bd 385 nptbins,ptmin,ptmax,700,-2.,5.);
386 fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #phi");
387 fhDeltaPhiUeLeftNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
388 outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
389
390 fhDeltaPhiUeRightNeutral = new TH2F
dd4871cd 391 ("DeltaPhiUeRightNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T Ueh^{0}} with neutral UE right side range of trigger particles",
123fc3bd 392 nptbins,ptmin,ptmax,700,-2.,5.);
393 fhDeltaPhiUeRightNeutral->SetYTitle("#Delta #phi");
394 fhDeltaPhiUeRightNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
395 outputContainer->Add(fhDeltaPhiUeRightNeutral) ;
396
397 fhPtImbalanceUeLeftNeutral =
dd4871cd 398 new TH2F("CorrelationUeNeutralLeft","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE left side of trigger",
123fc3bd 399 nptbins,ptmin,ptmax,1000,0.,2.);
400 fhPtImbalanceUeLeftNeutral->SetYTitle("z_{trigger Ueh^{0}}");
401 fhPtImbalanceUeLeftNeutral->SetXTitle("p_{T trigger}");
402 outputContainer->Add(fhPtImbalanceUeLeftNeutral) ;
403
404 fhPtImbalanceUeRightNeutral =
dd4871cd 405 new TH2F("CorrelationUeNeutralRight","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE right side of trigger",
123fc3bd 406 nptbins,ptmin,ptmax,1000,0.,2.);
407 fhPtImbalanceUeRightNeutral->SetYTitle("z_{trigger Ueh^{0}}");
408 fhPtImbalanceUeRightNeutral->SetXTitle("p_{T trigger}");
409 outputContainer->Add(fhPtImbalanceUeRightNeutral) ;
dd4871cd 410
411 fhPtHbpUeLeftNeutral =
412 new TH2F("HbpUeNeutralLeft","#xi = ln(1/x_{E}) with neutral UE left side of trigger",
413 nptbins,ptmin,ptmax,1000,0.,10.);
414 fhPtHbpUeLeftNeutral->SetYTitle("ln(1/x_{E})");
415 fhPtHbpUeLeftNeutral->SetXTitle("p_{T trigger}");
416 outputContainer->Add(fhPtHbpUeLeftNeutral) ;
417
418 fhPtHbpUeRightNeutral =
419 new TH2F("HbpUeNeutralRight","#xi = ln(1/x_{E}) with neutral UE right side of trigger",
420 nptbins,ptmin,ptmax,1000,0.,10.);
421 fhPtHbpUeRightNeutral->SetYTitle("ln(1/x_{E})");
422 fhPtHbpUeRightNeutral->SetXTitle("p_{T trigger}");
423 outputContainer->Add(fhPtHbpUeRightNeutral) ;
424
123fc3bd 425 }
426
477d6cee 427 //Keep neutral meson selection histograms if requiered
428 //Setting done in AliNeutralMesonSelection
429 if(GetNeutralMesonSelection()){
430 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
431 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
432 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
433 }
a3aebfff 434
477d6cee 435 }//Correlation with neutral hadrons
436
437 return outputContainer;
a3aebfff 438
1c5acb87 439}
440
477d6cee 441//____________________________________________________________________________
1c5acb87 442void AliAnaParticleHadronCorrelation::InitParameters()
443{
477d6cee 444
1c5acb87 445 //Initialize the parameters of the analysis.
a3aebfff 446 SetInputAODName("PWG4Particle");
591cc579 447 SetAODObjArrayName("Hadrons");
a3aebfff 448 AddToHistogramsName("AnaHadronCorr_");
449
9415d854 450 //Correlation with neutrals
451 //SetOutputAODClassName("AliAODPWG4Particle");
452 //SetOutputAODName("Pi0Correlated");
453
dde5a268 454 SetPtCutRange(0.,300);
1c5acb87 455 fDeltaPhiMinCut = 1.5 ;
456 fDeltaPhiMaxCut = 4.5 ;
6639984f 457 fSelectIsolated = kFALSE;
dde5a268 458 fMakeSeveralUE = kFALSE;
459 fUeDeltaPhiMinCut = 1. ;
460 fUeDeltaPhiMaxCut = 1.5 ;
1c5acb87 461}
462
463//__________________________________________________________________
464void AliAnaParticleHadronCorrelation::Print(const Option_t * opt) const
465{
466
467 //Print some relevant parameters set for the analysis
468 if(! opt)
469 return;
470
a3aebfff 471 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
472 AliAnaPartCorrBaseClass::Print(" ");
473
1c5acb87 474 printf("Phi trigger particle-Hadron < %3.2f\n", fDeltaPhiMaxCut) ;
475 printf("Phi trigger particle-Hadron > %3.2f\n", fDeltaPhiMinCut) ;
6639984f 476 printf("Isolated Trigger? %d\n", fSelectIsolated) ;
dde5a268 477 printf("Phi trigger particle-UeHadron < %3.2f\n", fUeDeltaPhiMaxCut) ;
478 printf("Phi trigger particle-UeHadron > %3.2f\n", fUeDeltaPhiMinCut) ;
479 printf("Several UE? %d\n", fMakeSeveralUE) ;
1c5acb87 480}
481
482//____________________________________________________________________________
483void AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD()
484{
477d6cee 485 //Particle-Hadron Correlation Analysis, fill AODs
486
487 if(!GetInputAODBranch()){
591cc579 488 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 489 abort();
490 }
9415d854 491
492 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
493 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Wrong type of AOD object, change AOD class name in input AOD: It should be <AliAODPWG4ParticleCorrelation> and not <%s> \n",GetInputAODBranch()->GetClass()->GetName());
494 abort();
495 }
496
477d6cee 497 if(GetDebug() > 1){
a3aebfff 498 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Begin hadron correlation analysis, fill AODs \n");
499 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
500 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetAODCTS()->GetEntriesFast());
501 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetAODEMCAL()->GetEntriesFast());
502 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n", GetAODPHOS()->GetEntriesFast());
477d6cee 503 }
504
505 //Loop on stored AOD particles, trigger
506 Int_t naod = GetInputAODBranch()->GetEntriesFast();
507 for(Int_t iaod = 0; iaod < naod ; iaod++){
508 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
9415d854 509
510 //Make correlation with charged hadrons
477d6cee 511 if(GetReader()->IsCTSSwitchedOn() )
a3aebfff 512 MakeChargedCorrelation(particle, GetAODCTS(),kFALSE);
477d6cee 513
514 //Make correlation with neutral pions
515 //Trigger particle in PHOS, correlation with EMCAL
516 if(particle->GetDetector()=="PHOS" && GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
9415d854 517 MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
477d6cee 518 //Trigger particle in EMCAL, correlation with PHOS
519 else if(particle->GetDetector()=="EMCAL" && GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
9415d854 520 MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
477d6cee 521 //Trigger particle in CTS, correlation with PHOS, EMCAL and CTS
522 else if(particle->GetDetector()=="CTS" ){
523 if(GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
9415d854 524 MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
477d6cee 525 if(GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
9415d854 526 MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
477d6cee 527 }
a3aebfff 528
529
477d6cee 530 }//Aod branch loop
531
a3aebfff 532 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
477d6cee 533
1c5acb87 534}
535
536//____________________________________________________________________________
537void AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms()
538{
477d6cee 539 //Particle-Hadron Correlation Analysis, fill histograms
540
541 if(!GetInputAODBranch()){
591cc579 542 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 543 abort();
544 }
9415d854 545
477d6cee 546 if(GetDebug() > 1){
a3aebfff 547 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Begin hadron correlation analysis, fill histograms \n");
548 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
477d6cee 549 }
9415d854 550
477d6cee 551 //Loop on stored AOD particles
552 Int_t naod = GetInputAODBranch()->GetEntriesFast();
9415d854 553 for(Int_t iaod = 0; iaod < naod ; iaod++){
554 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
555
556 //check if the particle is isolated or if we want to take the isolation into account
557 if(OnlyIsolated() && !particle->IsIsolated()) continue;
558
559 //Make correlation with charged hadrons
591cc579 560 TObjArray * reftracks = particle->GetObjArray(GetAODObjArrayName()+"Tracks");
9415d854 561 if(reftracks){
562 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Track Refs entries %d\n", iaod, reftracks->GetEntriesFast());
563 if(reftracks->GetEntriesFast() > 0) MakeChargedCorrelation(particle, reftracks,kTRUE);
564 }
565
566 //Make correlation with neutral pions
567 if(GetOutputAODBranch() && GetOutputAODBranch()->GetEntriesFast() > 0){
568 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Cluster Refs entries %d\n",iaod, GetOutputAODBranch()->GetEntriesFast());
569 MakeNeutralCorrelationFillHistograms(particle);
570 }
571
477d6cee 572 }//Aod branch loop
573
a3aebfff 574 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
477d6cee 575
1c5acb87 576}
577
578//____________________________________________________________________________
233e0df8 579void AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle, TObjArray* const pl, const Bool_t bFillHisto)
1c5acb87 580{
dde5a268 581 // Charged Hadron Correlation Analysis
582 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
583
584 Double_t ptTrig = aodParticle->Pt();
dd4871cd 585 Double_t pxTrig = aodParticle->Px();
586 Double_t pyTrig = aodParticle->Py();
587
dde5a268 588 Double_t phiTrig = aodParticle->Phi();
dd4871cd 589
dde5a268 590 Double_t pt = -100.;
dd4871cd 591 Double_t px = -100.;
592 Double_t py = -100.;
dde5a268 593 Double_t rat = -100.;
dd4871cd 594 Double_t xE = -100.;
595 Double_t cosi = -100.;
dde5a268 596 Double_t phi = -100. ;
597 Double_t eta = -100. ;
598 Double_t p[3];
dde5a268 599
591cc579 600 TObjArray * reftracks =0x0;
dde5a268 601 if(!bFillHisto)
591cc579 602 reftracks = new TObjArray;
dde5a268 603
604 //Track loop, select tracks with good pt, phi and fill AODs or histograms
605 for(Int_t ipr = 0;ipr < pl->GetEntries() ; ipr ++ ){
606 AliAODTrack * track = (AliAODTrack *) (pl->At(ipr)) ;
607 track->GetPxPyPz(p) ;
608 TLorentzVector mom(p[0],p[1],p[2],0);
609 pt = mom.Pt();
dd4871cd 610 px = mom.Px();
611 py = mom.Py();
dde5a268 612 eta = mom.Eta();
613 phi = mom.Phi() ;
614 if(phi < 0) phi+=TMath::TwoPi();
615 rat = pt/ptTrig ;
dd4871cd 616 xE = -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
617 if(xE <0.)xE =-xE;
618 cosi = TMath::Log(1/xE);
619 // printf("rat = %f, xE = %f, cosi =%f \n", rat, xE, cosi);
620 // printf("phi = %f \n", phi);
dde5a268 621
ff45398a 622 if(IsFiducialCutOn()){
623 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,"CTS") ;
dde5a268 624 if(! in ) continue ;
625 }
626
627 //Select only hadrons in pt range
628 if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
629
630 //Selection within angular range
631 Float_t deltaphi = phiTrig-phi;
632 if(deltaphi<-TMath::PiOver2()) deltaphi+=TMath::TwoPi();
633 if(deltaphi>3*TMath::PiOver2()) deltaphi-=TMath::TwoPi();
634
635 if(GetDebug() > 2)
636 printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Charged hadron: pt %f, phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
637 pt,phi, phiTrig,fDeltaPhiMinCut, deltaphi, fDeltaPhiMaxCut, GetMinPt());
638
639 if(bFillHisto){
640 // Fill Histograms
dd4871cd 641 fhEtaCharged->Fill(pt,eta);
642 fhPhiCharged->Fill(pt,phi);
dde5a268 643 fhDeltaEtaCharged->Fill(ptTrig,aodParticle->Eta()-eta);
644 fhDeltaPhiCharged->Fill(ptTrig, deltaphi);
645
646 if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
123fc3bd 647
648 //delta phi cut for correlation
649 if( (deltaphi > fDeltaPhiMinCut) && ( deltaphi < fDeltaPhiMaxCut) ) {
650 fhDeltaPhiChargedPt->Fill(pt,deltaphi);
651 fhPtImbalanceCharged->Fill(ptTrig,rat);
dd4871cd 652 fhPtHbpCharged->Fill(ptTrig,cosi);
123fc3bd 653 }
654 else {
dd4871cd 655 fhDeltaPhiUeChargedPt->Fill(pt,deltaphi);
656 fhPtImbalanceUeCharged->Fill(ptTrig,rat);
657 fhPtHbpUeCharged->Fill(ptTrig,cosi);
123fc3bd 658 }
dde5a268 659 //several UE calculation
660 if(fMakeSeveralUE){
661 if((deltaphi<-fUeDeltaPhiMinCut) && (deltaphi >-fUeDeltaPhiMaxCut)){
123fc3bd 662 fhDeltaPhiUeLeftCharged->Fill(pt,deltaphi);
dde5a268 663 fhPtImbalanceUeLeftCharged->Fill(ptTrig,rat);
dd4871cd 664 fhPtHbpUeLeftCharged->Fill(ptTrig,cosi);
dde5a268 665 }
666 if((deltaphi>fUeDeltaPhiMinCut) && (deltaphi <fUeDeltaPhiMaxCut)){
123fc3bd 667 fhDeltaPhiUeRightCharged->Fill(pt,deltaphi);
dde5a268 668 fhPtImbalanceUeRightCharged->Fill(ptTrig,rat);
dd4871cd 669 fhPtHbpUeRightCharged->Fill(ptTrig,cosi);
670
dde5a268 671 }
672 } //several UE calculation
673
dde5a268 674 }
675 else{
dde5a268 676 reftracks->Add(track);
dde5a268 677 }//aod particle loop
678 }// track loop
9415d854 679
dde5a268 680 //Fill AOD with reference tracks, if not filling histograms
681 if(!bFillHisto && reftracks->GetEntriesFast() > 0) {
591cc579 682 reftracks->SetName(GetAODObjArrayName()+"Tracks");
683 aodParticle->AddObjArray(reftracks);
dde5a268 684 }
9415d854 685
1c5acb87 686}
9415d854 687
1c5acb87 688//____________________________________________________________________________
233e0df8 689void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD(AliAODPWG4ParticleCorrelation* const aodParticle,TObjArray* const pl, TString detector)
1c5acb87 690{
9415d854 691 // Neutral Pion Correlation Analysis, find pi0, put them in new output aod, if correlation cuts passed
692 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Make trigger particle - neutral hadron correlation \n");
693
694 if(!NewOutputAOD()){
591cc579 695 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Output aod not created, set AOD class name and branch name in the configuration file, STOP! \n");
9415d854 696 abort();
697 }
1c5acb87 698
1c5acb87 699 Double_t phiTrig = aodParticle->Phi();
591cc579 700 Int_t tag = 0;
1c5acb87 701 TLorentzVector gammai;
702 TLorentzVector gammaj;
703
233e0df8 704 //Get vertex for photon momentum calculation
705 Double_t vertex[] = {0,0,0} ; //vertex
706 Double_t vertex2[] = {0,0,0} ; //vertex of second input aod
edd59991 707 if(GetReader()->GetDataType()!= AliCaloTrackReader::kMC)
233e0df8 708 {
709 GetReader()->GetVertex(vertex);
710 if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex2);
711 }
712
1c5acb87 713 //Cluster loop, select pairs with good pt, phi and fill AODs or histograms
9415d854 714 //Int_t iEvent= GetReader()->GetEventNumber() ;
715 Int_t nclus = pl->GetEntriesFast();
716 for(Int_t iclus = 0;iclus < nclus ; iclus ++ ){
1c5acb87 717 AliAODCaloCluster * calo = (AliAODCaloCluster *) (pl->At(iclus)) ;
9415d854 718
233e0df8 719 //Input from second AOD?
720 Int_t inputi = 0;
721 if (aodParticle->GetDetector() == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= iclus) inputi = 1 ;
722 else if(aodParticle->GetDetector() == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= iclus) inputi = 1;
723
ff45398a 724 //Cluster selection, not charged, with photon or pi0 id and in fiducial cut
9415d854 725 Int_t pdg=0;
233e0df8 726 if (inputi == 0 && !SelectCluster(calo, vertex, gammai, pdg)) continue ;
727 else if(inputi == 1 && !SelectCluster(calo, vertex2, gammai, pdg)) continue ;
12524a23 728
729 if(GetDebug() > 2)
9415d854 730 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Neutral cluster in %s: pt %f, phi %f, phi trigger %f. Cuts: delta phi min %2.2f, max %2.2f, pT min %2.2f \n",
731 detector.Data(), gammai.Pt(),gammai.Phi(),phiTrig,fDeltaPhiMinCut, fDeltaPhiMaxCut, GetMinPt());
732
733 //2 gamma overlapped, found with PID
734 if(pdg == AliCaloPID::kPi0){
477d6cee 735
9415d854 736 //Select only hadrons in pt range
737 if(gammai.Pt() < GetMinPt() || gammai.Pt() > GetMaxPt()) continue ;
1c5acb87 738
9415d854 739 //Selection within angular range
740 Float_t phi = gammai.Phi();
741 if(phi < 0) phi+=TMath::TwoPi();
dde5a268 742 //Float_t deltaphi = TMath::Abs(phiTrig-phi);
743 //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
477d6cee 744
9415d854 745 AliAODPWG4Particle pi0 = AliAODPWG4Particle(gammai);
746 //pi0.SetLabel(calo->GetLabel(0));
747 pi0.SetPdg(AliCaloPID::kPi0);
748 pi0.SetDetector(detector);
749
750 if(IsDataMC()){
233e0df8 751 pi0.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabel(0),GetReader(),inputi));
752 if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of candidate %d\n",pi0.GetTag());
753 }//Work with stack also
9415d854 754 //Set the indeces of the original caloclusters
755 pi0.SetCaloLabel(calo->GetID(),-1);
756 AddAODParticle(pi0);
757
758 if(GetDebug() > 2)
759 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Correlated with selected pi0 (pid): pt %f, phi %f\n",pi0.Pt(),pi0.Phi());
760
761 }// pdg = 111
762
763 //Make invariant mass analysis
764 else if(pdg == AliCaloPID::kPhoton){
765 //Search the photon companion in case it comes from a Pi0 decay
766 //Apply several cuts to select the good pair;
767 for(Int_t jclus = iclus+1; jclus < pl->GetEntries() ; jclus ++ ){
768 AliAODCaloCluster * calo2 = (AliAODCaloCluster *) (pl->At(jclus)) ;
477d6cee 769
233e0df8 770 //Input from second AOD?
771 Int_t inputj = 0;
772 if (aodParticle->GetDetector() == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= jclus) inputj = 1;
773 else if(aodParticle->GetDetector() == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= jclus) inputj = 1;
774
ff45398a 775 //Cluster selection, not charged with photon or pi0 id and in fiducial cut
9415d854 776 Int_t pdgj=0;
233e0df8 777 if (inputj == 0 && !SelectCluster(calo2, vertex, gammaj, pdgj)) continue ;
778 else if(inputj == 1 && !SelectCluster(calo2, vertex2, gammaj, pdgj)) continue ;
779
9415d854 780 if(!SelectCluster(calo2,vertex, gammaj, pdgj)) continue ;
477d6cee 781
9415d854 782 if(pdgj == AliCaloPID::kPhoton ){
783
784 if((gammai+gammaj).Pt() < GetMinPt() || (gammai+gammaj).Pt() > GetMaxPt()) continue ;
1c5acb87 785
9415d854 786 //Selection within angular range
787 Float_t phi = (gammai+gammaj).Phi();
788 if(phi < 0) phi+=TMath::TwoPi();
dde5a268 789 //Float_t deltaphi = TMath::Abs(phiTrig-phi);
790 //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
477d6cee 791
9415d854 792 //Select good pair (aperture and invariant mass)
793 if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
1c5acb87 794
9415d854 795 if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Neutral Hadron Correlation: AOD Selected gamma pair: pt %2.2f, phi %2.2f, eta %2.2f, M %2.3f\n",
796 (gammai+gammaj).Pt(),(gammai+gammaj).Phi(),(gammai+gammaj).Eta(), (gammai+gammaj).M());
1c5acb87 797
9415d854 798 TLorentzVector pi0mom = gammai+gammaj;
799 AliAODPWG4Particle pi0 = AliAODPWG4Particle(pi0mom);
800 //pi0.SetLabel(calo->GetLabel(0));
801 pi0.SetPdg(AliCaloPID::kPi0);
802 pi0.SetDetector(detector);
803 if(IsDataMC()){
804 //Check origin of the candidates
591cc579 805
806 Int_t label1 = calo->GetLabel(0);
807 Int_t label2 = calo2->GetLabel(0);
233e0df8 808 Int_t tag1 = GetMCAnalysisUtils()->CheckOrigin(label1, GetReader(), inputi);
809 Int_t tag2 = GetMCAnalysisUtils()->CheckOrigin(label2, GetReader(), inputj);
591cc579 810
9415d854 811 if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
591cc579 812 if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay) && GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0Decay)){
9415d854 813
591cc579 814 //Check if pi0 mother is the same
815 if(GetReader()->ReadStack()){
816 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
817 label1 = mother1->GetFirstMother();
818 //mother1 = GetMCStack()->Particle(label1);//pi0
819
820 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
821 label2 = mother2->GetFirstMother();
822 //mother2 = GetMCStack()->Particle(label2);//pi0
823 }
824 else if(GetReader()->ReadAODMCParticles()){
233e0df8 825 AliAODMCParticle * mother1 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(inputi))->At(label1);//photon in kine tree
591cc579 826 label1 = mother1->GetMother();
827 //mother1 = GetMCStack()->Particle(label1);//pi0
233e0df8 828 AliAODMCParticle * mother2 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(inputj))->At(label2);//photon in kine tree
591cc579 829 label2 = mother2->GetMother();
830 //mother2 = GetMCStack()->Particle(label2);//pi0
831 }
9415d854 832
591cc579 833 //printf("mother1 %d, mother2 %d\n",label1,label2);
834 if(label1 == label2)
835 GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCPi0);
477d6cee 836 }
591cc579 837 }//Work with mc information also
9415d854 838 pi0.SetTag(tag);
839 //Set the indeces of the original caloclusters
840 pi0.SetCaloLabel(calo->GetID(), calo2->GetID());
841 AddAODParticle(pi0);
842
843
844 }//Pair selected
845 }//if pair of gammas
846 }//2nd loop
847 }// if pdg = 22
1c5acb87 848 }//1st loop
849
9415d854 850 if(GetDebug() > 1)
851 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - End, %d pi0's found \n",GetOutputAODBranch()->GetEntriesFast());
852}
853
854//____________________________________________________________________________
233e0df8 855void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms(AliAODPWG4ParticleCorrelation * const aodParticle)
9415d854 856{
857 // Neutral Pion Correlation Analysis
858 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistogramS() - Make trigger particle - pi0 correlation, %d pi0's \n",GetOutputAODBranch()->GetEntriesFast());
859
860 Double_t pt = -100.;
dd4871cd 861 Double_t px = -100.;
862 Double_t py = -100.;
9415d854 863 Double_t rat = -100.;
864 Double_t phi = -100.;
865 Double_t eta = -100.;
dd4871cd 866 Double_t xE = -100.;
867 Double_t cosi = -100.;
868
9415d854 869 Double_t ptTrig = aodParticle->Pt();
870 Double_t phiTrig = aodParticle->Phi();
871 Double_t etaTrig = aodParticle->Eta();
dd4871cd 872 Double_t pxTrig = aodParticle->Px();
873 Double_t pyTrig = aodParticle->Py();
874
9415d854 875
876 if(!GetOutputAODBranch()){
877 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - No output pi0 in AOD branch with name < %s >,STOP \n",GetOutputAODName().Data());
878 abort();
a3aebfff 879 }
9415d854 880
881 //Loop on stored AOD pi0
882 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
883 if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - aod branch entries %d\n", naod);
884 for(Int_t iaod = 0; iaod < naod ; iaod++){
885 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
886 Int_t pdg = pi0->GetPdg();
887
888 if(pdg != AliCaloPID::kPi0) continue;
889 pt = pi0->Pt();
dd4871cd 890 px = pi0->Px();
891 py = pi0->Py();
9415d854 892 if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
893
894 //Selection within angular range
895 phi = pi0->Phi();
dde5a268 896 //Float_t deltaphi = TMath::Abs(phiTrig-phi);
897 //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
898 Float_t deltaphi = phiTrig-phi;
899 if(deltaphi<-TMath::PiOver2()) deltaphi+=TMath::TwoPi();
900 if(deltaphi>3*TMath::PiOver2()) deltaphi-=TMath::TwoPi();
901
9415d854 902 rat = pt/ptTrig ;
903 phi = pi0->Phi() ;
904 eta = pi0->Eta() ;
dd4871cd 905 xE = -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
906 if(xE <0.)xE =-xE;
907 cosi = TMath::Log(1/xE);
9415d854 908
dd4871cd 909 fhEtaNeutral->Fill(pt,eta);
910 fhPhiNeutral->Fill(pt,phi);
9415d854 911 fhDeltaEtaNeutral->Fill(ptTrig,etaTrig-eta);
912 fhDeltaPhiNeutral->Fill(ptTrig,deltaphi);
123fc3bd 913
914 //delta phi cut for correlation
915 if( (deltaphi > fDeltaPhiMinCut) && ( deltaphi < fDeltaPhiMaxCut) ) {
916 fhDeltaPhiNeutralPt->Fill(pt,deltaphi);
917 fhPtImbalanceNeutral->Fill(ptTrig,rat);
dd4871cd 918 fhPtHbpNeutral->Fill(ptTrig,cosi);
123fc3bd 919 }
920 else {
921 fhDeltaPhiUeNeutralPt->Fill(pt,deltaphi);
922 fhPtImbalanceUeNeutral->Fill(ptTrig,rat);
dd4871cd 923 fhPtHbpUeNeutral->Fill(ptTrig,cosi);
123fc3bd 924 }
925 //several UE calculation
926 if(fMakeSeveralUE){
927 if((deltaphi<-fUeDeltaPhiMinCut) && (deltaphi >-fUeDeltaPhiMaxCut)){
928 fhDeltaPhiUeLeftNeutral->Fill(pt,deltaphi);
929 fhPtImbalanceUeLeftNeutral->Fill(ptTrig,rat);
dd4871cd 930 fhPtHbpUeLeftNeutral->Fill(ptTrig,cosi);
dde5a268 931 }
123fc3bd 932 if((deltaphi>fUeDeltaPhiMinCut) && (deltaphi <fUeDeltaPhiMaxCut)){
933 fhDeltaPhiUeRightNeutral->Fill(pt,deltaphi);
934 fhPtImbalanceUeRightNeutral->Fill(ptTrig,rat);
dd4871cd 935 fhPtHbpUeRightNeutral->Fill(ptTrig,cosi);
123fc3bd 936 }
937 } //several UE calculation
938
dde5a268 939
9415d854 940 if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
941
942 }//loop
1c5acb87 943}
944
9415d854 945
1c5acb87 946//____________________________________________________________________________
947Bool_t AliAnaParticleHadronCorrelation::SelectCluster(AliAODCaloCluster * calo, Double_t *vertex, TLorentzVector & mom, Int_t & pdg) const {
9415d854 948 //Select cluster depending on its pid and acceptance selections
949
950 //Skip matched clusters with tracks
12524a23 951 if(calo->GetNTracksMatched() > 0) return kFALSE;
9415d854 952
12524a23 953 TString detector = "";
954 if(calo->IsPHOSCluster()) detector= "PHOS";
955 else if(calo->IsEMCALCluster()) detector= "EMCAL";
956
9415d854 957 //Check PID
958 calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line
959 pdg = AliCaloPID::kPhoton;
960 if(IsCaloPIDOn()){
961 //Get most probable PID, 2 options check PID weights (in MC this option is mandatory)
962 //or redo PID, recommended option for EMCal.
9415d854 963
964 if(!IsCaloPIDRecalculationOn() || GetReader()->GetDataType() == AliCaloTrackReader::kMC )
965 pdg = GetCaloPID()->GetPdg(detector,calo->PID(),mom.E());//PID with weights
966 else
967 pdg = GetCaloPID()->GetPdg(detector,mom,calo);//PID recalculated
968
969 if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - PDG of identified particle %d\n",pdg);
970
971 //If it does not pass pid, skip
12524a23 972 if(pdg != AliCaloPID::kPhoton && pdg != AliCaloPID::kPi0) {
9415d854 973 return kFALSE ;
974 }
12524a23 975 }//PID on
9415d854 976
977 //Check acceptance selection
ff45398a 978 if(IsFiducialCutOn()){
979 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,detector) ;
12524a23 980 if(! in ) return kFALSE ;
9415d854 981 }
982
983 if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
984
985 return kTRUE;
986
987}