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