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