AliPhysics  dab84fb (dab84fb)
AliAnaParticleHadronCorrelation.cxx
Go to the documentation of this file.
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 
16 // --- ROOT system ---
17 #include <TClass.h>
18 #include <TMath.h>
19 #include <TH2F.h>
20 #include <TDatabasePDG.h>
21 #include <TClonesArray.h>
22 #include <TList.h>
23 #include <TObjString.h>
24 
25 //---- ANALYSIS system ----
26 #include "AliNeutralMesonSelection.h"
27 #include "AliAnaParticleHadronCorrelation.h"
28 #include "AliCaloTrackReader.h"
29 #include "AliCaloTrackParticleCorrelation.h"
30 #include "AliFiducialCut.h"
31 #include "AliVTrack.h"
32 #include "AliVCluster.h"
33 #include "AliMCAnalysisUtils.h"
34 #include "AliMCEvent.h"
35 #include "AliVParticle.h"
36 #include "AliMixedEvent.h"
37 #include "AliAnalysisManager.h"
38 #include "AliInputEventHandler.h"
39 #include "AliEventplane.h"
40 #include "AliIsolationCut.h"
41 #include "AliESDEvent.h"
42 #include "AliAODEvent.h"
43 
44 // --- Detectors ---
45 #include "AliEMCALGeometry.h"
46 #include "AliPHOSGeoUtils.h"
47 
49 ClassImp(AliAnaParticleHadronCorrelation) ;
51 
52 
53 //___________________________________________________________________
55 //___________________________________________________________________
56 AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation():
57 AliAnaCaloTrackCorrBaseClass(),
58 fFillAODWithReferences(0), fCheckLeadingWithNeutralClusters(0),
59 fMaxAssocPt(1000.), fMinAssocPt(0.),
60 fDeltaPhiMaxCut(0.), fDeltaPhiMinCut(0.),
61 fSelectIsolated(0), fMakeSeveralUE(0),
62 fUeDeltaPhiMaxCut(0.), fUeDeltaPhiMinCut(0.),
63 fPi0AODBranchName(""), fAODNamepTInConeHisto(""), fNeutralCorr(0),
64 fPi0Trigger(0), fDecayTrigger(0),
65 fNDecayBits(0), fDecayBits(),
66 fNBkgBin(0),
67 fMakeAbsoluteLeading(0), fMakeNearSideLeading(0),
68 fLeadingTriggerIndex(-1), fHMPIDCorrelation(0), fFillBradHisto(0),
69 fNAssocPtBins(0), fAssocPtBinLimit(),
70 fCorrelVzBin(0),
71 fListMixTrackEvents(), fListMixCaloEvents(),
72 fUseMixStoredInReader(0), fFillNeutralEventMixPool(0),
73 fM02MaxCut(0), fM02MinCut(0),
74 fSelectLeadingHadronAngle(0), fFillLeadHadOppositeHisto(0),
75 fMinLeadHadPhi(0), fMaxLeadHadPhi(0),
76 fMinLeadHadPt(0), fMaxLeadHadPt(0),
77 fFillEtaGapsHisto(1), fFillMomImbalancePtAssocBinsHisto(0),
78 fFillInvMassHisto(0), fFillBkgBinsHisto(0),
79 fFillPerSMHistograms(0),
80 fFillTaggedDecayHistograms(0), fDecayTagsM02Cut(0),
81 fMCGenTypeMin(0), fMCGenTypeMax(0),
82 fTrackVector(), fMomentum(), fMomentumIM(),
83 fDecayMom1(), fDecayMom2(),
84 //Histograms
85 fhPtTriggerInput(0), fhPtTriggerSSCut(0),
86 fhPtTriggerIsoCut(0), fhPtTriggerFidCut(0),
87 fhPtTrigger(0), fhPtTriggerVtxBC0(0),
88 fhPtTriggerVzBin(0), fhPtTriggerBin(0),
89 fhPhiTrigger(0), fhEtaTrigger(0),
90 fhPtTriggerMC(),
91 fhPtDecayTrigger(), fhPtDecayTriggerMC(),
92 fhPtTriggerCentrality(0), fhPtTriggerEventPlane(0),
93 fhTriggerEventPlaneCentrality(0),
94 fhPtTriggerMixed(0), fhPtTriggerMixedVzBin(0), fhPtTriggerMixedBin(0),
95 fhPhiTriggerMixed(0), fhEtaTriggerMixed(0),
96 fhPtLeadingOppositeHadron(0), fhPtDiffPhiLeadingOppositeHadron(0), fhPtDiffEtaLeadingOppositeHadron(0),
97 fhPtNoLeadingOppositeHadron(0), fhEtaPhiNoLeadingOppositeHadron(0),
98 fhDeltaPhiDeltaEtaCharged(0),
99 fhPhiCharged(0), fhEtaCharged(0),
100 fhDeltaPhiCharged(0), fhDeltaEtaCharged(0),
101 fhDeltaPhiChargedPt(0), fhDeltaPhiUeChargedPt(0),
102 fhUePart(0),
103 fhXECharged(0), fhXECharged_Cone2(0), fhXEUeCharged(0),
104 fhXEUeChargedSmallCone(0), fhXEUeChargedMediumCone(0), fhXEUeChargedLargeCone(0),
105 fhXEPosCharged(0), fhXENegCharged(0),
106 fhPtHbpXECharged(0), fhPtHbpXECharged_Cone2(0), fhPtHbpXEUeCharged(0),
107 fhZTCharged(0), fhZTUeCharged(0),
108 fhZTPosCharged(0), fhZTNegCharged(0),
109 fhPtHbpZTCharged(0), fhPtHbpZTUeCharged(0),
110 fhXEChargedMC(), fhDeltaPhiChargedMC(),
111 fhXEUeChargedRightMC(), fhXEUeChargedLeftMC(),
112 fhDeltaPhiDeltaEtaChargedPtA3GeV(0),
113 fhDeltaPhiChargedPtA3GeV(0), fhDeltaEtaChargedPtA3GeV(0),
114 // Pile-Up
115 fhDeltaPhiChargedPileUp(), fhDeltaEtaChargedPileUp(),
116 fhDeltaPhiChargedPtA3GeVPileUp(), fhDeltaEtaChargedPtA3GeVPileUp(),
117 fhXEChargedPileUp(), fhXEUeChargedPileUp(),
118 fhZTChargedPileUp(), fhZTUeChargedPileUp(),
119 fhPtTrigChargedPileUp(),
120 fhDeltaPhiChargedOtherBC(), fhDeltaPhiChargedPtA3GeVOtherBC(),
121 fhXEChargedOtherBC(), fhXEUeChargedOtherBC(),
122 fhZTChargedOtherBC(), fhZTUeChargedOtherBC(),
123 fhPtTrigChargedOtherBC(),
124 fhDeltaPhiChargedBC0(), fhDeltaPhiChargedPtA3GeVBC0(),
125 fhXEChargedBC0(), fhXEUeChargedBC0(),
126 fhZTChargedBC0(), fhZTUeChargedBC0(),
127 fhPtTrigChargedBC0(),
128 fhDeltaPhiChargedVtxBC0(), fhDeltaPhiChargedPtA3GeVVtxBC0(),
129 fhXEChargedVtxBC0(), fhXEUeChargedVtxBC0(),
130 fhZTChargedVtxBC0(), fhZTUeChargedVtxBC0(),
131 fhPtTrigChargedVtxBC0(),
132 fhDeltaPhiUeLeftCharged(0),
133 fhDeltaPhiUeLeftUpCharged(0), fhDeltaPhiUeRightUpCharged(0),
134 fhDeltaPhiUeLeftDownCharged(0), fhDeltaPhiUeRightDownCharged(0),
135 fhXEUeLeftCharged(0),
136 fhXEUeLeftUpCharged(0), fhXEUeRightUpCharged(0),
137 fhXEUeLeftDownCharged(0), fhXEUeRightDownCharged(0),
138 fhPtHbpXEUeLeftCharged(0), fhZTUeLeftCharged(0),
139 fhPtHbpZTUeLeftCharged(0),
140 fhPtTrigPout(0), fhPtTrigCharged(0),
141 fhDeltaPhiChargedMult(0x0), fhDeltaEtaChargedMult(0x0),
142 fhXEMult(0x0), fhXEUeMult(0x0),
143 fhZTMult(0x0), fhZTUeMult(0x0),
144 fhAssocPtBkg(0), fhDeltaPhiDeltaEtaAssocPtBin(0),
145 fhDeltaPhiAssocPtBin(0),
146 fhDeltaPhiAssocPtBinDEta08(0), fhDeltaPhiAssocPtBinDEta0(0),
147 fhDeltaPhiAssocPtBinHMPID(0), fhDeltaPhiAssocPtBinHMPIDAcc(0),
148 fhDeltaPhiBradAssocPtBin(0), fhDeltaPhiBrad(0),
149 fhXEAssocPtBin(0), fhZTAssocPtBin(0),
150 fhXEVZ(0), fhZTVZ(0),
151 fhDeltaPhiDeltaEtaNeutral(0),
152 fhPhiNeutral(0), fhEtaNeutral(0),
153 fhDeltaPhiNeutral(0), fhDeltaEtaNeutral(0),
154 fhDeltaPhiNeutralPt(0), fhDeltaPhiUeNeutralPt(0),
155 fhXENeutral(0), fhXEUeNeutral(0),
156 fhPtHbpXENeutral(0), fhPtHbpXEUeNeutral(0),
157 fhZTNeutral(0), fhZTUeNeutral(0),
158 fhPtHbpZTNeutral(0), fhPtHbpZTUeNeutral(0),
159 fhDeltaPhiUeLeftNeutral(0), fhXEUeLeftNeutral(0),
160 fhPtHbpXEUeLeftNeutral(0), fhZTUeLeftNeutral(0),
161 fhPtHbpZTUeLeftNeutral(0), fhPtPi0DecayRatio(0),
162 fhDeltaPhiPi0DecayCharged(0), fhXEPi0DecayCharged(0), fhZTPi0DecayCharged(0),
163 fhDeltaPhiPi0DecayNeutral(0), fhXEPi0DecayNeutral(0), fhZTPi0DecayNeutral(0),
164 fhDeltaPhiDecayCharged(), fhXEDecayCharged(), fhZTDecayCharged(),
165 fhDeltaPhiDecayChargedAssocPtBin(),
166 fhMCPtTrigger(), fhMCPhiTrigger(), fhMCEtaTrigger(),
167 fhMCPtTriggerNotLeading(), fhMCPhiTriggerNotLeading(), fhMCEtaTriggerNotLeading(),
168 fhMCEtaCharged(), fhMCPhiCharged(),
169 fhMCDeltaEtaCharged(), fhMCDeltaPhiCharged(),
170 fhMCDeltaPhiDeltaEtaCharged(), fhMCDeltaPhiChargedPt(),
171 fhMCPtXECharged(), fhMCPtXEUeCharged(),
172 fhMCPtXEUeLeftCharged(),
173 fhMCPtHbpXECharged(), fhMCPtHbpXEUeCharged(),
174 fhMCPtHbpXEUeLeftCharged(),
175 fhMCUePart(),
176 fhMCPtZTCharged(), fhMCPtZTUeCharged(),
177 fhMCPtZTUeLeftCharged(),
178 fhMCPtHbpZTCharged(), fhMCPtHbpZTUeCharged(),
179 fhMCPtHbpZTUeLeftCharged(),
180 fhMCPtTrigPout(), fhMCPtAssocDeltaPhi(),
181 // Mixing
182 fhNEventsTrigger(0), fhNtracksMB(0), fhNclustersMB(0),
183 fhMixDeltaPhiCharged(0), fhMixDeltaPhiDeltaEtaCharged(0),
184 fhMixXECharged(0), fhMixXEUeCharged(0), fhMixHbpXECharged(0),
185 fhMixDeltaPhiChargedAssocPtBin(),
186 fhMixDeltaPhiChargedAssocPtBinDEta08(),
187 fhMixDeltaPhiChargedAssocPtBinDEta0(),
188 fhMixDeltaPhiDeltaEtaChargedAssocPtBin(),
189 fhEventBin(0), fhEventMixBin(0), fhEventMBBin(0),
190 fhMassPtTrigger(0), fhMCMassPtTrigger(),
191 fhPtLeadInConeBin(), fhPtSumInConeBin(),
192 fhPtLeadConeBinDecay(), fhSumPtConeBinDecay(),
193 fhPtLeadConeBinMC(), fhSumPtConeBinMC(),
194 fhTrackResolution(0), fhTrackResolutionUE(0),
195 fhPtTriggerPerSM(0)
196 {
197  InitParameters();
198 
199  for(Int_t i = 0; i < fgkNmcTypes; i++)
200  {
201  fhPtTriggerMC[i] = 0;
202  fhXEChargedMC[i] = 0;
203  fhDeltaPhiChargedMC[i] = 0;
204  fhMCMassPtTrigger [i] = 0 ;
205  fhXEUeChargedRightMC[i] = 0;
206  fhXEUeChargedLeftMC[i] = 0;
207 
208  for(Int_t ib = 0; ib < AliNeutralMesonSelection::fgkMaxNDecayBits; ib++) fhPtDecayTriggerMC[ib][i] = 0;
209  }
210 
211  for(Int_t ib = 0; ib < AliNeutralMesonSelection::fgkMaxNDecayBits; ib++) fhPtDecayTrigger[ib] = 0;
212 
213  for(Int_t i = 0; i < 7; i++)
214  {
215  fhPtTriggerPileUp [i] = 0 ;
216  fhDeltaPhiChargedPileUp [i] = 0 ; fhDeltaEtaChargedPileUp [i] = 0 ;
217  fhXEChargedPileUp [i] = 0 ; fhXEUeChargedPileUp [i] = 0 ;
218  fhZTChargedPileUp [i] = 0 ; fhZTUeChargedPileUp [i] = 0 ;
219  fhPtTrigChargedPileUp [i] = 0 ;
220  fhDeltaPhiChargedPtA3GeVPileUp[i] = 0 ; fhDeltaEtaChargedPtA3GeVPileUp[i] = 0 ;
221  }
222 
223  for(Int_t ism = 0; ism < 20; ism++)
224  {
225  fhXEChargedPerSM [ism] = 0 ;
226  fhXEUeChargedPerSM [ism] = 0 ;
227  fhDeltaPhiChargedPerSM [ism] = 0 ;
228  fhDeltaPhiChargedPtA3GeVPerSM[ism] = 0 ;
229  }
230 }
231 
232 //_________________________________________________________________
234 //_________________________________________________________________
235 AliAnaParticleHadronCorrelation::~AliAnaParticleHadronCorrelation()
236 {
237  if(DoOwnMix())
238  {
239  if(fListMixTrackEvents)
240  {
241  for(Int_t iz=0; iz < GetNZvertBin(); iz++)
242  {
243  for(Int_t ic=0; ic < GetNCentrBin(); ic++)
244  {
245  for(Int_t irp=0; irp<GetNRPBin(); irp++)
246  {
247  Int_t bin = GetEventMixBin(ic, iz, irp);
248  fListMixTrackEvents[bin]->Delete() ;
249  delete fListMixTrackEvents[bin] ;
250  }
251  }
252  }
253  }
254 
255  delete[] fListMixTrackEvents;
256 
257  if(fListMixCaloEvents)
258  {
259  for(Int_t iz=0; iz < GetNZvertBin(); iz++)
260  {
261  for(Int_t ic=0; ic < GetNCentrBin(); ic++)
262  {
263  for(Int_t irp=0; irp<GetNRPBin(); irp++)
264  {
265  Int_t bin = GetEventMixBin(ic, iz, irp);
266  fListMixCaloEvents[bin]->Delete() ;
267  delete fListMixCaloEvents[bin] ;
268  }
269  }
270  }
271  }
272 
273  delete[] fListMixCaloEvents;
274  }
275 }
276 
277 //____________________________________________________________________________________________________________________________________
279 //____________________________________________________________________________________________________________________________________
280 void AliAnaParticleHadronCorrelation::FillChargedAngularCorrelationHistograms(Float_t ptAssoc, Float_t ptTrig, Int_t bin,
281  Float_t phiAssoc, Float_t phiTrig, Float_t deltaPhi,
282  Float_t etaAssoc, Float_t etaTrig, Int_t sm,
283  Int_t decayTag, Float_t hmpidSignal, Int_t outTOF,
284  Int_t cen, Int_t mcTag)
285 {
286  Float_t deltaEta = etaTrig-etaAssoc;
287  Float_t deltaPhiOrg = phiTrig-phiAssoc;
288 
289  fhEtaCharged ->Fill(ptAssoc, etaAssoc, GetEventWeight());
290  fhPhiCharged ->Fill(ptAssoc, phiAssoc, GetEventWeight());
291  fhDeltaEtaCharged ->Fill(ptTrig , deltaEta, GetEventWeight());
292  fhDeltaPhiCharged ->Fill(ptTrig , deltaPhi, GetEventWeight());
293  fhDeltaPhiChargedPt->Fill(ptAssoc, deltaPhi, GetEventWeight());
294  fhDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta, GetEventWeight());
295 
296  if(ptAssoc > 3 )
297  {
298  fhDeltaEtaChargedPtA3GeV ->Fill(ptTrig , deltaEta, GetEventWeight());
299  fhDeltaPhiChargedPtA3GeV ->Fill(ptTrig , deltaPhi, GetEventWeight());
300  fhDeltaPhiDeltaEtaChargedPtA3GeV->Fill(deltaPhi, deltaEta, GetEventWeight());
301  }
302 
303  if( fFillPerSMHistograms )
304  {
305  fhDeltaPhiChargedPerSM[sm]->Fill(ptTrig, deltaPhi, GetEventWeight());
306  if(ptAssoc > 3 )
307  fhDeltaPhiChargedPtA3GeVPerSM[sm]->Fill(ptTrig, deltaPhi, GetEventWeight());
308  }
309 
310  // Pile up studies
311 
312  if(IsPileUpAnalysisOn())
313  {
314  if (outTOF==1)
315  {
316  fhDeltaPhiChargedOtherBC->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
317  if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVOtherBC->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
318  }
319  else if(outTOF==0)
320  {
321  fhDeltaPhiChargedBC0->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
322  if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVBC0->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
323  }
324 
325  Int_t vtxBC = GetReader()->GetVertexBC();
326  if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
327  {
328  fhDeltaPhiChargedVtxBC0->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
329  if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVVtxBC0->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
330  }
331 
332  if(GetReader()->IsPileUpFromSPD())
333  {
334  fhDeltaEtaChargedPileUp[0]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
335  fhDeltaPhiChargedPileUp[0]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
336  }
337  if(GetReader()->IsPileUpFromEMCal())
338  {
339  fhDeltaEtaChargedPileUp[1]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
340  fhDeltaPhiChargedPileUp[1]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
341  }
342  if(GetReader()->IsPileUpFromSPDOrEMCal())
343  {
344  fhDeltaEtaChargedPileUp[2]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
345  fhDeltaPhiChargedPileUp[2]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
346  }
347  if(GetReader()->IsPileUpFromSPDAndEMCal())
348  {
349  fhDeltaEtaChargedPileUp[3]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
350  fhDeltaPhiChargedPileUp[3]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
351  }
352  if(GetReader()->IsPileUpFromSPDAndNotEMCal())
353  {
354  fhDeltaEtaChargedPileUp[4]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
355  fhDeltaPhiChargedPileUp[4]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
356  }
357  if(GetReader()->IsPileUpFromEMCalAndNotSPD())
358  {
359  fhDeltaEtaChargedPileUp[5]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
360  fhDeltaPhiChargedPileUp[5]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
361  }
362  if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
363  {
364  fhDeltaEtaChargedPileUp[6]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
365  fhDeltaPhiChargedPileUp[6]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
366  }
367 
368  if(ptAssoc > 3 )
369  {
370  if(GetReader()->IsPileUpFromSPD())
371  {
372  fhDeltaEtaChargedPtA3GeVPileUp[0]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
373  fhDeltaPhiChargedPtA3GeVPileUp[0]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
374  }
375  if(GetReader()->IsPileUpFromEMCal())
376  {
377  fhDeltaEtaChargedPtA3GeVPileUp[1]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
378  fhDeltaPhiChargedPtA3GeVPileUp[1]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
379  }
380  if(GetReader()->IsPileUpFromSPDOrEMCal())
381  {
382  fhDeltaEtaChargedPtA3GeVPileUp[2]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
383  fhDeltaPhiChargedPtA3GeVPileUp[2]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
384  }
385  if(GetReader()->IsPileUpFromSPDAndEMCal())
386  {
387  fhDeltaEtaChargedPtA3GeVPileUp[3]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
388  fhDeltaPhiChargedPtA3GeVPileUp[3]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
389  }
390  if(GetReader()->IsPileUpFromSPDAndNotEMCal())
391  {
392  fhDeltaEtaChargedPtA3GeVPileUp[4]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
393  fhDeltaPhiChargedPtA3GeVPileUp[4]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
394  }
395  if(GetReader()->IsPileUpFromEMCalAndNotSPD())
396  {
397  fhDeltaEtaChargedPtA3GeVPileUp[5]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
398  fhDeltaPhiChargedPtA3GeVPileUp[5]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
399  }
400  if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
401  {
402  fhDeltaEtaChargedPtA3GeVPileUp[6]->Fill(ptTrig, deltaEta, GetEventWeight()) ;
403  fhDeltaPhiChargedPtA3GeVPileUp[6]->Fill(ptTrig, deltaPhi, GetEventWeight()) ;
404  }
405  }
406  }
407 
408  if(IsDataMC())
409  {
410  Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
411  fhDeltaPhiChargedMC[mcIndex]->Fill(ptTrig , deltaPhi, GetEventWeight());
412 
413  if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) )
414  {
415  // check index in GetMCTagIndexHistogram
416  if ( mcIndex == 2 ) fhDeltaPhiChargedMC[8]->Fill(ptTrig, deltaPhi, GetEventWeight()); // pi0 decay
417  else if( mcIndex == 4 ) fhDeltaPhiChargedMC[9]->Fill(ptTrig, deltaPhi, GetEventWeight()); // eta decay
418  }
419  }
420 
421  if(fDecayTrigger && decayTag > 0)
422  {
423  for(Int_t ibit = 0; ibit<fNDecayBits; ibit++)
424  {
425  if(GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[ibit]))
426  fhDeltaPhiDecayCharged[ibit]->Fill(ptTrig, deltaPhi, GetEventWeight());
427  }
428  }
429 
430  Double_t dphiBrad = -100;
431  if(fFillBradHisto)
432  {
433  dphiBrad = atan2(sin(deltaPhiOrg), cos(deltaPhiOrg))/TMath::Pi();//-1 to 1
434  if( TMath::Abs(dphiBrad) > 0.325 && TMath::Abs(dphiBrad) < 0.475 ) // Hardcoded values, BAD, FIXME
435  {
436  fhAssocPtBkg->Fill(ptTrig, ptAssoc, GetEventWeight());
437  }
438 
439  if( dphiBrad < -1./3 ) dphiBrad += 2;
440  fhDeltaPhiBrad->Fill(ptTrig, dphiBrad, GetEventWeight());
441  }
442 
443  // Fill histograms in bins of associated particle pT
444  if ( bin >= 0 )
445  {
446  fhDeltaPhiDeltaEtaAssocPtBin[bin]->Fill(deltaPhi,deltaEta, GetEventWeight());
447 
448  fhDeltaPhiAssocPtBin [bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
449 
450  if(fFillEtaGapsHisto)
451  {
452  if(TMath::Abs(deltaEta)> 0.8)
453  fhDeltaPhiAssocPtBinDEta08[bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
454 
455  if(TMath::Abs(deltaEta)< 0.01)
456  fhDeltaPhiAssocPtBinDEta0 [bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
457  }
458 
459  if (fFillBradHisto)
460  fhDeltaPhiBradAssocPtBin [bin]->Fill(ptTrig, dphiBrad, GetEventWeight());
461 
462  if(fDecayTrigger && decayTag > 0 && fNDecayBits > 0 &&
463  GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[0]))
464  {
465  fhDeltaPhiDecayChargedAssocPtBin[bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
466  }
467 
468  if(fHMPIDCorrelation)
469  {
470  if( hmpidSignal > 0 )
471  {
472  //printf("Track pt %f with HMPID signal %f \n",pt,hmpidSignal);
473  fhDeltaPhiAssocPtBinHMPID[bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
474  }
475 
476  if(phiAssoc > 5*TMath::DegToRad() && phiAssoc < 20*TMath::DegToRad())
477  {
478  //printf("Track pt %f in HMPID acceptance phi %f \n ",pt,phi*TMath::RadToDeg() );
479  fhDeltaPhiAssocPtBinHMPIDAcc[bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
480  }
481  }
482  }
483 
484  // Fill different multiplicity/centrality histogram
485  if(IsHighMultiplicityAnalysisOn() && cen >= 0 && cen < GetNCentrBin())
486  {
487  fhDeltaPhiChargedMult[cen]->Fill(ptTrig, deltaPhi, GetEventWeight());
488  fhDeltaEtaChargedMult[cen]->Fill(ptTrig, deltaEta, GetEventWeight());
489  }
490 }
491 
492 //___________________________________________________________________________________________________________________________________
494 //___________________________________________________________________________________________________________________________________
495 Bool_t AliAnaParticleHadronCorrelation::FillChargedMCCorrelationHistograms(Float_t mcAssocPt, Float_t mcAssocPhi, Float_t mcAssocEta,
496  Float_t mcTrigPt, Float_t mcTrigPhi, Float_t mcTrigEta,
497  Int_t histoIndex, Bool_t lostDecayPair)
498 {
499  Bool_t lead = kTRUE;
500 
501  // In case we requested the trigger to be a leading particle,
502  // check if this is true at the MC level.
503  // Not sure if it is correct to skip or not skip this.
504  // Absolute leading?
505  if( fMakeAbsoluteLeading && mcAssocPt > mcTrigPt ) lead = kFALSE; // skip event
506 
507  // Skip this event if near side associated particle pt larger than trigger
508  if( mcAssocPhi < 0 ) mcAssocPhi+=TMath::TwoPi();
509 
510  Float_t mcdeltaPhi= mcTrigPhi-mcAssocPhi;
511  if(mcdeltaPhi <= -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
512  if(mcdeltaPhi > 3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
513 
514  if( fMakeNearSideLeading)
515  {
516  if( mcAssocPt > mcTrigPt && mcdeltaPhi < TMath::PiOver2() ) lead = kFALSE; // skip event
517  }
518 
519  //
520  // Select only hadrons in pt range
521  if ( mcAssocPt < fMinAssocPt || mcAssocPt > fMaxAssocPt ) return lead ; // exclude but continue
522  if ( mcAssocPt < GetReader()->GetCTSPtMin()) return lead ;
523 
524 
525  //
526  // Remove trigger itself for correlation when use charged triggers
527  if(TMath::Abs(mcAssocPt -mcTrigPt ) < 1e-6 &&
528  mcdeltaPhi < 1e-6 &&
529  TMath::Abs(mcAssocEta-mcTrigEta) < 1e-6) return lead ; // exclude but continue
530 
531  Float_t mcxE =-mcAssocPt/mcTrigPt*TMath::Cos(mcdeltaPhi);// -(mcAssocPx*pxprim+mcAssocPy*pyprim)/(mcTrigPt*mcTrigPt);
532  Float_t mchbpXE =-100 ;
533  if(mcxE > 0 ) mchbpXE = TMath::Log(1./mcxE);
534 
535  Float_t mczT = mcAssocPt/mcTrigPt ;
536  Float_t mchbpZT =-100 ;
537  if(mczT > 0 ) mchbpZT = TMath::Log(1./mczT);
538 
539  Double_t mcpout = mcAssocPt*TMath::Sin(mcdeltaPhi) ;
540 
541  AliDebug(1,Form("Charged hadron: track Pt %f, track Phi %f, phi trigger %f. Cuts: delta phi %2.2f < %2.2f < %2.2f",
542  mcAssocPt,mcAssocPhi, mcTrigPhi,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut));
543 
544  // Fill Histograms
545  fhMCEtaCharged [histoIndex]->Fill(mcAssocPt, mcAssocEta, GetEventWeight());
546  fhMCPhiCharged [histoIndex]->Fill(mcAssocPt, mcAssocPhi, GetEventWeight());
547  fhMCDeltaPhiCharged[histoIndex]->Fill(mcTrigPt , mcdeltaPhi, GetEventWeight());
548  fhMCPtAssocDeltaPhi[histoIndex]->Fill(mcAssocPt, mcdeltaPhi, GetEventWeight());
549 
550  fhMCDeltaEtaCharged [histoIndex]->Fill(mcTrigPt , mcTrigEta-mcAssocEta, GetEventWeight());
551  fhMCDeltaPhiDeltaEtaCharged[histoIndex]->Fill(mcdeltaPhi, mcTrigEta-mcAssocEta, GetEventWeight());
552 
553  // Delta phi cut for correlation
554  if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) )
555  {
556  fhMCDeltaPhiChargedPt[histoIndex]->Fill(mcAssocPt, mcdeltaPhi, GetEventWeight());
557  fhMCPtXECharged [histoIndex]->Fill(mcTrigPt , mcxE , GetEventWeight());
558  fhMCPtHbpXECharged [histoIndex]->Fill(mcTrigPt , mchbpXE , GetEventWeight());
559  fhMCPtZTCharged [histoIndex]->Fill(mcTrigPt , mczT , GetEventWeight());
560  fhMCPtHbpZTCharged [histoIndex]->Fill(mcTrigPt , mchbpZT , GetEventWeight());
561  fhMCPtTrigPout [histoIndex]->Fill(mcTrigPt , mcpout , GetEventWeight());
562  }
563 
564  if(lostDecayPair)
565  {
566  // check index in GetMCTagIndexHistogram
567  if(histoIndex == 2 && 8 >= fMCGenTypeMin && 8 <= fMCGenTypeMax )
568  {
569  // pi0 decay
570  fhMCEtaCharged [8]->Fill(mcAssocPt, mcAssocEta, GetEventWeight());
571  fhMCPhiCharged [8]->Fill(mcAssocPt, mcAssocPhi, GetEventWeight());
572  fhMCDeltaPhiCharged[8]->Fill(mcTrigPt , mcdeltaPhi, GetEventWeight());
573  fhMCPtAssocDeltaPhi[8]->Fill(mcAssocPt, mcdeltaPhi, GetEventWeight());
574 
575  fhMCDeltaEtaCharged [8]->Fill(mcTrigPt , mcTrigEta-mcAssocEta, GetEventWeight());
576  fhMCDeltaPhiDeltaEtaCharged[8]->Fill(mcdeltaPhi, mcTrigEta-mcAssocEta, GetEventWeight());
577 
578  //delta phi cut for correlation
579  if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) )
580  {
581  fhMCDeltaPhiChargedPt[8]->Fill(mcAssocPt, mcdeltaPhi, GetEventWeight());
582  fhMCPtXECharged [8]->Fill(mcTrigPt , mcxE , GetEventWeight());
583  fhMCPtHbpXECharged [8]->Fill(mcTrigPt , mchbpXE , GetEventWeight());
584  fhMCPtZTCharged [8]->Fill(mcTrigPt , mczT , GetEventWeight());
585  fhMCPtHbpZTCharged [8]->Fill(mcTrigPt , mchbpZT , GetEventWeight());
586  fhMCPtTrigPout [8]->Fill(mcTrigPt , mcpout , GetEventWeight());
587  }
588  } // pi0 decay lost pair
589 
590  if(histoIndex == 4 && 9 >= fMCGenTypeMin && 9 <= fMCGenTypeMax )
591  {
592  // eta decay
593  fhMCEtaCharged [9]->Fill(mcAssocPt, mcAssocEta, GetEventWeight());
594  fhMCPhiCharged [9]->Fill(mcAssocPt, mcAssocPhi, GetEventWeight());
595  fhMCDeltaPhiCharged[9]->Fill(mcTrigPt , mcdeltaPhi, GetEventWeight());
596  fhMCPtAssocDeltaPhi[9]->Fill(mcAssocPt, mcdeltaPhi, GetEventWeight());
597 
598  fhMCDeltaEtaCharged [9]->Fill(mcTrigPt , mcTrigEta-mcAssocEta, GetEventWeight());
599  fhMCDeltaPhiDeltaEtaCharged[9]->Fill(mcdeltaPhi, mcTrigEta-mcAssocEta, GetEventWeight());
600 
601  // Delta phi cut for correlation
602  if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) )
603  {
604  fhMCDeltaPhiChargedPt[9]->Fill(mcAssocPt, mcdeltaPhi, GetEventWeight());
605  fhMCPtXECharged [9]->Fill(mcTrigPt , mcxE , GetEventWeight());
606  fhMCPtHbpXECharged [9]->Fill(mcTrigPt , mchbpXE , GetEventWeight());
607  fhMCPtZTCharged [9]->Fill(mcTrigPt , mczT , GetEventWeight());
608  fhMCPtHbpZTCharged [9]->Fill(mcTrigPt , mchbpZT , GetEventWeight());
609  fhMCPtTrigPout [9]->Fill(mcTrigPt , mcpout , GetEventWeight()) ;
610  }
611  } // eta decay lost pair
612  }
613  // Underlying event
614 
615  // Right
616  if ( (mcdeltaPhi > fUeDeltaPhiMinCut) && (mcdeltaPhi < fUeDeltaPhiMaxCut) )
617  {
618  //Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
619  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
620  Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
621  Double_t mcUezT = mcAssocPt/mcTrigPt;
622 
623  if(mcUexE < 0.)
624  AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
625  mcUexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
626 
627  fhMCPtXEUeCharged[histoIndex]->Fill(mcTrigPt, mcUexE, GetEventWeight());
628  if(mcUexE > 0) fhMCPtHbpXEUeCharged[histoIndex]->Fill(mcTrigPt, TMath::Log(1/mcUexE), GetEventWeight());
629 
630  fhMCPtZTUeCharged[histoIndex]->Fill(mcTrigPt, mcUezT, GetEventWeight());
631  if(mcUezT > 0) fhMCPtHbpZTUeCharged[histoIndex]->Fill(mcTrigPt, TMath::Log(1/mcUezT), GetEventWeight());
632 
633  fhMCUePart[histoIndex]->Fill(mcTrigPt, GetEventWeight());
634 
635  if(lostDecayPair)
636  {
637  // check index in GetMCTagIndexHistogram
638  if(histoIndex == 2 && 8 >= fMCGenTypeMin && 8 <= fMCGenTypeMax )
639  {
640  // pi0 decay
641  fhMCPtXEUeCharged[8]->Fill(mcTrigPt, mcUexE, GetEventWeight());
642  if(mcUexE > 0) fhMCPtHbpXEUeCharged[8]->Fill(mcTrigPt, TMath::Log(1/mcUexE), GetEventWeight());
643 
644  fhMCPtZTUeCharged[8]->Fill(mcTrigPt, mcUezT, GetEventWeight());
645  if(mcUezT > 0) fhMCPtHbpZTUeCharged[8]->Fill(mcTrigPt, TMath::Log(1/mcUezT), GetEventWeight());
646 
647  fhMCUePart[8]->Fill(mcTrigPt, GetEventWeight());
648  }
649  else if(histoIndex == 4 && 9 >= fMCGenTypeMin && 9 <= fMCGenTypeMax )
650  {
651  // eta decay
652  fhMCPtXEUeCharged[9]->Fill(mcTrigPt, mcUexE, GetEventWeight());
653  if(mcUexE > 0) fhMCPtHbpXEUeCharged[9]->Fill(mcTrigPt, TMath::Log(1/mcUexE), GetEventWeight());
654 
655  fhMCPtZTUeCharged[9]->Fill(mcTrigPt, mcUezT, GetEventWeight());
656  if(mcUezT > 0) fhMCPtHbpZTUeCharged[9]->Fill(mcTrigPt, TMath::Log(1/mcUezT), GetEventWeight());
657 
658  fhMCUePart[9]->Fill(mcTrigPt, GetEventWeight());
659  }
660  }
661  }
662 
663  if(fMakeSeveralUE)
664  {
665  // Left
666  if((mcdeltaPhi<-fUeDeltaPhiMinCut) || (mcdeltaPhi >2*fUeDeltaPhiMaxCut))
667  {
668  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
669  Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
670  Double_t mcUezT = mcAssocPt/mcTrigPt;
671 
672  if(mcUexE < 0.)
673  AliWarning(Form("Careful!!, negative xE %2.2f for left UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
674  mcUexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
675 
676  fhMCPtXEUeLeftCharged[histoIndex]->Fill(mcTrigPt, mcUexE, GetEventWeight());
677  if(mcUexE > 0) fhMCPtHbpXEUeLeftCharged[histoIndex]->Fill(mcTrigPt, TMath::Log(1/mcUexE), GetEventWeight());
678 
679  fhMCPtZTUeLeftCharged[histoIndex]->Fill(mcTrigPt, mcUezT, GetEventWeight());
680  if(mcUezT > 0) fhMCPtHbpZTUeLeftCharged[histoIndex]->Fill(mcTrigPt, TMath::Log(1/mcUezT), GetEventWeight());
681 
682  if(lostDecayPair)
683  {
684  // Check index in GetMCTagIndexHistogram
685  if(histoIndex == 2 && 8 >= fMCGenTypeMin && 8 <= fMCGenTypeMax )
686  {
687  // pi0 decay
688  fhMCPtXEUeLeftCharged[8]->Fill(mcTrigPt,mcUexE, GetEventWeight());
689  if(mcUexE > 0) fhMCPtHbpXEUeLeftCharged[8]->Fill(mcTrigPt, TMath::Log(1/mcUexE), GetEventWeight());
690 
691  fhMCPtZTUeLeftCharged[8]->Fill(mcTrigPt,mcUezT, GetEventWeight());
692  if(mcUezT > 0) fhMCPtHbpZTUeLeftCharged[8]->Fill(mcTrigPt, TMath::Log(1/mcUezT), GetEventWeight());
693  }
694  else if(histoIndex == 4 && 9 >= fMCGenTypeMin && 9 <= fMCGenTypeMax )
695  {
696  // eta decay
697  fhMCPtXEUeLeftCharged[9]->Fill(mcTrigPt, mcUexE, GetEventWeight());
698  if(mcUexE > 0) fhMCPtHbpXEUeLeftCharged[9]->Fill(mcTrigPt, TMath::Log(1/mcUexE), GetEventWeight());
699 
700  fhMCPtZTUeLeftCharged[9]->Fill(mcTrigPt, mcUezT, GetEventWeight());
701  if(mcUezT > 0) fhMCPtHbpZTUeLeftCharged[9]->Fill(mcTrigPt, TMath::Log(1/mcUezT), GetEventWeight());
702  }
703  }
704  }
705  }
706 
707  return lead;
708 }
709 
710 //______________________________________________________________________________________________________________
712 //______________________________________________________________________________________________________________
713 void AliAnaParticleHadronCorrelation::FillChargedMomentumImbalanceHistograms(Float_t ptTrig, Float_t ptAssoc,
714  Float_t deltaPhi, Int_t sm,
715  Int_t cen, Int_t charge,
716  Int_t bin, Int_t decayTag,
717  Int_t outTOF, Int_t mcTag)
718 
719 {
720  Float_t zT = ptAssoc/ptTrig ;
721  Float_t xE =-ptAssoc/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
722  Float_t pout = ptAssoc*TMath::Sin(deltaPhi) ;
723 
724  if(xE < 0.)
725  AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
726  xE,deltaPhi*TMath::RadToDeg(),TMath::Cos(deltaPhi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
727 
728  Float_t hbpXE = -100;
729  Float_t hbpZT = -100;
730 
731  if(xE > 0 ) hbpXE = TMath::Log(1./xE);
732  if(zT > 0 ) hbpZT = TMath::Log(1./zT);
733 
734  fhXECharged ->Fill(ptTrig, xE , GetEventWeight());
735  fhPtHbpXECharged ->Fill(ptTrig, hbpXE , GetEventWeight());
736  fhZTCharged ->Fill(ptTrig, zT , GetEventWeight());
737  fhPtHbpZTCharged ->Fill(ptTrig, hbpZT , GetEventWeight());
738  fhPtTrigPout ->Fill(ptTrig, pout , GetEventWeight());
739  fhPtTrigCharged ->Fill(ptTrig, ptAssoc, GetEventWeight());
740 
741  if ( fFillPerSMHistograms ) fhXEChargedPerSM[sm]->Fill(ptTrig, xE, GetEventWeight());
742 
743  if((deltaPhi > 5*TMath::Pi()/6.) && (deltaPhi < 7*TMath::Pi()/6.))
744  {
745  fhXECharged_Cone2 ->Fill(ptTrig, xE , GetEventWeight());
746  fhPtHbpXECharged_Cone2 ->Fill(ptTrig, hbpXE, GetEventWeight());
747  }
748 
749  // MC
750  if(IsDataMC())
751  {
752  Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
753  fhXEChargedMC[mcIndex]->Fill(ptTrig, xE, GetEventWeight());
754 
755  if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) )
756  {
757  // check index in GetMCTagIndexHistogram
758  if ( mcIndex == 2 ) fhXEChargedMC[8]->Fill(ptTrig, xE, GetEventWeight()); // pi0 decay
759  else if ( mcIndex == 4 ) fhXEChargedMC[9]->Fill(ptTrig, xE, GetEventWeight()); // eta decay
760  }
761  }
762 
763  // Pile-up studies
764  if(IsPileUpAnalysisOn())
765  {
766  if (outTOF==1)
767  {
768  fhXEChargedOtherBC ->Fill(ptTrig, xE , GetEventWeight());
769  fhZTChargedOtherBC ->Fill(ptTrig, zT , GetEventWeight());
770  fhPtTrigChargedOtherBC->Fill(ptTrig, ptAssoc, GetEventWeight());
771  }
772  else if(outTOF==0)
773  {
774  fhXEChargedBC0 ->Fill(ptTrig, xE , GetEventWeight());
775  fhZTChargedBC0 ->Fill(ptTrig, zT , GetEventWeight());
776  fhPtTrigChargedBC0->Fill(ptTrig, ptAssoc, GetEventWeight());
777  }
778 
779  Int_t vtxBC = GetReader()->GetVertexBC();
780  if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
781  {
782  fhXEChargedVtxBC0 ->Fill(ptTrig, xE , GetEventWeight());
783  fhZTChargedVtxBC0 ->Fill(ptTrig, zT , GetEventWeight());
784  fhPtTrigChargedVtxBC0->Fill(ptTrig, ptAssoc, GetEventWeight());
785  }
786 
787  if(GetReader()->IsPileUpFromSPD())
788  {
789  fhXEChargedPileUp [0]->Fill(ptTrig, xE , GetEventWeight());
790  fhZTChargedPileUp [0]->Fill(ptTrig, zT , GetEventWeight());
791  fhPtTrigChargedPileUp[0]->Fill(ptTrig, ptAssoc, GetEventWeight());
792  }
793  if(GetReader()->IsPileUpFromEMCal())
794  {
795  fhXEChargedPileUp [1]->Fill(ptTrig, xE , GetEventWeight());
796  fhZTChargedPileUp [1]->Fill(ptTrig, zT , GetEventWeight());
797  fhPtTrigChargedPileUp[1]->Fill(ptTrig, ptAssoc, GetEventWeight());
798  }
799  if(GetReader()->IsPileUpFromSPDOrEMCal())
800  {
801  fhXEChargedPileUp [2]->Fill(ptTrig, xE , GetEventWeight());
802  fhZTChargedPileUp [2]->Fill(ptTrig, zT , GetEventWeight());
803  fhPtTrigChargedPileUp[2]->Fill(ptTrig, ptAssoc, GetEventWeight());
804  }
805  if(GetReader()->IsPileUpFromSPDAndEMCal())
806  {
807  fhXEChargedPileUp [3]->Fill(ptTrig, xE , GetEventWeight());
808  fhZTChargedPileUp [3]->Fill(ptTrig, zT , GetEventWeight());
809  fhPtTrigChargedPileUp[3]->Fill(ptTrig, ptAssoc, GetEventWeight());
810  }
811  if(GetReader()->IsPileUpFromSPDAndNotEMCal())
812  {
813  fhXEChargedPileUp [4]->Fill(ptTrig, xE , GetEventWeight());
814  fhZTChargedPileUp [4]->Fill(ptTrig, zT , GetEventWeight());
815  fhPtTrigChargedPileUp[4]->Fill(ptTrig, ptAssoc, GetEventWeight());
816  }
817  if(GetReader()->IsPileUpFromEMCalAndNotSPD())
818  {
819  fhXEChargedPileUp [5]->Fill(ptTrig, xE , GetEventWeight());
820  fhZTChargedPileUp [5]->Fill(ptTrig, zT , GetEventWeight());
821  fhPtTrigChargedPileUp[5]->Fill(ptTrig, ptAssoc, GetEventWeight());
822  }
823  if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
824  {
825  fhXEChargedPileUp [6]->Fill(ptTrig, xE , GetEventWeight());
826  fhZTChargedPileUp [6]->Fill(ptTrig, zT , GetEventWeight());
827  fhPtTrigChargedPileUp[6]->Fill(ptTrig, ptAssoc, GetEventWeight());
828  }
829  }
830 
831  if(fDecayTrigger && decayTag > 0)
832  {
833  for(Int_t ibit = 0; ibit<fNDecayBits; ibit++)
834  {
835  if(GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[ibit]))
836  {
837  fhXEDecayCharged[ibit]->Fill(ptTrig, xE, GetEventWeight());
838  fhZTDecayCharged[ibit]->Fill(ptTrig, zT, GetEventWeight());
839  }
840  }
841  } // photon decay pi0/eta trigger
842 
843  if(bin >= 0 && fFillMomImbalancePtAssocBinsHisto)//away side
844  {
845  fhXEAssocPtBin[bin]->Fill(ptTrig, xE, GetEventWeight()) ;
846  fhZTAssocPtBin[bin]->Fill(ptTrig, zT, GetEventWeight()) ;
847  }
848 
849  if(fCorrelVzBin)
850  {
851  Int_t vz = GetEventVzBin();
852  fhXEVZ[vz]->Fill(ptTrig, xE, GetEventWeight()) ;
853  fhZTVZ[vz]->Fill(ptTrig, zT, GetEventWeight()) ;
854  }
855 
856  if(charge > 0)
857  {
858  fhXEPosCharged->Fill(ptTrig, xE, GetEventWeight()) ;
859  fhZTPosCharged->Fill(ptTrig, zT, GetEventWeight()) ;
860  }
861  else
862  {
863  fhXENegCharged->Fill(ptTrig, xE, GetEventWeight()) ;
864  fhZTNegCharged->Fill(ptTrig, zT, GetEventWeight()) ;
865  }
866 
867  // Fill different multiplicity/centrality histogram
868  if(IsHighMultiplicityAnalysisOn() && cen >= 0 && cen < GetNCentrBin())
869  {
870  fhXEMult[cen]->Fill(ptTrig, xE, GetEventWeight());
871  fhZTMult[cen]->Fill(ptTrig, zT, GetEventWeight());
872  } // multiplicity/centrality events selection
873 }
874 
875 //___________________________________________________________________________________________________________________
877 //___________________________________________________________________________________________________________________
878 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventHistograms(Float_t ptTrig, Float_t ptAssoc, Float_t deltaPhi,
879  Int_t sm, Int_t cen, Int_t outTOF, Int_t mcTag)
880 {
881  fhUePart->Fill(ptTrig, GetEventWeight());
882 
883  fhDeltaPhiUeChargedPt->Fill(ptAssoc, deltaPhi, GetEventWeight());
884 
885  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
886  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
887  Double_t uezT = ptAssoc/ptTrig;
888 
889  if(uexE < 0.)
890  AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
891  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
892 
893  fhXEUeCharged->Fill(ptTrig, uexE, GetEventWeight());
894 
895  if ( fFillPerSMHistograms ) fhXEUeChargedPerSM[sm]->Fill(ptTrig, uexE, GetEventWeight());
896 
897  if(deltaPhi > TMath::DegToRad()*80 && deltaPhi < TMath::DegToRad()*100)fhXEUeChargedSmallCone->Fill(ptTrig, uexE, GetEventWeight());
898  if(deltaPhi > TMath::DegToRad()*70 && deltaPhi < TMath::DegToRad()*110)fhXEUeChargedMediumCone->Fill(ptTrig, uexE, GetEventWeight());
899  if(deltaPhi > TMath::DegToRad()*60 && deltaPhi < TMath::DegToRad()*120)fhXEUeChargedLargeCone->Fill(ptTrig, uexE, GetEventWeight());
900 
901  if(uexE > 0) fhPtHbpXEUeCharged->Fill(ptTrig, TMath::Log(1/uexE), GetEventWeight());
902 
903  fhZTUeCharged->Fill(ptTrig, uezT, GetEventWeight());
904  if(uezT > 0) fhPtHbpZTUeCharged->Fill(ptTrig, TMath::Log(1/uezT), GetEventWeight());
905 
906  if(IsDataMC())
907  {
908  Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
909  fhXEUeChargedRightMC[mcIndex]->Fill(ptTrig, uexE, GetEventWeight());
910 
911  if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) )
912  {
913  // check index in GetMCTagIndexHistogram
914  if ( mcIndex == 2 ) fhXEUeChargedRightMC[8]->Fill(ptTrig, uexE, GetEventWeight()); // pi0 decay
915  else if ( mcIndex == 4 ) fhXEUeChargedRightMC[9]->Fill(ptTrig, uexE, GetEventWeight()); // eta decay
916  }
917  }
918 
919  // Pile-up studies
920 
921  if(IsPileUpAnalysisOn())
922  {
923  if (outTOF==1)
924  {
925  fhXEUeChargedOtherBC->Fill(ptTrig, uexE, GetEventWeight());
926  fhZTUeChargedOtherBC->Fill(ptTrig, uezT, GetEventWeight());
927  }
928  else if(outTOF==0)
929  {
930  fhXEUeChargedBC0->Fill(ptTrig, uexE, GetEventWeight());
931  fhZTUeChargedBC0->Fill(ptTrig, uezT, GetEventWeight());
932  }
933 
934  Int_t vtxBC = GetReader()->GetVertexBC();
935  if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
936  {
937  fhXEUeChargedVtxBC0->Fill(ptTrig, uexE, GetEventWeight());
938  fhZTUeChargedVtxBC0->Fill(ptTrig, uezT, GetEventWeight());
939  }
940 
941  if(GetReader()->IsPileUpFromSPD())
942  {
943  fhXEUeChargedPileUp[0]->Fill(ptTrig, uexE, GetEventWeight());
944  fhZTUeChargedPileUp[0]->Fill(ptTrig, uezT, GetEventWeight());
945  }
946  if(GetReader()->IsPileUpFromEMCal())
947  {
948  fhXEUeChargedPileUp[1]->Fill(ptTrig, uexE, GetEventWeight());
949  fhZTUeChargedPileUp[1]->Fill(ptTrig, uezT, GetEventWeight());
950  }
951  if(GetReader()->IsPileUpFromSPDOrEMCal())
952  {
953  fhXEUeChargedPileUp[2]->Fill(ptTrig, uexE, GetEventWeight());
954  fhZTUeChargedPileUp[2]->Fill(ptTrig, uezT, GetEventWeight());
955  }
956  if(GetReader()->IsPileUpFromSPDAndEMCal())
957  {
958  fhXEUeChargedPileUp[3]->Fill(ptTrig, uexE, GetEventWeight());
959  fhZTUeChargedPileUp[3]->Fill(ptTrig, uezT, GetEventWeight());
960  }
961  if(GetReader()->IsPileUpFromSPDAndNotEMCal())
962  {
963  fhXEUeChargedPileUp[4]->Fill(ptTrig, uexE, GetEventWeight());
964  fhZTUeChargedPileUp[4]->Fill(ptTrig, uezT, GetEventWeight());
965  }
966  if(GetReader()->IsPileUpFromEMCalAndNotSPD())
967  {
968  fhXEUeChargedPileUp[5]->Fill(ptTrig, uexE, GetEventWeight());
969  fhZTUeChargedPileUp[5]->Fill(ptTrig, uezT, GetEventWeight());
970  }
971  if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
972  {
973  fhXEUeChargedPileUp[6]->Fill(ptTrig, uexE, GetEventWeight());
974  fhZTUeChargedPileUp[6]->Fill(ptTrig, uezT, GetEventWeight());
975  }
976  }
977 
978  // Fill different multiplicity/centrality histogram
979  if(IsHighMultiplicityAnalysisOn() && cen >= 0 && cen < GetNCentrBin())
980  {
981  fhXEUeMult[cen]->Fill(ptTrig, uexE, GetEventWeight());
982  fhZTUeMult[cen]->Fill(ptTrig, uezT, GetEventWeight());
983  }
984 }
985 
986 //_______________________________________________________________________________________________
989 //_______________________________________________________________________________________________
990 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventSidesHistograms(Float_t ptTrig,
991  Float_t ptAssoc,
992  Float_t deltaPhi,
993  Int_t mcTag)
994 {
995  if((deltaPhi<-fUeDeltaPhiMinCut) || (deltaPhi >2*fUeDeltaPhiMaxCut))
996  {
997  fhDeltaPhiUeLeftCharged->Fill(ptAssoc, deltaPhi, GetEventWeight());
998 
999  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
1000  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
1001  Double_t uezT = ptAssoc/ptTrig;
1002 
1003  if(uexE < 0.)
1004  AliWarning(Form("Careful!!, negative xE %2.2f for left UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
1005  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
1006 
1007  fhXEUeLeftCharged->Fill(ptTrig, uexE, GetEventWeight());
1008  if(uexE > 0) fhPtHbpXEUeLeftCharged->Fill(ptTrig, TMath::Log(1/uexE), GetEventWeight());
1009 
1010  fhZTUeLeftCharged->Fill(ptTrig, uezT, GetEventWeight());
1011  if(uezT > 0) fhPtHbpZTUeLeftCharged->Fill(ptTrig, TMath::Log(1/uezT), GetEventWeight());
1012 
1013  fhDeltaPhiUeLeftCharged->Fill(ptAssoc, deltaPhi, GetEventWeight());
1014 
1015  if(IsDataMC())
1016  {
1017  Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
1018  fhXEUeChargedLeftMC[mcIndex]->Fill(ptTrig, uexE, GetEventWeight());
1019 
1020  if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) )
1021  {
1022  // check index in GetMCTagIndexHistogram
1023  if ( mcIndex == 2 ) fhXEUeChargedLeftMC[8]->Fill(ptTrig, uexE, GetEventWeight()); // pi0 decay
1024  else if ( mcIndex == 4 ) fhXEUeChargedLeftMC[9]->Fill(ptTrig, uexE, GetEventWeight()); // eta decay
1025  }
1026  }
1027  }
1028 
1029  if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-TMath::Pi()/2))
1030  {
1031  fhDeltaPhiUeLeftDownCharged->Fill(ptAssoc, deltaPhi, GetEventWeight());
1032 
1033  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
1034  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
1035 
1036  if(uexE < 0.)
1037  AliWarning(Form("Careful!!, negative xE %2.2f for left-down UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
1038  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
1039 
1040  fhXEUeLeftDownCharged->Fill(ptTrig, uexE, GetEventWeight());
1041  }
1042 
1043  if((deltaPhi>2*fUeDeltaPhiMaxCut) && (deltaPhi <3*TMath::Pi()/2))
1044  {
1045  fhDeltaPhiUeLeftUpCharged->Fill(ptAssoc, deltaPhi, GetEventWeight());
1046 
1047  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
1048  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
1049 
1050  if(uexE < 0.)
1051  AliWarning(Form("Careful!!, negative xE %2.2f for left-up UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
1052  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
1053 
1054  fhXEUeLeftUpCharged->Fill(ptTrig, uexE, GetEventWeight());
1055  }
1056 
1057  if((deltaPhi > TMath::Pi()/2) && (deltaPhi < fUeDeltaPhiMaxCut))
1058  {
1059  fhDeltaPhiUeRightUpCharged->Fill(ptAssoc, deltaPhi, GetEventWeight());
1060 
1061  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
1062  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
1063 
1064  if(uexE < 0.)
1065  AliWarning(Form("Careful!!, negative xE %2.2f for right-up UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
1066  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
1067 
1068  fhXEUeRightUpCharged->Fill(ptTrig, uexE, GetEventWeight());
1069  }
1070 
1071  if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < TMath::Pi()/2))
1072  {
1073  fhDeltaPhiUeRightDownCharged->Fill(ptAssoc, deltaPhi, GetEventWeight());
1074 
1075  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
1076  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
1077 
1078  if(uexE < 0.)
1079  AliWarning(Form("Careful!!, negative xE %2.2f for right-down UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
1080  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
1081 
1082  fhXEUeRightDownCharged->Fill(ptTrig, uexE, GetEventWeight());
1083  }
1084 }
1085 
1086 //_____________________________________________________________________________________________________________________________________
1088 //_____________________________________________________________________________________________________________________________________
1089 void AliAnaParticleHadronCorrelation::FillDecayPhotonCorrelationHistograms(Float_t ptAssoc, Float_t phiAssoc, Bool_t bChargedOrNeutral)
1090 {
1091  // Calculate the correlation parameters
1092  Float_t ptDecay1 = fDecayMom1.Pt();
1093  Float_t ptDecay2 = fDecayMom2.Pt();
1094 
1095  Float_t zTDecay1 = -100, zTDecay2 = -100;
1096  if(ptDecay1 > 0) zTDecay1 = ptAssoc/ptDecay1 ;
1097  if(ptDecay2 > 0) zTDecay2 = ptAssoc/ptDecay2 ;
1098 
1099  Float_t deltaPhiDecay1 = fDecayMom1.Phi()-phiAssoc;
1100  if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
1101  if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
1102 
1103  Float_t deltaPhiDecay2 = fDecayMom2.Phi()-phiAssoc;
1104  if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
1105  if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
1106 
1107  Float_t xEDecay1 =-zTDecay1*TMath::Cos(deltaPhiDecay1); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
1108  Float_t xEDecay2 =-zTDecay2*TMath::Cos(deltaPhiDecay2); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
1109 
1110  if(bChargedOrNeutral) // correlate with charges
1111  {
1112  fhDeltaPhiPi0DecayCharged->Fill(ptDecay1, deltaPhiDecay1, GetEventWeight());
1113  fhDeltaPhiPi0DecayCharged->Fill(ptDecay2, deltaPhiDecay2, GetEventWeight());
1114 
1115  AliDebug(2,Form("deltaPhoton1 = %f, deltaPhoton2 = %f", deltaPhiDecay1, deltaPhiDecay2));
1116 
1117  if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
1118  {
1119  fhZTPi0DecayCharged->Fill(ptDecay1, zTDecay1, GetEventWeight());
1120  fhXEPi0DecayCharged->Fill(ptDecay1, xEDecay1, GetEventWeight());
1121  }
1122  if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
1123  {
1124  fhZTPi0DecayCharged->Fill(ptDecay2, zTDecay2, GetEventWeight());
1125  fhXEPi0DecayCharged->Fill(ptDecay2, xEDecay2, GetEventWeight());
1126  }
1127  }
1128  else // correlate with neutrals
1129  {
1130  fhDeltaPhiPi0DecayNeutral->Fill(ptDecay1, deltaPhiDecay1, GetEventWeight());
1131  fhDeltaPhiPi0DecayNeutral->Fill(ptDecay2, deltaPhiDecay2, GetEventWeight());
1132 
1133  AliDebug(2,Form("deltaPhoton1 = %f, deltaPhoton2 = %f", deltaPhiDecay1, deltaPhiDecay2));
1134 
1135  if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
1136  {
1137  fhZTPi0DecayNeutral->Fill(ptDecay1, zTDecay1, GetEventWeight());
1138  fhXEPi0DecayNeutral->Fill(ptDecay1, xEDecay1, GetEventWeight());
1139  }
1140  if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
1141  {
1142  fhZTPi0DecayNeutral->Fill(ptDecay2, zTDecay2, GetEventWeight());
1143  fhXEPi0DecayNeutral->Fill(ptDecay2, xEDecay2, GetEventWeight());
1144  }
1145  }
1146 }
1147 
1148 //_____________________________________________________________________________________________________________________________
1151 //_____________________________________________________________________________________________________________________________
1152 void AliAnaParticleHadronCorrelation::FillNeutralUnderlyingEventSidesHistograms(Float_t ptTrig, Float_t ptAssoc,
1153  Float_t zT, Float_t hbpZT,
1154  Float_t deltaPhi)
1155 {
1156  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
1157 
1158  Float_t xE =-ptAssoc/ptTrig*TMath::Cos(randomphi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
1159  Float_t hbpXE = -100;
1160  if(xE > 0 ) hbpXE = TMath::Log(1./xE);
1161 
1162  if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut))
1163  {
1164  fhDeltaPhiUeLeftNeutral->Fill(ptAssoc, deltaPhi, GetEventWeight());
1165  fhXEUeLeftNeutral ->Fill(ptTrig , xE , GetEventWeight());
1166  fhPtHbpXEUeLeftNeutral ->Fill(ptTrig , hbpXE , GetEventWeight());
1167  fhZTUeLeftNeutral ->Fill(ptTrig , zT , GetEventWeight());
1168  fhPtHbpZTUeLeftNeutral ->Fill(ptTrig , hbpZT , GetEventWeight());
1169  }
1170 }
1171 
1172 //______________________________________________________
1174 //______________________________________________________
1175 void AliAnaParticleHadronCorrelation::FillEventMixPool()
1176 {
1177  if ( !DoOwnMix() ) return;
1178 
1179  FillChargedEventMixPool();
1180 
1181  // Do the cluster pool filling only if requested
1182  // or in case of isolation cut using clusters in the cone.
1183  Bool_t isoCase = OnlyIsolated() && (GetIsolationCut()->GetParticleTypeInCone() != AliIsolationCut::kOnlyCharged);
1184 
1185  if( !fFillNeutralEventMixPool && !isoCase) return;
1186 
1187  FillNeutralEventMixPool();
1188 }
1189 
1190 //_____________________________________________________________
1192 //_____________________________________________________________
1193 void AliAnaParticleHadronCorrelation::FillChargedEventMixPool()
1194 {
1195  if(fUseMixStoredInReader && GetReader()->GetLastTracksMixedEvent() == GetEventNumber())
1196  {
1197  //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
1198  return ; // pool filled previously for another trigger
1199  }
1200 
1201  AliAnalysisManager * manager = AliAnalysisManager::GetAnalysisManager();
1202  AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
1203 
1204  if(!inputHandler) return ;
1205 
1206  // Do mixing only with MB event (or the chosen mask), if not skip
1207  if( !(inputHandler->IsEventSelected( ) & GetReader()->GetMixEventTriggerMask()) ) return ;
1208 
1209  Int_t eventBin = GetEventMixBin();
1210 
1211  //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
1212  if(eventBin < 0) return;
1213 
1214  fhEventMBBin->Fill(eventBin, GetEventWeight());
1215 
1216  TObjArray * mixEventTracks = new TObjArray;
1217 
1218  if(fUseMixStoredInReader)
1219  {
1220  fListMixTrackEvents[eventBin] = GetReader()->GetListWithMixedEventsForTracks(eventBin);
1221  }
1222 
1223  if(!fListMixTrackEvents[eventBin]) fListMixTrackEvents[eventBin] = new TList();
1224 
1225  //printf("%s ***** Pool Event bin : %d - nTracks %d\n",GetInputAODName().Data(),eventBin, GetCTSTracks()->GetEntriesFast());
1226 
1227  TList * pool = fListMixTrackEvents[eventBin];
1228 
1229  for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
1230  {
1231  AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
1232 
1233  fTrackVector.SetXYZ(track->Px(),track->Py(),track->Pz());
1234  Float_t pt = fTrackVector.Pt();
1235 
1236  // Select only hadrons in pt range
1237  if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
1238 
1239  AliCaloTrackParticle * mixedTrack = new AliCaloTrackParticle(track->Px(),track->Py(),track->Pz(),0);
1240  mixedTrack->SetDetectorTag(kCTS);
1241  mixedTrack->SetChargedBit(track->Charge()>0);
1242  mixEventTracks->Add(mixedTrack);
1243  }
1244 
1245  fhNtracksMB->Fill(mixEventTracks->GetEntriesFast(), eventBin, GetEventWeight());
1246 
1247  // Set the event number where the last event was added, to avoid double pool filling
1248  GetReader()->SetLastTracksMixedEvent(GetEventNumber());
1249 
1250  //printf("Add event to pool with %d tracks \n ",mixEventTracks->GetEntries());
1251  pool->AddFirst(mixEventTracks);
1252  mixEventTracks = 0;
1253 
1254  //printf("Pool size %d, max %d\n",pool->GetSize(), GetNMaxEvMix());
1255 
1256  if(pool->GetSize() > GetNMaxEvMix())
1257  {
1258  // Remove last event
1259  TClonesArray * tmp = static_cast<TClonesArray*>(pool->Last()) ;
1260  pool->RemoveLast() ;
1261  delete tmp ;
1262  }
1263 }
1264 
1265 //_____________________________________________________________
1268 //_____________________________________________________________
1269 void AliAnaParticleHadronCorrelation::FillNeutralEventMixPool()
1270 {
1271  //printf("FillNeutralEventMixPool for %s\n",GetInputAODName().Data());
1272 
1273  if(fUseMixStoredInReader && GetReader()->GetLastCaloMixedEvent() == GetEventNumber())
1274  {
1275  //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
1276  return ; // pool filled previously for another trigger
1277  }
1278 
1279  TObjArray * pl = GetEMCALClusters();
1280  //if (GetAODObjArrayName.Contains("PHOS") )pl = GetPHOSClusters();
1281  //else pl = GetEMCALClusters();
1282 
1283  AliAnalysisManager * manager = AliAnalysisManager::GetAnalysisManager();
1284  AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
1285 
1286  if(!inputHandler) return ;
1287 
1288  // Do mixing only with MB event (or the chosen mask), if not skip
1289  if( !(inputHandler->IsEventSelected( ) & GetReader()->GetMixEventTriggerMask()) ) return ;
1290 
1291  Int_t eventBin = GetEventMixBin();
1292 
1293  // Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
1294  if(eventBin < 0) return;
1295 
1296  TObjArray * mixEventCalo = new TObjArray;
1297 
1298  if(fUseMixStoredInReader)
1299  {
1300  fListMixCaloEvents[eventBin] = GetReader()->GetListWithMixedEventsForCalo(eventBin);
1301  }
1302 
1303  if(!fListMixCaloEvents[eventBin]) fListMixCaloEvents[eventBin] = new TList();
1304 
1305  TList * poolCalo = fListMixCaloEvents[eventBin];
1306 
1307  for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ )
1308  {
1309  AliVCluster * calo = (AliVCluster *) (pl->At(ipr)) ;
1310 
1311  // Remove matched clusters
1312  if( IsTrackMatched( calo, GetReader()->GetInputEvent() ) ) continue ;
1313 
1314  // Cluster momentum calculation
1315  if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1316  {
1317  calo->GetMomentum(fMomentum,GetVertex(0)) ;
1318  }// Assume that come from vertex in straight line
1319  else
1320  {
1321  Double_t vertex[]={0,0,0};
1322  calo->GetMomentum(fMomentum,vertex) ;
1323  }
1324 
1325  Float_t pt = fMomentum.Pt();
1326 
1327  // Select only clusters in pt range
1328  if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
1329 
1330  AliCaloTrackParticle * mixedCalo = new AliCaloTrackParticle(fMomentum);
1331  mixedCalo->SetDetectorTag(kEMCAL);
1332  mixEventCalo->Add(mixedCalo);
1333  }
1334 
1335  fhNclustersMB->Fill(mixEventCalo->GetEntriesFast(), eventBin, GetEventWeight());
1336 
1337  // Set the event number where the last event was added, to avoid double pool filling
1338  GetReader()->SetLastCaloMixedEvent(GetEventNumber());
1339 
1340  //printf("Add event to pool with %d clusters \n ",mixEventCalo->GetEntries());
1341  poolCalo->AddFirst(mixEventCalo);
1342  mixEventCalo = 0;
1343 
1344  //printf("Pool size %d, max %d\n",poolCalo->GetSize(), GetNMaxEvMix());
1345 
1346  if(poolCalo->GetSize() > GetNMaxEvMix())
1347  {
1348  // Remove last event
1349  TClonesArray * tmp = static_cast<TClonesArray*>(poolCalo->Last()) ;
1350  poolCalo->RemoveLast() ;
1351  delete tmp ;
1352  }
1353 }
1354 
1355 //_________________________________________________________________________________________________________________
1358 //_________________________________________________________________________________________________________________
1359 Bool_t AliAnaParticleHadronCorrelation::FindLeadingOppositeHadronInWindow(AliCaloTrackParticleCorrelation * particle)
1360 {
1361  Float_t etaTrig = particle->Eta();
1362  Float_t ptTrig = particle->Pt();
1363  Float_t phiTrig = particle->Phi();
1364  if(phiTrig < 0 ) phiTrig+= TMath::TwoPi();
1365 
1366  Float_t ptLeadHad = 0 ;
1367  Float_t dphiLeadHad= -100 ;
1368  Float_t phiLeadHad = -100 ;
1369  Float_t etaLeadHad = -100 ;
1370  Int_t nTrack = 0;
1371 
1372  for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
1373  {
1374  AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
1375 
1376  fTrackVector.SetXYZ(track->Px(),track->Py(),track->Pz());
1377 
1378  Float_t pt = fTrackVector.Pt();
1379  Float_t phi = fTrackVector.Phi() ;
1380  if(phi < 0 ) phi+= TMath::TwoPi();
1381 
1382  Float_t deltaPhi = phiTrig-phi;
1383  //
1384  // Calculate deltaPhi shift so that for the particles on the opposite side
1385  // it is defined between 90 and 270 degrees
1386  // Shift [-360,-90] to [0, 270]
1387  // and [270,360] to [-90,0]
1388  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
1389  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
1390 
1391  if(pt > ptLeadHad && deltaPhi > TMath::PiOver2()) // in opposite hemisphere
1392  {
1393  ptLeadHad = pt ;
1394  phiLeadHad = phi;
1395  dphiLeadHad= deltaPhi;
1396  etaLeadHad = fTrackVector.Eta();
1397  nTrack++;
1398  }
1399  }// track loop
1400 
1401  if(fFillLeadHadOppositeHisto)
1402  {
1403  if(nTrack == 0)
1404  {
1405  fhPtNoLeadingOppositeHadron ->Fill(ptTrig , GetEventWeight());
1406  fhEtaPhiNoLeadingOppositeHadron->Fill(etaTrig, phiTrig, GetEventWeight());
1407  }
1408  else
1409  {
1410  fhPtLeadingOppositeHadron ->Fill(ptTrig, ptLeadHad, GetEventWeight());
1411  fhPtDiffPhiLeadingOppositeHadron->Fill(ptTrig, dphiLeadHad, GetEventWeight());
1412  fhPtDiffEtaLeadingOppositeHadron->Fill(ptTrig, etaLeadHad-etaTrig, GetEventWeight());
1413  }
1414  }
1415 
1416 
1417  AliDebug(1,Form("pT %2.2f, phi %2.2f, eta %2.2f, nTracks away %d, total tracks %d",
1418  ptLeadHad,phiLeadHad*TMath::RadToDeg(),etaLeadHad,nTrack, GetTrackMultiplicity()));
1419  AliDebug(1,Form("\t pT trig %2.2f, Dphi (trigger-hadron) %2.2f, Deta (trigger-hadron) %2.2f",
1420  ptTrig, dphiLeadHad*TMath::RadToDeg(), etaLeadHad-etaTrig));
1421  AliDebug(1,Form("\t cuts pT: min %2.2f, max %2.2f; DPhi: min %2.2f, max %2.2f",
1422  fMinLeadHadPt,fMaxLeadHadPt,fMinLeadHadPhi*TMath::RadToDeg(),fMaxLeadHadPhi*TMath::RadToDeg()));
1423 
1424 
1425  // reject the trigger if the leading hadron is not in the requested pt or phi window and
1426 
1427  if( nTrack == 0 ) return kFALSE; // No track found in opposite hemisphere
1428 
1429  if( ptLeadHad < fMinLeadHadPt || ptLeadHad > fMaxLeadHadPt ) return kFALSE;
1430 
1431  //printf("Accept leading hadron pT \n");
1432 
1433  if( dphiLeadHad < fMinLeadHadPhi || dphiLeadHad > fMaxLeadHadPhi ) return kFALSE;
1434 
1435  //printf("Accept leading hadron phi \n");
1436 
1437  return kTRUE ;
1438 }
1439 
1440 //____________________________________________________________
1442 //____________________________________________________________
1443 TObjString* AliAnaParticleHadronCorrelation::GetAnalysisCuts()
1444 {
1445  TString parList ; //this will be list of parameters used for this analysis.
1446  const Int_t buffersize = 560;
1447  char onePar[buffersize] ;
1448 
1449  snprintf(onePar,buffersize,"--- AliAnaPaticleHadronCorrelation ---:") ;
1450  parList+=onePar ;
1451  snprintf(onePar,buffersize," %3.2f < Pt associated < %3.2f; ", fMinAssocPt, fMaxAssocPt) ;
1452  parList+=onePar ;
1453  snprintf(onePar,buffersize," %3.2f < Phi trigger particle-Hadron < %3.2f; ", fDeltaPhiMinCut, fDeltaPhiMaxCut) ;
1454  parList+=onePar ;
1455  snprintf(onePar,buffersize," %3.2f < Phi trigger particle-UeHadron < %3.2f; ", fUeDeltaPhiMinCut, fUeDeltaPhiMaxCut) ;
1456  parList+=onePar ;
1457  snprintf(onePar,buffersize,"Isolated Trigger? %d;", fSelectIsolated) ;
1458  parList+=onePar ;
1459  snprintf(onePar,buffersize,"Several UE? %d;", fMakeSeveralUE) ;
1460  parList+=onePar ;
1461  snprintf(onePar,buffersize,"Name of AOD Pi0 Branch %s;", fPi0AODBranchName.Data());
1462  parList+=onePar ;
1463  snprintf(onePar,buffersize,"Do Decay-hadron correlation ? pi0 %d, decay %d;", fPi0Trigger, fDecayTrigger) ;
1464  parList+=onePar ;
1465  snprintf(onePar,buffersize,"Select absolute leading for cluster triggers ? %d or Near Side Leading %d;",
1466  fMakeAbsoluteLeading, fMakeNearSideLeading) ;
1467  parList+=onePar ;
1468  snprintf(onePar,buffersize,"Associated particle pt bins %d: ", fNAssocPtBins) ;
1469  parList+=onePar ;
1470  for (Int_t ibin = 0; ibin<fNAssocPtBins; ibin++) {
1471  snprintf(onePar,buffersize,"bin %d = [%2.1f,%2.1f];", ibin, fAssocPtBinLimit[ibin], fAssocPtBinLimit[ibin+1]) ;
1472  }
1473  parList+=onePar ;
1474 
1475  //Get parameters set in base class.
1476  parList += GetBaseParametersList() ;
1477 
1478  //Get parameters set in FiducialCut class (not available yet)
1479  //parlist += GetFidCut()->GetFidCutParametersList()
1480 
1481  return new TObjString(parList) ;
1482 }
1483 
1484 //________________________________________________________________
1487 //________________________________________________________________
1488 TList * AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
1489 {
1490  TList * outputContainer = new TList() ;
1491  outputContainer->SetName("CorrelationHistos") ;
1492 
1493  Int_t nptbins = GetHistogramRanges()->GetHistoPtBins(); Int_t nphibins = GetHistogramRanges()->GetHistoPhiBins(); Int_t netabins = GetHistogramRanges()->GetHistoEtaBins();
1494  Float_t ptmax = GetHistogramRanges()->GetHistoPtMax(); Float_t phimax = GetHistogramRanges()->GetHistoPhiMax(); Float_t etamax = GetHistogramRanges()->GetHistoEtaMax();
1495  Float_t ptmin = GetHistogramRanges()->GetHistoPtMin(); Float_t phimin = GetHistogramRanges()->GetHistoPhiMin(); Float_t etamin = GetHistogramRanges()->GetHistoEtaMin();
1496 
1497  Int_t ndeltaphibins = GetHistogramRanges()->GetHistoDeltaPhiBins(); Int_t ndeltaetabins = GetHistogramRanges()->GetHistoDeltaEtaBins();
1498  Float_t deltaphimax = GetHistogramRanges()->GetHistoDeltaPhiMax(); Float_t deltaetamax = GetHistogramRanges()->GetHistoDeltaEtaMax();
1499  Float_t deltaphimin = GetHistogramRanges()->GetHistoDeltaPhiMin(); Float_t deltaetamin = GetHistogramRanges()->GetHistoDeltaEtaMin();
1500 
1501  Int_t ntrbins = GetHistogramRanges()->GetHistoTrackMultiplicityBins(); Int_t nclbins = GetHistogramRanges()->GetHistoNClustersBins();
1502  Int_t trmax = GetHistogramRanges()->GetHistoTrackMultiplicityMax(); Int_t clmax = GetHistogramRanges()->GetHistoNClustersMax();
1503  Int_t trmin = GetHistogramRanges()->GetHistoTrackMultiplicityMin(); Int_t clmin = GetHistogramRanges()->GetHistoNClustersMin();
1504 
1505  Int_t nxeztbins = GetHistogramRanges()->GetHistoRatioBins(); Int_t nhbpbins = GetHistogramRanges()->GetHistoHBPBins();
1506  Float_t xeztmax = GetHistogramRanges()->GetHistoRatioMax(); Float_t hbpmax = GetHistogramRanges()->GetHistoHBPMax();
1507  Float_t xeztmin = GetHistogramRanges()->GetHistoRatioMin(); Float_t hbpmin = GetHistogramRanges()->GetHistoHBPMin();
1508 
1509  Int_t nMixBins = GetNCentrBin()*GetNZvertBin()*GetNRPBin();
1510 
1511  Int_t nmassbins = GetHistogramRanges()->GetHistoMassBins();
1512  Float_t massmin = GetHistogramRanges()->GetHistoMassMin();
1513  Float_t massmax = GetHistogramRanges()->GetHistoMassMax();
1514 
1515  TString nameMC[] = {"Photon","Pi0","Pi0Decay","Eta","EtaDecay","OtherDecay","Electron","Hadron","Pi0DecayLostPair","EtaDecayLostPair"};
1516  TString mcPartType[] = {"#gamma", "#pi^{0} (merged #gamma)", "#gamma_{#pi decay}","#eta (merged #gamma)" , "#gamma_{#eta decay}", "#gamma_{other decay}",
1517  "e^{#pm}" , "hadrons?" , "#gamma_{#pi decay} lost companion", "#gamma_{#eta decay} lost companion"} ;
1518  TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
1519 
1520  TString parTitle = Form("#it{R} = %2.2f",GetIsolationCut()->GetConeSize());
1521 
1522  // For vz dependent histograms, if option ON
1523  Int_t nz = 1 ;
1524  if(fCorrelVzBin) nz = GetNZvertBin();
1525  TString sz = "" ;
1526  TString tz = "" ;
1527 
1528  // Init the number of modules, set in the class AliCalorimeterUtils
1529  //
1530  InitCaloParameters(); // See AliCaloTrackCorrBaseClass
1531 
1532  // Fill histograms for neutral clusters in mixing?
1533  Bool_t isoCase = OnlyIsolated() && (GetIsolationCut()->GetParticleTypeInCone() != AliIsolationCut::kOnlyCharged);
1534  Bool_t neutralMix = fFillNeutralEventMixPool || isoCase ;
1535 
1536  if( GetReader()->GetDataType() == AliCaloTrackReader::kESD )
1537  {
1538  fhTrackResolution = new TH2F ("hTrackResolution","Track resolution: #sigma_{#it{p}_{T}} vs #it{p}_{T}, away side, ESDs",
1539  nptbins,ptmin,ptmax,600,0,0.3);
1540  fhTrackResolution->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1541  fhTrackResolution->SetYTitle("#sigma_{#it{p}_{T}} / #it{p}_{T} (GeV/#it{c})");
1542  outputContainer->Add(fhTrackResolution);
1543 
1544  fhTrackResolutionUE = new TH2F ("hTrackResolutionUE","Track resolution: #sigma_{#it{p}_{T}} vs #it{p}_{T}, UE, ESDs",
1545  nptbins,ptmin,ptmax,600,0,0.3);
1546  fhTrackResolutionUE->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1547  fhTrackResolutionUE->SetYTitle("#sigma_{#it{p}_{T}} / #it{p}_{T} (GeV/#it{c})");
1548  outputContainer->Add(fhTrackResolutionUE);
1549  }
1550 
1551  fhPtTriggerInput = new TH1F("hPtTriggerInput","Input trigger #it{p}_{T}", nptbins,ptmin,ptmax);
1552  fhPtTriggerInput->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1553  outputContainer->Add(fhPtTriggerInput);
1554 
1555  if( fM02MaxCut > 0 && fM02MinCut > 0 )
1556  {
1557  fhPtTriggerSSCut = new TH1F("hPtTriggerSSCut","Trigger #it{p}_{T} after #lambda^{2}_{0} cut", nptbins,ptmin,ptmax);
1558  fhPtTriggerSSCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1559  outputContainer->Add(fhPtTriggerSSCut);
1560  }
1561 
1562  if( OnlyIsolated() )
1563  {
1564  fhPtTriggerIsoCut = new TH1F("hPtTriggerIsoCut","Trigger #it{p}_{T} after isolation (and #lambda^{2}_{0} cut)", nptbins,ptmin,ptmax);
1565  fhPtTriggerIsoCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1566  outputContainer->Add(fhPtTriggerIsoCut);
1567  }
1568 
1569  fhPtTriggerFidCut = new TH1F("hPtTriggerFidCut","Trigger #it{p}_{T} after fiducial (isolation and #lambda^{2}_{0}) cut", nptbins,ptmin,ptmax);
1570  fhPtTriggerFidCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1571  outputContainer->Add(fhPtTriggerFidCut);
1572 
1573  fhPtTrigger = new TH1F("hPtTrigger","#it{p}_{T} distribution of trigger particles (after opposite hadron leading cut and rest)", nptbins,ptmin,ptmax);
1574  fhPtTrigger->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1575  outputContainer->Add(fhPtTrigger);
1576 
1577  if(fFillInvMassHisto)
1578  {
1579  fhMassPtTrigger = new TH2F("hMassPtTrigger","2 photons invariant mass vs p_{T}^{trig}",
1580  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
1581  fhMassPtTrigger->SetYTitle("M_{#gamma#gamma} (GeV/#it{c}^{2})");
1582  fhMassPtTrigger->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1583  outputContainer->Add(fhMassPtTrigger);
1584  }
1585 
1586  if(fFillBkgBinsHisto)
1587  {
1588 
1589  fhPtLeadInConeBin = new TH1F*[fNBkgBin] ;
1590 
1591  fhPtSumInConeBin = new TH1F*[fNBkgBin] ;
1592 
1593  fhPtLeadConeBinDecay = new TH1F*[fNBkgBin*fNDecayBits] ;
1594 
1595  fhSumPtConeBinDecay = new TH1F*[fNBkgBin*fNDecayBits] ;
1596 
1597  fhPtLeadConeBinMC = new TH1F*[fNBkgBin*fgkNmcTypes] ;
1598 
1599  fhSumPtConeBinMC = new TH1F*[fNBkgBin*fgkNmcTypes] ;
1600 
1601  for(Int_t ibin = 0; ibin < fNBkgBin; ibin++)
1602  {
1603  fhPtLeadInConeBin[ibin] = new TH1F
1604  (Form("hPtLeadCone_Bin%d",ibin),
1605  Form("cone %2.2f<#it{p}_{T}^{leading}<%2.2f GeV/#it{c}, %s",
1606  fBkgBinLimit[ibin],fBkgBinLimit[ibin+1], parTitle.Data()),nptbins,ptmin,ptmax);
1607  fhPtLeadInConeBin[ibin]->SetYTitle("d #it{N} / d #it{p}_{T}");
1608  fhPtLeadInConeBin[ibin]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1609  outputContainer->Add(fhPtLeadInConeBin[ibin]) ;
1610 
1611  fhPtSumInConeBin[ibin] = new TH1F
1612  (Form("hSumPtCone_Bin%d",ibin),
1613  Form("in cone %2.2f <#Sigma #it{p}_{T}< %2.2f GeV/#it{c}, %s",
1614  fBkgBinLimit[ibin],fBkgBinLimit[ibin+1], parTitle.Data()),nptbins,ptmin,ptmax);
1615  fhPtSumInConeBin[ibin]->SetYTitle("d #it{N} / d #it{p}_{T}");
1616  fhPtSumInConeBin[ibin]->SetXTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
1617  outputContainer->Add(fhPtSumInConeBin[ibin]) ;
1618 
1619  if(fFillTaggedDecayHistograms)
1620  {
1621  for(Int_t idecay = 0; idecay < fNDecayBits; idecay++)
1622  {
1623  Int_t bindecay = ibin+idecay*fNBkgBin;
1624 
1625  fhPtLeadConeBinDecay[bindecay] = new TH1F
1626  (Form("hPtLeadCone_Bin%d_DecayBit%d",ibin,fDecayBits[idecay]),
1627  Form("Decay bit %d, cone %2.2f<#it{p}_{T}^{leading}<%2.2f GeV/#it{c}, %s",
1628  fDecayBits[idecay],fBkgBinLimit[ibin],fBkgBinLimit[ibin+1], parTitle.Data()),nptbins,ptmin,ptmax);
1629  fhPtLeadConeBinDecay[bindecay]->SetYTitle("d #it{N} / d #it{p}_{T}");
1630  fhPtLeadConeBinDecay[bindecay]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1631  outputContainer->Add(fhPtLeadConeBinDecay[bindecay]) ;
1632 
1633  fhSumPtConeBinDecay[bindecay] = new TH1F
1634  (Form("hSumPtCone_Bin%d_DecayBit%d",ibin,fDecayBits[idecay]),
1635  Form("Decay bit %d, in cone %2.2f <#Sigma #it{p}_{T}< %2.2f GeV/#it{c}, %s",
1636  fDecayBits[idecay],fBkgBinLimit[ibin],fBkgBinLimit[ibin+1], parTitle.Data()),nptbins,ptmin,ptmax);
1637  fhSumPtConeBinDecay[bindecay]->SetYTitle("d #it{N} / d #it{p}_{T}");
1638  fhSumPtConeBinDecay[bindecay]->SetXTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
1639  outputContainer->Add(fhSumPtConeBinDecay[bindecay]) ;
1640  }
1641  }
1642 
1643  if(IsDataMC())
1644  {
1645  for(Int_t imc = 0; imc < fgkNmcTypes; imc++)
1646  {
1647  Int_t binmc = ibin+imc*fNBkgBin;
1648  fhPtLeadConeBinMC[binmc] = new TH1F
1649  (Form("hPtLeadCone_Bin%d_MC%s",ibin, nameMC[imc].Data()),
1650  Form("in cone %2.2f<#it{p}_{T}^{leading}<%2.2f GeV/#it{c}, MC %s, %s",
1651  fBkgBinLimit[ibin],fBkgBinLimit[ibin+1], mcPartType[imc].Data(), parTitle.Data()),nptbins,ptmin,ptmax);
1652  fhPtLeadConeBinMC[binmc]->SetYTitle("d #it{N} / d #it{p}_{T}");
1653  fhPtLeadConeBinMC[binmc]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
1654  outputContainer->Add(fhPtLeadConeBinMC[binmc]) ;
1655 
1656  fhSumPtConeBinMC[binmc] = new TH1F
1657  (Form("hSumPtCone_Bin%d_MC%s",ibin,nameMC[imc].Data()),
1658  Form("in cone %2.2f <#Sigma #it{p}_{T}< %2.2f GeV/#it{c}, MC %s, %s",
1659  fBkgBinLimit[ibin],fBkgBinLimit[ibin+1], mcPartType[imc].Data(), parTitle.Data()),nptbins,ptmin,ptmax);
1660  fhSumPtConeBinMC[binmc]->SetYTitle("d #it{N} / d #it{p}_{T}");
1661  fhSumPtConeBinMC[binmc]->SetXTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
1662  outputContainer->Add(fhSumPtConeBinMC[binmc]) ;
1663  } // MC particle loop
1664  }
1665  }
1666  }
1667 
1668  if(IsDataMC())
1669  {
1670  for(Int_t i=0; i < fgkNmcTypes; i++)
1671  {
1672  fhPtTriggerMC[i] = new TH1F(Form("hPtTrigger_MC%s",nameMC[i].Data()),
1673  Form("#it{p}_{T} distribution of trigger particles, trigger origin is %s",nameMC[i].Data()),
1674  nptbins,ptmin,ptmax);
1675  fhPtTriggerMC[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1676  outputContainer->Add(fhPtTriggerMC[i]);
1677 
1678  if(fFillInvMassHisto)
1679  {
1680  fhMCMassPtTrigger[i] = new TH2F(Form("hMassPtTrigger_MC%s",nameMC[i].Data()),
1681  Form("2 photons invariant mass, trigger origin is %s",nameMC[i].Data()),
1682  nptbins,ptmin,ptmax,nmassbins,massmin,massmax);
1683  fhMCMassPtTrigger[i]->SetYTitle("m_{#gamma#gamma}");
1684  fhMCMassPtTrigger[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1685  outputContainer->Add(fhMCMassPtTrigger[i]) ;
1686  }
1687  }
1688  }
1689 
1690  if(fDecayTrigger)
1691  {
1692  for(Int_t ibit = 0; ibit < fNDecayBits; ibit++)
1693  {
1694  fhPtDecayTrigger[ibit] = new TH1F(Form("hPtDecayTrigger_bit%d",fDecayBits[ibit]),
1695  Form("#it{p}_{T} distribution of trigger particles, decay Bit %d",fDecayBits[ibit]),
1696  nptbins,ptmin,ptmax);
1697  fhPtDecayTrigger[ibit]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1698  outputContainer->Add(fhPtDecayTrigger[ibit]);
1699 
1700  if(IsDataMC())
1701  {
1702  for(Int_t i=0; i < fgkNmcTypes; i++)
1703  {
1704  fhPtDecayTriggerMC[ibit][i] = new TH1F(Form("hPtDecayTrigger_bit%d_MC%s",fDecayBits[ibit], nameMC[i].Data()),
1705  Form("#it{p}_{T} distribution of trigger particles, decay Bit %d, trigger origin is %s",fDecayBits[ibit], nameMC[i].Data()),
1706  nptbins,ptmin,ptmax);
1707  fhPtDecayTriggerMC[ibit][i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1708  outputContainer->Add(fhPtDecayTriggerMC[ibit][i]);
1709  }
1710  }
1711  }
1712  }
1713 
1714  if(fCorrelVzBin)
1715  {
1716  fhPtTriggerVzBin = new TH2F("hPtTriggerVzBin","#it{p}_{T} distribution of trigger particles vs vz bin", nptbins,ptmin,ptmax,GetNZvertBin(),0,GetNZvertBin());
1717  fhPtTriggerVzBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1718  fhPtTriggerVzBin->SetYTitle("#it{v}_{#it{z}} bin");
1719  outputContainer->Add(fhPtTriggerVzBin);
1720  }
1721 
1722  fhPtTriggerBin = new TH2F ("hPtTriggerBin","#it{p}_{T} distribution of trigger particles", nptbins,ptmin,ptmax,nMixBins,0,nMixBins);
1723  fhPtTriggerBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1724  fhPtTriggerBin->SetYTitle("Bin");
1725  outputContainer->Add(fhPtTriggerBin);
1726 
1727  fhPhiTrigger = new TH2F ("hPhiTrigger","#varphi distribution of trigger Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
1728  fhPhiTrigger->SetYTitle("#varphi (rad)");
1729  outputContainer->Add(fhPhiTrigger);
1730 
1731  fhEtaTrigger = new TH2F ("hEtaTrigger","#eta distribution of trigger",nptbins,ptmin,ptmax, netabins,etamin,etamax);
1732  fhEtaTrigger->SetYTitle("#eta ");
1733  outputContainer->Add(fhEtaTrigger);
1734 
1735  if(IsHighMultiplicityAnalysisOn())
1736  {
1737  fhPtTriggerCentrality = new TH2F("hPtTriggerCentrality","Trigger particle #it{p}_{T} vs centrality",nptbins,ptmin,ptmax,100,0.,100) ;
1738  fhPtTriggerCentrality->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1739  fhPtTriggerCentrality->SetYTitle("Centrality (%)");
1740  outputContainer->Add(fhPtTriggerCentrality) ;
1741 
1742  fhPtTriggerEventPlane = new TH2F("hPtTriggerEventPlane","Trigger particle #it{p}_{T} vs event plane angle",nptbins,ptmin,ptmax, 100,0.,TMath::Pi()) ;
1743  fhPtTriggerEventPlane->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1744  fhPtTriggerEventPlane->SetXTitle("EP angle (rad)");
1745  outputContainer->Add(fhPtTriggerEventPlane) ;
1746 
1747  fhTriggerEventPlaneCentrality = new TH2F("hTriggerEventPlaneCentrality","Trigger particle centrality vs event plane angle",100,0.,100,100,0.,TMath::Pi()) ;
1748  fhTriggerEventPlaneCentrality->SetXTitle("Centrality (%)");
1749  fhTriggerEventPlaneCentrality->SetYTitle("EP angle (rad)");
1750  outputContainer->Add(fhTriggerEventPlaneCentrality) ;
1751  }
1752 
1753  // Leading hadron in oposite side
1754  if(fFillLeadHadOppositeHisto)
1755  {
1756  fhPtLeadingOppositeHadron = new TH2F("hPtTriggerPtLeadingOppositeHadron","Leading hadron opposite to trigger vs trigger #it{p}_{T}",
1757  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1758  fhPtLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1759  fhPtLeadingOppositeHadron->SetYTitle("#it{p}_{T}^{lead hadron} (GeV/#it{c})");
1760  outputContainer->Add(fhPtLeadingOppositeHadron);
1761 
1762  fhPtNoLeadingOppositeHadron = new TH1F("hPtTriggerNoLeadingOppositeHadron","No Leading hadron opposite to trigger #it{p}_{T}",
1763  nptbins,ptmin,ptmax);
1764  fhPtNoLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1765  outputContainer->Add(fhPtNoLeadingOppositeHadron);
1766 
1767  fhEtaPhiNoLeadingOppositeHadron = new TH2F("hEtaPhiTriggerNoLeadingOppositeHadron","No Leading hadron opposite to trigger #eta:#varphi",
1768  netabins,etamin,etamax,nphibins,phimin,phimax);
1769  fhEtaPhiNoLeadingOppositeHadron->SetXTitle("#eta");
1770  fhEtaPhiNoLeadingOppositeHadron->SetYTitle("#varphi");
1771  outputContainer->Add(fhEtaPhiNoLeadingOppositeHadron);
1772 
1773 
1774  fhPtDiffPhiLeadingOppositeHadron = new TH2F("hPtTriggerDiffPhiTriggerLeadingOppositeHadron","#varphi_{trigger}-#varphi_{leading opposite hadron} vs #it{p}_{T}^{trig}",
1775  nptbins,ptmin,ptmax,ndeltaphibins,deltaphimin,deltaphimax);
1776  fhPtDiffPhiLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1777  fhPtDiffPhiLeadingOppositeHadron->SetYTitle("#varphi_{trigger}-#varphi_{leading opposite hadron} (rad)");
1778  outputContainer->Add(fhPtDiffPhiLeadingOppositeHadron);
1779 
1780  fhPtDiffEtaLeadingOppositeHadron = new TH2F("hPtTriggerDiffEtaTriggerPhiLeadingOppositeHadron","#eta_{trigger}-#eta_{leading opposite hadron} vs #it{p}_{T}^{trig}",
1781  nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1782  fhPtDiffEtaLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1783  fhPtDiffEtaLeadingOppositeHadron->SetYTitle("#eta_{trigger}-#eta_{leading opposite hadron}");
1784  outputContainer->Add(fhPtDiffEtaLeadingOppositeHadron);
1785  }
1786 
1787  // Correlation with charged hadrons
1788 
1789  fhDeltaPhiDeltaEtaCharged = new TH2F
1790  ("hDeltaPhiDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #varphi_{trigger} - #varphi_{h^{#pm}}",
1791  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1792  fhDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #varphi (rad)");
1793  fhDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
1794 
1795  fhDeltaPhiDeltaEtaChargedPtA3GeV = new TH2F
1796  ("hDeltaPhiDeltaEtaChargedPtA3GeV","#eta_{trigger} - #eta_{h^{#pm}} vs #varphi_{trigger} - #varphi_{h^{#pm}, #it{p}_{TA}>3 GeV/#it{c}}",
1797  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1798  fhDeltaPhiDeltaEtaChargedPtA3GeV->SetXTitle("#Delta #varphi (rad)");
1799  fhDeltaPhiDeltaEtaChargedPtA3GeV->SetYTitle("#Delta #eta");
1800 
1801  fhPhiCharged = new TH2F
1802  ("hPhiCharged","#varphi_{h^{#pm}} vs #it{p}_{T #pm}",
1803  nptbins,ptmin,ptmax,180,0,TMath::TwoPi());
1804  fhPhiCharged->SetYTitle("#varphi_{h^{#pm}} (rad)");
1805  fhPhiCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
1806 
1807  fhEtaCharged = new TH2F
1808  ("hEtaCharged","#eta_{h^{#pm}} vs #it{p}_{T #pm}",
1809  nptbins,ptmin,ptmax,100,-1.,1.);
1810  fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
1811  fhEtaCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
1812 
1813  fhDeltaPhiCharged = new TH2F
1814  ("hDeltaPhiCharged","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}",
1815  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1816  fhDeltaPhiCharged->SetYTitle("#Delta #varphi (rad)");
1817  fhDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1818 
1819  fhDeltaPhiChargedPtA3GeV = new TH2F
1820  ("hDeltaPhiChargedPtA3GeV","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}",
1821  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1822  fhDeltaPhiChargedPtA3GeV->SetYTitle("#Delta #varphi (rad)");
1823  fhDeltaPhiChargedPtA3GeV->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1824 
1825 
1826  fhDeltaPhiChargedPt = new TH2F
1827  ("hDeltaPhiChargedPt","#varphi_{trigger} - #varphi_{#h^{#pm}} vs #it{p}_{T h^{#pm}}",
1828  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1829  fhDeltaPhiChargedPt->SetYTitle("#Delta #varphi (rad)");
1830  fhDeltaPhiChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1831 
1832  fhDeltaEtaCharged = new TH2F
1833  ("hDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}",
1834  nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1835  fhDeltaEtaCharged->SetYTitle("#Delta #eta");
1836  fhDeltaEtaCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1837 
1838  fhDeltaEtaChargedPtA3GeV = new TH2F
1839  ("hDeltaEtaChargedPtA3GeV","#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}",
1840  nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1841  fhDeltaEtaChargedPtA3GeV->SetYTitle("#Delta #eta");
1842  fhDeltaEtaChargedPtA3GeV->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1843 
1844  fhXECharged =
1845  new TH2F("hXECharged","#it{x}_{#it{E}} for charged tracks",
1846  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1847  fhXECharged->SetYTitle("#it{x}_{#it{E}}");
1848  fhXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1849 
1850  fhXECharged_Cone2 =
1851  new TH2F("hXECharged_Cone2","#it{x}_{#it{E}} for charged tracks in cone 2 (5#pi/6-7#pi/6)",
1852  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1853  fhXECharged_Cone2->SetYTitle("#it{x}_{#it{E}}");
1854  fhXECharged_Cone2->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1855 
1856  fhXEPosCharged =
1857  new TH2F("hXEPositiveCharged","#it{x}_{#it{E}} for positive charged tracks",
1858  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1859  fhXEPosCharged->SetYTitle("#it{x}_{#it{E}}");
1860  fhXEPosCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1861 
1862  fhXENegCharged =
1863  new TH2F("hXENegativeCharged","#it{x}_{#it{E}} for negative charged tracks",
1864  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1865  fhXENegCharged->SetYTitle("#it{x}_{#it{E}}");
1866  fhXENegCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1867 
1868  fhPtHbpXECharged =
1869  new TH2F("hHbpXECharged","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
1870  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1871  fhPtHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1872  fhPtHbpXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1873 
1874  fhPtHbpXECharged_Cone2 =
1875  new TH2F("hHbpXECharged_Cone2","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons in cone 2 (5#pi/6-7#pi/6)",
1876  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1877  fhPtHbpXECharged_Cone2->SetYTitle("ln(1/#it{x}_{#it{E}})");
1878  fhPtHbpXECharged_Cone2->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1879 
1880  fhZTCharged =
1881  new TH2F("hZTCharged","#it{z}_{T} for charged tracks",
1882  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1883  fhZTCharged->SetYTitle("#it{z}_{T}");
1884  fhZTCharged->SetXTitle("#it{p}_{T trigger}");
1885 
1886  fhZTPosCharged =
1887  new TH2F("hZTPositiveCharged","#it{z}_{T} for positive charged tracks",
1888  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1889  fhZTPosCharged->SetYTitle("#it{z}_{T}");
1890  fhZTPosCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1891 
1892  fhZTNegCharged =
1893  new TH2F("hZTNegativeCharged","#it{z}_{T} for negative charged tracks",
1894  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1895  fhZTNegCharged->SetYTitle("#it{z}_{T}");
1896  fhZTNegCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1897 
1898  fhPtHbpZTCharged =
1899  new TH2F("hHbpZTCharged","#xi = ln(1/#it{z}_{T}) with charged hadrons",
1900  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1901  fhPtHbpZTCharged->SetYTitle("ln(1/#it{z}_{T})");
1902  fhPtHbpZTCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1903 
1904  fhPtTrigPout =
1905  new TH2F("hPtTrigPout","Pout with triggers",
1906  nptbins,ptmin,ptmax,nptbins,-1.*ptmax/2.,ptmax/2.);
1907  fhPtTrigPout->SetYTitle("#it{p}_{out} (GeV/#it{c})");
1908  fhPtTrigPout->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1909 
1910  fhPtTrigCharged =
1911  new TH2F("hPtTrigCharged","trigger and charged tracks pt distribution",
1912  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1913  fhPtTrigCharged->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1914  fhPtTrigCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1915 
1916  outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
1917  outputContainer->Add(fhDeltaPhiDeltaEtaChargedPtA3GeV);
1918  outputContainer->Add(fhPhiCharged) ;
1919  outputContainer->Add(fhEtaCharged) ;
1920  outputContainer->Add(fhDeltaPhiCharged) ;
1921  outputContainer->Add(fhDeltaPhiChargedPtA3GeV) ;
1922  outputContainer->Add(fhDeltaEtaCharged) ;
1923  outputContainer->Add(fhDeltaEtaChargedPtA3GeV) ;
1924  outputContainer->Add(fhDeltaPhiChargedPt) ;
1925 
1926  outputContainer->Add(fhXECharged) ;
1927  outputContainer->Add(fhXECharged_Cone2) ;
1928 
1929  if(IsDataMC())
1930  {
1931  for(Int_t i=0; i < fgkNmcTypes; i++)
1932  {
1933 
1934  fhDeltaPhiChargedMC[i] = new TH2F(Form("hDeltaPhiCharged_MC%s",nameMC[i].Data()),
1935  Form("#Delta #varphi for charged tracks, trigger origin is %s",nameMC[i].Data()),
1936  nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
1937  fhDeltaPhiChargedMC[i]->SetYTitle("#it{x}_{#it{E}}");
1938  fhDeltaPhiChargedMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1939  outputContainer->Add(fhDeltaPhiChargedMC[i]) ;
1940 
1941  fhXEChargedMC[i] = new TH2F(Form("hXECharged_MC%s",nameMC[i].Data()),
1942  Form("#it{x}_{#it{E}} for charged tracks, trigger origin is %s",nameMC[i].Data()),
1943  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1944  fhXEChargedMC[i]->SetYTitle("#it{x}_{#it{E}}");
1945  fhXEChargedMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1946  outputContainer->Add(fhXEChargedMC[i]) ;
1947 
1948  fhXEUeChargedRightMC[i] = new TH2F(Form("hXEUeChargedRight_MC%s",nameMC[i].Data()),
1949  Form("#it{x}_{#it{E}} for charged tracks in right UE cone, trigger origin is %s",nameMC[i].Data()),
1950  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1951  fhXEUeChargedRightMC[i]->SetYTitle("#it{x}_{#it{E}}");
1952  fhXEUeChargedRightMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1953  outputContainer->Add(fhXEUeChargedRightMC[i]) ;
1954 
1955  fhXEUeChargedLeftMC[i] = new TH2F(Form("hXEUeChargedLeft_MC%s",nameMC[i].Data()),
1956  Form("#it{x}_{#it{E}} for charged tracks in left UE cone, trigger origin is %s",nameMC[i].Data()),
1957  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1958  fhXEUeChargedLeftMC[i]->SetYTitle("#it{x}_{#it{E}}");
1959  fhXEUeChargedLeftMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1960  outputContainer->Add(fhXEUeChargedLeftMC[i]) ;
1961  }
1962  }
1963 
1964  outputContainer->Add(fhXEPosCharged) ;
1965  outputContainer->Add(fhXENegCharged) ;
1966  outputContainer->Add(fhPtHbpXECharged) ;
1967  outputContainer->Add(fhPtHbpXECharged_Cone2) ;
1968 
1969  outputContainer->Add(fhZTCharged) ;
1970  outputContainer->Add(fhZTPosCharged) ;
1971  outputContainer->Add(fhZTNegCharged) ;
1972  outputContainer->Add(fhPtHbpZTCharged) ;
1973 
1974  outputContainer->Add(fhPtTrigPout) ;
1975  outputContainer->Add(fhPtTrigCharged) ;
1976 
1977  TString right = "";
1978  if(fMakeSeveralUE) right = "Right";
1979 
1980  fhUePart = new TH1F("hUePart","UE particles distribution vs pt trig",
1981  nptbins,ptmin,ptmax);
1982  fhUePart->SetYTitle("dNch");
1983  fhUePart->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1984 
1985  fhDeltaPhiUeChargedPt = new TH2F
1986  (Form("hDeltaPhiUe%sChargedPt",right.Data()),"#varphi_{trigger} - #varphi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}}",
1987  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1988  fhDeltaPhiUeChargedPt->SetYTitle("#Delta #varphi (rad)");
1989  fhDeltaPhiUeChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1990 
1991  fhXEUeCharged =
1992  new TH2F(Form("hXEUeCharged%s",right.Data()),"#it{x}_{#it{E}} for Underlying Event",
1993  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1994  fhXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
1995  fhXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1996 
1997  fhXEUeChargedSmallCone =
1998  new TH2F(Form("hXEUeChargedSmallCone%s",right.Data()),"#it{x}_{#it{E}} for Underlying Event in a cone [80,100] deg",
1999  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2000  fhXEUeChargedSmallCone->SetYTitle("#it{x}_{#it{E}}");
2001  fhXEUeChargedSmallCone->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2002 
2003  fhXEUeChargedMediumCone =
2004  new TH2F(Form("hXEUeChargedMediumCone%s",right.Data()),"#it{x}_{#it{E}} for Underlying Event in a cone [70,110] deg",
2005  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2006  fhXEUeChargedMediumCone->SetYTitle("#it{x}_{#it{E}}");
2007  fhXEUeChargedMediumCone->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2008 
2009  fhXEUeChargedLargeCone =
2010  new TH2F(Form("hXEUeChargedLargeCone%s",right.Data()),"#it{x}_{#it{E}} for Underlying Event in a cone [60,120] deg",
2011  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2012  fhXEUeChargedLargeCone->SetYTitle("#it{x}_{#it{E}}");
2013  fhXEUeChargedLargeCone->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2014 
2015  fhPtHbpXEUeCharged =
2016  new TH2F(Form("hHbpXEUeCharged%s",right.Data()),"#xi = ln(1/#it{x}_{#it{E}}) for Underlying Event",
2017  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2018  fhPtHbpXEUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2019  fhPtHbpXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2020 
2021  fhZTUeCharged =
2022  new TH2F(Form("hZTUeCharged%s",right.Data()),"#it{z}_{T} for Underlying Event",
2023  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2024  fhZTUeCharged->SetYTitle("#it{z}_{T}");
2025  fhZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2026 
2027  fhPtHbpZTUeCharged =
2028  new TH2F(Form("hHbpZTUeCharged%s",right.Data()),"#xi = ln(1/#it{z}_{T}) for Underlying Event",
2029  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2030  fhPtHbpZTUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2031  fhPtHbpZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2032 
2033  outputContainer->Add(fhUePart);
2034  outputContainer->Add(fhDeltaPhiUeChargedPt) ;
2035  outputContainer->Add(fhXEUeCharged) ;
2036  outputContainer->Add(fhXEUeChargedSmallCone) ;
2037  outputContainer->Add(fhXEUeChargedMediumCone) ;
2038  outputContainer->Add(fhXEUeChargedLargeCone) ;
2039  outputContainer->Add(fhPtHbpXEUeCharged) ;
2040  outputContainer->Add(fhZTUeCharged) ;
2041  outputContainer->Add(fhPtHbpZTUeCharged) ;
2042 
2043  if(fMakeSeveralUE)
2044  {
2045  fhDeltaPhiUeLeftCharged = new TH2F
2046  ("hDeltaPhiUeLeftChargedPt","#varphi_{trigger} - #varphi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left side range of trigger particles",
2047  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2048  fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #varphi (rad)");
2049  fhDeltaPhiUeLeftCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2050  outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
2051 
2052  fhDeltaPhiUeLeftUpCharged = new TH2F
2053  ("hDeltaPhiUeLeftUpChargedPt","#varphi_{trigger} - #varphi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left Up side range of trigger particles",
2054  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2055  fhDeltaPhiUeLeftUpCharged->SetYTitle("#Delta #varphi (rad)");
2056  fhDeltaPhiUeLeftUpCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2057  outputContainer->Add(fhDeltaPhiUeLeftUpCharged) ;
2058 
2059  fhDeltaPhiUeRightUpCharged = new TH2F
2060  ("hDeltaPhiUeRightUpChargedPt","#varphi_{trigger} - #varphi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right Up side range of trigger particles",
2061  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2062  fhDeltaPhiUeRightUpCharged->SetYTitle("#Delta #varphi (rad)");
2063  fhDeltaPhiUeRightUpCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2064  outputContainer->Add(fhDeltaPhiUeRightUpCharged) ;
2065 
2066  fhDeltaPhiUeLeftDownCharged = new TH2F
2067  ("hDeltaPhiUeLeftDownChargedPt","#varphi_{trigger} - #varphi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left Down side range of trigger particles",
2068  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2069  fhDeltaPhiUeLeftDownCharged->SetYTitle("#Delta #varphi (rad)");
2070  fhDeltaPhiUeLeftDownCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2071  outputContainer->Add(fhDeltaPhiUeLeftDownCharged) ;
2072 
2073  fhDeltaPhiUeRightDownCharged = new TH2F
2074  ("hDeltaPhiUeRightDownChargedPt","#varphi_{trigger} - #varphi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right Down side range of trigger particles",
2075  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2076  fhDeltaPhiUeRightDownCharged->SetYTitle("#Delta #varphi (rad)");
2077  fhDeltaPhiUeRightDownCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2078  outputContainer->Add(fhDeltaPhiUeRightDownCharged) ;
2079 
2080  fhXEUeLeftCharged =
2081  new TH2F("hXEUeChargedLeft","#it{x}_{#it{E}} with UE left side of trigger",
2082  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2083  fhXEUeLeftCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
2084  fhXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2085  outputContainer->Add(fhXEUeLeftCharged) ;
2086 
2087  fhXEUeLeftUpCharged =
2088  new TH2F("hXEUeChargedLeftUp","#it{x}_{#it{E}} with UE left Up side of trigger",
2089  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2090  fhXEUeLeftUpCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
2091  fhXEUeLeftUpCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2092  outputContainer->Add(fhXEUeLeftUpCharged) ;
2093 
2094  fhXEUeRightUpCharged =
2095  new TH2F("hXEUeChargedRightUp","#it{x}_{#it{E} h^{#pm}} with UE right Up side of trigger",
2096  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2097  fhXEUeRightUpCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
2098  fhXEUeRightUpCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2099  outputContainer->Add(fhXEUeRightUpCharged) ;
2100 
2101  fhXEUeLeftDownCharged =
2102  new TH2F("hXEUeChargedLeftDown","#it{x}_{#it{E}} with UE left Down side of trigger",
2103  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2104  fhXEUeLeftDownCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
2105  fhXEUeLeftDownCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2106  outputContainer->Add(fhXEUeLeftDownCharged) ;
2107 
2108  fhXEUeRightDownCharged =
2109  new TH2F("hXEUeChargedRightDown","#it{x}_{#it{E} h^{#pm}} with UE right Down side of trigger",
2110  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2111  fhXEUeRightDownCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
2112  fhXEUeRightDownCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2113  outputContainer->Add(fhXEUeRightDownCharged) ;
2114 
2115  fhPtHbpXEUeLeftCharged =
2116  new TH2F("hHbpXEUeChargedLeft","#xi = ln(1/#it{x}_{#it{E}}) with charged UE left side of trigger",
2117  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2118  fhPtHbpXEUeLeftCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2119  fhPtHbpXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2120  outputContainer->Add(fhPtHbpXEUeLeftCharged) ;
2121 
2122  fhZTUeLeftCharged =
2123  new TH2F("hZTUeChargedLeft","#it{z}_{trigger h^{#pm}} = #it{p}_{T Ueh^{#pm}} / #it{p}_{T trigger} with UE left side of trigger",
2124  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2125  fhZTUeLeftCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
2126  fhZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2127  outputContainer->Add(fhZTUeLeftCharged) ;
2128 
2129  fhPtHbpZTUeLeftCharged =
2130  new TH2F("hHbpZTUeChargedLeft","#xi = ln(1/#it{z}_{T}) with charged UE left side of trigger",
2131  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2132  fhPtHbpZTUeLeftCharged->SetYTitle("ln(1/#it{z}_{T})");
2133  fhPtHbpZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2134  outputContainer->Add(fhPtHbpZTUeLeftCharged) ;
2135  }
2136 
2137  if(IsPileUpAnalysisOn())
2138  {
2139  fhDeltaPhiChargedOtherBC = new TH2F
2140  ("hDeltaPhiChargedOtherBC","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, track BC!=0",
2141  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2142  fhDeltaPhiChargedOtherBC->SetYTitle("#Delta #varphi (rad)");
2143  fhDeltaPhiChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2144 
2145  fhDeltaPhiChargedPtA3GeVOtherBC = new TH2F
2146  ("hDeltaPhiChargedPtA3GeVOtherBC","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC!=0",
2147  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2148  fhDeltaPhiChargedPtA3GeVOtherBC->SetYTitle("#Delta #varphi (rad)");
2149  fhDeltaPhiChargedPtA3GeVOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2150 
2151  fhPtTrigChargedOtherBC =
2152  new TH2F("hPtTrigChargedOtherBC","trigger and charged tracks pt distribution, track BC!=0",
2153  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
2154  fhPtTrigChargedOtherBC->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2155  fhPtTrigChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2156 
2157  fhXEChargedOtherBC =
2158  new TH2F("hXEChargedOtherBC","#it{x}_{#it{E}} for charged tracks, track BC!=0",
2159  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2160  fhXEChargedOtherBC->SetYTitle("#it{x}_{#it{E}}");
2161  fhXEChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2162 
2163  fhXEUeChargedOtherBC =
2164  new TH2F("hXEUeChargedOtherBC","#it{x}_{#it{E}} for Underlying Event, track BC!=0",
2165  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2166  fhXEUeChargedOtherBC->SetYTitle("#it{x}_{#it{E}}");
2167  fhXEUeChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2168 
2169  fhZTChargedOtherBC =
2170  new TH2F("hZTChargedOtherBC","#it{z}_{T} for charged tracks, track BC!=0",
2171  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2172  fhZTChargedOtherBC->SetYTitle("#it{z}_{T}");
2173  fhZTChargedOtherBC->SetXTitle("#it{p}_{T trigger}");
2174 
2175  fhZTUeChargedOtherBC =
2176  new TH2F("hZTUeChargedOtherBC","#it{z}_{T} for Underlying Event, track BC!=0",
2177  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2178  fhZTUeChargedOtherBC->SetYTitle("#it{z}_{T}");
2179  fhZTUeChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2180 
2181  outputContainer->Add(fhDeltaPhiChargedOtherBC) ;
2182  outputContainer->Add(fhDeltaPhiChargedPtA3GeVOtherBC) ;
2183  outputContainer->Add(fhXEChargedOtherBC) ;
2184  outputContainer->Add(fhXEUeChargedOtherBC) ;
2185  outputContainer->Add(fhZTChargedOtherBC) ;
2186  outputContainer->Add(fhZTUeChargedOtherBC) ;
2187  outputContainer->Add(fhPtTrigChargedOtherBC) ;
2188 
2189  fhDeltaPhiChargedBC0 = new TH2F
2190  ("hDeltaPhiChargedBC0","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, track BC==0",
2191  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2192  fhDeltaPhiChargedBC0->SetYTitle("#Delta #varphi (rad)");
2193  fhDeltaPhiChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2194 
2195  fhDeltaPhiChargedPtA3GeVBC0 = new TH2F
2196  ("hDeltaPhiChargedPtA3GeVBC0","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC==0",
2197  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2198  fhDeltaPhiChargedPtA3GeVBC0->SetYTitle("#Delta #varphi (rad)");
2199  fhDeltaPhiChargedPtA3GeVBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2200 
2201  fhPtTrigChargedBC0 =
2202  new TH2F("hPtTrigChargedBC0","trigger and charged tracks pt distribution, track BC==0",
2203  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
2204  fhPtTrigChargedBC0->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2205  fhPtTrigChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2206 
2207  fhXEChargedBC0 =
2208  new TH2F("hXEChargedBC0","#it{x}_{#it{E}} for charged tracks, track BC==0",
2209  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2210  fhXEChargedBC0->SetYTitle("#it{x}_{#it{E}}");
2211  fhXEChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2212 
2213  fhXEUeChargedBC0 =
2214  new TH2F("hXEUeChargedBC0","#it{x}_{#it{E}} for Underlying Event, track BC==0",
2215  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2216  fhXEUeChargedBC0->SetYTitle("#it{x}_{#it{E}}");
2217  fhXEUeChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2218 
2219  fhZTChargedBC0 =
2220  new TH2F("hZTChargedBC0","#it{z}_{T} for charged tracks, track BC==0",
2221  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2222  fhZTChargedBC0->SetYTitle("#it{z}_{T}");
2223  fhZTChargedBC0->SetXTitle("#it{p}_{T trigger}");
2224 
2225  fhZTUeChargedBC0 =
2226  new TH2F("hZTUeChargedBC0","#it{z}_{T} for Underlying Event, track BC==0",
2227  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2228  fhZTUeChargedBC0->SetYTitle("#it{z}_{T}");
2229  fhZTUeChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2230 
2231  outputContainer->Add(fhDeltaPhiChargedBC0) ;
2232  outputContainer->Add(fhDeltaPhiChargedPtA3GeVBC0) ;
2233  outputContainer->Add(fhXEChargedBC0) ;
2234  outputContainer->Add(fhXEUeChargedBC0) ;
2235  outputContainer->Add(fhZTChargedBC0) ;
2236  outputContainer->Add(fhZTUeChargedBC0) ;
2237  outputContainer->Add(fhPtTrigChargedBC0) ;
2238 
2239  fhPtTriggerVtxBC0 = new TH1F("hPtTriggerVtxBC0","#it{p}_{T} distribution of trigger particles", nptbins,ptmin,ptmax);
2240  fhPtTriggerVtxBC0->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2241 
2242  fhDeltaPhiChargedVtxBC0 = new TH2F
2243  ("hDeltaPhiChargedVtxBC0","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, track BC==0",
2244  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2245  fhDeltaPhiChargedVtxBC0->SetYTitle("#Delta #varphi (rad)");
2246  fhDeltaPhiChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2247 
2248  fhDeltaPhiChargedPtA3GeVVtxBC0 = new TH2F
2249  ("hDeltaPhiChargedPtA3GeVVtxBC0","#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC==0",
2250  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2251  fhDeltaPhiChargedPtA3GeVVtxBC0->SetYTitle("#Delta #varphi (rad)");
2252  fhDeltaPhiChargedPtA3GeVVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2253 
2254  fhPtTrigChargedVtxBC0 =
2255  new TH2F("hPtTrigChargedVtxBC0","trigger and charged tracks pt distribution, track BC==0",
2256  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
2257  fhPtTrigChargedVtxBC0->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2258  fhPtTrigChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2259 
2260  fhXEChargedVtxBC0 =
2261  new TH2F("hXEChargedVtxBC0","#it{x}_{#it{E}} for charged tracks, track BC==0",
2262  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2263  fhXEChargedVtxBC0->SetYTitle("#it{x}_{#it{E}}");
2264  fhXEChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2265 
2266  fhXEUeChargedVtxBC0 =
2267  new TH2F("hXEUeChargedVtxBC0","#it{x}_{#it{E}} for Underlying Event, track BC==0",
2268  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2269  fhXEUeChargedVtxBC0->SetYTitle("#it{x}_{#it{E}}");
2270  fhXEUeChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2271 
2272  fhZTChargedVtxBC0 =
2273  new TH2F("hZTChargedVtxBC0","#it{z}_{T} for charged tracks, track BC==0",
2274  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2275  fhZTChargedVtxBC0->SetYTitle("#it{z}_{T}");
2276  fhZTChargedVtxBC0->SetXTitle("#it{p}_{T trigger}");
2277 
2278  fhZTUeChargedVtxBC0 =
2279  new TH2F("hZTUeChargedVtxBC0","#it{z}_{T} for Underlying Event, track BC==0",
2280  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2281  fhZTUeChargedVtxBC0->SetYTitle("#it{z}_{T}");
2282  fhZTUeChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2283 
2284  outputContainer->Add(fhPtTriggerVtxBC0);
2285  outputContainer->Add(fhDeltaPhiChargedVtxBC0) ;
2286  outputContainer->Add(fhDeltaPhiChargedPtA3GeVVtxBC0) ;
2287  outputContainer->Add(fhXEChargedVtxBC0) ;
2288  outputContainer->Add(fhXEUeChargedVtxBC0) ;
2289  outputContainer->Add(fhZTChargedVtxBC0) ;
2290  outputContainer->Add(fhZTUeChargedVtxBC0) ;
2291  outputContainer->Add(fhPtTrigChargedVtxBC0) ;
2292 
2293  for(Int_t i = 0 ; i < 7 ; i++)
2294  {
2295  fhPtTriggerPileUp[i] = new TH1F(Form("hPtTriggerPileUp%s",pileUpName[i].Data()),
2296  Form("#it{p}_{T} distribution of trigger particles, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
2297  fhPtTriggerPileUp[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2298  outputContainer->Add(fhPtTriggerPileUp[i]);
2299 
2300  fhDeltaPhiChargedPileUp[i] = new TH2F(Form("hDeltaPhiChargedPileUp%s",pileUpName[i].Data()),
2301  Form("#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, %s Pile-Up event",pileUpName[i].Data()),
2302  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2303  fhDeltaPhiChargedPileUp[i]->SetYTitle("#Delta #varphi (rad)");
2304  fhDeltaPhiChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2305  outputContainer->Add(fhDeltaPhiChargedPileUp[i]) ;
2306 
2307  fhDeltaPhiChargedPtA3GeVPileUp[i] = new TH2F(Form("hDeltaPhiChargedPtA3GeVPileUp%s",pileUpName[i].Data()),
2308  Form("#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, %s Pile-Up event",pileUpName[i].Data()),
2309  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2310  fhDeltaPhiChargedPtA3GeVPileUp[i]->SetYTitle("#Delta #varphi (rad)");
2311  fhDeltaPhiChargedPtA3GeVPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2312  outputContainer->Add(fhDeltaPhiChargedPtA3GeVPileUp[i]) ;
2313 
2314  fhDeltaEtaChargedPileUp[i] = new TH2F(Form("hDeltaEtaChargedPileUp%s",pileUpName[i].Data()),
2315  Form("#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, %s Pile-Up event",pileUpName[i].Data()),
2316  nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
2317  fhDeltaEtaChargedPileUp[i]->SetYTitle("#Delta #eta");
2318  fhDeltaEtaChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2319  outputContainer->Add(fhDeltaEtaChargedPileUp[i]) ;
2320 
2321  fhDeltaEtaChargedPtA3GeVPileUp[i] = new TH2F(Form("hDeltaEtaChargedPtA3GeVPileUp%s",pileUpName[i].Data()),
2322  Form("#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, %s Pile-Up event",pileUpName[i].Data()),
2323  nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
2324  fhDeltaEtaChargedPtA3GeVPileUp[i]->SetYTitle("#Delta #eta");
2325  fhDeltaEtaChargedPtA3GeVPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2326  outputContainer->Add(fhDeltaEtaChargedPtA3GeVPileUp[i]) ;
2327 
2328  fhXEChargedPileUp[i] = new TH2F(Form("hXEChargedPileUp%s",pileUpName[i].Data()),
2329  Form("#it{x}_{#it{E}} for charged tracks, %s Pile-Up event",pileUpName[i].Data()),
2330  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2331  fhXEChargedPileUp[i]->SetYTitle("#it{x}_{#it{E}}");
2332  fhXEChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2333  outputContainer->Add(fhXEChargedPileUp[i]) ;
2334 
2335  fhXEUeChargedPileUp[i] = new TH2F(Form("hXEUeChargedPileUp%s",pileUpName[i].Data()),
2336  Form("#it{x}_{#it{E}} for Underlying Event, %s Pile-Up event",pileUpName[i].Data()),
2337  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2338  fhXEUeChargedPileUp[i]->SetYTitle("#it{x}_{#it{E}}");
2339  fhXEUeChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2340  outputContainer->Add(fhXEUeChargedPileUp[i]) ;
2341 
2342  fhZTChargedPileUp[i] = new TH2F(Form("hZTChargedPileUp%s",pileUpName[i].Data()),
2343  Form("#it{z}_{T} for charged tracks, %s Pile-Up event",pileUpName[i].Data()),
2344  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2345  fhZTChargedPileUp[i]->SetYTitle("#it{z}_{T}");
2346  fhZTChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2347  outputContainer->Add(fhZTChargedPileUp[i]) ;
2348 
2349  fhZTUeChargedPileUp[i] = new TH2F(Form("hZTUeChargedPileUp%s",pileUpName[i].Data()),
2350  Form("#it{z}_{T} for Underlying Event, %s Pile-Up event",pileUpName[i].Data()),
2351  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2352  fhZTUeChargedPileUp[i]->SetYTitle("#it{z}_{T}");
2353  fhZTUeChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2354  outputContainer->Add(fhZTUeChargedPileUp[i]) ;
2355 
2356  fhPtTrigChargedPileUp[i] = new TH2F(Form("hPtTrigChargedPileUp%s",pileUpName[i].Data()),
2357  Form("trigger and charged tracks pt distribution, %s Pile-Up event",pileUpName[i].Data()),
2358  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
2359  fhPtTrigChargedPileUp[i]->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2360  fhPtTrigChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2361  outputContainer->Add(fhPtTrigChargedPileUp[i]) ;
2362  }
2363  }
2364 
2365  if(IsHighMultiplicityAnalysisOn())
2366  {
2367  Int_t nMultiBins = GetNCentrBin();
2368  fhDeltaPhiChargedMult = new TH2F*[nMultiBins] ;
2369  fhDeltaEtaChargedMult = new TH2F*[nMultiBins] ;
2370  fhXEMult = new TH2F*[nMultiBins] ;
2371  fhXEUeMult = new TH2F*[nMultiBins] ;
2372  fhZTMult = new TH2F*[nMultiBins] ;
2373  fhZTUeMult = new TH2F*[nMultiBins] ;
2374 
2375  for(Int_t im=0; im<nMultiBins; im++)
2376  {
2377  fhDeltaPhiChargedMult[im] = new TH2F
2378  (Form("hDeltaPhiCharged_Mult%d",im),Form("#Delta #varphi charged Mult bin %d",im), nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2379  fhDeltaPhiChargedMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2380  fhDeltaPhiChargedMult[im]->SetYTitle("#Delta #varphi (rad)");
2381 
2382  fhDeltaEtaChargedMult[im] = new TH2F
2383  (Form("hDeltaEtaCharged_Mult%d",im),Form("#Delta #eta charged Mult bin %d",im), nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);
2384  fhDeltaEtaChargedMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2385  fhDeltaEtaChargedMult[im]->SetYTitle("#Delta #eta");
2386 
2387  fhXEMult[im] = new TH2F
2388  (Form("hXECharged_Mult%d",im),Form("#it{x}_{E} charged Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2389  fhXEMult[im]->SetYTitle("#it{x}_{E}");
2390  fhXEMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2391 
2392  fhXEUeMult[im] = new TH2F
2393  (Form("hXEUeCharged_Mult%d",im),Form("#it{x}_{E} UE charged Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2394  fhXEUeMult[im]->SetYTitle("#it{x}_{E}");
2395  fhXEUeMult[im]->SetXTitle("#it{p}_{T trigger}(GeV/#it{c})");
2396 
2397  fhZTMult[im] = new TH2F
2398  (Form("hZTCharged_Mult%d",im),Form("#it{z}_{T} charged Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2399  fhZTMult[im]->SetYTitle("#it{z}_{T}");
2400  fhZTMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2401 
2402  fhZTUeMult[im] = new TH2F
2403  (Form("hZTUeCharged_Mult%d",im),Form("#it{z}_{T} UE charged Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2404  fhZTUeMult[im]->SetYTitle("#it{z}_{T}");
2405  fhZTUeMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2406 
2407  outputContainer->Add(fhDeltaPhiChargedMult[im]) ;
2408  outputContainer->Add(fhDeltaEtaChargedMult[im]) ;
2409  outputContainer->Add(fhXEMult [im]);
2410  outputContainer->Add(fhXEUeMult[im]);
2411  outputContainer->Add(fhZTMult [im]);
2412  outputContainer->Add(fhZTUeMult[im]);
2413  }
2414  }
2415 
2416  if(fFillBradHisto)
2417  {
2418  fhAssocPtBkg = new TH2F("hAssocPtBkg", " Trigger #it{p}_{T} vs associated hadron #it{p}_{T} from background",
2419  nptbins, ptmin, ptmax,nptbins,ptmin,ptmax);
2420  fhAssocPtBkg->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2421  fhAssocPtBkg->SetYTitle("#it{p}_{T associated} (GeV/#it{c})");
2422  outputContainer->Add(fhAssocPtBkg) ;
2423 
2424  fhDeltaPhiBrad = new TH2F("hDeltaPhiBrad","atan2(sin(#Delta #varphi), cos(#Delta #varphi))/#pi vs #it{p}_{T trigger} ",
2425  nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
2426  fhDeltaPhiBrad->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2427  fhDeltaPhiBrad->SetYTitle("atan2(sin(#Delta #varphi), cos(#Delta #varphi))/#pi");
2428  outputContainer->Add(fhDeltaPhiBrad) ;
2429  }
2430 
2431  fhDeltaPhiDeltaEtaAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2432  fhDeltaPhiAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2433  fhDeltaPhiAssocPtBinDEta08 = new TH2F*[fNAssocPtBins*nz];
2434  fhDeltaPhiAssocPtBinDEta0 = new TH2F*[fNAssocPtBins*nz];
2435  if(fFillMomImbalancePtAssocBinsHisto)
2436  {
2437  fhXEAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2438  fhZTAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2439  }
2440 
2441  if(fCorrelVzBin)
2442  {
2443  fhXEVZ = new TH2F*[nz];
2444  fhZTVZ = new TH2F*[nz];
2445  }
2446 
2447  if(fFillBradHisto)
2448  fhDeltaPhiBradAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2449 
2450 
2451  fhDeltaPhiAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2452  if(fFillEtaGapsHisto)fhDeltaPhiAssocPtBinDEta08 = new TH2F*[fNAssocPtBins*nz];
2453  if(fDecayTrigger) fhDeltaPhiDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2454 
2455 
2456  if(fHMPIDCorrelation)
2457  {
2458  fhDeltaPhiAssocPtBinHMPID = new TH2F*[fNAssocPtBins*nz];
2459  fhDeltaPhiAssocPtBinHMPIDAcc= new TH2F*[fNAssocPtBins*nz];
2460  }
2461 
2462  for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2463  {
2464  for(Int_t z = 0 ; z < nz ; z++)
2465  {
2466  Int_t bin = i*nz+z;
2467 
2468  if(fCorrelVzBin)
2469  {
2470  sz = Form("_vz%d",z);
2471  tz = Form(", #it{v}_{#it{z}} bin %d",z);
2472  }
2473 
2474  //printf("iAssoc %d, Vz %d, bin %d - sz %s, tz %s \n",i,z,bin,sz.Data(),tz.Data());
2475 
2476  fhDeltaPhiDeltaEtaAssocPtBin[bin] = new TH2F(Form("hDeltaPhiDeltaEtaPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2477  Form("#Delta #varphi vs #Delta #eta for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2478  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
2479  fhDeltaPhiDeltaEtaAssocPtBin[bin]->SetXTitle("#Delta #varphi (rad)");
2480  fhDeltaPhiDeltaEtaAssocPtBin[bin]->SetYTitle("#Delta #eta");
2481 
2482  fhDeltaPhiAssocPtBin[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2483  Form("#Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2484  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2485  fhDeltaPhiAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2486  fhDeltaPhiAssocPtBin[bin]->SetYTitle("#Delta #varphi (rad)");
2487 
2488  outputContainer->Add(fhDeltaPhiDeltaEtaAssocPtBin[bin]) ;
2489  outputContainer->Add(fhDeltaPhiAssocPtBin[bin]) ;
2490 
2491  if(fFillEtaGapsHisto)
2492  {
2493  fhDeltaPhiAssocPtBinDEta08[bin] = new TH2F(Form("hDeltaPhiDeltaEta0.8PtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2494  Form("#Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta > 0.8", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2495  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2496  fhDeltaPhiAssocPtBinDEta08[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2497  fhDeltaPhiAssocPtBinDEta08[bin]->SetYTitle("#Delta #varphi (rad)");
2498 
2499  fhDeltaPhiAssocPtBinDEta0[bin] = new TH2F(Form("hDeltaPhiDeltaEta0PtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2500  Form("#Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta = 0.", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2501  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2502  fhDeltaPhiAssocPtBinDEta0[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2503  fhDeltaPhiAssocPtBinDEta0[bin]->SetYTitle("#Delta #varphi (rad)");
2504 
2505  outputContainer->Add(fhDeltaPhiAssocPtBinDEta08[bin]) ;
2506  outputContainer->Add(fhDeltaPhiAssocPtBinDEta0[bin]) ;
2507  }
2508 
2509  if(fDecayTrigger)
2510  {
2511  fhDeltaPhiDecayChargedAssocPtBin[bin] = new TH2F(Form("hDeltaPhiPtDecayChargedAssocPt%2.1f_%2.1f%s_bit%d", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data(),fDecayBits[0]),
2512  Form("#Delta #varphi vs #it{p}_{T trigger} tagged as decay for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, Bit %d", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data(),fDecayBits[0]),
2513  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2514  fhDeltaPhiDecayChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2515  fhDeltaPhiDecayChargedAssocPtBin[bin]->SetYTitle("#Delta #varphi (rad)");
2516 
2517  outputContainer->Add(fhDeltaPhiDecayChargedAssocPtBin[bin]) ;
2518  }
2519 
2520  if(fFillBradHisto)
2521  {
2522  fhDeltaPhiBradAssocPtBin[bin] = new TH2F(Form("hDeltaPhiBradPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2523  Form("atan2(sin(#Delta #varphi), cos(#Delta #varphi))/#pi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2524  nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
2525  fhDeltaPhiBradAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2526  fhDeltaPhiBradAssocPtBin[bin]->SetYTitle("atan2(sin(#Delta #varphi), cos(#Delta #varphi))/#pi");
2527  outputContainer->Add(fhDeltaPhiBradAssocPtBin[bin]) ;
2528  }
2529 
2530  if(fHMPIDCorrelation)
2531  {
2532  fhDeltaPhiAssocPtBinHMPID[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%sHMPID", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2533  Form("#Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, with track having HMPID signal", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2534  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2535  fhDeltaPhiAssocPtBinHMPID[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})" );
2536  fhDeltaPhiAssocPtBinHMPID[bin]->SetYTitle("#Delta #varphi (rad)");
2537 
2538  fhDeltaPhiAssocPtBinHMPIDAcc[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%sHMPIDAcc", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2539  Form("#Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, with track within 5<phi<20 deg", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2540  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2541  fhDeltaPhiAssocPtBinHMPIDAcc[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2542  fhDeltaPhiAssocPtBinHMPIDAcc[bin]->SetYTitle("#Delta #varphi (rad)");
2543 
2544  outputContainer->Add(fhDeltaPhiAssocPtBinHMPID [bin]) ;
2545  outputContainer->Add(fhDeltaPhiAssocPtBinHMPIDAcc[bin]) ;
2546  }
2547  }
2548  }
2549 
2550  if(fFillMomImbalancePtAssocBinsHisto)
2551  {
2552  for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2553  {
2554  fhXEAssocPtBin[i] = new TH2F(Form("hXEAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2555  Form("#it{x}_{#it{E}} vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2556  nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2557  fhXEAssocPtBin[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2558  fhXEAssocPtBin[i]->SetYTitle("#it{x}_{#it{E}}");
2559 
2560  fhZTAssocPtBin[i] = new TH2F(Form("hZTAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2561  Form("#it{z}_{T} vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2562  nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2563  fhZTAssocPtBin[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2564  fhZTAssocPtBin[i]->SetYTitle("#it{z}_{T}");
2565 
2566 
2567  outputContainer->Add(fhXEAssocPtBin[i]);
2568  outputContainer->Add(fhZTAssocPtBin[i]);
2569  }
2570  }
2571 
2572  if(fCorrelVzBin)
2573  {
2574  for(Int_t z = 0 ; z < nz ; z++)
2575  {
2576  sz = Form("_vz%d",z);
2577  tz = Form(", #it{v}_{#it{z}} bin %d",z);
2578 
2579  fhXEVZ[z] = new TH2F(Form("hXE%s", sz.Data()),
2580  Form("#it{x}_{#it{E}} vs #it{p}_{T trigger}%s", tz.Data()),
2581  nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2582  fhXEVZ[z]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2583  fhXEVZ[z]->SetYTitle("#it{x}_{#it{E}}");
2584 
2585  fhZTVZ[z] = new TH2F(Form("hZT%s",sz.Data()),
2586  Form("#it{z}_{T} vs #it{p}_{T trigger}%s", tz.Data()),
2587  nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2588  fhZTVZ[z]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2589  fhZTVZ[z]->SetYTitle("#it{z}_{T}");
2590 
2591  outputContainer->Add(fhXEVZ[z]);
2592  outputContainer->Add(fhZTVZ[z]);
2593  }
2594  }
2595 
2596  if(fPi0Trigger)
2597  {
2598  fhPtPi0DecayRatio = new TH2F
2599  ("hPtPi0DecayRatio","#it{p}_{T} of #pi^{0} and the ratio of pt for two decay",
2600  nptbins,ptmin,ptmax, 100,0.,2.);
2601  fhPtPi0DecayRatio->SetXTitle("#it{p}_{T}^{#pi^{0}} (GeV/#it{c})");
2602  fhPtPi0DecayRatio->SetYTitle("#it{p}_{T}^{Decay}/#it{p}_{T}^{#pi^{0}}");
2603  outputContainer->Add(fhPtPi0DecayRatio) ;
2604 
2605  fhDeltaPhiPi0DecayCharged = new TH2F
2606  ("hDeltaPhiPi0DecayCharged","#varphi_{Decay} - #varphi_{h^{#pm}} vs #it{p}_{T Decay}",
2607  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2608  fhDeltaPhiPi0DecayCharged->SetYTitle("#Delta #varphi (rad)");
2609  fhDeltaPhiPi0DecayCharged->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2610 
2611  fhXEPi0DecayCharged =
2612  new TH2F("hXEPi0DecayCharged","#it{x}_{#it{E}} Decay",
2613  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2614  fhXEPi0DecayCharged->SetYTitle("#it{x}_{#it{E}}");
2615  fhXEPi0DecayCharged->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2616 
2617  fhZTPi0DecayCharged =
2618  new TH2F("hZTPi0DecayCharged","#it{z}_{trigger h^{#pm}} = #it{p}_{T h^{#pm}} / #it{p}_{T Decay}",
2619  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2620  fhZTPi0DecayCharged->SetYTitle("#it{z}_{decay h^{#pm}}");
2621  fhZTPi0DecayCharged->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2622 
2623  outputContainer->Add(fhDeltaPhiPi0DecayCharged) ;
2624  outputContainer->Add(fhXEPi0DecayCharged) ;
2625  outputContainer->Add(fhZTPi0DecayCharged) ;
2626  }
2627 
2628  if(fDecayTrigger)
2629  {
2630  for(Int_t ibit = 0; ibit< fNDecayBits; ibit++)
2631  {
2632  fhDeltaPhiDecayCharged[ibit] = new TH2F
2633  (Form("hDeltaPhiDecayCharged_bit%d",fDecayBits[ibit]),
2634  Form("#varphi_{Decay} - #varphi_{h^{#pm}} vs #it{p}_{T Decay}, Bit %d",fDecayBits[ibit]),
2635  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2636  fhDeltaPhiDecayCharged[ibit]->SetYTitle("#Delta #varphi (rad)");
2637  fhDeltaPhiDecayCharged[ibit]->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2638 
2639  fhXEDecayCharged[ibit] =
2640  new TH2F(Form("hXEDecayCharged_bit%d",fDecayBits[ibit]),
2641  Form("#it{x}_{#it{E}} Decay, Bit %d",fDecayBits[ibit]),
2642  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2643  fhXEDecayCharged[ibit]->SetYTitle("#it{x}_{#it{E}}");
2644  fhXEDecayCharged[ibit]->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2645 
2646  fhZTDecayCharged[ibit] =
2647  new TH2F(Form("hZTDecayCharged_bit%d",fDecayBits[ibit]),
2648  Form("#it{z}_{trigger h^{#pm}} = #it{p}_{T h^{#pm}} / #it{p}_{T Decay}, Bit %d",fDecayBits[ibit]),
2649  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2650  fhZTDecayCharged[ibit]->SetYTitle("#it{z}_{decay h^{#pm}}");
2651  fhZTDecayCharged[ibit]->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2652 
2653  outputContainer->Add(fhDeltaPhiDecayCharged[ibit]) ;
2654  outputContainer->Add(fhXEDecayCharged[ibit]) ;
2655  outputContainer->Add(fhZTDecayCharged[ibit]) ;
2656  }
2657  }
2658 
2659  //Correlation with neutral hadrons
2660  if(fNeutralCorr)
2661  {
2662  fhDeltaPhiDeltaEtaNeutral = new TH2F
2663  ("hDeltaPhiDeltaEtaNeutral","#varphi_{trigger} - #varphi_{h^{0}} vs #eta_{trigger} - #eta_{h^{0}}",
2664  ndeltaphibins ,deltaphimin,deltaphimax, ndeltaetabins ,deltaetamin,deltaetamax);
2665  fhDeltaPhiDeltaEtaNeutral->SetXTitle("#Delta #varphi (rad)");
2666  fhDeltaPhiDeltaEtaNeutral->SetYTitle("#Delta #eta");
2667 
2668  fhPhiNeutral = new TH2F
2669  ("hPhiNeutral","#varphi_{#pi^{0}} vs #it{p}_{T #pi^{0}}",
2670  nptbins,ptmin,ptmax,180,0,TMath::TwoPi());
2671  fhPhiNeutral->SetYTitle("#varphi_{#pi^{0}} (rad)");
2672  fhPhiNeutral->SetXTitle("#it{p}_{T #pi^{0}} (GeV/#it{c})");
2673 
2674  fhEtaNeutral = new TH2F
2675  ("hEtaNeutral","#eta_{#pi^{0}} vs #it{p}_{T #pi^{0}}",
2676  nptbins,ptmin,ptmax,200,-1.,1.);
2677  fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
2678  fhEtaNeutral->SetXTitle("#it{p}_{T #pi^{0}} (GeV/#it{c})");
2679 
2680  fhDeltaPhiNeutral = new TH2F
2681  ("hDeltaPhiNeutral","#varphi_{trigger} - #varphi_{#pi^{0}} vs #it{p}_{T trigger}",
2682  nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2683  fhDeltaPhiNeutral->SetYTitle("#Delta #varphi (rad)");
2684  fhDeltaPhiNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2685 
2686  fhDeltaPhiNeutralPt = new TH2F
2687  ("hDeltaPhiNeutralPt","#varphi_{trigger} - #varphi_{#pi^{0}} vs #it{p}_{T #pi^{0}}}",
2688  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2689  fhDeltaPhiNeutralPt->SetYTitle("#Delta #varphi (rad)");
2690  fhDeltaPhiNeutralPt->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2691 
2692  fhDeltaEtaNeutral = new TH2F
2693  ("hDeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs #it{p}_{T trigger}",
2694  nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);
2695  fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
2696  fhDeltaEtaNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2697 
2698  fhXENeutral =
2699  new TH2F("hXENeutral","#it{x}_{#it{E}} for #pi^{0} associated",
2700  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2701  fhXENeutral->SetYTitle("#it{x}_{#it{E}}");
2702  fhXENeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2703 
2704  fhPtHbpXENeutral =
2705  new TH2F("hHbpXENeutral","#xi = ln(1/#it{x}_{#it{E}})for #pi^{0} associated",
2706  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2707  fhPtHbpXENeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2708  fhPtHbpXENeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2709 
2710  fhZTNeutral =
2711  new TH2F("hZTNeutral","#it{z}_{trigger #pi} = #it{p}_{T #pi^{0}} / #it{p}_{T trigger} for #pi^{0} associated",
2712  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2713  fhZTNeutral->SetYTitle("#it{z}_{trigger #pi^{0}}");
2714  fhZTNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2715 
2716  fhPtHbpZTNeutral =
2717  new TH2F("hHbpZTNeutral","#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2718  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2719  fhPtHbpZTNeutral->SetYTitle("ln(1/#it{z}_{T})");
2720  fhPtHbpZTNeutral->SetXTitle("#it{p}_{T trigger}");
2721 
2722  outputContainer->Add(fhDeltaPhiDeltaEtaNeutral);
2723  outputContainer->Add(fhPhiNeutral) ;
2724  outputContainer->Add(fhEtaNeutral) ;
2725  outputContainer->Add(fhDeltaPhiNeutral) ;
2726  outputContainer->Add(fhDeltaPhiNeutralPt) ;
2727  outputContainer->Add(fhDeltaEtaNeutral) ;
2728  outputContainer->Add(fhXENeutral) ;
2729  outputContainer->Add(fhPtHbpXENeutral) ;
2730  outputContainer->Add(fhZTNeutral) ;
2731  outputContainer->Add(fhPtHbpZTNeutral) ;
2732 
2733  fhDeltaPhiUeNeutralPt = new TH2F
2734  (Form("hDeltaPhiUe%sNeutralPt",right.Data()),"#varphi_{trigger} - #varphi_{#pi^{0}} vs #it{p}_{T #pi^{0}}}",
2735  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2736  fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #varphi (rad)");
2737  fhDeltaPhiUeNeutralPt->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2738 
2739  fhXEUeNeutral =
2740  new TH2F(Form("hXEUeNeutral%s",right.Data()),"#it{x}_{#it{E}} for #pi^{0} associated",
2741  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2742  fhXEUeNeutral->SetYTitle("#it{x}_{#it{E}}");
2743  fhXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2744 
2745  fhPtHbpXEUeNeutral =
2746  new TH2F(Form("hHbpXEUeNeutral%s",right.Data()),"#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2747  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2748  fhPtHbpXEUeNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2749  fhPtHbpXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2750 
2751  fhZTUeNeutral =
2752  new TH2F(Form("hZTUeNeutral%s",right.Data()),"#it{z}_{trigger #pi} = #it{p}_{T #pi^{0}} / #it{p}_{T trigger} for #pi^{0} associated",
2753  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2754  fhZTUeNeutral->SetYTitle("#it{z}_{trigger #pi^{0}}");
2755  fhZTUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2756 
2757  fhPtHbpZTUeNeutral =
2758  new TH2F(Form("hHbpZTUeNeutral%s",right.Data()),"#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2759  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2760  fhPtHbpXEUeNeutral->SetYTitle("ln(1/#it{z}_{T})");
2761  fhPtHbpXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2762 
2763  outputContainer->Add(fhDeltaPhiUeNeutralPt) ;
2764  outputContainer->Add(fhXEUeNeutral) ;
2765  outputContainer->Add(fhPtHbpXEUeNeutral) ;
2766  outputContainer->Add(fhZTUeNeutral) ;
2767  outputContainer->Add(fhPtHbpZTUeNeutral) ;
2768 
2769  if(fMakeSeveralUE)
2770  {
2771  fhDeltaPhiUeLeftNeutral = new TH2F
2772  ("hDeltaPhiUeLeftNeutralPt","#varphi_{trigger} - #varphi_{#Ueh^{0}} vs #it{p}_{T h^{0}} with neutral UE left side range of trigger particles",
2773  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2774  fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #varphi (rad)");
2775  fhDeltaPhiUeLeftNeutral->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2776  outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
2777 
2778  fhXEUeLeftNeutral =
2779  new TH2F("hXEUeNeutralLeft","#it{x}_{#it{E}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE left side of trigger",
2780  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2781  fhXEUeLeftNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2782  fhXEUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2783  outputContainer->Add(fhXEUeLeftNeutral) ;
2784 
2785  fhPtHbpXEUeLeftNeutral =
2786  new TH2F("hHbpXEUeNeutralLeft","#xi = ln(1/#it{x}_{#it{E}}) with neutral UE left side of trigger",
2787  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2788  fhPtHbpXEUeLeftNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2789  fhPtHbpXEUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2790  outputContainer->Add(fhPtHbpXEUeLeftNeutral) ;
2791 
2792  fhZTUeLeftNeutral =
2793  new TH2F("hZTUeNeutralLeft","#it{z}_{trigger h^{0}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE left side of trigger",
2794  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2795  fhZTUeLeftNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2796  fhZTUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2797  outputContainer->Add(fhZTUeLeftNeutral) ;
2798 
2799  fhPtHbpZTUeLeftNeutral =
2800  new TH2F("hHbpZTUeNeutralLeft","#xi = ln(1/#it{z}_{T}) with neutral UE left side of trigger",
2801  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2802  fhPtHbpZTUeLeftNeutral->SetYTitle("ln(1/#it{z}_{T})");
2803  fhPtHbpZTUeLeftNeutral->SetXTitle("#it{p}_{T trigger}");
2804  outputContainer->Add(fhPtHbpZTUeLeftNeutral) ;
2805  }
2806 
2807  if(fPi0Trigger)
2808  {
2809  fhDeltaPhiPi0DecayNeutral = new TH2F
2810  ("hDeltaPhiPi0DecayNeutral","#varphi_{Decay} - #varphi_{h^{0}} vs #it{p}_{T Decay}",
2811  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2812  fhDeltaPhiPi0DecayNeutral->SetYTitle("#Delta #varphi (rad)");
2813  fhDeltaPhiPi0DecayNeutral->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2814 
2815  fhXEPi0DecayNeutral =
2816  new TH2F("hXEPi0DecayNeutral","#it{x}_{#it{E}} for decay trigger",
2817  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2818  fhXEPi0DecayNeutral->SetYTitle("#it{x}_{#it{E}}");
2819  fhXEPi0DecayNeutral->SetXTitle("#it{p}_{T decay}");
2820 
2821  fhZTPi0DecayNeutral =
2822  new TH2F("hZTPi0DecayNeutral","#it{z}_{trigger h^{0}} = #it{p}_{T h^{0}} / #it{p}_{T Decay}",
2823  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2824  fhZTPi0DecayNeutral->SetYTitle("#it{z}_{h^{0}}");
2825  fhZTPi0DecayNeutral->SetXTitle("#it{p}_{T decay}");
2826 
2827  outputContainer->Add(fhDeltaPhiPi0DecayNeutral) ;
2828  outputContainer->Add(fhXEPi0DecayNeutral) ;
2829  outputContainer->Add(fhZTPi0DecayNeutral) ;
2830  }
2831  }//Correlation with neutral hadrons
2832 
2833  // If data is MC, fill more histograms, depending on origin
2834  if(IsDataMC())
2835  {
2836  for(Int_t i= fMCGenTypeMin; i <= fMCGenTypeMax; i++)
2837  {
2838  fhMCPtTrigger[i] = new TH1F (Form("hMCPtTrigger_%s",nameMC[i].Data()),
2839  Form("MC %s: generated trigger #it{p}_{T}",nameMC[i].Data()),
2840  nptbins,ptmin,ptmax);
2841  fhMCPtTrigger[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2842 
2843  fhMCPhiTrigger[i] = new TH2F (Form("hMCPhiTrigger_%s",nameMC[i].Data()),
2844  Form("MC %s: generated trigger #varphi",nameMC[i].Data()),
2845  nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2846  fhMCPhiTrigger[i]->SetYTitle("#varphi (rad)");
2847  fhMCPhiTrigger[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2848 
2849  fhMCEtaTrigger[i] = new TH2F (Form("hMCEtaTrigger_%s",nameMC[i].Data()),
2850  Form("MC %s: generated trigger #eta",nameMC[i].Data()),
2851  nptbins,ptmin,ptmax, netabins,etamin,etamax);
2852  fhMCEtaTrigger[i]->SetYTitle("#eta");
2853  fhMCEtaTrigger[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2854 
2855  if(fMakeAbsoluteLeading || fMakeNearSideLeading)
2856  {
2857  fhMCPtTriggerNotLeading[i] = new TH1F (Form("hMCPtTriggerNotLeading_%s",nameMC[i].Data()),
2858  Form("MC %s: generated trigger #it{p}_{T}, when not leading of primaries",nameMC[i].Data()),
2859  nptbins,ptmin,ptmax);
2860  fhMCPtTriggerNotLeading[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2861 
2862  fhMCPhiTriggerNotLeading[i] = new TH2F (Form("hMCPhiTriggerNotLeading_%s",nameMC[i].Data()),
2863  Form("MC %s: generated trigger #varphi, when not leading of primaries",nameMC[i].Data()),
2864  nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2865  fhMCPhiTriggerNotLeading[i]->SetYTitle("#varphi (rad)");
2866  fhMCPhiTriggerNotLeading[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2867 
2868 
2869  fhMCEtaTriggerNotLeading[i] = new TH2F (Form("hMCEtaTriggerNotLeading_%s",nameMC[i].Data()),
2870  Form("MC %s: generated triogger #eta, when not leading of primaries",nameMC[i].Data()),
2871  nptbins,ptmin,ptmax, netabins,etamin,etamax);
2872  fhMCEtaTriggerNotLeading[i]->SetYTitle("#eta ");
2873  fhMCEtaTriggerNotLeading[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2874  }
2875 
2876  fhMCEtaCharged[i] = new TH2F (Form("hMCEtaCharged_%s",nameMC[i].Data()),
2877  Form("MC %s: #eta_{h^{#pm}} vs #it{p}_{T #pm}",nameMC[i].Data()),
2878  nptbins,ptmin,ptmax,100,-1.,1.);
2879  fhMCEtaCharged[i]->SetYTitle("#eta_{h^{#pm}} (rad)");
2880  fhMCEtaCharged[i]->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
2881 
2882  fhMCPhiCharged[i] = new TH2F(Form("hMCPhiCharged_%s",nameMC[i].Data()),
2883  Form("MC %s: phi_{h^{#pm}} vs #it{p}_{T #pm}",nameMC[i].Data()),
2884  nptbins,ptmin,ptmax,180,0,TMath::TwoPi());
2885  fhMCPhiCharged[i]->SetYTitle("MC #varphi_{h^{#pm}} (rad)");
2886  fhMCPhiCharged[i]->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
2887 
2888  fhMCDeltaPhiDeltaEtaCharged[i] = new TH2F (Form("hMCDeltaPhiDeltaEtaCharged_%s",nameMC[i].Data()),
2889  Form("MC %s: phi_{trigger} - #varphi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",nameMC[i].Data()),
2890  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
2891  fhMCDeltaPhiDeltaEtaCharged[i]->SetXTitle("#Delta #varphi (rad)");
2892  fhMCDeltaPhiDeltaEtaCharged[i]->SetYTitle("#Delta #eta");
2893 
2894  fhMCDeltaEtaCharged[i] = new TH2F (Form("hMCDeltaEtaCharged_%s",nameMC[i].Data()),
2895  Form("MC %s: #eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger} and #it{p}_{T assoc}",nameMC[i].Data()),
2896  nptbins,ptmin,ptmax,ndeltaetabins ,deltaetamin,deltaetamax);
2897  fhMCDeltaEtaCharged[i]->SetYTitle("#Delta #eta");
2898  fhMCDeltaEtaCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2899 
2900  fhMCDeltaPhiCharged[i] = new TH2F (Form("hMCDeltaPhiCharged_%s",nameMC[i].Data()),
2901  Form("MC %s: #varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}",nameMC[i].Data()),
2902  nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2903  fhMCDeltaPhiCharged[i]->SetYTitle("#Delta #varphi (rad)");
2904  fhMCDeltaPhiCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2905 
2906  fhMCDeltaPhiChargedPt[i] = new TH2F (Form("hMCDeltaPhiChargedPt_%s",nameMC[i].Data()),
2907  Form("MC %s: #varphi_{trigger} - #varphi_{#h^{#pm}} vs #it{p}_{T h^{#pm}}",nameMC[i].Data()),
2908  nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2909  fhMCDeltaPhiChargedPt[i]->SetYTitle("#Delta #varphi (rad)");
2910  fhMCDeltaPhiChargedPt[i]->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2911 
2912  fhMCPtXECharged[i] = new TH2F (Form("hMCPtXECharged_%s",nameMC[i].Data()),
2913  Form("MC %s: #it{x}_{#it{E}} with charged hadrons",nameMC[i].Data()),
2914  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2915  fhMCPtXECharged[i]->SetYTitle("#it{x}_{#it{E}}");
2916  fhMCPtXECharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2917 
2918  fhMCPtHbpXECharged[i] = new TH2F(Form("hMCHbpXECharged_%s",nameMC[i].Data()),
2919  Form("MC %s: #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",nameMC[i].Data()),
2920  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2921  fhMCPtHbpXECharged[i]->SetYTitle("ln(1/#it{x}_{#it{E}})");
2922  fhMCPtHbpXECharged[i]->SetXTitle("#it{p}_{T trigger}");
2923 
2924  fhMCPtZTCharged[i] = new TH2F(Form("hMCPtZTCharged_%s",nameMC[i].Data()),
2925  Form("MC %s: #it{z}_{T} with charged hadrons",nameMC[i].Data()),
2926  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2927  fhMCPtZTCharged[i]->SetYTitle("#it{z}_{T}");
2928  fhMCPtZTCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2929 
2930  fhMCPtHbpZTCharged[i] = new TH2F(Form("hMCHbpZTCharged_%s",nameMC[i].Data()),
2931  Form("MC %s: #xi = ln(1/#it{z}_{T}) with charged hadrons",nameMC[i].Data()),
2932  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2933  fhMCPtHbpZTCharged[i]->SetYTitle("ln(1/#it{z}_{T})");
2934  fhMCPtHbpZTCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2935 
2936  fhMCPtTrigPout[i] = new TH2F(Form("hMCPtTrigPout_%s",nameMC[i].Data()),
2937  Form("MC %s: #it{p}_{out} with triggers",nameMC[i].Data()),
2938  nptbins,ptmin,ptmax,nptbins,-1.*ptmax/2.,ptmax/2.);
2939  fhMCPtTrigPout[i]->SetYTitle("#it{p}_{out} (GeV/#it{c})");
2940  fhMCPtTrigPout[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2941 
2942  fhMCPtAssocDeltaPhi[i] = new TH2F(Form("hMCPtAssocDeltaPhi_%s",nameMC[i].Data()),
2943  Form("MC %s: #Delta #varphi with associated charged hadrons",nameMC[i].Data()),
2944  nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2945  fhMCPtAssocDeltaPhi[i]->SetYTitle("#Delta #varphi (rad)");
2946  fhMCPtAssocDeltaPhi[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2947 
2948  outputContainer->Add(fhMCPtTrigger[i]);
2949  outputContainer->Add(fhMCPhiTrigger[i]);
2950  outputContainer->Add(fhMCEtaTrigger[i]);
2951 
2952  if(fMakeAbsoluteLeading || fMakeNearSideLeading)
2953  {
2954  outputContainer->Add(fhMCPtTriggerNotLeading[i]);
2955  outputContainer->Add(fhMCPhiTriggerNotLeading[i]);
2956  outputContainer->Add(fhMCEtaTriggerNotLeading[i]);
2957  }
2958 
2959  outputContainer->Add(fhMCDeltaPhiDeltaEtaCharged[i]);
2960  outputContainer->Add(fhMCPhiCharged[i]) ;
2961  outputContainer->Add(fhMCEtaCharged[i]) ;
2962  outputContainer->Add(fhMCDeltaEtaCharged[i]) ;
2963  outputContainer->Add(fhMCDeltaPhiCharged[i]) ;
2964 
2965  outputContainer->Add(fhMCDeltaPhiChargedPt[i]) ;
2966  outputContainer->Add(fhMCPtXECharged[i]) ;
2967  outputContainer->Add(fhMCPtZTCharged[i]) ;
2968  outputContainer->Add(fhMCPtHbpXECharged[i]) ;
2969  outputContainer->Add(fhMCPtHbpZTCharged[i]) ;
2970  outputContainer->Add(fhMCPtTrigPout[i]) ;
2971  outputContainer->Add(fhMCPtAssocDeltaPhi[i]) ;
2972 
2973  // Underlying event
2974 
2975  fhMCUePart[i] =
2976  new TH1F(Form("hMCUePart_%s",nameMC[i].Data()),
2977  Form("MC %s: UE particles distribution vs #it{p}_{T trigger}",nameMC[i].Data()),
2978  nptbins,ptmin,ptmax);
2979  fhMCUePart[i]->SetYTitle("#it{dN}^{ch}");
2980  fhMCUePart[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2981 
2982  fhMCPtXEUeCharged[i] =
2983  new TH2F(Form("hMCPtXEUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2984  Form("MC %s: #it{x}_{#it{E}} with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2985  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2986  fhMCPtXEUeCharged[i]->SetYTitle("#it{x}_{#it{E}}");
2987  fhMCPtXEUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2988 
2989  fhMCPtHbpXEUeCharged[i] =
2990  new TH2F(Form("hMCPtHbpXEUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2991  Form("MC %s: #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2992  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2993  fhMCPtHbpXEUeCharged[i]->SetYTitle("ln(1/#it{x}_{#it{E}})");
2994  fhMCPtHbpXEUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2995 
2996  fhMCPtZTUeCharged[i] =
2997  new TH2F(Form("hMCPtZTUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2998  Form("MC %s: #it{z}_{T} with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2999  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3000  fhMCPtZTUeCharged[i]->SetYTitle("#it{z}_{T}");
3001  fhMCPtZTUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3002 
3003  fhMCPtHbpZTUeCharged[i] =
3004  new TH2F(Form("hMCPtHbpZTUeCharged%s_%s",right.Data(),nameMC[i].Data()),
3005  Form("MC %s: #xi = ln(1/#it{z}_{T}) with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
3006  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
3007  fhMCPtHbpZTUeCharged[i]->SetYTitle("ln(1/#it{z}_{T})");
3008  fhMCPtHbpZTUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3009 
3010  outputContainer->Add(fhMCUePart[i]);
3011  outputContainer->Add(fhMCPtXEUeCharged[i]) ;
3012  outputContainer->Add(fhMCPtZTUeCharged[i]) ;
3013  outputContainer->Add(fhMCPtHbpZTUeCharged[i]);
3014  outputContainer->Add(fhMCPtHbpXEUeCharged[i]);
3015 
3016  if(fMakeSeveralUE)
3017  {
3018  fhMCPtXEUeLeftCharged[i] = new TH2F(Form("hMCPtXEUeChargedLeft_%s",nameMC[i].Data()),
3019  Form("MC %s: #it{x}_{#it{E}} with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
3020  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3021  fhMCPtXEUeLeftCharged[i]->SetYTitle("#it{x}_{#it{E}}");
3022  fhMCPtXEUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3023 
3024  fhMCPtHbpXEUeLeftCharged[i] = new TH2F(Form("hMCPtHbpXEUeChargedLeft_%s",nameMC[i].Data()),
3025  Form("MC %s: #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
3026  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
3027  fhMCPtHbpXEUeLeftCharged[i]->SetYTitle("ln(1/#it{x}_{#it{E}})");
3028  fhMCPtHbpXEUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3029 
3030  fhMCPtZTUeLeftCharged[i] = new TH2F(Form("hMCPtZTUeChargedLeft_%s",nameMC[i].Data()),
3031  Form("MC %s: #it{z}_{T} with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
3032  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3033  fhMCPtZTUeLeftCharged[i]->SetYTitle("#it{z}_{T}");
3034  fhMCPtZTUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3035 
3036  fhMCPtHbpZTUeLeftCharged[i] = new TH2F(Form("hMCPtHbpZTUeChargedLeft_%s",nameMC[i].Data()),
3037  Form("MC %s: #xi = ln(1/#it{z}_{T}) with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
3038  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
3039  fhMCPtHbpZTUeLeftCharged[i]->SetYTitle("ln(1/#it{z}_{T})");
3040  fhMCPtHbpZTUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3041 
3042  outputContainer->Add(fhMCPtXEUeLeftCharged[i]) ;
3043  outputContainer->Add(fhMCPtZTUeLeftCharged[i]) ;
3044  outputContainer->Add(fhMCPtHbpXEUeLeftCharged[i]);
3045  outputContainer->Add(fhMCPtHbpZTUeLeftCharged[i]) ;
3046 
3047  }
3048  }
3049  } // For MC histogram
3050 
3051  if(DoOwnMix())
3052  {
3053  // Create event containers
3054 
3055  if(!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForTracksExists()))
3056  {
3057  Int_t nvz = GetNZvertBin();
3058  Int_t nrp = GetNRPBin();
3059  Int_t nce = GetNCentrBin();
3060 
3061  fListMixTrackEvents= new TList*[nvz*nrp*nce] ;
3062 
3063  for( Int_t ice = 0 ; ice < nce ; ice++ )
3064  {
3065  for( Int_t ivz = 0 ; ivz < nvz ; ivz++ )
3066  {
3067  for( Int_t irp = 0 ; irp < nrp ; irp++ )
3068  {
3069  Int_t bin = GetEventMixBin(ice,ivz,irp); //ic*nvz*nrp+iz*nrp+irp;
3070 
3071  //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
3072  // ic,iz, irp, bin);
3073 
3074  fListMixTrackEvents[bin] = new TList() ;
3075  fListMixTrackEvents[bin]->SetOwner(kFALSE);
3076  }
3077  }
3078  }
3079  }
3080 
3081  fhPtTriggerMixed = new TH1F ("hPtTriggerMixed","#it{p}_{T} distribution of trigger particles, used for mixing", nptbins,ptmin,ptmax);
3082  fhPtTriggerMixed->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
3083 
3084  if(fCorrelVzBin)
3085  {
3086  fhPtTriggerMixedVzBin = new TH2F ("hPtTriggerMixedVzBin","#it{p}_{T} distribution of trigger particles, used for mixing", nptbins,ptmin,ptmax,GetNZvertBin(),0,GetNZvertBin());
3087  fhPtTriggerMixedVzBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
3088  fhPtTriggerMixedVzBin->SetYTitle("#it{v}_{#it{z}} bin");
3089  outputContainer->Add(fhPtTriggerMixedVzBin);
3090  }
3091 
3092  fhPtTriggerMixedBin = new TH2F ("hPtTriggerMixedBin","#it{p}_{T} distribution of trigger particles vs mixing bin", nptbins,ptmin,ptmax,nMixBins,0,nMixBins);
3093  fhPtTriggerMixedBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
3094  fhPtTriggerMixedBin->SetYTitle("Bin");
3095 
3096  fhPhiTriggerMixed = new TH2F ("hPhiTriggerMixed","#varphi distribution of trigger Particles, used for mixing",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
3097  fhPhiTriggerMixed->SetYTitle("#varphi (rad)");
3098 
3099  fhEtaTriggerMixed = new TH2F ("hEtaTriggerMixed","#eta distribution of trigger, used for mixing",nptbins,ptmin,ptmax, netabins,etamin,etamax);
3100  fhEtaTriggerMixed->SetYTitle("#eta ");
3101 
3102  outputContainer->Add(fhPtTriggerMixed);
3103  outputContainer->Add(fhPtTriggerMixedBin);
3104  outputContainer->Add(fhPhiTriggerMixed);
3105  outputContainer->Add(fhEtaTriggerMixed);
3106 
3107  // Fill the cluster pool only in isolation analysis or if requested
3108  if( neutralMix && (!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForCaloExists())))
3109  {
3110  Int_t nvz = GetNZvertBin();
3111  Int_t nrp = GetNRPBin();
3112  Int_t nce = GetNCentrBin();
3113 
3114  fListMixCaloEvents= new TList*[nvz*nrp*nce] ;
3115 
3116  for( Int_t ice = 0 ; ice < nce ; ice++ )
3117  {
3118  for( Int_t ivz = 0 ; ivz < nvz ; ivz++ )
3119  {
3120  for( Int_t irp = 0 ; irp < nrp ; irp++ )
3121  {
3122  Int_t bin = GetEventMixBin(ice,ivz,irp); //ic*nvz*nrp+iz*nrp+irp;
3123 
3124  //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
3125  // ic,iz, irp, bin);
3126 
3127  fListMixCaloEvents[bin] = new TList() ;
3128  fListMixCaloEvents[bin]->SetOwner(kFALSE);
3129  }
3130  }
3131  }
3132  }
3133 
3134  // Init the list in the reader if not done previously
3135  if(fUseMixStoredInReader)
3136  {
3137  if( !GetReader()->ListWithMixedEventsForTracksExists() )
3138  GetReader()->SetListWithMixedEventsForTracks(fListMixTrackEvents);
3139 
3140  if( !GetReader()->ListWithMixedEventsForCaloExists() )
3141  GetReader()->SetListWithMixedEventsForCalo (fListMixCaloEvents );
3142  }
3143 
3144  fhEventBin=new TH1I("hEventBin","Number of triggers per bin(cen,vz,rp)",
3145  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
3146  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
3147  fhEventBin->SetXTitle("event bin");
3148  outputContainer->Add(fhEventBin) ;
3149 
3150  fhEventMixBin=new TH1I("hEventMixBin","Number of triggers mixed per event bin(cen,vz,rp)",
3151  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
3152  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
3153  fhEventMixBin->SetXTitle("event bin");
3154  outputContainer->Add(fhEventMixBin) ;
3155 
3156  fhEventMBBin=new TH1I("hEventMBBin","Number of min bias events per bin(cen,vz,rp)",
3157  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
3158  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
3159  fhEventMBBin->SetXTitle("event bin");
3160  outputContainer->Add(fhEventMBBin) ;
3161 
3162  fhNtracksMB=new TH2F("hNtracksMBEvent","Number of filtered tracks in MB event per event bin",ntrbins,trmin,trmax,
3163  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
3164  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
3165  fhNtracksMB->SetYTitle("event bin");
3166  fhNtracksMB->SetXTitle("#it{N}_{track}");
3167  outputContainer->Add(fhNtracksMB);
3168 
3169  if( neutralMix )
3170  {
3171  fhNclustersMB=new TH2F("hNclustersMBEvent","Number of filtered clusters in MB events per event bin",nclbins,clmin,clmax,
3172  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
3173  GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
3174  fhNclustersMB->SetYTitle("event bin");
3175  fhNclustersMB->SetXTitle("#it{N}_{cluster}");
3176  outputContainer->Add(fhNclustersMB);
3177  }
3178 
3179  fhMixDeltaPhiCharged = new TH2F
3180  ("hMixDeltaPhiCharged","Mixed event : #varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}",
3181  nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
3182  fhMixDeltaPhiCharged->SetYTitle("#Delta #varphi (rad)");
3183  fhMixDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3184  outputContainer->Add(fhMixDeltaPhiCharged);
3185 
3186  fhMixDeltaPhiDeltaEtaCharged = new TH2F
3187  ("hMixDeltaPhiDeltaEtaCharged","Mixed event : #varphi_{trigger} - #varphi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
3188  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
3189  fhMixDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #varphi (rad)");
3190  fhMixDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
3191  outputContainer->Add(fhMixDeltaPhiDeltaEtaCharged);
3192 
3193  fhMixXECharged =
3194  new TH2F("hMixXECharged","Mixed event : #it{x}_{#it{E}} for charged tracks",
3195  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3196  fhMixXECharged->SetYTitle("#it{x}_{#it{E}}");
3197  fhMixXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3198  outputContainer->Add(fhMixXECharged);
3199 
3200  fhMixXEUeCharged =
3201  new TH2F("hMixXEUeCharged","Mixed event : #it{x}_{#it{E}} for charged tracks in Ue region",
3202  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3203  fhMixXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
3204  fhMixXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3205  outputContainer->Add(fhMixXEUeCharged);
3206 
3207  fhMixHbpXECharged =
3208  new TH2F("hMixHbpXECharged","mixed event : #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
3209  nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
3210  fhMixHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
3211  fhMixHbpXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3212  outputContainer->Add(fhMixHbpXECharged);
3213 
3214  fhMixDeltaPhiChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
3215  fhMixDeltaPhiChargedAssocPtBinDEta08 = new TH2F*[fNAssocPtBins*nz];
3216  fhMixDeltaPhiChargedAssocPtBinDEta0 = new TH2F*[fNAssocPtBins*nz];
3217  fhMixDeltaPhiDeltaEtaChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
3218 
3219  for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
3220  {
3221  for(Int_t z = 0 ; z < nz ; z++)
3222  {
3223  Int_t bin = i*nz+z;
3224 
3225  if(fCorrelVzBin)
3226  {
3227  sz = Form("_vz%d",z);
3228  tz = Form(", #it{v}_{#it{z}} bin %d",z);
3229  }
3230 
3231  //printf("MIX : iAssoc %d, Vz %d, bin %d - sz %s, tz %s \n",i,z,bin,sz.Data(),tz.Data());
3232 
3233  fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin] = new TH2F(Form("hMixDeltaPhiDeltaEtaChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
3234  Form("Mixed event #Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
3235  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
3236  fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->SetXTitle("#Delta #varphi (rad)");
3237  fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->SetYTitle("#Delta #eta");
3238 
3239  outputContainer->Add(fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]);
3240 
3241  fhMixDeltaPhiChargedAssocPtBin[bin] = new TH2F(Form("hMixDeltaPhiChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
3242  Form("Mixed event #Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
3243  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
3244  fhMixDeltaPhiChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3245  fhMixDeltaPhiChargedAssocPtBin[bin]->SetYTitle("#Delta #varphi (rad)");
3246 
3247  outputContainer->Add(fhMixDeltaPhiChargedAssocPtBin[bin]);
3248 
3249  if(fFillEtaGapsHisto)
3250  {
3251  fhMixDeltaPhiChargedAssocPtBinDEta08[bin] = new TH2F(Form("hMixDeltaPhiDeltaEta0.8ChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
3252  Form("Mixed event #Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta > 0.8", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
3253  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
3254  fhMixDeltaPhiChargedAssocPtBinDEta08[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3255  fhMixDeltaPhiChargedAssocPtBinDEta08[bin]->SetYTitle("#Delta #varphi (rad)");
3256 
3257  fhMixDeltaPhiChargedAssocPtBinDEta0[bin] = new TH2F(Form("hMixDeltaPhiDeltaEta0ChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
3258  Form("Mixed event #Delta #varphi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta = 0", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
3259  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
3260  fhMixDeltaPhiChargedAssocPtBinDEta0[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3261  fhMixDeltaPhiChargedAssocPtBinDEta0[bin]->SetYTitle("#Delta #varphi (rad)");
3262 
3263  outputContainer->Add(fhMixDeltaPhiChargedAssocPtBinDEta08[bin]);
3264  outputContainer->Add(fhMixDeltaPhiChargedAssocPtBinDEta0[bin]);
3265  }
3266  }
3267  }
3268  }
3269 
3270  if ( fFillPerSMHistograms )
3271  {
3272  Int_t totalSM = fLastModule-fFirstModule+1;
3273 
3274  fhPtTriggerPerSM = new TH2F("hPtTriggerPerSM","Selected triggers #it{p}_{T} and super-module number",
3275  nptbins,ptmin,ptmax,
3276  totalSM,fFirstModule-0.5,fLastModule+0.5);
3277  fhPtTriggerPerSM->SetYTitle("SuperModule ");
3278  fhPtTriggerPerSM->SetXTitle("#it{p}_{T} (GeV/#it{c})");
3279  outputContainer->Add(fhPtTriggerPerSM) ;
3280 
3281  for(Int_t ism = 0; ism < fNModules; ism++)
3282  {
3283  if ( ism < fFirstModule || ism > fLastModule ) continue;
3284 
3285  fhXEChargedPerSM[ism] = new TH2F
3286  (Form("hXECharged_SM%d",ism),Form("#it{x}_{#it{E}} for charged tracks, SM %d",ism),
3287  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3288  fhXEChargedPerSM[ism]->SetYTitle("#it{x}_{#it{E}}");
3289  fhXEChargedPerSM[ism]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3290  outputContainer->Add(fhXEChargedPerSM[ism]) ;
3291 
3292  fhXEUeChargedPerSM[ism] = new TH2F
3293  (Form("hXEUeCharged%s_SM%d",right.Data(),ism),Form("#it{x}_{#it{E}} for Underlying event, SM %d",ism),
3294  nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
3295  fhXEUeChargedPerSM[ism]->SetYTitle("#it{x}_{#it{E}}");
3296  fhXEUeChargedPerSM[ism]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3297  outputContainer->Add(fhXEUeChargedPerSM[ism]) ;
3298 
3299  fhDeltaPhiChargedPerSM[ism] = new TH2F
3300  (Form("hDeltaPhiCharged_SM%d",ism),
3301  Form("#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, SM %d",ism),
3302  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
3303  fhDeltaPhiChargedPerSM[ism]->SetYTitle("#Delta #varphi (rad)");
3304  fhDeltaPhiChargedPerSM[ism]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3305  outputContainer->Add(fhDeltaPhiChargedPerSM[ism]) ;
3306 
3307  fhDeltaPhiChargedPtA3GeVPerSM[ism] = new TH2F
3308  (Form("hDeltaPhiChargedPtA3GeV_SM%d",ism),
3309  Form("#varphi_{trigger} - #varphi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, SM %d",ism),
3310  nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
3311  fhDeltaPhiChargedPtA3GeVPerSM[ism]->SetYTitle("#Delta #varphi (rad)");
3312  fhDeltaPhiChargedPtA3GeVPerSM[ism]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
3313  outputContainer->Add(fhDeltaPhiChargedPtA3GeVPerSM[ism]) ;
3314  }
3315  } // Per SM
3316 
3317  return outputContainer;
3318 }
3319 
3320 //_____________________________________________________________________________________________________________________
3328 //_____________________________________________________________________________________________________________________
3329 Bool_t AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum(Int_t indexPhoton1, Int_t indexPhoton2, Int_t idetector)
3330 {
3331  if(indexPhoton1!=-1 || indexPhoton2!=-1) return kFALSE;
3332 
3333  AliDebug(1,Form("indexPhoton1 = %d, indexPhoton2 = %d", indexPhoton1, indexPhoton2));
3334 
3335  TObjArray * clusters = 0x0 ;
3336  if(idetector==kEMCAL) clusters = GetEMCALClusters() ;
3337  else clusters = GetPHOSClusters() ;
3338 
3339  for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++)
3340  {
3341  AliVCluster * photon = (AliVCluster*) (clusters->At(iclus));
3342 
3343  if(photon->GetID()==indexPhoton1) photon->GetMomentum(fDecayMom1,GetVertex(0)) ;
3344  if(photon->GetID()==indexPhoton2) photon->GetMomentum(fDecayMom2,GetVertex(0)) ;
3345 
3346  AliDebug(1,Form("Photon1 = %f, Photon2 = %f", fDecayMom1.Pt(), fDecayMom2.Pt()));
3347  } //cluster loop
3348 
3349  return kTRUE;
3350 }
3351 
3352 //________________________________________________________________________
3355 //________________________________________________________________________
3356 Int_t AliAnaParticleHadronCorrelation::GetMCTagHistogramIndex(Int_t mcTag)
3357 {
3358  if ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt) ||
3359  GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation) ||
3360  GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCISR) ) return 0;
3361  else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0) ) return 1;
3362  else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay) ) return 2;
3363  else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEta) ) return 3;
3364  else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay) ) return 4;
3365  else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton) ) return 5; // other decays
3366  else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron) ) return 6;
3367  else return 7;
3368 }
3369 
3370 //_________________________________________
3372 //_________________________________________
3373 void AliAnaParticleHadronCorrelation::Init()
3374 {
3375  if(!GetReader()->IsCTSSwitchedOn())
3376  AliFatal("STOP!: You want to use CTS tracks in analysis but not read!! \n!!Check the configuration file!!");
3377 }
3378 
3379 //____________________________________________________
3381 //____________________________________________________
3382 void AliAnaParticleHadronCorrelation::InitParameters()
3383 {
3384  SetInputAODName("Particle");
3385  SetAODObjArrayName("Hadrons");
3386  AddToHistogramsName("AnaHadronCorr_");
3387 
3388  SetPtCutRange(0.,300);
3389  fDeltaPhiMinCut = TMath::DegToRad()*120.;
3390  fDeltaPhiMaxCut = TMath::DegToRad()*240. ;
3391  fSelectIsolated = kFALSE;
3392  fMakeSeveralUE = kFALSE;
3393  fUeDeltaPhiMinCut = TMath::DegToRad()*60.;
3394  fUeDeltaPhiMaxCut = TMath::DegToRad()*120 ;
3395 
3396  fNeutralCorr = kFALSE ;
3397  fPi0Trigger = kFALSE ;
3398  fDecayTrigger = kFALSE ;
3399  fHMPIDCorrelation = kFALSE ;
3400 
3401  fMakeAbsoluteLeading = kTRUE;
3402  fMakeNearSideLeading = kFALSE;
3403 
3404  fNAssocPtBins = 9 ;
3405  fAssocPtBinLimit[0] = 0.2 ;
3406  fAssocPtBinLimit[1] = 0.5 ;
3407  fAssocPtBinLimit[2] = 1.0 ;
3408  fAssocPtBinLimit[3] = 2.0 ;
3409  fAssocPtBinLimit[4] = 3.0 ;
3410  fAssocPtBinLimit[5] = 4.0 ;
3411  fAssocPtBinLimit[6] = 5.0 ;
3412  fAssocPtBinLimit[7] = 6.0 ;
3413  fAssocPtBinLimit[8] = 7.0 ;
3414  fAssocPtBinLimit[9] = 8.0 ;
3415  fAssocPtBinLimit[10] = 9.0 ;
3416  fAssocPtBinLimit[11] = 10.0 ;
3417  fAssocPtBinLimit[12] = 12.0 ;
3418  fAssocPtBinLimit[13] = 14.0 ;
3419  fAssocPtBinLimit[14] = 16.0 ;
3420  fAssocPtBinLimit[15] = 20.0 ;
3421  fAssocPtBinLimit[16] = 30.0 ;
3422  fAssocPtBinLimit[17] = 40.0 ;
3423  fAssocPtBinLimit[18] = 50.0 ;
3424  fAssocPtBinLimit[19] = 200.0 ;
3425 
3426  fUseMixStoredInReader = kTRUE;
3427 
3428  fM02MinCut = -1 ;
3429  fM02MaxCut = -1 ;
3430 
3431  fDecayTagsM02Cut = 0.27;
3432 
3433  fSelectLeadingHadronAngle = kFALSE;
3434  fFillLeadHadOppositeHisto = kFALSE;
3435  fMinLeadHadPhi = 150*TMath::DegToRad();
3436  fMaxLeadHadPhi = 210*TMath::DegToRad();
3437 
3438  fMinLeadHadPt = 1;
3439  fMaxLeadHadPt = 100;
3440 
3441  fMCGenTypeMin = 0;
3442  fMCGenTypeMax = fgkNmcTypes-1;
3443 
3444  fNDecayBits = 1;
3445  fDecayBits[0] = AliNeutralMesonSelection::kPi0;
3446  fDecayBits[1] = AliNeutralMesonSelection::kEta;
3447  fDecayBits[2] = AliNeutralMesonSelection::kPi0RightSide;
3448  fDecayBits[3] = AliNeutralMesonSelection::kEtaRightSide;
3449  fDecayBits[4] = AliNeutralMesonSelection::kEtaLeftSide ;
3450  fDecayBits[5] = AliNeutralMesonSelection::kEtaBothSides;
3451  fDecayBits[6] = AliNeutralMesonSelection::kPi0LeftSide ; // Leave it last since likely not used
3452  fDecayBits[7] = AliNeutralMesonSelection::kPi0BothSides; // Leave it last since likely not used
3453 
3454  fNBkgBin = 11;
3455  fBkgBinLimit[ 0] = 00.0; fBkgBinLimit[ 1] = 00.2; fBkgBinLimit[ 2] = 00.3; fBkgBinLimit[ 3] = 00.4; fBkgBinLimit[ 4] = 00.5;
3456  fBkgBinLimit[ 5] = 01.0; fBkgBinLimit[ 6] = 01.5; fBkgBinLimit[ 7] = 02.0; fBkgBinLimit[ 8] = 03.0; fBkgBinLimit[ 9] = 05.0;
3457  fBkgBinLimit[10] = 10.0; fBkgBinLimit[11] = 100.;
3458  for(Int_t ibin = fNBkgBin+1; ibin < 20; ibin++)
3459  {
3460  fBkgBinLimit[ibin] = 00.0;
3461  }
3462 
3463 }
3464 
3465 //________________________________________________________________________________________________________
3470 //________________________________________________________________________________________________________
3471 void AliAnaParticleHadronCorrelation::InvMassHisto(AliCaloTrackParticleCorrelation * trigger, Int_t mcIndex)
3472 {
3473  Int_t idTrig = trigger->GetCaloLabel(0);
3474  Float_t ptTrig = trigger->Pt();
3475  Float_t tofTrig = trigger->GetTime();
3476 
3477  fMomentum = *(trigger->GetMomentum());
3478 
3479  // Loop on second cluster to be associated to trigger particle to have the invariant mass
3480  Int_t nphoton = GetInputAODBranch()->GetEntriesFast();
3481 
3482  for(Int_t iphoton = 0; iphoton < nphoton; iphoton++)
3483  {
3484  AliCaloTrackParticleCorrelation * photon1 = (AliCaloTrackParticleCorrelation*) (GetInputAODBranch()->At(iphoton));
3485 
3486  if(idTrig == photon1->GetCaloLabel(0)) continue; //Do not the invariant mass for the same particle
3487 
3488  fMomentumIM = *(photon1->GetMomentum());
3489 
3490  // Select secondary photon with proper invariant mass
3491  if(fM02MaxCut > 0 && fM02MinCut > 0) //clID1 > 0 && clID2 < 0 &&
3492  {
3493  Float_t m02 = photon1->GetM02();
3494 
3495  if(m02 > fM02MaxCut || m02 < fM02MinCut) continue ;
3496  }
3497 
3498  // Select clusters with good time window difference
3499  Double_t tdiff = tofTrig - photon1->GetTime();
3500  if(TMath::Abs(tdiff) > GetPairTimeCut()) continue;
3501 
3502  // Add momenta of trigger and possible secondary decay
3503  fMomentumIM += fMomentum;
3504 
3505  Double_t mass = fMomentumIM.M();
3506  // Double_t epair = fMomentumIM.E();
3507  // Float_t ptpair = fMomentumIM.Pt();
3508 
3509  // Mass of all pairs
3510  fhMassPtTrigger->Fill(ptTrig, mass, GetEventWeight());
3511  if ( IsDataMC() ) fhMCMassPtTrigger[mcIndex]->Fill(ptTrig, mass, GetEventWeight());
3512  }
3513 }
3514 
3515 //_________________________________________________________________________
3518 //_________________________________________________________________________
3519 Bool_t AliAnaParticleHadronCorrelation::IsTriggerTheEventLeadingParticle()
3520 {
3521  Double_t ptTrig = GetMinPt();
3522  Double_t phiTrig = 0 ;
3523  fLeadingTriggerIndex =-1 ;
3524  Int_t index =-1 ;
3525  AliCaloTrackParticleCorrelation* pLeading = 0;
3526 
3527  // Loop on stored AOD particles, find leading trigger on the selected list, with at least min pT.
3528 
3529  for(Int_t iaod = 0; iaod < GetInputAODBranch()->GetEntriesFast() ; iaod++)
3530  {
3531  AliCaloTrackParticleCorrelation* particle = (AliCaloTrackParticleCorrelation*) (GetInputAODBranch()->At(iaod));
3532  particle->SetLeadingParticle(kFALSE); // set it later
3533 
3534  // Vertex cut in case of mixing
3535  Int_t check = CheckMixedEventVertex(particle->GetCaloLabel(0), particle->GetTrackLabel(0));
3536  if(check == 0) continue;
3537  if(check == -1) return kFALSE; // not sure if it is correct.
3538 
3539  // Find the leading particles with highest momentum
3540  if (particle->Pt() > ptTrig)
3541  {
3542  ptTrig = particle->Pt() ;
3543  phiTrig = particle->Phi();
3544  index = iaod ;
3545  pLeading = particle ;
3546  }
3547  }// Finish search of leading trigger particle on the AOD branch.
3548 
3549  if(index < 0) return kFALSE;
3550 
3551  //printf("AOD leading pT %2.2f, ID %d\n", pLeading->Pt(),pLeading->GetCaloLabel(0));
3552 
3553  if(phiTrig < 0 ) phiTrig += TMath::TwoPi();
3554 
3555  // Compare if it is the leading of all tracks
3556 
3557  for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
3558  {
3559  AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
3560 
3561  // In case of isolation of single tracks or conversion photon (2 tracks) or pi0 (4 tracks),
3562  // do not count the candidate or the daughters of the candidate
3563  // in the isolation conte
3564  if ( pLeading->GetDetectorTag() == kCTS ) // make sure conversions are tagged as kCTS!!!
3565  {
3566  Int_t trackID = GetReader()->GetTrackID(track) ; // needed instead of track->GetID() since AOD needs some manipulations
3567  Bool_t contained = kFALSE;
3568 
3569  for(Int_t i = 0; i < 4; i++)
3570  {
3571  if( trackID == pLeading->GetTrackLabel(i) ) contained = kTRUE;
3572  }
3573 
3574  if ( contained ) continue ;
3575  }
3576 
3577  fTrackVector.SetXYZ(track->Px(),track->Py(),track->Pz());
3578  Float_t pt = fTrackVector.Pt();
3579  Float_t phi = fTrackVector.Phi() ;
3580  if(phi < 0) phi+=TMath::TwoPi();
3581 
3582  //jump out this event if near side associated particle pt larger than trigger
3583  if (fMakeNearSideLeading)
3584  {
3585  Float_t deltaPhi = phiTrig-phi;
3586  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
3587  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
3588 
3589  if(pt > ptTrig && deltaPhi < TMath::PiOver2()) return kFALSE;
3590  }
3591  //jump out this event if there is any other particle with pt larger than trigger
3592  else
3593  {
3594  if(pt > ptTrig) return kFALSE ;
3595  }
3596  }// track loop
3597 
3598  // Compare if it is leading of all calorimeter clusters
3599 
3600  if(fCheckLeadingWithNeutralClusters)
3601  {
3602  // Select the calorimeter cluster list
3603  TObjArray * nePl = 0x0;
3604  if (pLeading->GetDetectorTag() == kPHOS )
3605  nePl = GetPHOSClusters();
3606  else
3607  nePl = GetEMCALClusters();
3608 
3609  if(!nePl) return kTRUE; // Do the selection just with the tracks if no calorimeter is available.
3610 
3611  for(Int_t ipr = 0;ipr < nePl->GetEntriesFast() ; ipr ++ )
3612  {
3613  AliVCluster * cluster = (AliVCluster *) (nePl->At(ipr)) ;
3614 
3615  if(cluster->GetID() == pLeading->GetCaloLabel(0) || cluster->GetID() == pLeading->GetCaloLabel(1) ) continue ;
3616 
3617  cluster->GetMomentum(fMomentum,GetVertex(0));
3618 
3619  Float_t pt = fMomentum.Pt();
3620  Float_t phi = fMomentum.Phi() ;
3621  if(phi < 0) phi+=TMath::TwoPi();
3622 
3623  if(IsTrackMatched(cluster,GetReader()->GetInputEvent())) continue ; // avoid charged clusters, already covered by tracks, or cluster merging with track.
3624 
3625  //jump out this event if near side associated particle pt larger than trigger
3626  // not really needed for calorimeter, unless DCal is included
3627  if (fMakeNearSideLeading)
3628  {
3629  Float_t deltaPhi = phiTrig-phi;
3630  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
3631  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
3632 
3633  if(pt > ptTrig && deltaPhi < TMath::PiOver2()) return kFALSE ;
3634  }
3635  //jump out this event if there is any other particle with pt larger than trigger
3636  else
3637  {
3638  if(pt > ptTrig) return kFALSE ;
3639  }
3640  }// cluster loop
3641  } // check neutral clusters
3642 
3643  fLeadingTriggerIndex = index ;
3644  pLeading->SetLeadingParticle(kTRUE);
3645 
3646  AliDebug(1,Form("\t particle AOD with index %d is leading with pT %2.2f", fLeadingTriggerIndex, pLeading->Pt()));
3647 
3648  return kTRUE;
3649 }
3650 
3651 //_________________________________________________________________
3665 //_________________________________________________________________
3666 void AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms()
3667 {
3668  if(!GetInputAODBranch())
3669  {
3670  AliFatal(Form("No input particles in AOD with name branch < %s >, STOP",GetInputAODName().Data()));
3671  return ; // coverity
3672  }
3673 
3674  Int_t naod = GetInputAODBranch()->GetEntriesFast();
3675  if( naod == 0 )
3676  {
3677  AliDebug(1,"No particle AOD found!");
3678  return ; // no trigger particles found.
3679  }
3680 
3681  AliDebug(1,Form("Begin hadron correlation analysis, fill histograms"));
3682  AliDebug(1,Form("n particle branch aod entries %d", naod));
3683  AliDebug(1,Form("In CTS aod entries %d",GetCTSTracks()->GetEntriesFast()));
3684 
3685  //------------------------------------------------------
3686  // Find leading trigger if analysis request only leading,
3687  // if there is no leading trigger, then skip the event
3688 
3689  Int_t iaod = 0 ;
3690  if( fMakeAbsoluteLeading || fMakeNearSideLeading )
3691  {
3692  Bool_t leading = IsTriggerTheEventLeadingParticle();
3693 
3694  AliDebug(1,Form("AOD Leading trigger? %d, with index %d",leading,fLeadingTriggerIndex));
3695 
3696  if(!leading)
3697  {
3698  AliDebug(1,"Leading was requested and not found");
3699  return ;
3700  }
3701  else
3702  {
3703  // Select only the leading in the trigger AOD loop
3704  naod = fLeadingTriggerIndex+1 ;
3705  iaod = fLeadingTriggerIndex ;
3706  }
3707  }
3708 
3709  //------------------------------------------------------
3710  // Get event multiplicity and bins
3711 
3712  Float_t cen = GetEventCentrality();
3713  Float_t ep = GetEventPlaneAngle();
3714  if(IsHighMultiplicityAnalysisOn()) fhTriggerEventPlaneCentrality->Fill(cen, ep, GetEventWeight());
3715 
3716  Int_t mixEventBin = GetEventMixBin();
3717  Int_t vzbin = GetEventVzBin();
3718 
3719  //------------------------------------------------------
3720  // Loop on trigger AOD
3721 
3722  for( iaod = 0; iaod < naod; iaod++ )
3723  {
3724  AliCaloTrackParticleCorrelation* particle = (AliCaloTrackParticleCorrelation*) (GetInputAODBranch()->At(iaod));
3725 
3726  //
3727  // Trigger particle selection criteria:
3728  //
3729  Float_t pt = particle->Pt();
3730 
3731  if(pt < GetMinPt() || pt > GetMaxPt() ) continue ;
3732 
3733  fhPtTriggerInput->Fill(pt, GetEventWeight());
3734 
3735  //
3736  // check if it was a calorimeter cluster
3737  // and if the shower shape cut was requested apply it.
3738  // Not needed if already done at the particle identification level,
3739  // but for isolation studies, it is preferred not to remove so we do it here
3740  //
3741 
3742  AliDebug(1,Form("%s Trigger : min %f, max %f, det %d",
3743  GetInputAODName().Data(),fM02MinCut,fM02MaxCut,particle->GetDetectorTag()));
3744 
3745  if(fM02MaxCut > 0 && fM02MinCut > 0) //clID1 > 0 && clID2 < 0 &&
3746  {
3747  // Int_t iclus = -1;
3748  // TObjArray* clusters = 0x0;
3749  // if (particle->GetDetectorTag() == kEMCAL) clusters = GetEMCALClusters();
3750  // else if(particle->GetDetectorTag() == kPHOS ) clusters = GetPHOSClusters();
3751  //
3752  // if(clusters)
3753  // {
3754  // AliVCluster *cluster = FindCluster(clusters,clID1,iclus);
3755  // Float_t m02 = cluster->GetM02();
3756  // if(m02 > fM02MaxCut || m02 < fM02MinCut) continue ;
3757  // }
3758 
3759  Float_t m02 = particle->GetM02();
3760 
3761  if(m02 > fM02MaxCut || m02 < fM02MinCut) continue ;
3762 
3763  fhPtTriggerSSCut->Fill(pt, GetEventWeight());
3764 
3765  AliDebug(1,"Pass the shower shape cut");
3766  }
3767 
3768  //
3769  // Check if the particle is isolated or if we want to take the isolation into account
3770  // This bool is set in AliAnaParticleIsolation
3771  //
3772  if(OnlyIsolated())
3773  {
3774  if( !particle->IsIsolated() ) continue;
3775 
3776  fhPtTriggerIsoCut->Fill(pt, GetEventWeight());
3777 
3778  AliDebug(1,"Pass the isolation cut");
3779  }
3780 
3781  //
3782  // Check if trigger is in fiducial region
3783  //
3784  if(IsFiducialCutOn())
3785  {
3786  Bool_t in = GetFiducialCut()->IsInFiducialCut(particle->Eta(),particle->Phi(),particle->GetDetectorTag()) ;
3787 
3788  if(! in ) continue ;
3789 
3790  AliDebug(1,"Pass the fiducial cut");
3791  }
3792 
3793  fhPtTriggerFidCut->Fill(pt, GetEventWeight());
3794 
3795  //---------------------------------------
3796  // Make correlation
3797 
3798  // Find the leading hadron in the opposite hemisphere to the trigger
3799  // and accept the trigger if leading is in defined window.
3800  Bool_t okLeadHad = kTRUE;
3801  if(fSelectLeadingHadronAngle || fFillLeadHadOppositeHisto)
3802  {
3803  okLeadHad = FindLeadingOppositeHadronInWindow(particle);
3804  if(!okLeadHad && fSelectLeadingHadronAngle) continue;
3805  }
3806 
3807  //
3808  // Charged particles correlation
3809  //
3810  MakeChargedCorrelation(particle);
3811 
3812  // MC
3813  Int_t mcIndex = -1;
3814  Int_t mcTag = particle->GetTag();
3815  Bool_t lostDecayPair = kFALSE;
3816  if(IsDataMC())
3817  {
3818  mcIndex = GetMCTagHistogramIndex(mcTag);
3819  lostDecayPair = GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost);
3820  MakeMCChargedCorrelation(particle->GetLabel(), mcIndex,lostDecayPair);
3821  }
3822 
3823  // Do own mixed event with charged,
3824  // add event and remove previous or fill the mixed histograms
3825  if(DoOwnMix())
3826  MakeChargedMixCorrelation(particle);
3827 
3828  //
3829  // Neutral particles correlation
3830  //
3831  if(fNeutralCorr)
3832  MakeNeutralCorrelation(particle);
3833 
3834  //----------------------------------------------------------------
3835  // Fill trigger pT related histograms if not absolute leading
3836 
3837  //
3838  // pT of the trigger, vs trigger origin if MC
3839  //
3840  fhPtTrigger->Fill(pt, GetEventWeight());
3841 
3842  if ( fFillPerSMHistograms ) fhPtTriggerPerSM->Fill(pt, particle->GetSModNumber(),GetEventWeight());
3843 
3844 
3845  if(IsDataMC() && mcIndex >=0 && mcIndex < fgkNmcTypes)
3846  {
3847  fhPtTriggerMC[mcIndex]->Fill(pt, GetEventWeight());
3848  if( lostDecayPair )
3849  {
3850  // check index in GetMCTagIndexHistogram
3851  if ( mcIndex == 2 ) fhPtTriggerMC[8]->Fill(pt, GetEventWeight()); // pi0 decay
3852  else if( mcIndex == 4 ) fhPtTriggerMC[9]->Fill(pt, GetEventWeight()); // eta decay
3853  }
3854  }
3855 
3856  //
3857  // pT lead cone, pT sum cone background bins distribution
3858  //
3859  if(fFillBkgBinsHisto)
3860  {
3861  Float_t m02 = particle->GetM02();
3862  Float_t pTLeadTrackInCone = 0;
3863  Float_t pTSumTrackInCone = 0;
3864  Float_t pTLeadClusterInCone = 0;
3865  Float_t pTSumClusterInCone = 0;
3866 
3867  pTLeadTrackInCone = particle->GetChargedLeadPtInCone();
3868  pTLeadClusterInCone = particle->GetNeutralLeadPtInCone();
3869 
3870  pTSumTrackInCone = particle->GetChargedPtSumInCone();
3871  pTSumClusterInCone = particle->GetNeutralPtSumInCone();
3872 
3873  Float_t pTLeadInCone = pTLeadTrackInCone;
3874  if(pTLeadClusterInCone > pTLeadInCone) pTLeadInCone = pTLeadClusterInCone;
3875  Float_t pTSumInCone = pTSumTrackInCone + pTSumClusterInCone;
3876 
3877  // Get the background bin for this cone and trigger
3878  Int_t pTSumBin = -1;
3879  Int_t pTLeadBin = -1;
3880 
3881  for(Int_t ibin = 0; ibin < fNBkgBin; ibin++)
3882  {
3883  if( pTSumInCone >= fBkgBinLimit[ibin] && pTSumInCone < fBkgBinLimit[ibin+1]) pTSumBin = ibin;
3884  if( pTLeadInCone >= fBkgBinLimit[ibin] && pTLeadInCone < fBkgBinLimit[ibin+1]) pTLeadBin = ibin;
3885  }
3886  if(pTSumBin > 0)
3887  {
3888  fhPtSumInConeBin[pTSumBin]->Fill(pt, GetEventWeight());
3889  }
3890 
3891  if(pTLeadBin > 0)
3892  {
3893  fhPtLeadInConeBin[pTLeadBin]->Fill(pt, GetEventWeight());
3894  }
3895 
3896  // Check if it was a decay
3897  if( fFillTaggedDecayHistograms && m02 < fDecayTagsM02Cut )
3898  {
3899  Int_t decayTag = particle->DecayTag();
3900  if(decayTag < 0) decayTag = 0;
3901 
3902  for(Int_t ibit = 0; ibit < fNDecayBits; ibit++)
3903  {
3904  if(GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[ibit]))
3905  {
3906  Int_t pTLeadBinDecay = pTLeadBin+ibit*fNBkgBin;
3907  Int_t pTSumBinDecay = pTSumBin+ibit*fNBkgBin;
3908  if( pTLeadBin >=0 ) fhPtLeadConeBinDecay[pTLeadBinDecay]->Fill(pt, GetEventWeight());
3909  if( pTSumBin >=0 ) fhSumPtConeBinDecay [pTSumBinDecay] ->Fill(pt, GetEventWeight());
3910  }
3911  }
3912  }
3913 
3914  if(IsDataMC())
3915  {
3916  Int_t pTLeadBinMC = pTLeadBin+mcIndex*fNBkgBin;
3917  Int_t pTSumBinMC = pTSumBin+mcIndex*fNBkgBin;
3918 
3919  if( pTLeadBin >=0 )
3920  {
3921  fhPtLeadConeBinMC[pTLeadBinMC]->Fill(pt, GetEventWeight());
3922  }
3923 
3924  if( pTSumBin >=0 )
3925  {
3926  fhSumPtConeBinMC [pTSumBinMC]->Fill(pt, GetEventWeight());
3927  }
3928 
3929  if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPhoton))
3930  {
3931  pTLeadBinMC = pTLeadBin+kmcPhoton*fNBkgBin;
3932  pTSumBinMC = pTSumBin+kmcPhoton*fNBkgBin;
3933  if( pTLeadBin >=0 )
3934  {
3935  fhPtLeadConeBinMC[pTLeadBinMC]->Fill(pt, GetEventWeight());
3936  }
3937 
3938  if( pTSumBin >=0 )
3939  {
3940  fhSumPtConeBinMC [ pTSumBinMC]->Fill(pt, GetEventWeight());
3941  }
3942  }
3943 
3944  // Check if decay and if pair is lost
3945  if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) )
3946  {
3947  if ( mcIndex == kmcPi0Decay )
3948  {
3949  pTLeadBinMC = pTLeadBin+kmcPi0DecayLostPair*fNBkgBin;
3950  pTSumBinMC = pTSumBin+kmcPi0DecayLostPair*fNBkgBin;
3951  }
3952  else if(mcIndex == kmcEtaDecay)
3953  {
3954  pTLeadBinMC = pTLeadBin+kmcEtaDecayLostPair*fNBkgBin;
3955  pTSumBinMC = pTSumBin+kmcEtaDecayLostPair*fNBkgBin;
3956  }
3957  else
3958  AliFatal(Form("Lost decay Bit assigned to bad case, mcIndex %d",mcIndex));
3959 
3960  if( pTLeadBin >=0 )
3961  {
3962  fhPtLeadConeBinMC[pTLeadBinMC]->Fill(pt, GetEventWeight());
3963  }
3964 
3965  if( pTSumBin >=0 )
3966  {
3967  fhSumPtConeBinMC [ pTSumBinMC]->Fill(pt, GetEventWeight());
3968  }
3969 
3970  }
3971 
3972  }
3973  }
3974 
3975  // Invariant mass of trigger particle
3976  if(fFillInvMassHisto) InvMassHisto(particle,mcIndex);
3977 
3978  if(fDecayTrigger)
3979  {
3980  Int_t decayTag = particle->DecayTag();
3981  if(decayTag < 0) decayTag = 0;
3982 
3983  for(Int_t ibit = 0; ibit<fNDecayBits; ibit++)
3984  {
3985  if(GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[ibit]))
3986  {
3987  fhPtDecayTrigger[ibit]->Fill(pt, GetEventWeight());
3988 
3989  if(IsDataMC() && mcIndex >=0 && mcIndex < fgkNmcTypes)
3990  {
3991  fhPtDecayTriggerMC[ibit][mcIndex]->Fill(pt, GetEventWeight());
3992  if( lostDecayPair )
3993  {
3994  // check index in GetMCTagIndexHistogram
3995  if( mcIndex == 2 ) fhPtDecayTriggerMC[ibit][8]->Fill(pt, GetEventWeight()); // pi0 decay
3996  if( mcIndex == 4 ) fhPtDecayTriggerMC[ibit][9]->Fill(pt, GetEventWeight()); // eta decay
3997  }
3998  }
3999  }// check bit
4000  }// bit loop
4001  }
4002 
4003  //
4004  // Acceptance of the trigger
4005  //
4006  Float_t phi = particle->Phi();
4007  if( phi < 0 ) phi+=TMath::TwoPi();
4008  fhPhiTrigger->Fill(pt, phi, GetEventWeight());
4009 
4010  fhEtaTrigger->Fill(pt, particle->Eta(), GetEventWeight());
4011  //printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Trigger particle : pt %f, eta %f, phi %f\n",particle->Pt(),particle->Eta(),phi);
4012 
4013  //----------------------------------
4014  // Trigger particle pT vs event bins
4015 
4016  fhPtTriggerBin->Fill(pt, mixEventBin, GetEventWeight());
4017  if(fCorrelVzBin)
4018  fhPtTriggerVzBin->Fill(pt, vzbin, GetEventWeight());
4019 
4020  if(IsHighMultiplicityAnalysisOn())
4021  {
4022  fhPtTriggerCentrality->Fill(pt, cen, GetEventWeight());
4023  fhPtTriggerEventPlane->Fill(pt, ep , GetEventWeight());
4024  }
4025 
4026  //----------------------------------
4027  // Trigger particle pT vs pile-up
4028 
4029  if(IsPileUpAnalysisOn())
4030  {
4031  Int_t vtxBC = GetReader()->GetVertexBC();
4032  if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA) fhPtTriggerVtxBC0->Fill(pt, GetEventWeight());
4033 
4034  if(GetReader()->IsPileUpFromSPD()) fhPtTriggerPileUp[0]->Fill(pt, GetEventWeight());
4035  if(GetReader()->IsPileUpFromEMCal()) fhPtTriggerPileUp[1]->Fill(pt, GetEventWeight());
4036  if(GetReader()->IsPileUpFromSPDOrEMCal()) fhPtTriggerPileUp[2]->Fill(pt, GetEventWeight());
4037  if(GetReader()->IsPileUpFromSPDAndEMCal()) fhPtTriggerPileUp[3]->Fill(pt, GetEventWeight());
4038  if(GetReader()->IsPileUpFromSPDAndNotEMCal()) fhPtTriggerPileUp[4]->Fill(pt, GetEventWeight());
4039  if(GetReader()->IsPileUpFromEMCalAndNotSPD()) fhPtTriggerPileUp[5]->Fill(pt, GetEventWeight());
4040  if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtTriggerPileUp[6]->Fill(pt, GetEventWeight());
4041  }
4042  } // AOD trigger loop
4043 
4044  // Re-init for next event
4045  fLeadingTriggerIndex = -1;
4046 
4047  AliDebug(1,"End fill histograms");
4048 }
4049 
4050 //_______________________________________________________________________________________________________
4056 //_______________________________________________________________________________________________________
4057 void AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliCaloTrackParticleCorrelation *aodParticle)
4058 {
4059  AliDebug(1,"Make trigger particle - charged hadron correlation");
4060 
4061  Float_t phiTrig = aodParticle->Phi();
4062  Float_t etaTrig = aodParticle->Eta();
4063  Float_t ptTrig = aodParticle->Pt();
4064  Int_t mcTag = aodParticle->GetTag();
4065  Double_t bz = GetReader()->GetInputEvent()->GetMagneticField();
4066  Int_t sm = aodParticle->GetSModNumber();
4067 
4068  Int_t decayTag = 0;
4069  if(fDecayTrigger)
4070  {
4071  //decay = aodParticle->IsTagged();
4072  decayTag = aodParticle->DecayTag();
4073  if(decayTag < 0) decayTag = 0;
4074  // printf("Correlation: pT %2.2f, BTag %d, Tagged %d\n",ptTrig, decayTag, aodParticle->IsTagged());
4075  // printf("\t check bit Pi0 %d, Eta %d, Pi0Side %d, EtaSide %d\n",
4076  // GetNeutralMesonSelection()->CheckDecayBit(decayTag,AliNeutralMesonSelection::kPi0),
4077  // GetNeutralMesonSelection()->CheckDecayBit(decayTag,AliNeutralMesonSelection::kEta),
4078  // GetNeutralMesonSelection()->CheckDecayBit(decayTag,AliNeutralMesonSelection::kPi0Side),
4079  // GetNeutralMesonSelection()->CheckDecayBit(decayTag,AliNeutralMesonSelection::kEtaSide));
4080  }
4081 
4082  Float_t pt = -100. ;
4083  Float_t phi = -100. ;
4084  Float_t eta = -100. ;
4085  Float_t deltaPhi = -100. ;
4086 
4087  TObjArray * reftracks = 0x0;
4088  Int_t nrefs = 0;
4089 
4090  // Mixed event settings
4091  Int_t evtIndex11 = -1 ; // cluster trigger or pi0 trigger
4092  Int_t evtIndex12 = -1 ; // pi0 trigger
4093  Int_t evtIndex13 = -1 ; // charged trigger
4094 
4095  if (GetMixedEvent())
4096  {
4097  evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
4098  evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
4099  evtIndex13 = GetMixedEvent()->EventIndex(aodParticle->GetTrackLabel(0)) ;
4100  }
4101 
4102  // Track multiplicity or cent bin
4103  Int_t cenbin = 0;
4104  if(IsHighMultiplicityAnalysisOn()) cenbin = GetEventCentralityBin();
4105 
4106  //
4107  // In case of pi0/eta trigger, we may want to check their decay correlation,
4108  // get their decay children
4109  //
4110 
4111  Bool_t decayFound = kFALSE;
4112  if( fPi0Trigger )
4113  {
4114  decayFound = GetDecayPhotonMomentum(aodParticle->GetCaloLabel(0),aodParticle->GetCaloLabel(1),aodParticle->GetDetectorTag());
4115  if(decayFound)
4116  {
4117  fhPtPi0DecayRatio->Fill(ptTrig, fDecayMom1.Pt()/ptTrig, GetEventWeight());
4118  fhPtPi0DecayRatio->Fill(ptTrig, fDecayMom2.Pt()/ptTrig, GetEventWeight());
4119  }
4120  }
4121 
4122  //-----------------------------------------------------------------------
4123  // Track loop, select tracks with good pt, phi and fill AODs or histograms
4124  //-----------------------------------------------------------------------
4125 
4126  for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
4127  {
4128  AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
4129 
4130  fTrackVector.SetXYZ(track->Px(),track->Py(),track->Pz());
4131  pt = fTrackVector.Pt();
4132  eta = fTrackVector.Eta();
4133  phi = fTrackVector.Phi() ;
4134  if(phi < 0) phi+=TMath::TwoPi();
4135 
4136  //Select only hadrons in pt range
4137  if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
4138 
4139  // In case of isolation of single tracks or conversion photon (2 tracks) or pi0 (4 tracks),
4140  // do not count the candidate or the daughters of the candidate
4141  // in the isolation conte
4142  if ( aodParticle->GetDetectorTag() == kCTS ) // make sure conversions are tagged as kCTS!!!
4143  {
4144  Int_t trackID = GetReader()->GetTrackID(track) ; // needed instead of track->GetID() since AOD needs some manipulations
4145  Bool_t contained = kFALSE;
4146 
4147  for(Int_t i = 0; i < 4; i++)
4148  {
4149  if( trackID == aodParticle->GetTrackLabel(i) ) contained = kTRUE;
4150  }
4151 
4152  if ( contained ) continue ;
4153  }
4154 
4155  //Only for mixed event frame
4156  Int_t evtIndex2 = 0 ;
4157  if (GetMixedEvent())
4158  {
4159  evtIndex2 = GetMixedEvent()->EventIndex(GetReader()->GetTrackID(track)) ;
4160  if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 || evtIndex13 == evtIndex2 ) // photon and track from different events
4161  continue ;
4162  //vertex cut
4163  if (TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut())
4164  continue;
4165  }
4166 
4167  AliDebug(2,Form("Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f",pt,phi,eta));
4168 
4169  // ------------------------------
4170  // Track type bin or bits setting
4171  //
4172 
4173  //
4174  // * Set the pt associated bin for the defined bins *
4175  //
4176  Int_t assocBin = -1;
4177  for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
4178  {
4179  if(pt > fAssocPtBinLimit[i] && pt < fAssocPtBinLimit[i+1]) assocBin= i;
4180  }
4181 
4182  //
4183  // * Assign to the histogram array a bin corresponding
4184  // to a combination of pTa and vz bins *
4185  //
4186  Int_t nz = 1;
4187  Int_t vz = 0;
4188 
4189  if(fCorrelVzBin)
4190  {
4191  nz = GetNZvertBin();
4192  vz = GetEventVzBin();
4193  }
4194 
4195  Int_t bin = assocBin*nz+vz;
4196 
4197  //printf("assoc Bin = %d, vZ bin = %d, bin = %d \n", assocBin,GetEventVzBin(),bin);
4198 
4199  //
4200  // * Get the status of the TOF bit *
4201  //
4202  ULong_t status = track->GetStatus();
4203  Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
4204  //Double32_t tof = track->GetTOFsignal()*1e-3;
4205  Int_t trackBC = track->GetTOFBunchCrossing(bz);
4206 
4207  Int_t outTOF = -1;
4208  if (okTOF && trackBC!=0) outTOF = 1;
4209  else if(okTOF && trackBC==0) outTOF = 0;
4210 
4211  //----------------
4212  // Fill Histograms
4213 
4214  //
4215  // Azimuthal Angle histograms
4216  //
4217 
4218  deltaPhi = phiTrig-phi;
4219 
4220  //
4221  // Calculate deltaPhi shift so that for the particles on the opposite side
4222  // it is defined between 90 and 270 degrees
4223  // Shift [-360,-90] to [0, 270]
4224  // and [270,360] to [-90,0]
4225  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
4226  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
4227 
4228  FillChargedAngularCorrelationHistograms(pt, ptTrig, bin, phi, phiTrig, deltaPhi,
4229  eta, etaTrig, sm, decayTag, track->GetHMPIDsignal(),
4230  outTOF, cenbin, mcTag);
4231 
4232  //
4233  // Imbalance zT/xE/pOut histograms
4234  //
4235 
4236  AliESDtrack * esdTrack = dynamic_cast<AliESDtrack*>(track);
4237 
4238  //
4239  // Delta phi cut for momentum imbalance correlation
4240  //
4241  if ( (deltaPhi > fDeltaPhiMinCut) && (deltaPhi < fDeltaPhiMaxCut) )
4242  {
4243  // Check track resolution
4244  if ( esdTrack )
4245  {
4246  fhTrackResolution->Fill(pt, TMath::Sqrt(esdTrack->GetCovariance()[14])*pt,
4247  GetEventWeight());
4248  }
4249 
4250  FillChargedMomentumImbalanceHistograms(ptTrig, pt, deltaPhi, sm, cenbin, track->Charge(),
4251  assocBin, decayTag, outTOF, mcTag);
4252  }
4253 
4254  //
4255  // Underlying event, right side, default case
4256  //
4257  if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )
4258  {
4259  // Check track resolution
4260  if ( esdTrack )
4261  {
4262  fhTrackResolutionUE->Fill(pt, TMath::Sqrt(esdTrack->GetCovariance()[14])*pt,
4263  GetEventWeight());
4264  }
4265 
4266  FillChargedUnderlyingEventHistograms(ptTrig, pt, deltaPhi, sm, cenbin, outTOF,mcTag);
4267  }
4268 
4269  //
4270  // Several UE calculation, in different perpendicular regions, up to 6:
4271  // left, right, upper-left, lower left, upper-right, lower-right
4272  //
4273  if(fMakeSeveralUE)
4274  FillChargedUnderlyingEventSidesHistograms(ptTrig,pt,deltaPhi,mcTag);
4275 
4276  //
4277  if(fPi0Trigger && decayFound)
4278  FillDecayPhotonCorrelationHistograms(pt, phi, kTRUE) ;
4279 
4280  //
4281  // Add track reference to array
4282  //
4283  if(fFillAODWithReferences)
4284  {
4285  nrefs++;
4286  if(nrefs==1)
4287  {
4288  reftracks = new TObjArray(0);
4289  TString trackname = Form("%sTracks", GetAODObjArrayName().Data());
4290  reftracks->SetName(trackname.Data());
4291  reftracks->SetOwner(kFALSE);
4292  }
4293 
4294  reftracks->Add(track);
4295  }// reference track to AOD
4296  }// track loop
4297 
4298  //Fill AOD with reference tracks, if not filling histograms
4299  if(fFillAODWithReferences && reftracks)
4300  {
4301  aodParticle->AddObjArray(reftracks);
4302  }
4303 }
4304 
4305 //_________________________________________________________________________________________________________
4307 //_________________________________________________________________________________________________________
4308 void AliAnaParticleHadronCorrelation::MakeChargedMixCorrelation(AliCaloTrackParticleCorrelation *aodParticle)
4309 {
4310  AliDebug(1,Form("Make trigger particle - charged hadron mixed event correlation"));
4311 
4312  if(GetMixedEvent()) return; // This is not the mixed event from general mixing frame
4313 
4314  // Get the event with similar caracteristics
4315  //printf("MakeChargedMixCorrelation for %s\n",GetInputAODName().Data());
4316 
4317  AliAnalysisManager * manager = AliAnalysisManager::GetAnalysisManager();
4318 
4319  AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
4320 
4321  if(!inputHandler) return;
4322 
4323  if(!(inputHandler->IsEventSelected( ) & GetReader()->GetEventTriggerMask())) return;
4324 
4325  // Get the pool, check if it exits
4326  Int_t eventBin = GetEventMixBin();
4327 
4328  //Check that the bin exists, if not (bad determination of RP, ntrality or vz bin) do nothing
4329  if(eventBin < 0) return;
4330 
4331  fhEventBin->Fill(eventBin, GetEventWeight());
4332 
4333  // Get neutral clusters pool?
4334  Bool_t isoCase = OnlyIsolated() && (GetIsolationCut()->GetParticleTypeInCone() != AliIsolationCut::kOnlyCharged);
4335  Bool_t neutralMix = fFillNeutralEventMixPool || isoCase ;
4336 
4337  TList * pool = 0;
4338  TList * poolCalo = 0;
4339  if(fUseMixStoredInReader)
4340  {
4341  pool = GetReader()->GetListWithMixedEventsForTracks(eventBin);
4342  if(neutralMix) poolCalo = GetReader()->GetListWithMixedEventsForCalo (eventBin);
4343  }
4344  else
4345  {
4346  pool = fListMixTrackEvents[eventBin];
4347  if(neutralMix) poolCalo = fListMixCaloEvents [eventBin];
4348  }
4349 
4350  if(!pool) return ;
4351 
4352  if( neutralMix && !poolCalo )
4353  AliWarning("Careful, cluster pool not available");
4354 
4355  Double_t ptTrig = aodParticle->Pt();
4356  Double_t etaTrig = aodParticle->Eta();
4357  Double_t phiTrig = aodParticle->Phi();
4358  if(phiTrig < 0.) phiTrig+=TMath::TwoPi();
4359 
4360  AliDebug(1,Form("Pool bin %d size %d, trigger trigger pt=%f, phi=%f, eta=%f",
4361  eventBin,pool->GetSize(), ptTrig,phiTrig,etaTrig));
4362 
4363  Double_t ptAssoc = -999.;
4364  Double_t phiAssoc = -999.;
4365  Double_t etaAssoc = -999.;
4366  Double_t deltaPhi = -999.;
4367  Double_t deltaEta = -999.;
4368  Double_t xE = -999.;
4369 
4370  // Start from first event in pool except if in this same event the pool was filled
4371  Int_t ev0 = 0;
4372  if(GetReader()->GetLastTracksMixedEvent() == GetEventNumber()) ev0 = 1;
4373 
4374  for(Int_t ev=ev0; ev < pool->GetSize(); ev++)
4375  {
4376  //
4377  // Recover the lists of tracks or clusters
4378  //
4379  TObjArray* bgTracks = static_cast<TObjArray*>(pool->At(ev));
4380  TObjArray* bgCalo = 0;
4381 
4382  // Recover the clusters list if requested
4383  if( neutralMix && poolCalo )
4384  {
4385  if(pool->GetSize()!=poolCalo->GetSize())
4386  AliWarning("Different size of calo and track pools");
4387 
4388  bgCalo = static_cast<TObjArray*>(poolCalo->At(ev));
4389 
4390  if(!bgCalo) AliDebug(1,Form("Event %d in calo pool not available?",ev));
4391  }
4392 
4393  //
4394  // Isolate the trigger in the mixed event with mixed tracks and clusters
4395  //
4396  if( OnlyIsolated() )
4397  {
4398  Int_t n=0, nfrac = 0;
4399  Bool_t isolated = kFALSE;
4400  Float_t coneptsum = 0, coneptlead = 0;
4401  GetIsolationCut()->MakeIsolationCut(bgTracks,bgCalo,
4402  GetReader(), GetCaloPID(),
4403  kFALSE, aodParticle, "",
4404  n,nfrac,coneptsum,coneptlead,isolated);
4405 
4406  //printf("AliAnaParticleHadronCorrelation::MakeChargedMixCorrelation() - Isolated? %d - cone %f, ptthres %f",
4407  // isolated,GetIsolationCut()->GetConeSize(),GetIsolationCut()->GetPtThreshold());
4408  //if(bgTracks)printf(" - n track %d", bgTracks->GetEntriesFast());
4409  //printf("\n");
4410 
4411  if(!isolated) continue ;
4412  }
4413 
4414  //
4415  // Check if the trigger is leading of mixed event
4416  //
4417  Int_t nTracks=bgTracks->GetEntriesFast();
4418 
4419  if(fMakeNearSideLeading || fMakeAbsoluteLeading)
4420  {
4421  Bool_t leading = kTRUE;
4422  for(Int_t jlead = 0;jlead < nTracks; jlead++ )
4423  {
4424  AliCaloTrackParticle *track = (AliCaloTrackParticle*) bgTracks->At(jlead) ;
4425 
4426  ptAssoc = track->Pt();
4427  phiAssoc = track->Phi() ;
4428  if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
4429 
4430  if (fMakeNearSideLeading)
4431  {
4432  deltaPhi = phiTrig-phiAssoc;
4433  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
4434  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
4435 
4436  if(ptAssoc > ptTrig && deltaPhi < TMath::PiOver2())
4437  {
4438  leading = kFALSE;
4439  break;
4440  }
4441  }
4442  //jump out this event if there is any other particle with pt larger than trigger
4443  else if(fMakeAbsoluteLeading)
4444  {
4445  if(ptAssoc > ptTrig)
4446  {
4447  leading = kFALSE;
4448  break;
4449  }
4450  }
4451  }
4452 
4453  if( !neutralMix && fCheckLeadingWithNeutralClusters )
4454  AliWarning("Leading of clusters requested but no clusters in mixed event");
4455 
4456  if(neutralMix && fCheckLeadingWithNeutralClusters && bgCalo)
4457  {
4458  Int_t nClusters=bgCalo->GetEntriesFast();
4459  for(Int_t jlead = 0;jlead <nClusters; jlead++ )
4460  {
4461  AliCaloTrackParticle *cluster= (AliCaloTrackParticle*) bgCalo->At(jlead) ;
4462 
4463  ptAssoc = cluster->Pt();
4464  phiAssoc = cluster->Phi() ;
4465  if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
4466 
4467  if (fMakeNearSideLeading)
4468  {
4469  deltaPhi = phiTrig-phiAssoc;
4470  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
4471  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
4472 
4473  if(ptAssoc > ptTrig && deltaPhi < TMath::PiOver2())
4474  {
4475  leading = kFALSE;
4476  break;
4477  }
4478  }
4479  //jump out this event if there is any other particle with pt larger than trigger
4480  else if(fMakeAbsoluteLeading)
4481  {
4482  if(ptAssoc > ptTrig)
4483  {
4484  leading = kFALSE;
4485  break;
4486  }
4487  }
4488  }
4489  }
4490 
4491  if(!leading) continue; // not leading, check the next event in pool
4492  }
4493 
4494  //
4495  // Fill histograms for selected triggers
4496  //
4497 
4498  fhEventMixBin->Fill(eventBin, GetEventWeight());
4499 
4500  //printf("\t Read Pool event %d, nTracks %d\n",ev,nTracks);
4501 
4502  fhPtTriggerMixed ->Fill(ptTrig, GetEventWeight());
4503  fhPhiTriggerMixed ->Fill(ptTrig, phiTrig, GetEventWeight());
4504  fhEtaTriggerMixed ->Fill(ptTrig, etaTrig, GetEventWeight());
4505  fhPtTriggerMixedBin->Fill(ptTrig, eventBin, GetEventWeight());
4506 
4507  if ( fCorrelVzBin ) fhPtTriggerMixedVzBin->Fill(ptTrig, GetEventVzBin(), GetEventWeight());
4508 
4509  //
4510  // Correlation histograms
4511  //
4512  for(Int_t j1 = 0;j1 <nTracks; j1++ )
4513  {
4514  AliCaloTrackParticle *track = (AliCaloTrackParticle*) bgTracks->At(j1) ;
4515 
4516  if(!track) continue;
4517 
4518  ptAssoc = track->Pt();
4519  etaAssoc = track->Eta();
4520  phiAssoc = track->Phi() ;
4521  if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
4522 
4523  deltaPhi = phiTrig-phiAssoc;
4524  if(deltaPhi < -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
4525  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
4526  deltaEta = etaTrig-etaAssoc;
4527 
4528  AliDebug(1,Form("deltaPhi= %f, deltaEta=%f",deltaPhi, deltaEta));
4529 
4530  // Angular correlation
4531  fhMixDeltaPhiCharged ->Fill(ptTrig , deltaPhi, GetEventWeight());
4532  fhMixDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta, GetEventWeight());
4533 
4534  //
4535  // Momentum imbalance
4536  //
4537  if ( (deltaPhi > fDeltaPhiMinCut) && (deltaPhi < fDeltaPhiMaxCut) )
4538  {
4539  xE = -ptAssoc/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
4540 
4541  if( xE < 0.)
4542  AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
4543  xE,deltaPhi*TMath::RadToDeg(),TMath::Cos(deltaPhi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
4544 
4545  fhMixXECharged->Fill(ptTrig, xE, GetEventWeight());
4546  if ( xE > 0 ) fhMixHbpXECharged->Fill(ptTrig, TMath::Log(1./xE), GetEventWeight());
4547  }
4548 
4549  //
4550  // Underlying event momentum imbalance
4551  //
4552  if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )
4553  {
4554  //Underlying event region
4555  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
4556  Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
4557 
4558  if(uexE < 0.)
4559  AliWarning(Form("Careful!!, negative xE %2.2f for left UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
4560  uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
4561 
4562  fhMixXEUeCharged->Fill(ptTrig, uexE, GetEventWeight());
4563  }
4564 
4565  // Set the pt associated bin for the defined bins
4566  Int_t assocBin = -1;
4567  for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
4568  {
4569  if(ptAssoc > fAssocPtBinLimit[i] && ptAssoc < fAssocPtBinLimit[i+1]) assocBin= i;
4570  }
4571 
4572  //
4573  // Assign to the histogram array a bin corresponding to a combination of pTa and vz bins
4574  //
4575  Int_t nz = 1;
4576  Int_t vz = 0;
4577 
4578  if(fCorrelVzBin)
4579  {
4580  nz = GetNZvertBin();
4581  vz = GetEventVzBin();
4582  }
4583 
4584  Int_t bin = assocBin*nz+vz;
4585 
4586  if(bin < 0) continue ; // this pt bin was not considered
4587 
4588  fhMixDeltaPhiChargedAssocPtBin [bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
4589  fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->Fill(deltaPhi, deltaEta, GetEventWeight());
4590 
4591  if(fFillEtaGapsHisto)
4592  {
4593  if(TMath::Abs(deltaEta) > 0.8)
4594  fhMixDeltaPhiChargedAssocPtBinDEta08 [bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
4595  if(TMath::Abs(deltaEta) < 0.01)
4596  fhMixDeltaPhiChargedAssocPtBinDEta0 [bin]->Fill(ptTrig, deltaPhi, GetEventWeight());
4597  }
4598 
4599  } // track loop
4600  } // mixed event loop
4601 }
4602 
4603 //_______________________________________________________________________________________________________
4606 //_______________________________________________________________________________________________________
4607 void AliAnaParticleHadronCorrelation::MakeNeutralCorrelation(AliCaloTrackParticleCorrelation * aodParticle)
4608 {
4609  TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); // For the future, foresee more possible pi0 lists
4610  if(!pi0list) return ;
4611 
4612  Int_t npi0 = pi0list->GetEntriesFast();
4613  if(npi0 == 0) return ;
4614 
4615  AliDebug(1,Form("Particle - pi0 correlation, %d pi0's",npi0));
4616 
4617  Int_t evtIndex11 = 0 ;
4618  Int_t evtIndex12 = 0 ;
4619  if (GetMixedEvent())
4620  {
4621  evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
4622  evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
4623  }
4624 
4625  Float_t pt = -100. ;
4626  Float_t zT = -100. ;
4627  Float_t phi = -100. ;
4628  Float_t eta = -100. ;
4629  Float_t xE = -100. ;
4630  Float_t hbpXE= -100. ;
4631  Float_t hbpZT= -100. ;
4632 
4633  Float_t ptTrig = aodParticle->Pt();
4634  Float_t phiTrig = aodParticle->Phi();
4635  Float_t etaTrig = aodParticle->Eta();
4636  Float_t deltaPhi= -100. ;
4637  Float_t deltaEta= -100. ;
4638 
4639  // In case of pi0/eta trigger, we may want to check their decay correlation,
4640  // get their decay children
4641 
4642  Bool_t decayFound = kFALSE;
4643  if(fPi0Trigger) decayFound = GetDecayPhotonMomentum(aodParticle->GetCaloLabel(0),aodParticle->GetCaloLabel(1),aodParticle->GetDetectorTag());
4644 
4645  TObjArray * refpi0 = 0x0;
4646  Int_t nrefs = 0;
4647 
4648  // Loop on stored AOD pi0
4649 
4650  for(Int_t iaod = 0; iaod < npi0 ; iaod++)
4651  {
4652  AliCaloTrackParticle* pi0 = (AliCaloTrackParticle*) (pi0list->At(iaod));
4653 
4654  Int_t evtIndex2 = 0 ;
4655  Int_t evtIndex3 = 0 ;
4656  if (GetMixedEvent())
4657  {
4658  evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(0)) ;
4659  evtIndex3 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(1)) ;
4660 
4661  if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 ||
4662  evtIndex11 == evtIndex3 || evtIndex12 == evtIndex3) // trigger and pi0 are not from different events
4663  continue ;
4664  }
4665 
4666  pt = pi0->Pt();
4667 
4668  if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
4669 
4670  // Remove trigger itself for correlation when use charged triggers
4671  if(aodParticle->GetCaloLabel(0) >= 0 &&
4672  (pi0->GetCaloLabel(0) == aodParticle->GetCaloLabel(0) || pi0->GetCaloLabel(1) == aodParticle->GetCaloLabel(0))) continue ;
4673 
4674  if( aodParticle->GetCaloLabel(1) >= 0 &&
4675  (pi0->GetCaloLabel(0) == aodParticle->GetCaloLabel(1) || pi0->GetCaloLabel(1) == aodParticle->GetCaloLabel(1))) continue ;
4676 
4677  //
4678  // Angular correlations
4679  //
4680  phi = pi0->Phi() ;
4681  eta = pi0->Eta() ;
4682  deltaEta = etaTrig-eta;
4683  deltaPhi = phiTrig-phi;
4684  if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
4685  if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
4686 
4687  fhEtaNeutral ->Fill(pt , eta , GetEventWeight());
4688  fhPhiNeutral ->Fill(pt , phi , GetEventWeight());
4689  fhDeltaEtaNeutral->Fill(ptTrig, deltaEta, GetEventWeight());
4690  fhDeltaPhiNeutral->Fill(ptTrig, deltaPhi, GetEventWeight());
4691 
4692  if ( pt > 2 ) fhDeltaPhiDeltaEtaNeutral->Fill(deltaPhi, deltaEta, GetEventWeight());
4693 
4694  //
4695  // Momentum imbalance
4696  //
4697  zT = pt/ptTrig ;
4698 
4699  hbpZT = -100;
4700  hbpXE = -100;
4701 
4702  if(zT > 0 ) hbpZT = TMath::Log(1./zT);
4703 
4704  // Delta phi cut for correlation
4705  if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) )
4706  {
4707  xE =-pt/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
4708 
4709  if( xE < 0. )
4710  AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
4711  xE,deltaPhi*TMath::RadToDeg(),TMath::Cos(deltaPhi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
4712 
4713  if( xE > 0 ) hbpXE = TMath::Log(1./xE);
4714 
4715  fhDeltaPhiNeutralPt->Fill(pt , deltaPhi, GetEventWeight());
4716  fhXENeutral ->Fill(ptTrig, xE , GetEventWeight());
4717  fhPtHbpXENeutral ->Fill(ptTrig, hbpXE , GetEventWeight());
4718  fhZTNeutral ->Fill(ptTrig, zT , GetEventWeight());
4719  fhPtHbpZTNeutral ->Fill(ptTrig, hbpZT , GetEventWeight());
4720  }
4721  else if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )
4722  {
4723  // Randomize angle for xE calculation
4724  Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
4725 
4726  xE = -(pt/ptTrig)*TMath::Cos(randomphi);
4727  if(xE > 0 ) hbpXE = TMath::Log(1./xE);
4728 
4729  fhDeltaPhiUeNeutralPt->Fill(pt , deltaPhi, GetEventWeight());
4730  fhZTUeNeutral ->Fill(ptTrig, zT , GetEventWeight());
4731  fhPtHbpZTUeNeutral ->Fill(ptTrig, hbpZT , GetEventWeight());
4732  fhXEUeNeutral ->Fill(ptTrig, xE , GetEventWeight());
4733  fhPtHbpXEUeNeutral ->Fill(ptTrig, hbpXE , GetEventWeight());
4734  }
4735 
4736  // Several UE calculation, not sure it is useful
4737  // with partical calorimter acceptance
4738  if(fMakeSeveralUE) FillNeutralUnderlyingEventSidesHistograms(ptTrig,pt,zT,hbpZT,deltaPhi);
4739 
4740  //
4741  // Decay photon correlations
4742  //
4743  if(fPi0Trigger && decayFound)
4744  FillDecayPhotonCorrelationHistograms(pt, phi, kFALSE) ;
4745 
4746  if(fFillAODWithReferences)
4747  {
4748  nrefs++;
4749  if(nrefs==1)
4750  {
4751  refpi0 = new TObjArray(0);
4752  refpi0->SetName(GetAODObjArrayName()+"Pi0s");
4753  refpi0->SetOwner(kFALSE);
4754  }
4755  refpi0->Add(pi0);
4756  } // put references in trigger AOD
4757 
4758  AliDebug(1,Form("Selected pi0: pt %2.2f, phi %2.2f, eta %2.2f",pt,phi,eta));
4759 
4760  } // loop
4761 
4762  // Fill AOD with reference tracks, if not filling histograms
4763  if(fFillAODWithReferences && refpi0)
4764  {
4765  aodParticle->AddObjArray(refpi0);
4766  }
4767 }
4768 
4769 //__________________________________________________________________________________________________________________
4771 //__________________________________________________________________________________________________________________
4772 void AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation(Int_t label, Int_t histoIndex, Bool_t lostDecayPair)
4773 {
4774  AliDebug(1,"Make trigger particle - charged hadron correlation in AOD MC level");
4775 
4776  if( label < 0 )
4777  {
4778  AliDebug(1,Form(" *** bad label ***: label %d", label));
4779  return;
4780  }
4781 
4782  if ( !GetMC() )
4783  {
4784  AliFatal("Stack not available, is the MC handler called? STOP");
4785  return;
4786  }
4787 
4788  //Int_t nTracks = GetMC()->GetNumberOfTracks() ;
4789  Int_t nTracks = GetMC()->GetNumberOfPrimaries();
4790 
4791  if( label >= nTracks )
4792  {
4793  if(GetDebug() > 2)
4794  AliInfo(Form(" *** large label ***: label %d, n prim tracks %d", label,nTracks));
4795 
4796  return;
4797  }
4798