]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/PartCorrDep/AliAnaParticleHadronCorrelation.cxx
Bug corrected in QA checker
[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}
78219bac 77/*
1c5acb87 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}
78219bac 168*/
1c5acb87 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
5a2dbc3c 178 Int_t nptbins = GetHistoPtBins();
179 Int_t nphibins = GetHistoPhiBins();
180 Int_t netabins = GetHistoEtaBins();
477d6cee 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}",
4df35693 204 nptbins,ptmin,ptmax,140,-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}}",
4df35693 210 nptbins,ptmin,ptmax,140,-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}}",
4df35693 216 nptbins,ptmin,ptmax,140,-2.,5.);
dde5a268 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}",
4df35693 228 nptbins,ptmin,ptmax,200,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}",
4df35693 234 nptbins,ptmin,ptmax,200,0.,2.);
dde5a268 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",
4df35693 240 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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",
4df35693 246 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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",
4df35693 263 nptbins,ptmin,ptmax,140,-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",
4df35693 270 nptbins,ptmin,ptmax,140,-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",
4df35693 277 nptbins,ptmin,ptmax,200,0.,2.);
dde5a268 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",
4df35693 284 nptbins,ptmin,ptmax,200,0.,2.);
dde5a268 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",
4df35693 291 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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",
4df35693 298 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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}}}",
4df35693 329 nptbins,ptmin,ptmax,140,-2.,5.);
477d6cee 330 fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi");
de8a210f 331 fhDeltaPhiNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
123fc3bd 332
333 fhDeltaPhiUeNeutralPt = new TH2F
334 ("DeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
4df35693 335 nptbins,ptmin,ptmax,140,-2.,5.);
123fc3bd 336 fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi");
de8a210f 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}",
4df35693 347 nptbins,ptmin,ptmax,200,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}",
4df35693 353 nptbins,ptmin,ptmax,200,0.,2.);
123fc3bd 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",
4df35693 359 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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",
4df35693 365 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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",
4df35693 385 nptbins,ptmin,ptmax,140,-2.,5.);
123fc3bd 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",
4df35693 392 nptbins,ptmin,ptmax,140,-2.,5.);
123fc3bd 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",
4df35693 399 nptbins,ptmin,ptmax,140,0.,2.);
123fc3bd 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",
4df35693 406 nptbins,ptmin,ptmax,200,0.,2.);
123fc3bd 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",
4df35693 413 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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",
4df35693 420 nptbins,ptmin,ptmax,200,0.,10.);
dd4871cd 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)) ;
a14fd526 433 delete nmsHistos;
477d6cee 434 }
a3aebfff 435
477d6cee 436 }//Correlation with neutral hadrons
437
438 return outputContainer;
a3aebfff 439
1c5acb87 440}
441
477d6cee 442//____________________________________________________________________________
1c5acb87 443void AliAnaParticleHadronCorrelation::InitParameters()
444{
477d6cee 445
1c5acb87 446 //Initialize the parameters of the analysis.
a3aebfff 447 SetInputAODName("PWG4Particle");
591cc579 448 SetAODObjArrayName("Hadrons");
a3aebfff 449 AddToHistogramsName("AnaHadronCorr_");
450
9415d854 451 //Correlation with neutrals
452 //SetOutputAODClassName("AliAODPWG4Particle");
453 //SetOutputAODName("Pi0Correlated");
454
dde5a268 455 SetPtCutRange(0.,300);
1c5acb87 456 fDeltaPhiMinCut = 1.5 ;
457 fDeltaPhiMaxCut = 4.5 ;
6639984f 458 fSelectIsolated = kFALSE;
dde5a268 459 fMakeSeveralUE = kFALSE;
460 fUeDeltaPhiMinCut = 1. ;
461 fUeDeltaPhiMaxCut = 1.5 ;
1c5acb87 462}
463
464//__________________________________________________________________
465void AliAnaParticleHadronCorrelation::Print(const Option_t * opt) const
466{
467
468 //Print some relevant parameters set for the analysis
469 if(! opt)
470 return;
471
a3aebfff 472 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
473 AliAnaPartCorrBaseClass::Print(" ");
474
1c5acb87 475 printf("Phi trigger particle-Hadron < %3.2f\n", fDeltaPhiMaxCut) ;
476 printf("Phi trigger particle-Hadron > %3.2f\n", fDeltaPhiMinCut) ;
6639984f 477 printf("Isolated Trigger? %d\n", fSelectIsolated) ;
dde5a268 478 printf("Phi trigger particle-UeHadron < %3.2f\n", fUeDeltaPhiMaxCut) ;
479 printf("Phi trigger particle-UeHadron > %3.2f\n", fUeDeltaPhiMinCut) ;
480 printf("Several UE? %d\n", fMakeSeveralUE) ;
1c5acb87 481}
482
483//____________________________________________________________________________
484void AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD()
485{
477d6cee 486 //Particle-Hadron Correlation Analysis, fill AODs
487
488 if(!GetInputAODBranch()){
591cc579 489 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 490 abort();
491 }
9415d854 492
493 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
494 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());
495 abort();
496 }
497
477d6cee 498 if(GetDebug() > 1){
a3aebfff 499 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Begin hadron correlation analysis, fill AODs \n");
500 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
501 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetAODCTS()->GetEntriesFast());
502 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetAODEMCAL()->GetEntriesFast());
503 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n", GetAODPHOS()->GetEntriesFast());
477d6cee 504 }
505
506 //Loop on stored AOD particles, trigger
507 Int_t naod = GetInputAODBranch()->GetEntriesFast();
508 for(Int_t iaod = 0; iaod < naod ; iaod++){
509 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
9415d854 510
511 //Make correlation with charged hadrons
477d6cee 512 if(GetReader()->IsCTSSwitchedOn() )
a3aebfff 513 MakeChargedCorrelation(particle, GetAODCTS(),kFALSE);
477d6cee 514
515 //Make correlation with neutral pions
516 //Trigger particle in PHOS, correlation with EMCAL
517 if(particle->GetDetector()=="PHOS" && GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
9415d854 518 MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
477d6cee 519 //Trigger particle in EMCAL, correlation with PHOS
520 else if(particle->GetDetector()=="EMCAL" && GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
9415d854 521 MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
477d6cee 522 //Trigger particle in CTS, correlation with PHOS, EMCAL and CTS
523 else if(particle->GetDetector()=="CTS" ){
524 if(GetReader()->IsPHOSSwitchedOn() && GetAODPHOS()->GetEntriesFast() > 0)
9415d854 525 MakeNeutralCorrelationFillAOD(particle, GetAODPHOS(),"PHOS");
477d6cee 526 if(GetReader()->IsEMCALSwitchedOn() && GetAODEMCAL()->GetEntriesFast() > 0)
9415d854 527 MakeNeutralCorrelationFillAOD(particle, GetAODEMCAL(),"EMCAL");
477d6cee 528 }
a3aebfff 529
530
477d6cee 531 }//Aod branch loop
532
a3aebfff 533 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
477d6cee 534
1c5acb87 535}
536
537//____________________________________________________________________________
538void AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms()
539{
477d6cee 540 //Particle-Hadron Correlation Analysis, fill histograms
541
542 if(!GetInputAODBranch()){
591cc579 543 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
477d6cee 544 abort();
545 }
9415d854 546
477d6cee 547 if(GetDebug() > 1){
a3aebfff 548 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Begin hadron correlation analysis, fill histograms \n");
549 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
477d6cee 550 }
9415d854 551
477d6cee 552 //Loop on stored AOD particles
553 Int_t naod = GetInputAODBranch()->GetEntriesFast();
9415d854 554 for(Int_t iaod = 0; iaod < naod ; iaod++){
555 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
556
557 //check if the particle is isolated or if we want to take the isolation into account
558 if(OnlyIsolated() && !particle->IsIsolated()) continue;
559
560 //Make correlation with charged hadrons
591cc579 561 TObjArray * reftracks = particle->GetObjArray(GetAODObjArrayName()+"Tracks");
9415d854 562 if(reftracks){
563 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Track Refs entries %d\n", iaod, reftracks->GetEntriesFast());
564 if(reftracks->GetEntriesFast() > 0) MakeChargedCorrelation(particle, reftracks,kTRUE);
565 }
a70a8477 566
9415d854 567 //Make correlation with neutral pions
568 if(GetOutputAODBranch() && GetOutputAODBranch()->GetEntriesFast() > 0){
569 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Particle %d, In Cluster Refs entries %d\n",iaod, GetOutputAODBranch()->GetEntriesFast());
570 MakeNeutralCorrelationFillHistograms(particle);
571 }
572
477d6cee 573 }//Aod branch loop
574
a3aebfff 575 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
477d6cee 576
1c5acb87 577}
578
579//____________________________________________________________________________
233e0df8 580void AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle, TObjArray* const pl, const Bool_t bFillHisto)
1c5acb87 581{
dde5a268 582 // Charged Hadron Correlation Analysis
583 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
584
585 Double_t ptTrig = aodParticle->Pt();
dd4871cd 586 Double_t pxTrig = aodParticle->Px();
587 Double_t pyTrig = aodParticle->Py();
588
dde5a268 589 Double_t phiTrig = aodParticle->Phi();
dd4871cd 590
dde5a268 591 Double_t pt = -100.;
dd4871cd 592 Double_t px = -100.;
593 Double_t py = -100.;
dde5a268 594 Double_t rat = -100.;
9ccc5a0b 595 Double_t xE = -100.;
596 Double_t cosi = -100.;
dde5a268 597 Double_t phi = -100. ;
598 Double_t eta = -100. ;
599 Double_t p[3];
dde5a268 600
a33b76ed 601 TObjArray * reftracks =0x0;
602 Int_t nrefs = 0;
603
dde5a268 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 ;
9ccc5a0b 616 xE = -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
dd4871cd 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{
a33b76ed 676 nrefs++;
a70a8477 677 if(nrefs==1){
9ccc5a0b 678 reftracks = new TObjArray(0);
679 reftracks->SetName(GetAODObjArrayName()+"Tracks");
680 reftracks->SetOwner(kFALSE);
a70a8477 681 }
dde5a268 682 reftracks->Add(track);
dde5a268 683 }//aod particle loop
684 }// track loop
9415d854 685
dde5a268 686 //Fill AOD with reference tracks, if not filling histograms
a70a8477 687 if(!bFillHisto && reftracks) {
591cc579 688 aodParticle->AddObjArray(reftracks);
dde5a268 689 }
a70a8477 690
691 //delete reftracks;
9415d854 692
1c5acb87 693}
9415d854 694
1c5acb87 695//____________________________________________________________________________
233e0df8 696void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD(AliAODPWG4ParticleCorrelation* const aodParticle,TObjArray* const pl, TString detector)
1c5acb87 697{
9415d854 698 // Neutral Pion Correlation Analysis, find pi0, put them in new output aod, if correlation cuts passed
699 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Make trigger particle - neutral hadron correlation \n");
700
701 if(!NewOutputAOD()){
591cc579 702 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Output aod not created, set AOD class name and branch name in the configuration file, STOP! \n");
9415d854 703 abort();
704 }
1c5acb87 705
1c5acb87 706 Double_t phiTrig = aodParticle->Phi();
591cc579 707 Int_t tag = 0;
1c5acb87 708 TLorentzVector gammai;
709 TLorentzVector gammaj;
710
233e0df8 711 //Get vertex for photon momentum calculation
712 Double_t vertex[] = {0,0,0} ; //vertex
713 Double_t vertex2[] = {0,0,0} ; //vertex of second input aod
edd59991 714 if(GetReader()->GetDataType()!= AliCaloTrackReader::kMC)
233e0df8 715 {
716 GetReader()->GetVertex(vertex);
717 if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex2);
718 }
719
1c5acb87 720 //Cluster loop, select pairs with good pt, phi and fill AODs or histograms
9415d854 721 //Int_t iEvent= GetReader()->GetEventNumber() ;
722 Int_t nclus = pl->GetEntriesFast();
723 for(Int_t iclus = 0;iclus < nclus ; iclus ++ ){
1c5acb87 724 AliAODCaloCluster * calo = (AliAODCaloCluster *) (pl->At(iclus)) ;
9415d854 725
233e0df8 726 //Input from second AOD?
727 Int_t inputi = 0;
728 if (aodParticle->GetDetector() == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= iclus) inputi = 1 ;
729 else if(aodParticle->GetDetector() == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= iclus) inputi = 1;
730
ff45398a 731 //Cluster selection, not charged, with photon or pi0 id and in fiducial cut
9415d854 732 Int_t pdg=0;
233e0df8 733 if (inputi == 0 && !SelectCluster(calo, vertex, gammai, pdg)) continue ;
734 else if(inputi == 1 && !SelectCluster(calo, vertex2, gammai, pdg)) continue ;
12524a23 735
736 if(GetDebug() > 2)
9415d854 737 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",
738 detector.Data(), gammai.Pt(),gammai.Phi(),phiTrig,fDeltaPhiMinCut, fDeltaPhiMaxCut, GetMinPt());
739
740 //2 gamma overlapped, found with PID
741 if(pdg == AliCaloPID::kPi0){
477d6cee 742
9415d854 743 //Select only hadrons in pt range
744 if(gammai.Pt() < GetMinPt() || gammai.Pt() > GetMaxPt()) continue ;
1c5acb87 745
9415d854 746 //Selection within angular range
747 Float_t phi = gammai.Phi();
748 if(phi < 0) phi+=TMath::TwoPi();
dde5a268 749 //Float_t deltaphi = TMath::Abs(phiTrig-phi);
750 //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
477d6cee 751
9415d854 752 AliAODPWG4Particle pi0 = AliAODPWG4Particle(gammai);
753 //pi0.SetLabel(calo->GetLabel(0));
754 pi0.SetPdg(AliCaloPID::kPi0);
755 pi0.SetDetector(detector);
756
757 if(IsDataMC()){
233e0df8 758 pi0.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabel(0),GetReader(),inputi));
759 if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of candidate %d\n",pi0.GetTag());
760 }//Work with stack also
9415d854 761 //Set the indeces of the original caloclusters
762 pi0.SetCaloLabel(calo->GetID(),-1);
763 AddAODParticle(pi0);
764
765 if(GetDebug() > 2)
766 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Correlated with selected pi0 (pid): pt %f, phi %f\n",pi0.Pt(),pi0.Phi());
767
768 }// pdg = 111
769
770 //Make invariant mass analysis
771 else if(pdg == AliCaloPID::kPhoton){
772 //Search the photon companion in case it comes from a Pi0 decay
773 //Apply several cuts to select the good pair;
774 for(Int_t jclus = iclus+1; jclus < pl->GetEntries() ; jclus ++ ){
775 AliAODCaloCluster * calo2 = (AliAODCaloCluster *) (pl->At(jclus)) ;
477d6cee 776
233e0df8 777 //Input from second AOD?
778 Int_t inputj = 0;
779 if (aodParticle->GetDetector() == "EMCAL" && GetReader()->GetAODEMCALNormalInputEntries() <= jclus) inputj = 1;
780 else if(aodParticle->GetDetector() == "PHOS" && GetReader()->GetAODPHOSNormalInputEntries() <= jclus) inputj = 1;
781
ff45398a 782 //Cluster selection, not charged with photon or pi0 id and in fiducial cut
9415d854 783 Int_t pdgj=0;
233e0df8 784 if (inputj == 0 && !SelectCluster(calo2, vertex, gammaj, pdgj)) continue ;
785 else if(inputj == 1 && !SelectCluster(calo2, vertex2, gammaj, pdgj)) continue ;
786
9415d854 787 if(!SelectCluster(calo2,vertex, gammaj, pdgj)) continue ;
477d6cee 788
9415d854 789 if(pdgj == AliCaloPID::kPhoton ){
790
791 if((gammai+gammaj).Pt() < GetMinPt() || (gammai+gammaj).Pt() > GetMaxPt()) continue ;
1c5acb87 792
9415d854 793 //Selection within angular range
794 Float_t phi = (gammai+gammaj).Phi();
795 if(phi < 0) phi+=TMath::TwoPi();
dde5a268 796 //Float_t deltaphi = TMath::Abs(phiTrig-phi);
797 //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
477d6cee 798
9415d854 799 //Select good pair (aperture and invariant mass)
800 if(GetNeutralMesonSelection()->SelectPair(gammai, gammaj)){
1c5acb87 801
9415d854 802 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",
803 (gammai+gammaj).Pt(),(gammai+gammaj).Phi(),(gammai+gammaj).Eta(), (gammai+gammaj).M());
1c5acb87 804
9415d854 805 TLorentzVector pi0mom = gammai+gammaj;
806 AliAODPWG4Particle pi0 = AliAODPWG4Particle(pi0mom);
807 //pi0.SetLabel(calo->GetLabel(0));
808 pi0.SetPdg(AliCaloPID::kPi0);
809 pi0.SetDetector(detector);
810 if(IsDataMC()){
811 //Check origin of the candidates
591cc579 812
813 Int_t label1 = calo->GetLabel(0);
814 Int_t label2 = calo2->GetLabel(0);
233e0df8 815 Int_t tag1 = GetMCAnalysisUtils()->CheckOrigin(label1, GetReader(), inputi);
816 Int_t tag2 = GetMCAnalysisUtils()->CheckOrigin(label2, GetReader(), inputj);
591cc579 817
9415d854 818 if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - Origin of: photon1 %d; photon2 %d \n",tag1, tag2);
591cc579 819 if(GetMCAnalysisUtils()->CheckTagBit(tag1,AliMCAnalysisUtils::kMCPi0Decay) && GetMCAnalysisUtils()->CheckTagBit(tag2,AliMCAnalysisUtils::kMCPi0Decay)){
9415d854 820
591cc579 821 //Check if pi0 mother is the same
822 if(GetReader()->ReadStack()){
823 TParticle * mother1 = GetMCStack()->Particle(label1);//photon in kine tree
824 label1 = mother1->GetFirstMother();
825 //mother1 = GetMCStack()->Particle(label1);//pi0
826
827 TParticle * mother2 = GetMCStack()->Particle(label2);//photon in kine tree
828 label2 = mother2->GetFirstMother();
829 //mother2 = GetMCStack()->Particle(label2);//pi0
830 }
831 else if(GetReader()->ReadAODMCParticles()){
233e0df8 832 AliAODMCParticle * mother1 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(inputi))->At(label1);//photon in kine tree
591cc579 833 label1 = mother1->GetMother();
834 //mother1 = GetMCStack()->Particle(label1);//pi0
233e0df8 835 AliAODMCParticle * mother2 = (AliAODMCParticle *) (GetReader()->GetAODMCParticles(inputj))->At(label2);//photon in kine tree
591cc579 836 label2 = mother2->GetMother();
837 //mother2 = GetMCStack()->Particle(label2);//pi0
838 }
9415d854 839
591cc579 840 //printf("mother1 %d, mother2 %d\n",label1,label2);
841 if(label1 == label2)
842 GetMCAnalysisUtils()->SetTagBit(tag,AliMCAnalysisUtils::kMCPi0);
477d6cee 843 }
591cc579 844 }//Work with mc information also
9415d854 845 pi0.SetTag(tag);
846 //Set the indeces of the original caloclusters
847 pi0.SetCaloLabel(calo->GetID(), calo2->GetID());
848 AddAODParticle(pi0);
849
850
851 }//Pair selected
852 }//if pair of gammas
853 }//2nd loop
854 }// if pdg = 22
1c5acb87 855 }//1st loop
856
9415d854 857 if(GetDebug() > 1)
858 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillAOD() - End, %d pi0's found \n",GetOutputAODBranch()->GetEntriesFast());
859}
860
861//____________________________________________________________________________
233e0df8 862void AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms(AliAODPWG4ParticleCorrelation * const aodParticle)
9415d854 863{
864 // Neutral Pion Correlation Analysis
865 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistogramS() - Make trigger particle - pi0 correlation, %d pi0's \n",GetOutputAODBranch()->GetEntriesFast());
866
867 Double_t pt = -100.;
dd4871cd 868 Double_t px = -100.;
869 Double_t py = -100.;
9415d854 870 Double_t rat = -100.;
871 Double_t phi = -100.;
872 Double_t eta = -100.;
dd4871cd 873 Double_t xE = -100.;
874 Double_t cosi = -100.;
875
9415d854 876 Double_t ptTrig = aodParticle->Pt();
877 Double_t phiTrig = aodParticle->Phi();
878 Double_t etaTrig = aodParticle->Eta();
dd4871cd 879 Double_t pxTrig = aodParticle->Px();
880 Double_t pyTrig = aodParticle->Py();
881
9415d854 882
883 if(!GetOutputAODBranch()){
884 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - No output pi0 in AOD branch with name < %s >,STOP \n",GetOutputAODName().Data());
885 abort();
a3aebfff 886 }
9415d854 887
888 //Loop on stored AOD pi0
889 Int_t naod = GetOutputAODBranch()->GetEntriesFast();
890 if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - aod branch entries %d\n", naod);
891 for(Int_t iaod = 0; iaod < naod ; iaod++){
892 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
893 Int_t pdg = pi0->GetPdg();
894
895 if(pdg != AliCaloPID::kPi0) continue;
896 pt = pi0->Pt();
dd4871cd 897 px = pi0->Px();
898 py = pi0->Py();
9415d854 899 if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
900
901 //Selection within angular range
902 phi = pi0->Phi();
dde5a268 903 //Float_t deltaphi = TMath::Abs(phiTrig-phi);
904 //if( (deltaphi < fDeltaPhiMinCut) || ( deltaphi > fDeltaPhiMaxCut) ) continue ;
905 Float_t deltaphi = phiTrig-phi;
906 if(deltaphi<-TMath::PiOver2()) deltaphi+=TMath::TwoPi();
907 if(deltaphi>3*TMath::PiOver2()) deltaphi-=TMath::TwoPi();
908
9415d854 909 rat = pt/ptTrig ;
910 phi = pi0->Phi() ;
911 eta = pi0->Eta() ;
dd4871cd 912 xE = -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
913 if(xE <0.)xE =-xE;
914 cosi = TMath::Log(1/xE);
9415d854 915
dd4871cd 916 fhEtaNeutral->Fill(pt,eta);
917 fhPhiNeutral->Fill(pt,phi);
9415d854 918 fhDeltaEtaNeutral->Fill(ptTrig,etaTrig-eta);
919 fhDeltaPhiNeutral->Fill(ptTrig,deltaphi);
123fc3bd 920
921 //delta phi cut for correlation
922 if( (deltaphi > fDeltaPhiMinCut) && ( deltaphi < fDeltaPhiMaxCut) ) {
923 fhDeltaPhiNeutralPt->Fill(pt,deltaphi);
924 fhPtImbalanceNeutral->Fill(ptTrig,rat);
dd4871cd 925 fhPtHbpNeutral->Fill(ptTrig,cosi);
123fc3bd 926 }
927 else {
928 fhDeltaPhiUeNeutralPt->Fill(pt,deltaphi);
929 fhPtImbalanceUeNeutral->Fill(ptTrig,rat);
dd4871cd 930 fhPtHbpUeNeutral->Fill(ptTrig,cosi);
123fc3bd 931 }
932 //several UE calculation
933 if(fMakeSeveralUE){
934 if((deltaphi<-fUeDeltaPhiMinCut) && (deltaphi >-fUeDeltaPhiMaxCut)){
935 fhDeltaPhiUeLeftNeutral->Fill(pt,deltaphi);
936 fhPtImbalanceUeLeftNeutral->Fill(ptTrig,rat);
dd4871cd 937 fhPtHbpUeLeftNeutral->Fill(ptTrig,cosi);
dde5a268 938 }
123fc3bd 939 if((deltaphi>fUeDeltaPhiMinCut) && (deltaphi <fUeDeltaPhiMaxCut)){
940 fhDeltaPhiUeRightNeutral->Fill(pt,deltaphi);
941 fhPtImbalanceUeRightNeutral->Fill(ptTrig,rat);
dd4871cd 942 fhPtHbpUeRightNeutral->Fill(ptTrig,cosi);
123fc3bd 943 }
944 } //several UE calculation
945
dde5a268 946
9415d854 947 if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
948
949 }//loop
1c5acb87 950}
951
9415d854 952
1c5acb87 953//____________________________________________________________________________
954Bool_t AliAnaParticleHadronCorrelation::SelectCluster(AliAODCaloCluster * calo, Double_t *vertex, TLorentzVector & mom, Int_t & pdg) const {
9415d854 955 //Select cluster depending on its pid and acceptance selections
956
957 //Skip matched clusters with tracks
12524a23 958 if(calo->GetNTracksMatched() > 0) return kFALSE;
9415d854 959
12524a23 960 TString detector = "";
961 if(calo->IsPHOSCluster()) detector= "PHOS";
962 else if(calo->IsEMCALCluster()) detector= "EMCAL";
963
9415d854 964 //Check PID
965 calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line
966 pdg = AliCaloPID::kPhoton;
967 if(IsCaloPIDOn()){
968 //Get most probable PID, 2 options check PID weights (in MC this option is mandatory)
969 //or redo PID, recommended option for EMCal.
9415d854 970
971 if(!IsCaloPIDRecalculationOn() || GetReader()->GetDataType() == AliCaloTrackReader::kMC )
972 pdg = GetCaloPID()->GetPdg(detector,calo->PID(),mom.E());//PID with weights
973 else
974 pdg = GetCaloPID()->GetPdg(detector,mom,calo);//PID recalculated
975
976 if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - PDG of identified particle %d\n",pdg);
977
978 //If it does not pass pid, skip
12524a23 979 if(pdg != AliCaloPID::kPhoton && pdg != AliCaloPID::kPi0) {
9415d854 980 return kFALSE ;
981 }
12524a23 982 }//PID on
9415d854 983
984 //Check acceptance selection
ff45398a 985 if(IsFiducialCutOn()){
986 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,detector) ;
12524a23 987 if(! in ) return kFALSE ;
9415d854 988 }
989
990 if(GetDebug() > 5) printf("AliAnaParticleHadronCorrelation::SelectCluster() - Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg);
991
992 return kTRUE;
993
994}