Remove depricated function,
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaParticleHadronCorrelation.cxx
CommitLineData
e09cf5ef 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)
21
22// Modified by Yaxian Mao:
23// 1. add the UE subtraction for corrlation study
24// 2. change the correlation variable
25// 3. Only use leading particle(cluster/track) as trigger for correlation (2010/07/02)
26// 4. Make decay photon-hadron correlations where decay contribute pi0 mass (2010/09/09)
27// 5. fill the pout to extract kt at the end, also to study charge asymmetry(2010/10/06)
045396c8 28// 6. Add the possibility for event selection analysis based on vertex and multiplicity bins (10/10/2010)
29// 7. change the way of delta phi cut for UE study due to memory issue (reduce histograms)
30// 8. Add the possibility to request the absolute leading particle at the near side or not, set trigger bins, general clean-up (08/2011)
e09cf5ef 31//////////////////////////////////////////////////////////////////////////////
32
33
34// --- ROOT system ---
35//#include "TClonesArray.h"
36#include "TClass.h"
37#include "TMath.h"
38#include "TH3D.h"
045396c8 39#include "TDatabasePDG.h"
e09cf5ef 40
41//---- ANALYSIS system ----
42#include "AliNeutralMesonSelection.h"
43#include "AliAnaParticleHadronCorrelation.h"
44#include "AliCaloTrackReader.h"
45#include "AliCaloPID.h"
46#include "AliAODPWG4ParticleCorrelation.h"
47#include "AliFiducialCut.h"
48#include "AliVTrack.h"
49#include "AliVCluster.h"
50#include "AliMCAnalysisUtils.h"
51#include "TParticle.h"
52#include "AliStack.h"
53#include "AliAODMCParticle.h"
54#include "AliMixedEvent.h"
55
56ClassImp(AliAnaParticleHadronCorrelation)
57
58
59//____________________________________________________________________________
60 AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation():
61 AliAnaPartCorrBaseClass(),
62 fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.), fSelectIsolated(0),
63 fMakeSeveralUE(0),
64 fUeDeltaPhiMaxCut(0.), fUeDeltaPhiMinCut(0.),
65 fPi0AODBranchName(""),fNeutralCorr(0), fPi0Trigger(0),
045396c8 66 fMakeAbsoluteLeading(0), fLeadingTriggerIndex(-1),
67 fNTriggerPtBins(0),
e09cf5ef 68 fhPtLeading(0),fhPhiLeading(0),fhEtaLeading(0),
69 fhDeltaPhiDeltaEtaCharged(0),
70 fhPhiCharged(0), fhEtaCharged(0),
71 fhDeltaPhiCharged(0),
72 fhDeltaEtaCharged(0),
73 fhDeltaPhiChargedPt(0),
74 fhDeltaPhiUeChargedPt(0),
75 fhPtImbalanceCharged(0),
76 fhPtImbalanceUeCharged(0),
77 fhPtImbalancePosCharged(0),fhPtImbalanceNegCharged(0),
78 fhPtHbpCharged(0), fhPtHbpUeCharged(0),
79 fhDeltaPhiUeLeftCharged(0),fhDeltaPhiUeRightCharged(0),
80 fhPtImbalanceUeLeftCharged(0),fhPtImbalanceUeRightCharged(0),
81 fhPtHbpUeLeftCharged(0),fhPtHbpUeRightCharged(0),
82 fhPtTrigPout(0), fhPtAssocDeltaPhi(0),
83 fhPtTrigCharged(0),
84 fhTrigDeltaPhiCharged(0x0), fhTrigDeltaEtaCharged(0x0),fhTrigCorr(0x0),fhTrigUeCorr(0x0),
85 fhDeltaPhiDeltaEtaNeutral(0),
86 fhPhiNeutral(0), fhEtaNeutral(0),
87 fhDeltaPhiNeutral(0), fhDeltaEtaNeutral(0),
88 fhDeltaPhiNeutralPt(0),fhDeltaPhiUeNeutralPt(0),
89 fhPtImbalanceNeutral(0),fhPtImbalanceUeNeutral(0),
90 fhPtHbpNeutral(0), fhPtHbpUeNeutral(0),
91 fhDeltaPhiUeLeftNeutral(0),fhDeltaPhiUeRightNeutral(0),
92 fhPtImbalanceUeLeftNeutral(0),fhPtImbalanceUeRightNeutral(0),
93 fhPtHbpUeLeftNeutral(0),fhPtHbpUeRightNeutral(0),
94 fhPtPi0DecayRatio(0),
95 fhDeltaPhiDecayCharged(0),
96 fhPtImbalanceDecayCharged(0),
97 fhDeltaPhiDecayNeutral(0),
045396c8 98 fhPtImbalanceDecayNeutral(0),
99 fh2phiLeadingParticle(0x0),
100 fhMCLeadingCount(0),
101 fhMCEtaCharged(0),
102 fhMCPhiCharged(0),
103 fhMCDeltaEtaCharged(0),
104 fhMCDeltaPhiCharged(0x0),
105 fhMCDeltaPhiDeltaEtaCharged(0),
106 fhMCDeltaPhiChargedPt(0),
107 fhMCPtImbalanceCharged(0),
108 fhMCPtHbpCharged(0),
109 fhMCPtTrigPout(0),
110 fhMCPtAssocDeltaPhi(0)
e09cf5ef 111{
112 //Default Ctor
113
114 //Initialize parameters
115 InitParameters();
116}
117
118//________________________________________________________________________
119TList * AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
120{
121
122 // Create histograms to be saved in output file and
123 // store them in fOutputContainer
124 TList * outputContainer = new TList() ;
125 outputContainer->SetName("CorrelationHistos") ;
126
127 Int_t nptbins = GetHistoPtBins();
128 Int_t nphibins = GetHistoPhiBins();
129 Int_t netabins = GetHistoEtaBins();
130 Float_t ptmax = GetHistoPtMax();
131 Float_t phimax = GetHistoPhiMax();
132 Float_t etamax = GetHistoEtaMax();
133 Float_t ptmin = GetHistoPtMin();
134 Float_t phimin = GetHistoPhiMin();
135 Float_t etamin = GetHistoEtaMin();
e09cf5ef 136
e09cf5ef 137 fhPtLeading = new TH1F ("hPtLeading","p_T distribution of leading particles", nptbins,ptmin,ptmax);
138 fhPtLeading->SetXTitle("p_{T}^{trig} (GeV/c)");
139
140 fhPhiLeading = new TH2F ("hPhiLeading","#phi distribution of leading Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
141 fhPhiLeading->SetYTitle("#phi (rad)");
142
143 fhEtaLeading = new TH2F ("hEtaLeading","#eta distribution of leading",nptbins,ptmin,ptmax, netabins,etamin,etamax);
144 fhEtaLeading->SetYTitle("#eta ");
045396c8 145
e09cf5ef 146 outputContainer->Add(fhPtLeading);
147 outputContainer->Add(fhPhiLeading);
148 outputContainer->Add(fhEtaLeading);
149
150 //Correlation with charged hadrons
151 if(GetReader()->IsCTSSwitchedOn()) {
152 fhDeltaPhiDeltaEtaCharged = new TH2F
153 ("DeltaPhiDeltaEtaCharged","#phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
154 140,-2.,5.,200,-2,2);
155 fhDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
156 fhDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
157
158 fhPhiCharged = new TH2F
159 ("PhiCharged","#phi_{h^{#pm}} vs p_{T #pm}",
160 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
161 fhPhiCharged->SetYTitle("#phi_{h^{#pm}} (rad)");
162 fhPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
163
164 fhEtaCharged = new TH2F
165 ("EtaCharged","#eta_{h^{#pm}} vs p_{T #pm}",
166 nptbins,ptmin,ptmax,netabins,etamin,etamax);
167 fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
168 fhEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
169
170 fhDeltaPhiCharged = new TH2F
171 ("DeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
172 nptbins,ptmin,ptmax,140,-2.,5.);
173 fhDeltaPhiCharged->SetYTitle("#Delta #phi");
174 fhDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
175
176 fhDeltaPhiChargedPt = new TH2F
177 ("DeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
178 nptbins,ptmin,ptmax,140,-2.,5.);
179 fhDeltaPhiChargedPt->SetYTitle("#Delta #phi");
180 fhDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
181
182 fhDeltaPhiUeChargedPt = new TH2F
183 ("DeltaPhiUeChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}}",
184 nptbins,ptmin,ptmax,140,-2.,5.);
185 fhDeltaPhiUeChargedPt->SetYTitle("#Delta #phi");
186 fhDeltaPhiUeChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
187
188 fhDeltaEtaCharged = new TH2F
189 ("DeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger}",
190 nptbins,ptmin,ptmax,200,-2,2);
191 fhDeltaEtaCharged->SetYTitle("#Delta #eta");
192 fhDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
193
194 fhPtImbalanceCharged =
195 new TH2F("CorrelationCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T trigger}",
196 nptbins,ptmin,ptmax,200,0.,2.);
197 fhPtImbalanceCharged->SetYTitle("z_{trigger h^{#pm}}");
198 fhPtImbalanceCharged->SetXTitle("p_{T trigger}");
199
200 fhPtImbalanceUeCharged =
201 new TH2F("CorrelationUeCharged","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger}",
202 nptbins,ptmin,ptmax,200,0.,2.);
203 fhPtImbalanceUeCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
204 fhPtImbalanceUeCharged->SetXTitle("p_{T trigger}");
205
206 fhPtImbalancePosCharged =
207 new TH2F("CorrelationPositiveCharged","z_{trigger h^{+}} = p_{T h^{+}} / p_{T trigger}",
208 nptbins,ptmin,ptmax,200,0.,2.);
209 fhPtImbalancePosCharged->SetYTitle("z_{trigger h^{+}}");
210 fhPtImbalancePosCharged->SetXTitle("p_{T trigger}");
211
212 fhPtImbalanceNegCharged =
213 new TH2F("CorrelationNegativeCharged","z_{trigger h^{-}} = p_{T h^{-}} / p_{T trigger}",
214 nptbins,ptmin,ptmax,200,0.,2.);
215 fhPtImbalanceNegCharged->SetYTitle("z_{trigger h^{-}}");
216 fhPtImbalanceNegCharged->SetXTitle("p_{T trigger}");
217
218 fhPtHbpCharged =
219 new TH2F("HbpCharged","#xi = ln(1/x_{E}) with charged hadrons",
220 nptbins,ptmin,ptmax,200,0.,10.);
221 fhPtHbpCharged->SetYTitle("ln(1/x_{E})");
222 fhPtHbpCharged->SetXTitle("p_{T trigger}");
223
224 fhPtHbpUeCharged =
225 new TH2F("HbpUeCharged","#xi = ln(1/x_{E}) with charged hadrons",
226 nptbins,ptmin,ptmax,200,0.,10.);
227 fhPtHbpUeCharged->SetYTitle("ln(1/x_{E})");
228 fhPtHbpUeCharged->SetXTitle("p_{T trigger}");
229
230 fhPtTrigPout =
231 new TH2F("PtTrigPout","Pout with triggers",
232 nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax);
233 fhPtTrigPout->SetYTitle("p_{out} (GeV/c)");
234 fhPtTrigPout->SetXTitle("p_{T trigger} (GeV/c)");
235
236 fhPtAssocDeltaPhi =
237 new TH2F("fhPtAssocDeltaPhi"," charged hadrons vs. delta phi",
238 nptbins,ptmin,ptmax,140,-2.,5.);
239 fhPtAssocDeltaPhi->SetXTitle("p_{T h^{#pm}} (GeV/c)");
240 fhPtAssocDeltaPhi->SetYTitle("#Delta #phi (GeV/c)");
241
242 fhPtTrigCharged =
243 new TH2F("PtTrigCharged","trgger and charged tracks pt distribution",
244 nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
245 fhPtTrigCharged->SetYTitle("p_{T h^{#pm}} (GeV/c)");
246 fhPtTrigCharged->SetXTitle("p_{T trigger} (GeV/c)");
247
248 outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
249 outputContainer->Add(fhPhiCharged) ;
250 outputContainer->Add(fhEtaCharged) ;
251 outputContainer->Add(fhDeltaPhiCharged) ;
252 outputContainer->Add(fhDeltaEtaCharged) ;
253 outputContainer->Add(fhDeltaPhiChargedPt) ;
254 outputContainer->Add(fhDeltaPhiUeChargedPt) ;
255 outputContainer->Add(fhPtImbalanceCharged) ;
256 outputContainer->Add(fhPtImbalancePosCharged) ;
257 outputContainer->Add(fhPtImbalanceNegCharged) ;
258 outputContainer->Add(fhPtImbalanceUeCharged) ;
259 outputContainer->Add(fhPtHbpCharged) ;
260 outputContainer->Add(fhPtHbpUeCharged) ;
261 outputContainer->Add(fhPtTrigPout) ;
262 outputContainer->Add(fhPtAssocDeltaPhi) ;
263 outputContainer->Add(fhPtTrigCharged) ;
264
265 if(DoEventSelect()){
266 Int_t nMultiBins = GetMultiBin();
267 fhTrigDeltaPhiCharged = new TH2F*[nMultiBins] ;
268 fhTrigDeltaEtaCharged = new TH2F*[nMultiBins] ;
269 fhTrigCorr = new TH2F*[nMultiBins];
270 fhTrigUeCorr = new TH2F*[nMultiBins];
271 for(Int_t im=0; im<nMultiBins; im++){
272 fhTrigDeltaPhiCharged[im] = new TH2F
273 (Form("fhTrigDeltaPhiCharged_%d",im),Form("fhTrigDeltaPhiCharged_%d",im), nptbins,ptmin,ptmax, 140,-2.,5.);
274 fhTrigDeltaPhiCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
275 fhTrigDeltaPhiCharged[im]->SetYTitle("#Delta #phi");
276 fhTrigDeltaEtaCharged[im] = new TH2F
277 (Form("fhTrigDeltaEtaCharged_%d",im),Form("fhTrigDeltaEtaCharged_%d",im), nptbins,ptmin,ptmax, 200,-2,2);
278 fhTrigDeltaEtaCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
279 fhTrigDeltaEtaCharged[im]->SetYTitle("#Delta #eta");
280 fhTrigCorr[im] = new TH2F
281 (Form("fhTrigPtCorr_%d",im),Form("fhTrigPtCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
282 fhTrigCorr[im]->SetYTitle("z_{trigger h^{#pm}}");
283 fhTrigCorr[im]->SetXTitle("p_{T trigger}");
284 fhTrigUeCorr[im] = new TH2F
285 (Form("fhTrigPtUeCorr_%d",im),Form("fhTrigPtUeCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
286 fhTrigUeCorr[im]->SetYTitle("z_{trigger h^{#pm}}");
287 fhTrigUeCorr[im]->SetXTitle("p_{T trigger}");
045396c8 288
e09cf5ef 289 outputContainer->Add(fhTrigDeltaPhiCharged[im]) ;
290 outputContainer->Add(fhTrigDeltaEtaCharged[im]) ;
291 outputContainer->Add(fhTrigCorr[im]);
292 outputContainer->Add(fhTrigUeCorr[im]);
293
045396c8 294 }
295 }
296
297 for(Int_t i = 0 ; i < fNTriggerPtBins ; i++){
298 fhTrigPt[i] = new TH1F(Form("fhTrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
299 Form("fhTrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
300 20, fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]);
301 fhDphiTrigPtAssocPt[i] = new TH1F(Form("fhDphiTrigPtAssocPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
302 Form("fhDphiTrigPtAssocPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
303 288, -1.0/3.0, 5.0/3.0);
304 fhAssocPtBkg[i] = new TH1F(Form("fhAssocPtBkg%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
305 Form("fhAssocPtBkg%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
306 600, 0., 30.0);
307 fhXE[i] = new TH1F(Form("fhXETrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
308 Form("fhXETrigPt%1.f_%1.f", fTriggerPtBinLimit[i], fTriggerPtBinLimit[i+1]),
309 50, 0.0, 2.0);
310 outputContainer->Add(fhTrigPt[i]) ;
311 outputContainer->Add(fhDphiTrigPtAssocPt[i]) ;
312 outputContainer->Add(fhAssocPtBkg[i]);
313 outputContainer->Add(fhXE[i]);
e09cf5ef 314 }
045396c8 315
e09cf5ef 316 if(fPi0Trigger){
317 fhPtPi0DecayRatio = new TH2F
318 ("hPtPi0DecayRatio","p_T of #pi^{0} and the ratio of pt for two decay",
319 nptbins,ptmin,ptmax, 100,0.,2.);
320 fhPtPi0DecayRatio->SetXTitle("p_{T}^{#pi^{0}} (GeV/c)");
321 fhPtPi0DecayRatio->SetYTitle("p_{T}^{Decay}/p_{T}^{#pi^{0}}");
322
323 fhDeltaPhiDecayCharged = new TH2F
324 ("DeltaPhiDecayCharged","#phi_{Decay} - #phi_{h^{#pm}} vs p_{T Decay}",
325 nptbins,ptmin,ptmax,140,-2.,5.);
326 fhDeltaPhiDecayCharged->SetYTitle("#Delta #phi");
327 fhDeltaPhiDecayCharged->SetXTitle("p_{T Decay} (GeV/c)");
328
329 fhPtImbalanceDecayCharged =
330 new TH2F("CorrelationDecayCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T Decay}",
331 nptbins,ptmin,ptmax,200,0.,2.);
332 fhPtImbalanceDecayCharged->SetYTitle("z_{decay h^{#pm}}");
333 fhPtImbalanceDecayCharged->SetXTitle("p_{T decay}");
334
335 outputContainer->Add(fhPtPi0DecayRatio) ;
336 outputContainer->Add(fhDeltaPhiDecayCharged) ;
337 outputContainer->Add(fhPtImbalanceDecayCharged) ;
338 }
339
340
341 if(fMakeSeveralUE){
342 fhDeltaPhiUeLeftCharged = new TH2F
343 ("DeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE left side range of trigger particles",
344 nptbins,ptmin,ptmax,140,-2.,5.);
345 fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #phi");
346 fhDeltaPhiUeLeftCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
347 outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
348
349 fhDeltaPhiUeRightCharged = new TH2F
350 ("DeltaPhiUeRightChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE right side range of trigger particles",
351 nptbins,ptmin,ptmax,140,-2.,5.);
352 fhDeltaPhiUeRightCharged->SetYTitle("#Delta #phi");
353 fhDeltaPhiUeRightCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
354 outputContainer->Add(fhDeltaPhiUeRightCharged) ;
355
356 fhPtImbalanceUeLeftCharged =
357 new TH2F("CorrelationUeChargedLeft","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE left side of trigger",
358 nptbins,ptmin,ptmax,200,0.,2.);
359 fhPtImbalanceUeLeftCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
360 fhPtImbalanceUeLeftCharged->SetXTitle("p_{T trigger}");
361 outputContainer->Add(fhPtImbalanceUeLeftCharged) ;
362
363 fhPtImbalanceUeRightCharged =
364 new TH2F("CorrelationUeChargedRight","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE right side of trigger",
365 nptbins,ptmin,ptmax,200,0.,2.);
366 fhPtImbalanceUeRightCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
367 fhPtImbalanceUeRightCharged->SetXTitle("p_{T trigger}");
368 outputContainer->Add(fhPtImbalanceUeRightCharged) ;
369
370 fhPtHbpUeLeftCharged =
371 new TH2F("HbpUeChargedLeft","#xi = ln(1/x_{E}) with charged UE left side of trigger",
372 nptbins,ptmin,ptmax,200,0.,10.);
373 fhPtHbpUeLeftCharged->SetYTitle("ln(1/x_{E})");
374 fhPtHbpUeLeftCharged->SetXTitle("p_{T trigger}");
375 outputContainer->Add(fhPtHbpUeLeftCharged) ;
376
377 fhPtHbpUeRightCharged =
378 new TH2F("HbpUeChargedRight","#xi = ln(1/x_{E}) with charged UE right side of trigger",
379 nptbins,ptmin,ptmax,200,0.,10.);
380 fhPtHbpUeRightCharged->SetYTitle("ln(1/x_{E})");
381 fhPtHbpUeRightCharged->SetXTitle("p_{T trigger}");
382 outputContainer->Add(fhPtHbpUeRightCharged) ;
383
384 }
385 } //Correlation with charged hadrons
386
387 //Correlation with neutral hadrons
388 if(fNeutralCorr){
389
390 fhDeltaPhiDeltaEtaNeutral = new TH2F
391 ("DeltaPhiDeltaEtaNeutral","#phi_{trigger} - #phi_{h^{0}} vs #eta_{trigger} - #eta_{h^{0}}",
392 140,-2.,5.,200,-2,2);
393 fhDeltaPhiDeltaEtaNeutral->SetXTitle("#Delta #phi");
394 fhDeltaPhiDeltaEtaNeutral->SetYTitle("#Delta #eta");
395
396 fhPhiNeutral = new TH2F
397 ("PhiNeutral","#phi_{#pi^{0}} vs p_{T #pi^{0}}",
398 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
399 fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
400 fhPhiNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
401
402 fhEtaNeutral = new TH2F
403 ("EtaNeutral","#eta_{#pi^{0}} vs p_{T #pi^{0}}",
404 nptbins,ptmin,ptmax,netabins,etamin,etamax);
405 fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
406 fhEtaNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
407
408 fhDeltaPhiNeutral = new TH2F
409 ("DeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T trigger}",
410 nptbins,ptmin,ptmax,nphibins,phimin,phimax);
411 fhDeltaPhiNeutral->SetYTitle("#Delta #phi");
412 fhDeltaPhiNeutral->SetXTitle("p_{T trigger} (GeV/c)");
413
414 fhDeltaPhiNeutralPt = new TH2F
415 ("DeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
416 nptbins,ptmin,ptmax,140,-2.,5.);
417 fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi");
418 fhDeltaPhiNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
419
420 fhDeltaPhiUeNeutralPt = new TH2F
421 ("DeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
422 nptbins,ptmin,ptmax,140,-2.,5.);
423 fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi");
424 fhDeltaPhiUeNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
425
426 fhDeltaEtaNeutral = new TH2F
427 ("DeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs p_{T trigger}",
428 nptbins,ptmin,ptmax,200,-2,2);
429 fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
430 fhDeltaEtaNeutral->SetXTitle("p_{T trigger} (GeV/c)");
431
432 fhPtImbalanceNeutral =
433 new TH2F("CorrelationNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
434 nptbins,ptmin,ptmax,200,0.,2.);
435 fhPtImbalanceNeutral->SetYTitle("z_{trigger #pi^{0}}");
436 fhPtImbalanceNeutral->SetXTitle("p_{T trigger}");
437
438 fhPtImbalanceUeNeutral =
439 new TH2F("CorrelationUeNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger}",
440 nptbins,ptmin,ptmax,200,0.,2.);
441 fhPtImbalanceUeNeutral->SetYTitle("z_{trigger #pi^{0}}");
442 fhPtImbalanceUeNeutral->SetXTitle("p_{T trigger}");
443
444 fhPtHbpNeutral =
445 new TH2F("HbpNeutral","#xi = ln(1/x_{E}) with neutral particles",
446 nptbins,ptmin,ptmax,200,0.,10.);
447 fhPtHbpNeutral->SetYTitle("ln(1/x_{E})");
448 fhPtHbpNeutral->SetXTitle("p_{T trigger}");
449
450 fhPtHbpUeNeutral =
451 new TH2F("HbpUeNeutral","#xi = ln(1/x_{E}) with neutral particles",
452 nptbins,ptmin,ptmax,200,0.,10.);
453 fhPtHbpUeNeutral->SetYTitle("ln(1/x_{E})");
454 fhPtHbpUeNeutral->SetXTitle("p_{T trigger}");
455
456
457 outputContainer->Add(fhDeltaPhiDeltaEtaNeutral);
458 outputContainer->Add(fhPhiNeutral) ;
459 outputContainer->Add(fhEtaNeutral) ;
460 outputContainer->Add(fhDeltaPhiNeutral) ;
461 outputContainer->Add(fhDeltaPhiNeutralPt) ;
462 outputContainer->Add(fhDeltaPhiUeNeutralPt) ;
463 outputContainer->Add(fhDeltaEtaNeutral) ;
464 outputContainer->Add(fhPtImbalanceNeutral) ;
465 outputContainer->Add(fhPtImbalanceUeNeutral) ;
466 outputContainer->Add(fhPtHbpNeutral) ;
467 outputContainer->Add(fhPtHbpUeNeutral) ;
468
469 if(fPi0Trigger){
470 fhDeltaPhiDecayNeutral = new TH2F
471 ("DeltaPhiDecayNeutral","#phi_{Decay} - #phi_{h^{0}} vs p_{T Decay}",
472 nptbins,ptmin,ptmax,140,-2.,5.);
473 fhDeltaPhiDecayNeutral->SetYTitle("#Delta #phi");
474 fhDeltaPhiDecayNeutral->SetXTitle("p_{T Decay} (GeV/c)");
475
476 fhPtImbalanceDecayNeutral =
477 new TH2F("CorrelationDecayNeutral","z_{trigger h^{0}} = p_{T h^{0}} / p_{T Decay}",
478 nptbins,ptmin,ptmax,200,0.,2.);
479 fhPtImbalanceDecayNeutral->SetYTitle("z_{decay h^{0}}");
480 fhPtImbalanceDecayNeutral->SetXTitle("p_{T decay}");
481
482 outputContainer->Add(fhDeltaPhiDecayNeutral) ;
483 outputContainer->Add(fhPtImbalanceDecayNeutral) ;
484 }
485
486
487 if(fMakeSeveralUE){
488 fhDeltaPhiUeLeftNeutral = new TH2F
489 ("DeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T h^{0}} with neutral UE left side range of trigger particles",
490 nptbins,ptmin,ptmax,140,-2.,5.);
491 fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #phi");
492 fhDeltaPhiUeLeftNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
493 outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
494
495 fhDeltaPhiUeRightNeutral = new TH2F
496 ("DeltaPhiUeRightNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T Ueh^{0}} with neutral UE right side range of trigger particles",
497 nptbins,ptmin,ptmax,140,-2.,5.);
498 fhDeltaPhiUeRightNeutral->SetYTitle("#Delta #phi");
499 fhDeltaPhiUeRightNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
500 outputContainer->Add(fhDeltaPhiUeRightNeutral) ;
501
502 fhPtImbalanceUeLeftNeutral =
503 new TH2F("CorrelationUeNeutralLeft","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE left side of trigger",
504 nptbins,ptmin,ptmax,140,0.,2.);
505 fhPtImbalanceUeLeftNeutral->SetYTitle("z_{trigger Ueh^{0}}");
506 fhPtImbalanceUeLeftNeutral->SetXTitle("p_{T trigger}");
507 outputContainer->Add(fhPtImbalanceUeLeftNeutral) ;
508
509 fhPtImbalanceUeRightNeutral =
510 new TH2F("CorrelationUeNeutralRight","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE right side of trigger",
511 nptbins,ptmin,ptmax,200,0.,2.);
512 fhPtImbalanceUeRightNeutral->SetYTitle("z_{trigger Ueh^{0}}");
513 fhPtImbalanceUeRightNeutral->SetXTitle("p_{T trigger}");
514 outputContainer->Add(fhPtImbalanceUeRightNeutral) ;
515
516 fhPtHbpUeLeftNeutral =
517 new TH2F("HbpUeNeutralLeft","#xi = ln(1/x_{E}) with neutral UE left side of trigger",
518 nptbins,ptmin,ptmax,200,0.,10.);
519 fhPtHbpUeLeftNeutral->SetYTitle("ln(1/x_{E})");
520 fhPtHbpUeLeftNeutral->SetXTitle("p_{T trigger}");
521 outputContainer->Add(fhPtHbpUeLeftNeutral) ;
522
523 fhPtHbpUeRightNeutral =
524 new TH2F("HbpUeNeutralRight","#xi = ln(1/x_{E}) with neutral UE right side of trigger",
525 nptbins,ptmin,ptmax,200,0.,10.);
526 fhPtHbpUeRightNeutral->SetYTitle("ln(1/x_{E})");
527 fhPtHbpUeRightNeutral->SetXTitle("p_{T trigger}");
528 outputContainer->Add(fhPtHbpUeRightNeutral) ;
529
530 }
531
045396c8 532 //if data is MC, fill more histograms
533 if(IsDataMC()){
534 fh2phiLeadingParticle=new TH2F("fh2phiLeadingParticle","#phi resolustion for trigger particles",nptbins,ptmin,ptmax,100,-1,1);
535 fh2phiLeadingParticle->GetXaxis()->SetTitle("p_{T gen Leading} (GeV/c)");
536 fh2phiLeadingParticle->GetYaxis()->SetTitle("(#phi_{rec}-#phi_{gen})/#phi_{gen}");
537
538 fhMCLeadingCount=new TH1F("MCLeadingTriggerCount","MCLeadingTriggerCount",nptbins,ptmin,ptmax);
539 fhMCLeadingCount->SetXTitle("p_{T trig}");
540
541 fhMCEtaCharged = new TH2F
542 ("MCEtaCharged","MC #eta_{h^{#pm}} vs p_{T #pm}",
543 nptbins,ptmin,ptmax,netabins,etamin,etamax);
544 fhMCEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
545 fhMCEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
546
547 fhMCPhiCharged = new TH2F
548 ("MCPhiCharged","#MC phi_{h^{#pm}} vs p_{T #pm}",
549 200,ptmin,ptmax,nphibins,phimin,phimax);
550 fhMCPhiCharged->SetYTitle("MC #phi_{h^{#pm}} (rad)");
551 fhMCPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
552
553 fhMCDeltaPhiDeltaEtaCharged = new TH2F
554 ("MCDeltaPhiDeltaEtaCharged","#MC phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
555 140,-2.,5.,200,-2,2);
556 fhMCDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
557 fhMCDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
558
559 fhMCDeltaEtaCharged = new TH2F
560 ("MCDeltaEtaCharged","MC #eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger} and p_{T assoc}",
561 nptbins,ptmin,ptmax,200,-2,2);
562 fhMCDeltaEtaCharged->SetYTitle("#Delta #eta");
563 fhMCDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
564
565 fhMCDeltaPhiCharged = new TH2F
566 ("MCDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
567 nptbins,ptmin,ptmax,140,-2.,5.);
568 fhMCDeltaPhiCharged->SetYTitle("#Delta #phi");
569 fhMCDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
570
571 fhMCDeltaPhiChargedPt = new TH2F
572 ("MCDeltaPhiChargedPt","MC #phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
573 nptbins,ptmin,ptmax,140,-2.,5.);
574 fhMCDeltaPhiChargedPt->SetYTitle("#Delta #phi");
575 fhMCDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
576
577 fhMCPtImbalanceCharged =
578 new TH2F("MCCorrelationCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T trigger}",
579 nptbins,ptmin,ptmax,200,0.,2.);
580 fhMCPtImbalanceCharged->SetYTitle("z_{trigger h^{#pm}}");
581 fhMCPtImbalanceCharged->SetXTitle("p_{T trigger}");
582
583 fhMCPtHbpCharged =
584 new TH2F("MCHbpCharged","MC #xi = ln(1/x_{E}) with charged hadrons",
585 nptbins,ptmin,ptmax,200,0.,10.);
586 fhMCPtHbpCharged->SetYTitle("ln(1/x_{E})");
587 fhMCPtHbpCharged->SetXTitle("p_{T trigger}");
588
589 fhMCPtTrigPout =
590 new TH2F("MCPtTrigPout","AOD MC Pout with triggers",
591 nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax);
592 fhMCPtTrigPout->SetYTitle("p_{out} (GeV/c)");
593 fhMCPtTrigPout->SetXTitle("p_{T trigger} (GeV/c)");
594
595 fhMCPtAssocDeltaPhi =
596 new TH2F("fhMCPtAssocDeltaPhi","AOD MC delta phi with associated charged hadrons",
597 nptbins,ptmin,ptmax,140,-2.,5.);
598 fhMCPtAssocDeltaPhi->SetYTitle("#Delta #phi");
599 fhMCPtAssocDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)");
600
601 outputContainer->Add(fh2phiLeadingParticle);
602 outputContainer->Add(fhMCLeadingCount);
603 outputContainer->Add(fhMCDeltaPhiDeltaEtaCharged);
604 outputContainer->Add(fhMCPhiCharged) ;
605 outputContainer->Add(fhMCEtaCharged) ;
606 outputContainer->Add(fhMCDeltaEtaCharged) ;
607 outputContainer->Add(fhMCDeltaPhiCharged) ;
608
609 outputContainer->Add(fhMCDeltaPhiChargedPt) ;
610 outputContainer->Add(fhMCPtImbalanceCharged) ;
611 outputContainer->Add(fhMCPtHbpCharged) ;
612 outputContainer->Add(fhMCPtTrigPout) ;
613 outputContainer->Add(fhMCPtAssocDeltaPhi) ;
614 } //for MC histogram
615
e09cf5ef 616 //Keep neutral meson selection histograms if requiered
617 //Setting done in AliNeutralMesonSelection
618 if(GetNeutralMesonSelection()){
619 TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ;
620 if(GetNeutralMesonSelection()->AreNeutralMesonSelectionHistosKept())
621 for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ;
622 delete nmsHistos;
623 }
624
625 }//Correlation with neutral hadrons
626
627 return outputContainer;
628
629}
630
631//____________________________________________________________________________
632void AliAnaParticleHadronCorrelation::InitParameters()
633{
634
635 //Initialize the parameters of the analysis.
636 SetInputAODName("PWG4Particle");
637 SetAODObjArrayName("Hadrons");
638 AddToHistogramsName("AnaHadronCorr_");
045396c8 639
e09cf5ef 640 SetPtCutRange(0.,300);
045396c8 641 fDeltaPhiMinCut = 1.5 ;
642 fDeltaPhiMaxCut = 4.5 ;
643 fSelectIsolated = kFALSE;
644 fMakeSeveralUE = kFALSE;
645 fUeDeltaPhiMinCut = 1. ;
646 fUeDeltaPhiMaxCut = 1.5 ;
647 fNeutralCorr = kFALSE ;
648 fPi0Trigger = kFALSE ;
649 fMakeAbsoluteLeading = kTRUE;
650
651 fNTriggerPtBins = 6 ;
652 fTriggerPtBinLimit[0] = 4.;
653 fTriggerPtBinLimit[1] = 5.;
654 fTriggerPtBinLimit[2] = 6.;
655 fTriggerPtBinLimit[3] = 7.;
656 fTriggerPtBinLimit[4] = 8.;
657 fTriggerPtBinLimit[5] =10.;
658 fTriggerPtBinLimit[6] =12.;
659 fTriggerPtBinLimit[7] =15.;
660 fTriggerPtBinLimit[8] =20.;
661 fTriggerPtBinLimit[9] =25.;
662 fTriggerPtBinLimit[6] =100.;
663
e09cf5ef 664}
665
666//__________________________________________________________________
667void AliAnaParticleHadronCorrelation::Print(const Option_t * opt) const
668{
669
670 //Print some relevant parameters set for the analysis
671 if(! opt)
672 return;
673
674 printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
675 AliAnaPartCorrBaseClass::Print(" ");
676
677 printf("Phi trigger particle-Hadron < %3.2f\n", fDeltaPhiMaxCut) ;
678 printf("Phi trigger particle-Hadron > %3.2f\n", fDeltaPhiMinCut) ;
679 printf("Isolated Trigger? %d\n", fSelectIsolated) ;
680 printf("Phi trigger particle-UeHadron < %3.2f\n", fUeDeltaPhiMaxCut) ;
681 printf("Phi trigger particle-UeHadron > %3.2f\n", fUeDeltaPhiMinCut) ;
682 printf("Several UE? %d\n", fMakeSeveralUE) ;
683 printf("Name of AOD Pi0 Branch %s \n",fPi0AODBranchName.Data());
684 printf("Do Decay-hadron correlation ? %d\n", fPi0Trigger) ;
045396c8 685 printf("Select absolute leading for cluster triggers ? %d\n", fMakeAbsoluteLeading) ;
686 printf("Trigger pt bins %d\n", fNTriggerPtBins) ;
687 for (Int_t ibin = 0; ibin<fNTriggerPtBins; ibin++) {
688 printf("\t bin %d = [%2.1f,%2.1f]\n", ibin, fTriggerPtBinLimit[ibin], fTriggerPtBinLimit[ibin+1]) ;
689 }
e09cf5ef 690
e09cf5ef 691}
692
693//__________________________________________________________________
694TObjString* AliAnaParticleHadronCorrelation::GetAnalysisCuts()
695{
696 //Save parameters used for analysis
697 TString parList ; //this will be list of parameters used for this analysis.
045396c8 698 const Int_t buffersize = 560;
e09cf5ef 699 char onePar[buffersize] ;
045396c8 700
e09cf5ef 701 snprintf(onePar,buffersize,"--- AliAnaPaticleHadronCorrelation ---\n") ;
702 parList+=onePar ;
703 snprintf(onePar,buffersize,"Phi trigger particle-Hadron < %3.2f ", fDeltaPhiMaxCut) ;
704 parList+=onePar ;
705 snprintf(onePar,buffersize,"Phi trigger particle-Hadron > %3.2f ", fDeltaPhiMinCut) ;
706 parList+=onePar ;
707 snprintf(onePar,buffersize,"Isolated Trigger? %d\n", fSelectIsolated) ;
708 parList+=onePar ;
709 snprintf(onePar,buffersize,"Phi trigger particle-UeHadron < %3.2f ", fUeDeltaPhiMaxCut) ;
710 parList+=onePar ;
711 snprintf(onePar,buffersize,"Phi trigger particle-UeHadron > %3.2f ", fUeDeltaPhiMinCut) ;
712 parList+=onePar ;
713 snprintf(onePar,buffersize,"Several UE? %d\n", fMakeSeveralUE) ;
714 parList+=onePar ;
715 snprintf(onePar,buffersize,"Name of AOD Pi0 Branch %s ",fPi0AODBranchName.Data());
716 parList+=onePar ;
717 snprintf(onePar,buffersize,"Do Decay-hadron correlation ? %d", fPi0Trigger) ;
718 parList+=onePar ;
045396c8 719 snprintf(onePar,buffersize,"Select absolute leading for cluster triggers ? %d\n", fMakeAbsoluteLeading) ;
720 parList+=onePar ;
721 snprintf(onePar,buffersize,"Trigger pt bins %d: ", fNTriggerPtBins) ;
722 parList+=onePar ;
723 for (Int_t ibin = 0; ibin<fNTriggerPtBins; ibin++) {
724 snprintf(onePar,buffersize,"bin %d = [%2.1f,%2.1f];", ibin, fTriggerPtBinLimit[ibin], fTriggerPtBinLimit[ibin+1]) ;
725 }
726 parList+=onePar ;
727
e09cf5ef 728 //Get parameters set in base class.
729 parList += GetBaseParametersList() ;
730
731 //Get parameters set in PID class.
732 //parList += GetCaloPID()->GetPIDParametersList() ;
733
734 //Get parameters set in FiducialCut class (not available yet)
735 //parlist += GetFidCut()->GetFidCutParametersList()
736
737 return new TObjString(parList) ;
045396c8 738
e09cf5ef 739}
740
741
742//____________________________________________________________________________
743void AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD()
744{
745 //Particle-Hadron Correlation Analysis, fill AODs
746
747 if(!GetInputAODBranch()){
748 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
749 abort();
750 }
751
752 if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation")){
753 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());
754 abort();
755 }
756
757 if(GetDebug() > 1){
758 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Begin hadron correlation analysis, fill AODs \n");
759 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
760 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n", GetCTSTracks()->GetEntriesFast());
761 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetEMCALClusters()->GetEntriesFast());
762 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n", GetPHOSClusters()->GetEntriesFast());
763 }
764
045396c8 765 //Get the vertex and check it is not too large in z
766 Double_t v[3] = {0,0,0}; //vertex ;
767 GetReader()->GetVertex(v);
768 if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;
769
770 //Loop on stored AOD particles, find leading trigger
771 Double_t ptTrig = GetMinPt() ;
772 fLeadingTriggerIndex = -1 ;
773 Int_t naod = GetInputAODBranch()->GetEntriesFast() ;
e09cf5ef 774 for(Int_t iaod = 0; iaod < naod ; iaod++){
775 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
045396c8 776
777 //vertex cut in case of mixing
778 if (GetMixedEvent()) {
779 Int_t evt=-1;
780 Int_t id =-1;
781 if (particle->GetCaloLabel(0) >= 0 ){
782 id=particle->GetCaloLabel(0);
783 if(id >= 0 )evt= GetMixedEvent()-> EventIndexForCaloCluster(id) ;
784 }
785 else if(particle->GetTrackLabel(0) >= 0 ){
786 id=particle->GetTrackLabel(0);
787 if(id >= 0 )evt= GetMixedEvent()->EventIndex(id) ;
788 }
789 else continue;
790
791 if (TMath::Abs(GetVertex(evt)[2]) > GetZvertexCut())
792 return ;
e09cf5ef 793 }
045396c8 794
795 // find the leading particles with highest momentum
796 if (particle->Pt() > ptTrig) {
797 ptTrig = particle->Pt() ;
798 fLeadingTriggerIndex = iaod ;
799 }
800 }// finish search of leading trigger particle
e09cf5ef 801
802 //Do correlation with leading particle
045396c8 803 if(fLeadingTriggerIndex >= 0){
e09cf5ef 804
045396c8 805 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
806
807 //check if the particle is isolated or if we want to take the isolation into account
808 if(OnlyIsolated() && !particle->IsIsolated()) return;
809
e09cf5ef 810 //Make correlation with charged hadrons
045396c8 811 Bool_t okcharged = kTRUE;
812 Bool_t okneutral = kTRUE;
e09cf5ef 813 if(GetReader()->IsCTSSwitchedOn() )
045396c8 814 okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kFALSE);
e09cf5ef 815
816 TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
817 if(fNeutralCorr && pi0list && pi0list->GetEntriesFast() > 0)
045396c8 818 okneutral = MakeNeutralCorrelation(particle, pi0list,kFALSE);
e09cf5ef 819
820 }//Correlate leading
821
822 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
823
824}
825
826//____________________________________________________________________________
827void AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms()
828{
829 //Particle-Hadron Correlation Analysis, fill histograms
830
831 if(!GetInputAODBranch()){
832 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
833 abort();
834 }
835
836 if(GetDebug() > 1){
837 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Begin hadron correlation analysis, fill histograms \n");
838 printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
839 }
840
841 //Get the vertex and check it is not too large in z
842 Double_t v[3] = {0,0,0}; //vertex ;
843 GetReader()->GetVertex(v);
844 if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;
845
846 //Loop on stored AOD particles, find leading
045396c8 847 Double_t ptTrig = GetMinPt() ;
848 if(fLeadingTriggerIndex < 0){//Search leading if not done before
849 Int_t naod = GetInputAODBranch()->GetEntriesFast() ;
850 for(Int_t iaod = 0; iaod < naod ; iaod++){ //loop on input trigger AOD file
851 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
852 //vertex cut in case of mixing
853 if (GetMixedEvent()) {
854 Int_t evt=-1;
855 Int_t id =-1;
856 if (particle->GetCaloLabel(0) >= 0 ){
857 id=particle->GetCaloLabel(0);
858 if(id >= 0 )evt= GetMixedEvent()-> EventIndexForCaloCluster(id) ;
859 }
860 else if(particle->GetTrackLabel(0) >= 0 ){
861 id=particle->GetTrackLabel(0);
862 if(id >= 0 )evt= GetMixedEvent()->EventIndex(id) ;
863 }
864 else continue;
865
866 if (TMath::Abs(GetVertex(evt)[2]) > GetZvertexCut())
867 return ;
e09cf5ef 868 }
045396c8 869 //check if the particle is isolated or if we want to take the isolation into account
870 if(OnlyIsolated() && !particle->IsIsolated()) continue;
e09cf5ef 871
045396c8 872 //check if inside the vertex cut
873 //find the leading particles with highest momentum
874 if (particle->Pt() > ptTrig) {
875 ptTrig = particle->Pt() ;
876 fLeadingTriggerIndex = iaod ;
e09cf5ef 877 }
045396c8 878 }//finish search of leading trigger particle
879 }//Search leading if not done before
880
881 if(fLeadingTriggerIndex >= 0 ){ //using trigger particle to do correlations
882 AliAODPWG4ParticleCorrelation* particle = (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
883
884 //check if the particle is isolated or if we want to take the isolation into account
885 if(OnlyIsolated() && !particle->IsIsolated()) return;
e09cf5ef 886
887 //Make correlation with charged hadrons
045396c8 888 Bool_t okcharged = kTRUE;
889 Bool_t okneutral = kTRUE;
890 if(GetReader()->IsCTSSwitchedOn() ){
891 okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kTRUE);
a5f1e836 892 if(IsDataMC()){
045396c8 893 MakeMCChargedCorrelation(particle);
894 }
a5f1e836 895 }
896
e09cf5ef 897 TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
a5f1e836 898 if(fNeutralCorr && pi0list){
899 if(pi0list->GetEntriesFast() > 0)
900 okneutral = MakeNeutralCorrelation(particle, pi0list,kTRUE);
901 }
e09cf5ef 902
045396c8 903 // Fill leading particle histogram if correlation went well and
904 // no problem was found, like not absolute leading, or bad vertex in mixing.
905 if(okcharged && okneutral){
906 fhPtLeading->Fill(particle->Pt());
907 Float_t phi = particle->Phi();
908 if(phi<0)phi+=TMath::TwoPi();
909 fhPhiLeading->Fill(particle->Pt(), phi);
910 fhEtaLeading->Fill(particle->Pt(), particle->Eta());
911 }//ok charged && neutral
e09cf5ef 912 }//Aod branch loop
913
045396c8 914 //Reinit for next event
915 fLeadingTriggerIndex = -1;
e09cf5ef 916
045396c8 917 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
e09cf5ef 918}
919
920//____________________________________________________________________________
045396c8 921Bool_t AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle,
922 TObjArray* pl, const Bool_t bFillHisto)
e09cf5ef 923{
924 // Charged Hadron Correlation Analysis
925 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
926
045396c8 927 Int_t evtIndex11 = -1 ; //cluster trigger or pi0 trigger
928 Int_t evtIndex12 = -1 ; // pi0 trigger
929 Int_t evtIndex13 = -1 ; // charged trigger
e09cf5ef 930 Int_t indexPhoton1 = -1 ;
931 Int_t indexPhoton2 = -1 ;
045396c8 932
933 Double_t v[3] = {0,0,0}; //vertex ;
e09cf5ef 934 GetReader()->GetVertex(v);
e09cf5ef 935
936 if (GetMixedEvent()) {
937 evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
938 evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
939 evtIndex13 = GetMixedEvent()->EventIndex(aodParticle->GetTrackLabel(0)) ;
940 }
941
e09cf5ef 942 Double_t phiTrig = aodParticle->Phi();
045396c8 943 Double_t etaTrig = aodParticle->Eta();
944 Double_t ptTrig = aodParticle->Pt();
945 Int_t trBin = -1;
946 for(Int_t i = 0 ; i < fNTriggerPtBins ; i++){
947 if(ptTrig > fTriggerPtBinLimit[i] && ptTrig < fTriggerPtBinLimit[i+1]) trBin= i;
948 }
949 if(trBin==-1) return kFALSE;
950 if(OnlyIsolated() && aodParticle->IsIsolated()){
951 fhTrigPt[trBin]->Fill(ptTrig);
952 }
953
954 Double_t pt = -100. ;
955 Double_t px = -100. ;
956 Double_t py = -100. ;
957 Double_t rat = -100. ;
958 Double_t xE = -100. ;
959 Double_t cosi = -100. ;
960 Double_t phi = -100. ;
961 Double_t eta = -100. ;
962 Double_t pout = -100. ;
e09cf5ef 963
045396c8 964 Double_t ptDecay1 = 0. ;
965 Double_t pxDecay1 = 0. ;
966 Double_t pyDecay1 = 0. ;
967 Double_t phiDecay1 = 0. ;
968 Double_t ptDecay2 = 0. ;
969 Double_t pxDecay2 = 0. ;
970 Double_t pyDecay2 = 0. ;
971 Double_t phiDecay2 = 0. ;
972
973 Double_t ratDecay1 = -100. ;
974 Double_t ratDecay2 = -100. ;
975 Double_t deltaPhi = -100. ;
976 Double_t deltaPhiOrg = -100. ;
977 Double_t deltaPhiDecay1 = -100. ;
978 Double_t deltaPhiDecay2 = -100. ;
979
980 TVector3 p3;
e09cf5ef 981 TLorentzVector photonMom ;
045396c8 982 TObjArray * clusters = 0x0 ;
983 TObjArray * reftracks = 0x0;
984 Int_t nrefs = 0;
985 Int_t nTracks = GetCTSTracks()->GetEntriesFast() ;
e09cf5ef 986
987 if(fPi0Trigger){
988 indexPhoton1 = aodParticle->GetCaloLabel (0);
989 indexPhoton2 = aodParticle->GetCaloLabel (1);
990 if(GetDebug() > 1)printf("indexPhoton1 = %d, indexPhoton2 = %d \n", indexPhoton1, indexPhoton2);
991
992 if(indexPhoton1!=-1 && indexPhoton2!=-1){
993 if(aodParticle->GetDetector()=="EMCAL") clusters = GetEMCALClusters() ;
045396c8 994 else clusters = GetPHOSClusters() ;
e09cf5ef 995 for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++){
996 AliVCluster * photon = (AliVCluster*) (clusters->At(iclus));
997 photon->GetMomentum(photonMom,GetVertex(0)) ;
998 if(photon->GetID()==indexPhoton1) {
999 ptDecay1 = photonMom.Pt();
1000 pxDecay1 = photonMom.Px();
1001 pyDecay1 = photonMom.Py();
1002 phiDecay1 = photonMom.Phi();
1003 if(ptTrig && bFillHisto) fhPtPi0DecayRatio->Fill(ptTrig, ptDecay1/ptTrig);
1004 }
1005 if(photon->GetID()==indexPhoton2) {
1006 ptDecay2 = photonMom.Pt();
1007 pxDecay2 = photonMom.Px();
1008 pyDecay2 = photonMom.Py();
1009 phiDecay2 = photonMom.Phi();
1010 if(ptTrig && bFillHisto) fhPtPi0DecayRatio->Fill(ptTrig, ptDecay2/ptTrig);
1011 }
1012 if(GetDebug() > 1)printf("Photon1 = %f, Photon2 = %f \n", ptDecay1, ptDecay2);
1013 } //cluster loop
1014 } //index of decay photons found
e09cf5ef 1015 } //make decay-hadron correlation
1016
1017 //Track loop, select tracks with good pt, phi and fill AODs or histograms
e09cf5ef 1018 for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ ){
1019 AliVTrack * track = (AliVTrack *) (pl->At(ipr)) ;
045396c8 1020
e09cf5ef 1021 Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
1022 p3.SetXYZ(mom[0],mom[1],mom[2]);
1023 pt = p3.Pt();
1024 px = p3.Px();
1025 py = p3.Py();
1026 eta = p3.Eta();
1027 phi = p3.Phi() ;
1028 if(phi < 0) phi+=TMath::TwoPi();
045396c8 1029
e09cf5ef 1030 //Select only hadrons in pt range
1031 if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
1032 //remove trigger itself for correlation when use charged triggers
1033 if( track->GetID() == aodParticle->GetTrackLabel(0) || track->GetID() == aodParticle->GetTrackLabel(1) ||
045396c8 1034 track->GetID() == aodParticle->GetTrackLabel(2) || track->GetID() == aodParticle->GetTrackLabel(3) )
1035 continue ;
1036
e09cf5ef 1037 if(IsFiducialCutOn()){
1038 Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,"CTS") ;
1039 if(! in ) continue ;
1040 }
045396c8 1041
1042 //jump out this event if near side associated particle pt larger than trigger
1043 if (fMakeAbsoluteLeading){
1044 if(pt > ptTrig && TMath::Abs(phi-phiTrig)<TMath::PiOver2()) return kFALSE;
1045 }
1046
1047 //Only for mixed event
1048 Int_t evtIndex2 = 0 ;
1049 if (GetMixedEvent()) {
1050 evtIndex2 = GetMixedEvent()->EventIndex(track->GetID()) ;
1051 if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 || evtIndex13 == evtIndex2 ) // photon and track from different events
1052 continue ;
1053 //vertex cut
1054 if (TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut())
1055 return kFALSE;
1056 }
1057
1058 if(fPi0Trigger){
e09cf5ef 1059 if(indexPhoton1!=-1 && indexPhoton2!=-1){
1060 if(ptDecay1) ratDecay1 = pt/ptDecay1 ;
1061 if(ptDecay2) ratDecay2 = pt/ptDecay2 ;
045396c8 1062 deltaPhiDecay1 = phiDecay1-phi;
1063 deltaPhiDecay2 = phiDecay2-phi;
1064 if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
1065 if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
1066 if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
1067 if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
e09cf5ef 1068 }
1069 } //do decay-hadron correlation
045396c8 1070
e09cf5ef 1071 //Selection within angular range
045396c8 1072 deltaPhi = phiTrig-phi;
1073 deltaPhiOrg = deltaPhi;
1074 if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
1075 if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
1076
1077 pout = pt*TMath::Sin(deltaPhi) ;
1078 rat = pt/ptTrig ;
1079 xE =-pt/ptTrig*TMath::Cos(deltaPhi);
c514f286 1080 if(xE > 0 )cosi = TMath::Log(1/xE);
1081 else {
1082 printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - xE=%f negative or null, check!, pT %f, ptTrig %f\n", xE,pt,ptTrig);
1083 return kFALSE;
1084 }
e09cf5ef 1085
1086 if(GetDebug() > 2)
1087 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",
045396c8 1088 pt,phi, phiTrig,fDeltaPhiMinCut, deltaPhi, fDeltaPhiMaxCut, GetMinPt());
e09cf5ef 1089
045396c8 1090 // Fill Histograms
e09cf5ef 1091 if(bFillHisto){
045396c8 1092
1093 if(OnlyIsolated() && aodParticle->IsIsolated()){ //FILIP, if only wantisolated histograms
1094 if(TMath::Cos(deltaPhi)<0){//away side
1095 fhXE[trBin]->Fill(xE) ;
1096 }
1097
1098 //Hardcoded values, BAD, FIXME
1099 Double_t dphiBrad = atan2(sin(deltaPhiOrg), cos(deltaPhiOrg))/TMath::Pi();//-1 to 1
1100 if(TMath::Abs(dphiBrad)>0.325 && TMath::Abs(dphiBrad)<0.475){
1101 fhAssocPtBkg[trBin]->Fill(pt);
1102 }
1103 if(dphiBrad<-1./3) dphiBrad += 2;
1104 fhDphiTrigPtAssocPt[trBin]->Fill(dphiBrad);
1105 }
1106
e09cf5ef 1107 fhEtaCharged->Fill(pt,eta);
1108 fhPhiCharged->Fill(pt,phi);
1109 fhDeltaEtaCharged->Fill(ptTrig,aodParticle->Eta()-eta);
045396c8 1110 fhDeltaPhiCharged->Fill(ptTrig, deltaPhi);
1111 fhDeltaPhiDeltaEtaCharged->Fill(deltaPhi,aodParticle->Eta()-eta);
1112 fhPtAssocDeltaPhi->Fill(pt, deltaPhi);
e09cf5ef 1113
1114 if(GetDebug() > 2 ) printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
1115 //fill different multiplicity histogram
1116 if(DoEventSelect()){
1117 for(Int_t im=0; im<GetMultiBin(); im++){
1118 if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1)){
045396c8 1119 fhTrigDeltaPhiCharged[im]->Fill(ptTrig,deltaPhi);
e09cf5ef 1120 fhTrigDeltaEtaCharged[im]->Fill(ptTrig,aodParticle->Eta()-eta);
045396c8 1121 }
e09cf5ef 1122 }
1123 }
1124 //delta phi cut for correlation
045396c8 1125 if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) ) {
1126 fhDeltaPhiChargedPt->Fill(pt,deltaPhi);
e09cf5ef 1127 fhPtImbalanceCharged->Fill(ptTrig,xE);
1128 fhPtHbpCharged->Fill(ptTrig,cosi);
1129 fhPtTrigPout->Fill(ptTrig, pout) ;
1130 fhPtTrigCharged->Fill(ptTrig, pt) ;
1131 if(track->Charge()>0) fhPtImbalancePosCharged->Fill(ptTrig,xE) ;
1132 else fhPtImbalanceNegCharged->Fill(ptTrig,xE) ;
1133 //fill different multiplicity histogram
1134 if(DoEventSelect()){
1135 for(Int_t im=0; im<GetMultiBin(); im++){
1136 if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
1137 fhTrigCorr[im]->Fill(ptTrig,xE);
1138 }
1139 } //multiplicity events selection
1140 } //delta phi cut for correlation
045396c8 1141 else if ((deltaPhi > fUeDeltaPhiMinCut) && ( deltaPhi < fUeDeltaPhiMaxCut)) { //UE study
1142 fhDeltaPhiUeChargedPt->Fill(pt,deltaPhi);
1143 Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
1144 Double_t uexE = -(pt/ptTrig)*TMath::Cos(randomphi);
1145 if(uexE < 0.) uexE = -uexE;
1146 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - xe = %f, uexE = %f \n", xE, uexE);
1147 fhPtImbalanceUeCharged->Fill(ptTrig,uexE);
c514f286 1148 if(uexE>0)fhPtHbpUeCharged->Fill(ptTrig,TMath::Log(1/uexE));
e09cf5ef 1149 if(DoEventSelect()){
1150 for(Int_t im=0; im<GetMultiBin(); im++){
1151 if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
1152 fhTrigUeCorr[im]->Fill(ptTrig,xE);
1153 }
1154 } //multiplicity events selection
1155
1156 } //UE study
1157
1158 if(fPi0Trigger){
1159 if(indexPhoton1!=-1 && indexPhoton2!=-1){
045396c8 1160 fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaPhiDecay1);
1161 fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaPhiDecay2);
1162 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
1163 if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
e09cf5ef 1164 fhPtImbalanceDecayCharged->Fill(ptDecay1,ratDecay1);
045396c8 1165 if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
e09cf5ef 1166 fhPtImbalanceDecayCharged->Fill(ptDecay2,ratDecay2);
045396c8 1167 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - ratPhoton1 = %f, ratPhoton2 = %f \n", pt/ptDecay1, pt/ptDecay2);
e09cf5ef 1168 } //index of decay photons found
1169 } //make decay-hadron correlation
1170
1171 //several UE calculation
1172 if(fMakeSeveralUE){
045396c8 1173 if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut)){
1174 fhDeltaPhiUeLeftCharged->Fill(pt,deltaPhi);
e09cf5ef 1175 fhPtImbalanceUeLeftCharged->Fill(ptTrig,rat);
1176 fhPtHbpUeLeftCharged->Fill(ptTrig,cosi);
1177 }
045396c8 1178 if((deltaPhi>fUeDeltaPhiMinCut) && (deltaPhi <fUeDeltaPhiMaxCut)){
1179 fhDeltaPhiUeRightCharged->Fill(pt,deltaPhi);
e09cf5ef 1180 fhPtImbalanceUeRightCharged->Fill(ptTrig,rat);
1181 fhPtHbpUeRightCharged->Fill(ptTrig,cosi);
1182
1183 }
1184 } //several UE calculation
1185
045396c8 1186 //Fill leading particle histogram
1187 fhPtLeading->Fill(ptTrig);
1188 if(phiTrig<0)phiTrig+=TMath::TwoPi();
1189 fhPhiLeading->Fill(ptTrig, phiTrig);
1190 fhEtaLeading->Fill(ptTrig, etaTrig);
1191
e09cf5ef 1192 } //Fill histogram
1193 else{
1194 nrefs++;
1195 if(nrefs==1){
1196 reftracks = new TObjArray(0);
045396c8 1197 TString trackname = Form("%s+Tracks", GetAODObjArrayName().Data());
1198 reftracks->SetName(trackname.Data());
e09cf5ef 1199 reftracks->SetOwner(kFALSE);
1200 }
1201 reftracks->Add(track);
1202 }//aod particle loop
1203 }// track loop
1204
1205 //Fill AOD with reference tracks, if not filling histograms
1206 if(!bFillHisto && reftracks) {
1207 aodParticle->AddObjArray(reftracks);
1208 }
1209
045396c8 1210 return kTRUE;
e09cf5ef 1211
1212}
1213
e09cf5ef 1214//____________________________________________________________________________
045396c8 1215Bool_t AliAnaParticleHadronCorrelation::MakeNeutralCorrelation(AliAODPWG4ParticleCorrelation * const aodParticle,
1216 TObjArray* pi0list, const Bool_t bFillHisto)
e09cf5ef 1217{
1218 // Neutral Pion Correlation Analysis
1219 if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Make trigger particle - pi0 correlation, %d pi0's \n",pi0list->GetEntriesFast());
1220
1221 Int_t evtIndex11 = 0 ;
1222 Int_t evtIndex12 = 0 ;
1223 if (GetMixedEvent()) {
1224 evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
1225 evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
1226 }
1227
045396c8 1228 Double_t pt = -100. ;
1229 Double_t px = -100. ;
1230 Double_t py = -100. ;
1231 Double_t rat = -100. ;
1232 Double_t phi = -100. ;
1233 Double_t eta = -100. ;
1234 Double_t xE = -100. ;
1235 Double_t cosi = -100. ;
e09cf5ef 1236
1237 Double_t ptTrig = aodParticle->Pt();
1238 Double_t phiTrig = aodParticle->Phi();
1239 Double_t etaTrig = aodParticle->Eta();
1240 Double_t pxTrig = aodParticle->Px();
1241 Double_t pyTrig = aodParticle->Py();
1242
045396c8 1243 Int_t indexPhoton1 =-1 ;
1244 Int_t indexPhoton2 =-1 ;
1245 Double_t ptDecay1 = 0. ;
e09cf5ef 1246 Double_t pxDecay1 = 0. ;
1247 Double_t pyDecay1 = 0. ;
045396c8 1248 Double_t phiDecay1 = 0. ;
1249 Double_t ptDecay2 = 0. ;
e09cf5ef 1250 Double_t pxDecay2 = 0. ;
1251 Double_t pyDecay2 = 0. ;
045396c8 1252 Double_t phiDecay2 = 0. ;
1253
1254 Double_t ratDecay1 = -100. ;
1255 Double_t ratDecay2 = -100. ;
1256 Double_t deltaPhi = -100. ;
1257 Double_t deltaPhiDecay1 = -100. ;
1258 Double_t deltaPhiDecay2 = -100. ;
e09cf5ef 1259
e09cf5ef 1260 TObjArray * clusters = 0x0 ;
045396c8 1261 TLorentzVector photonMom ;
1262
e09cf5ef 1263 if(fPi0Trigger){
1264 indexPhoton1 = aodParticle->GetCaloLabel (0);
1265 indexPhoton2 = aodParticle->GetCaloLabel (1);
045396c8 1266 if(GetDebug() > 1)
1267 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - indexPhoton1 = %d, indexPhoton2 = %d \n", indexPhoton1, indexPhoton2);
e09cf5ef 1268
1269 if(indexPhoton1!=-1 && indexPhoton2!=-1){
1270 if(aodParticle->GetDetector()=="EMCAL") clusters = GetEMCALClusters() ;
1271 else clusters = GetPHOSClusters() ;
1272 for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++){
1273 AliVCluster * photon = (AliVCluster*) (clusters->At(iclus));
1274 photon->GetMomentum(photonMom,GetVertex(0)) ;
1275 if(photon->GetID()==indexPhoton1) {
1276 ptDecay1 = photonMom.Pt();
1277 pxDecay1 = photonMom.Px();
1278 pyDecay1 = photonMom.Py();
1279 phiDecay1 = photonMom.Phi();
1280 }
1281 if(photon->GetID()==indexPhoton2) {
1282 ptDecay2 = photonMom.Pt();
1283 pxDecay2 = photonMom.Px();
1284 pyDecay2 = photonMom.Py();
1285 phiDecay2 = photonMom.Phi();
1286 }
045396c8 1287 if(GetDebug() > 1)
1288 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - Photon1 = %f, Photon2 = %f \n", ptDecay1, ptDecay2);
e09cf5ef 1289 } //photonAOD loop
1290 } //index of decay photons found
1291 if(ptTrig && bFillHisto) fhPtPi0DecayRatio->Fill(ptTrig, ptDecay1/ptTrig, ptDecay2/ptTrig);
1292 } //make decay-hadron correlation
1293
1294 TObjArray * refpi0 =0x0;
1295 Int_t nrefs = 0;
1296
1297 //Loop on stored AOD pi0
1298 Int_t naod = pi0list->GetEntriesFast();
045396c8 1299 if(GetDebug() > 0)
1300 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() - aod branch entries %d\n", naod);
e09cf5ef 1301 for(Int_t iaod = 0; iaod < naod ; iaod++){
1302 AliAODPWG4Particle* pi0 = (AliAODPWG4Particle*) (pi0list->At(iaod));
1303
1304 Int_t evtIndex2 = 0 ;
1305 Int_t evtIndex3 = 0 ;
1306 if (GetMixedEvent()) {
1307 evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(0)) ;
1308 evtIndex3 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(1)) ;
1309
045396c8 1310 if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 ||
1311 evtIndex11 == evtIndex3 || evtIndex12 == evtIndex3) // trigger and pi0 are not from different events
e09cf5ef 1312 continue ;
1313 }
1314
1315 //Int_t pdg = pi0->GetPdg();
1316 //if(pdg != AliCaloPID::kPi0) continue;
1317
1318 pt = pi0->Pt();
1319 px = pi0->Px();
1320 py = pi0->Py();
1321 if(pt < GetMinPt() || pt > GetMaxPt()) continue ;
1322 //jumped out this event if near side associated partile pt larger than trigger
1323 if(pt > ptTrig && TMath::Abs(phi-phiTrig)<TMath::PiOver2()) break ;
045396c8 1324
e09cf5ef 1325 //Selection within angular range
1326 phi = pi0->Phi();
045396c8 1327 //Float_t deltaPhi = TMath::Abs(phiTrig-phi);
1328 //if( (deltaPhi < fDeltaPhiMinCut) || ( deltaPhi > fDeltaPhiMaxCut) ) continue ;
e09cf5ef 1329
1330 if(bFillHisto){
1331
045396c8 1332 deltaPhi = phiTrig-phi;
1333 if(deltaPhi<-TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
1334 if(deltaPhi>3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
e09cf5ef 1335
1336 rat = pt/ptTrig ;
1337 phi = pi0->Phi() ;
1338 eta = pi0->Eta() ;
1339 xE = -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
1340 if(xE <0.)xE =-xE;
1341 cosi = TMath::Log(1/xE);
1342
1343 if(fPi0Trigger){
1344 if(indexPhoton1!=-1 && indexPhoton2!=-1){
1345 if(ptDecay1) ratDecay1 = pt/ptDecay1 ;
1346 if(ptDecay2) ratDecay2 = pt/ptDecay2 ;
045396c8 1347 deltaPhiDecay1 = phiDecay1-phi;
1348 deltaPhiDecay2 = phiDecay2-phi;
1349 if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
1350 if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
1351 if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
1352 if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
1353 fhDeltaPhiDecayNeutral->Fill(ptDecay1, deltaPhiDecay1);
1354 fhDeltaPhiDecayNeutral->Fill(ptDecay2, deltaPhiDecay2);
1355 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
1356 if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
e09cf5ef 1357 fhPtImbalanceDecayNeutral->Fill(ptDecay1,ratDecay1);
045396c8 1358 if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
e09cf5ef 1359 fhPtImbalanceDecayNeutral->Fill(ptDecay2,ratDecay2);
045396c8 1360 if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - ratPhoton1 = %f, ratPhoton2 = %f \n", pt/ptDecay1, pt/ptDecay2);
e09cf5ef 1361 }
1362 } //do decay-hadron correlation
1363
1364 fhEtaNeutral->Fill(pt,eta);
1365 fhPhiNeutral->Fill(pt,phi);
1366 fhDeltaEtaNeutral->Fill(ptTrig,etaTrig-eta);
045396c8 1367 fhDeltaPhiNeutral->Fill(ptTrig,deltaPhi);
1368 fhDeltaPhiDeltaEtaNeutral->Fill(deltaPhi,etaTrig-eta);
e09cf5ef 1369
1370 //delta phi cut for correlation
045396c8 1371 if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) ) {
1372 fhDeltaPhiNeutralPt->Fill(pt,deltaPhi);
e09cf5ef 1373 fhPtImbalanceNeutral->Fill(ptTrig,rat);
1374 fhPtHbpNeutral->Fill(ptTrig,cosi);
1375 }
1376 else {
045396c8 1377 fhDeltaPhiUeNeutralPt->Fill(pt,deltaPhi);
e09cf5ef 1378 fhPtImbalanceUeNeutral->Fill(ptTrig,rat);
1379 fhPtHbpUeNeutral->Fill(ptTrig,cosi);
1380 }
1381 //several UE calculation
1382 if(fMakeSeveralUE){
045396c8 1383 if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut)){
1384 fhDeltaPhiUeLeftNeutral->Fill(pt,deltaPhi);
e09cf5ef 1385 fhPtImbalanceUeLeftNeutral->Fill(ptTrig,rat);
1386 fhPtHbpUeLeftNeutral->Fill(ptTrig,cosi);
1387 }
045396c8 1388 if((deltaPhi>fUeDeltaPhiMinCut) && (deltaPhi <fUeDeltaPhiMaxCut)){
1389 fhDeltaPhiUeRightNeutral->Fill(pt,deltaPhi);
e09cf5ef 1390 fhPtImbalanceUeRightNeutral->Fill(ptTrig,rat);
1391 fhPtHbpUeRightNeutral->Fill(ptTrig,cosi);
1392 }
1393 } //several UE calculation
1394 }
1395 else{
1396 nrefs++;
1397 if(nrefs==1){
1398 refpi0 = new TObjArray(0);
1399 refpi0->SetName(GetAODObjArrayName()+"Pi0s");
1400 refpi0->SetOwner(kFALSE);
1401 }
1402 refpi0->Add(pi0);
1403 }//put references in trigger AOD
045396c8 1404
1405 if(GetDebug() > 2 )
1406 printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
1407
1408 }//loop
1409
1410 return kTRUE;
1411}
1412
1413//____________________________________________________________________________
1414void AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle)
1415{
1416 // Charged Hadron Correlation Analysis with MC information
1417 if(GetDebug()>1)
1418 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Make trigger particle - charged hadron correlation in AOD MC level\n");
1419
1420 AliStack * stack = 0x0 ;
1421 TParticle * primary = 0x0 ;
1422 TClonesArray * mcparticles0 = 0x0 ;
1423 TClonesArray * mcparticles = 0x0 ;
1424 AliAODMCParticle * aodprimary = 0x0 ;
1425
1426 Double_t eprim = 0 ;
1427 Double_t ptprim = 0 ;
1428 Double_t phiprim = 0 ;
1429 Double_t etaprim = 0 ;
1430 Double_t pxprim = 0 ;
1431 Double_t pyprim = 0 ;
1432 Double_t pzprim = 0 ;
1433 Int_t nTracks = 0 ;
1434 Int_t iParticle = 0 ;
1435 Double_t charge = 0.;
1436
1437 Double_t mcrat =-100 ;
1438 Double_t mcxE =-100 ;
1439 Double_t mccosi =-100 ;
1440
1441 //Track loop, select tracks with good pt, phi and fill AODs or histograms
1442 //Int_t currentIndex = -1 ;
1443 Double_t mcTrackPt = 0 ;
1444 Double_t mcTrackPhi = 0 ;
1445 Double_t mcTrackEta = 0 ;
1446 Double_t mcTrackPx = 0 ;
1447 Double_t mcTrackPy = 0 ;
1448 Double_t mcTrackPz = 0 ;
1449
1450 if(GetReader()->ReadStack()){
1451 nTracks = GetMCStack()->GetNtrack() ;
1452 }
1453 else nTracks = GetReader()->GetAODMCParticles()->GetEntriesFast() ;
1454 //Int_t trackIndex[nTracks];
1455
1456 Int_t label= aodParticle->GetLabel();
1457 if(label<0){
1458 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** bad label ***: label %d \n", label);
1459 return;
1460 }
1461
1462 if(GetReader()->ReadStack()){
1463 stack = GetMCStack() ;
1464 if(!stack) {
1465 printf(" AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation- Stack not available, is the MC handler called? STOP\n");
1466 abort();
1467 }
1468
1469 nTracks=stack->GetNprimary();
1470 if(label >= stack->GetNtrack()) {
1471 if(GetDebug() > 2) printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack());
1472 return ;
1473 }
1474 primary = stack->Particle(label);
1475 if(!primary){
1476 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** no primary ***: label %d \n", label);
1477 return;
1478 }
1479
1480 eprim = primary->Energy();
1481 ptprim = primary->Pt();
1482 phiprim = primary->Phi();
1483 etaprim = primary->Eta();
1484 pxprim = primary->Px();
1485 pyprim = primary->Py();
1486 pzprim = primary->Pz();
1487
1488 if(primary){
e09cf5ef 1489
045396c8 1490 for (iParticle = 0 ; iParticle < nTracks ; iParticle++) {
1491 TParticle * particle = stack->Particle(iParticle);
1492 TLorentzVector momentum;
1493 //keep only final state particles
1494 if(particle->GetStatusCode()!=1) continue ;
1495 Int_t pdg = particle->GetPdgCode();
1496 charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
1497 particle->Momentum(momentum);
1498
1499 //---------- Charged particles ----------------------
1500 if(charge != 0){
1501 //Particles in CTS acceptance
1502 Bool_t inCTS = GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
1503 if(TMath::Abs(pdg) == 11 && stack->Particle(particle->GetFirstMother())->GetPdgCode()==22) continue ;
1504 if(inCTS&&momentum.Pt() >GetMinPt())
1505 {
1506 mcTrackPt = particle->Pt();
1507 mcTrackPhi = particle->Phi();
1508 mcTrackEta = particle->Eta();
1509 mcTrackPx = particle->Px();
1510 mcTrackPy = particle->Py();
1511 mcTrackPz = particle->Pz();
1512 if(mcTrackPhi < 0) mcTrackPhi+=TMath::TwoPi();
1513 //Select only hadrons in pt range
1514 if(mcTrackPt < GetMinPt() || mcTrackPt > GetMaxPt()) continue ;
1515 //remove trigger itself for correlation when use charged triggers
1516 if(label==iParticle && mcTrackPt==ptprim && mcTrackPhi==phiprim && mcTrackEta==etaprim)
1517 continue ;
1518 //jumped out this event if near side associated partile pt larger than trigger
1519 if( mcTrackPt> ptprim && TMath::Abs(mcTrackPhi-phiprim)<TMath::PiOver2())
1520 return ;
1521
1522 mcrat = mcTrackPt/ptprim ;
1523 mcxE = -(mcTrackPx*pxprim+mcTrackPy*pyprim)/(ptprim*ptprim);
1524 if(mcxE <0.) mcxE =-mcxE;
1525 mccosi = TMath::Log(1/mcxE);
1526 // printf("rat = %f, xE = %f, cosi =%f \n", rat, xE, cosi);
1527 // printf("phi = %f \n", phi);
1528
1529 //Selection within angular range
1530 Double_t mcdeltaPhi = phiprim-mcTrackPhi;
1531 if( mcdeltaPhi< -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
1532 if( mcdeltaPhi>3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
1533 Double_t mcpout = mcTrackPt*TMath::Sin(mcdeltaPhi) ;
1534 if(GetDebug()>0 )
1535 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Charged hadron: track Pt %f, track Phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
1536 mcTrackPt,mcTrackPhi, phiprim,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut, GetMinPt());
1537 // Fill Histograms
1538 fhMCEtaCharged->Fill(mcTrackPt,mcTrackEta);
1539 fhMCPhiCharged->Fill(mcTrackPt,mcTrackPhi);
1540 fhMCDeltaEtaCharged->Fill(ptprim,etaprim-mcTrackEta);
1541 fhMCDeltaPhiCharged->Fill(ptprim,mcdeltaPhi);
1542 fhMCPtAssocDeltaPhi->Fill(mcTrackPt, mcdeltaPhi);
1543 // fhDeltaPhiCharged->Fill(ptTrig, deltaPhi);
1544 fhMCDeltaPhiDeltaEtaCharged->Fill(mcdeltaPhi,etaprim-mcTrackEta);
1545
1546 //delta phi cut for correlation
1547 if( (mcdeltaPhi > fDeltaPhiMinCut) && ( mcdeltaPhi < fDeltaPhiMaxCut) ) {
1548 fhMCDeltaPhiChargedPt->Fill(mcTrackPt,mcdeltaPhi);
1549 fhMCPtImbalanceCharged->Fill(ptprim,mcxE);
1550 fhMCPtHbpCharged->Fill(ptprim,mccosi);
1551 fhMCPtTrigPout->Fill(ptprim, mcpout) ;
1552 }//delta phi cut for correlation
1553 } //tracks after cuts
1554 }//Charged
1555 } //track loop
1556 } //when the leading particles could trace back to MC
1557 } //ESD MC
1558 else if(GetReader()->ReadAODMCParticles()){
1559 //Get the list of MC particles
1560 mcparticles0 = GetReader()->GetAODMCParticles(0);
1561 if(!mcparticles0) return;
1562 if(label >=mcparticles0->GetEntriesFast()){
1563 if(GetDebug() > 2)
1564 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** large label ***: label %d, n tracks %d \n", label,mcparticles0->GetEntriesFast());
1565 return;
1566 }
1567 //Get the particle
1568 aodprimary = (AliAODMCParticle*) mcparticles0->At(label);
1569 if(!aodprimary) {
1570 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** no AOD primary ***: label %d \n", label);
1571 return;
1572 }
1573
1574 ptprim = aodprimary->Pt();
1575 phiprim = aodprimary->Phi();
1576 etaprim =aodprimary->Eta();
1577 pxprim =aodprimary->Px();
1578 pyprim =aodprimary->Py();
1579 pzprim =aodprimary->Pz();
1580 if(aodprimary){
1581 mcparticles= GetReader()->GetAODMCParticles();
1582 for (Int_t i=0; i<nTracks;i++) {
1583 AliAODMCParticle *part = (AliAODMCParticle*) mcparticles->At(i);
1584 if (!part->IsPhysicalPrimary()) continue;
1585 Int_t pdg = part->GetPdgCode();
1586 charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
1587 TLorentzVector momentum(part->Px(),part->Py(),part->Pz(),part->E());
1588 if(charge != 0){
1589 if(part->Pt()> GetReader()->GetCTSPtMin()){
1590 //Particles in CTS acceptance
1591 Bool_t inCTS = GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
1592 Int_t indexmother=part->GetMother();
1593 if(indexmother>-1)
1594 {
1595 Int_t mPdg = ((AliAODMCParticle*) mcparticles->At(indexmother)) ->GetPdgCode();
1596 if(TMath::Abs(pdg) == 11 && mPdg == 22) continue;
1597 }
1598 if(inCTS&&momentum.Pt() >GetMinPt())
1599 {
1600 mcTrackPt=part->Pt();
1601 mcTrackPhi=part->Phi();
1602 mcTrackEta=part->Eta();
1603 mcTrackPx=part->Px();
1604 mcTrackPy=part->Py();
1605 mcTrackPz=part->Pz();
1606 if(mcTrackPhi < 0) mcTrackPhi+=TMath::TwoPi();
1607 //Select only hadrons in pt range
1608 if(mcTrackPt < GetMinPt() || mcTrackPt > GetMaxPt()) continue ;
1609 //remove trigger itself for correlation when use charged triggers
1610 if(label==i && mcTrackPt==ptprim && mcTrackPhi==phiprim && mcTrackEta==etaprim)
1611 continue ;
1612 //jumped out this event if near side associated partile pt larger than trigger
1613 if( mcTrackPt> ptprim && TMath::Abs(mcTrackPhi-phiprim)<TMath::PiOver2())
1614 return ;
1615
1616 mcrat = mcTrackPt/ptprim ;
1617 mcxE = -(mcTrackPx*pxprim+mcTrackPy*pyprim)/(ptprim*ptprim);
1618 if(mcxE <0.)mcxE =-mcxE;
1619 mccosi = TMath::Log(1/mcxE);
1620
1621 //Selection within angular range
1622 Double_t mcdeltaPhi = phiprim-mcTrackPhi;
1623 if( mcdeltaPhi< -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
1624 if( mcdeltaPhi>3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
1625 Double_t mcpout = mcTrackPt*TMath::Sin(mcdeltaPhi) ;
1626 if(GetDebug()>0)
1627 printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Charged hadron: track Pt %f, track Phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f, pT min %2.2f \n",
1628 mcTrackPt,mcTrackPhi, phiprim,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut, GetMinPt());
1629 // Fill Histograms
1630 fhMCEtaCharged->Fill(mcTrackPt,mcTrackEta);
1631 fhMCPhiCharged->Fill(mcTrackPt,mcTrackPhi);
1632 fhMCDeltaEtaCharged->Fill(ptprim,etaprim-mcTrackEta);
1633 fhMCDeltaPhiCharged->Fill(ptprim,mcdeltaPhi);
1634 fhMCPtAssocDeltaPhi->Fill(mcTrackPt, mcdeltaPhi);
1635 // fhDeltaPhiCharged->Fill(ptTrig, deltaPhi);
1636 fhMCDeltaPhiDeltaEtaCharged->Fill(mcdeltaPhi,etaprim-mcTrackEta);
1637
1638 //delta phi cut for correlation
1639 if( (mcdeltaPhi > fDeltaPhiMinCut) && ( mcdeltaPhi < fDeltaPhiMaxCut) ) {
1640 fhMCDeltaPhiChargedPt->Fill(mcTrackPt,mcdeltaPhi);
1641 fhMCPtImbalanceCharged->Fill(ptprim,mcxE);
1642 fhMCPtHbpCharged->Fill(ptprim,mccosi);
1643 fhMCPtTrigPout->Fill(ptprim, mcpout) ;
1644 }//delta phi cut for correlation
1645
1646 } //tracks after cuts
1647
1648 } //with minimum pt cut
1649 } //only charged particles
1650 } //MC particle loop
1651 } //when the leading particles could trace back to MC
1652 }// AOD MC
e09cf5ef 1653}
045396c8 1654
1655//____________________________________________________________________________
1656void AliAnaParticleHadronCorrelation::SetNTriggerPtBins(Int_t n)
1657{
1658 // Set number of bins
e09cf5ef 1659
045396c8 1660 if(n < 10 && n > 0)
1661 {
1662 fNTriggerPtBins = n ;
1663 }
1664 else
1665 {
1666 printf("n = larger than 9 or too small, set to 9 \n");
1667 fNTriggerPtBins = 9;
1668 }
1669}
e09cf5ef 1670
1671//____________________________________________________________________________
045396c8 1672void AliAnaParticleHadronCorrelation::SetTriggerPtBinLimit(Int_t ibin, Float_t pt)
1673{
1674 // Set the list of limits for the trigger pt bins
1675
1676 if(ibin <= fNTriggerPtBins || ibin >= 0)
1677 {
1678 fTriggerPtBinLimit[ibin] = pt ;
1679 }
1680 else {
1681 printf("AliAnaParticleHadronCorrelation::SetTriggerPtBinLimit() - bin number too large %d > %d or small, nothing done\n", ibin, fNTriggerPtBins) ;
1682
1683 }
1684}
1685