change all printf's to AliDebug/AliInfo/AliWarning
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaParticleHadronCorrelation.cxx
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 //_________________________________________________________________________
17 // Class for the analysis of particle - hadron correlations
18 // Particle (for example direct gamma) must be found
19 // in a previous analysis
20 //
21 //-- Author: Gustavo Conesa (LNF-INFN) (LPSC-IN2P3-CNRS)
22 //           Yaxian Mao (LPSC-IN2P3-CNRS) and (CNWU) first usable implementation.
23 //           Xiangrong Zhu (CNWU), implementtion of own mixing.
24 //////////////////////////////////////////////////////////////////////////////
25
26
27 // --- ROOT system ---
28 //#include "TClonesArray.h"
29 #include <TClass.h>
30 #include <TMath.h>
31 #include <TH2F.h>
32 #include <TDatabasePDG.h>
33
34 //---- ANALYSIS system ----
35 #include "AliNeutralMesonSelection.h" 
36 #include "AliAnaParticleHadronCorrelation.h" 
37 #include "AliCaloTrackReader.h"
38 #include "AliAODPWG4ParticleCorrelation.h"
39 #include "AliFiducialCut.h"
40 #include "AliVTrack.h"
41 #include "AliVCluster.h"
42 #include "AliMCAnalysisUtils.h"
43 #include "TParticle.h"
44 #include "AliStack.h"
45 #include "AliAODMCParticle.h"
46 #include "AliMixedEvent.h"
47 #include "AliAnalysisManager.h"
48 #include "AliInputEventHandler.h"
49 #include "AliEventplane.h"
50
51 ClassImp(AliAnaParticleHadronCorrelation)
52
53
54 //___________________________________________________________________
55   AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation(): 
56     AliAnaCaloTrackCorrBaseClass(),
57     fFillAODWithReferences(0),      fCheckLeadingWithNeutralClusters(0),
58     fMaxAssocPt(1000.),             fMinAssocPt(0.),
59     fDeltaPhiMaxCut(0.),            fDeltaPhiMinCut(0.),   
60     fSelectIsolated(0),             fMakeSeveralUE(0),              
61     fUeDeltaPhiMaxCut(0.),          fUeDeltaPhiMinCut(0.), 
62     fPi0AODBranchName(""),          fNeutralCorr(0),       
63     fPi0Trigger(0),                 fDecayTrigger(0),
64     fNDecayBits(0),                 fDecayBits(),
65     fMakeAbsoluteLeading(0),        fMakeNearSideLeading(0),      
66     fLeadingTriggerIndex(-1),       fHMPIDCorrelation(0),  fFillBradHisto(0),
67     fNAssocPtBins(0),               fAssocPtBinLimit(),
68     fCorrelVzBin(0),
69     fListMixTrackEvents(),          fListMixCaloEvents(),
70     fUseMixStoredInReader(0),       fFillNeutralEventMixPool(0),
71     fM02MaxCut(0),                  fM02MinCut(0),  
72     fSelectLeadingHadronAngle(0),   fFillLeadHadOppositeHisto(0),
73     fMinLeadHadPhi(0),              fMaxLeadHadPhi(0),
74     fMinLeadHadPt(0),               fMaxLeadHadPt(0),
75     fFillEtaGapsHisto(1),           fFillMomImbalancePtAssocBinsHisto(0),
76     fMCGenTypeMin(0),               fMCGenTypeMax(0),
77     fTrackVector(),                 fMomentum(),
78     fDecayMom1(),                   fDecayMom2(),
79     //Histograms
80     fhPtTriggerInput(0),            fhPtTriggerSSCut(0),
81     fhPtTriggerIsoCut(0),           fhPtTriggerFidCut(0),
82     fhPtTrigger(0),                 fhPtTriggerVtxBC0(0),
83     fhPtTriggerVzBin(0),            fhPtTriggerBin(0),                 
84     fhPhiTrigger(0),                fhEtaTrigger(0),   
85     fhPtTriggerMC(),
86     fhPtDecayTrigger(),             fhPtDecayTriggerMC(),
87     fhPtTriggerCentrality(0),       fhPtTriggerEventPlane(0), 
88     fhTriggerEventPlaneCentrality(0),
89     fhPtTriggerMixed(0),            fhPtTriggerMixedVzBin(0), fhPtTriggerMixedBin(0),              
90     fhPhiTriggerMixed(0),           fhEtaTriggerMixed(0),
91     fhPtLeadingOppositeHadron(0),   fhPtDiffPhiLeadingOppositeHadron(0), fhPtDiffEtaLeadingOppositeHadron(0),
92     fhPtNoLeadingOppositeHadron(0), fhEtaPhiNoLeadingOppositeHadron(0),
93     fhDeltaPhiDeltaEtaCharged(0),
94     fhPhiCharged(0),                fhEtaCharged(0), 
95     fhDeltaPhiCharged(0),           fhDeltaEtaCharged(0), 
96     fhDeltaPhiChargedPt(0),         fhDeltaPhiUeChargedPt(0), 
97     fhUePart(0),
98     fhXECharged(0),                 fhXECharged_Cone2(0),      fhXEUeCharged(0),
99     fhXEPosCharged(0),              fhXENegCharged(0),
100     fhPtHbpXECharged(0),            fhPtHbpXECharged_Cone2(0), fhPtHbpXEUeCharged(0),
101     fhZTCharged(0),                 fhZTUeCharged(0),
102     fhZTPosCharged(0),              fhZTNegCharged(0),
103     fhPtHbpZTCharged(0),            fhPtHbpZTUeCharged(0),
104     fhXEChargedMC(),                fhDeltaPhiChargedMC(),  
105     fhDeltaPhiDeltaEtaChargedPtA3GeV(0),
106     fhDeltaPhiChargedPtA3GeV(0),    fhDeltaEtaChargedPtA3GeV(0),
107     //Pile-Up
108     fhDeltaPhiChargedPileUp(),      fhDeltaEtaChargedPileUp(),
109     fhDeltaPhiChargedPtA3GeVPileUp(), fhDeltaEtaChargedPtA3GeVPileUp(),
110     fhXEChargedPileUp(),            fhXEUeChargedPileUp(),
111     fhZTChargedPileUp(),            fhZTUeChargedPileUp(), 
112     fhPtTrigChargedPileUp(),
113     fhDeltaPhiChargedOtherBC(),     fhDeltaPhiChargedPtA3GeVOtherBC(),
114     fhXEChargedOtherBC(),           fhXEUeChargedOtherBC(),
115     fhZTChargedOtherBC(),           fhZTUeChargedOtherBC(),
116     fhPtTrigChargedOtherBC(),
117     fhDeltaPhiChargedBC0(),         fhDeltaPhiChargedPtA3GeVBC0(),
118     fhXEChargedBC0(),               fhXEUeChargedBC0(),
119     fhZTChargedBC0(),               fhZTUeChargedBC0(),
120     fhPtTrigChargedBC0(),
121     fhDeltaPhiChargedVtxBC0(),      fhDeltaPhiChargedPtA3GeVVtxBC0(),
122     fhXEChargedVtxBC0(),            fhXEUeChargedVtxBC0(),
123     fhZTChargedVtxBC0(),            fhZTUeChargedVtxBC0(),
124     fhPtTrigChargedVtxBC0(),
125     fhDeltaPhiUeLeftCharged(0),
126     fhDeltaPhiUeLeftUpCharged(0),   fhDeltaPhiUeRightUpCharged(0),
127     fhDeltaPhiUeLeftDownCharged(0), fhDeltaPhiUeRightDownCharged(0),
128     fhXEUeLeftCharged(0),
129     fhXEUeLeftUpCharged(0),         fhXEUeRightUpCharged(0),
130     fhXEUeLeftDownCharged(0),       fhXEUeRightDownCharged(0),
131     fhPtHbpXEUeLeftCharged(0),      fhZTUeLeftCharged(0),
132     fhPtHbpZTUeLeftCharged(0),
133     fhPtTrigPout(0),                fhPtTrigCharged(0),
134     fhDeltaPhiChargedMult(0x0),     fhDeltaEtaChargedMult(0x0),
135     fhXEMult(0x0),                  fhXEUeMult(0x0),
136     fhZTMult(0x0),                  fhZTUeMult(0x0),
137     fhAssocPtBkg(0),                fhDeltaPhiDeltaEtaAssocPtBin(0),
138     fhDeltaPhiAssocPtBin(0),        
139     fhDeltaPhiAssocPtBinDEta08(0),  fhDeltaPhiAssocPtBinDEta0(0),
140     fhDeltaPhiAssocPtBinHMPID(0),   fhDeltaPhiAssocPtBinHMPIDAcc(0),
141     fhDeltaPhiBradAssocPtBin(0),    fhDeltaPhiBrad(0),
142     fhXEAssocPtBin(0),              fhZTAssocPtBin(0),
143     fhXEVZ(0),                      fhZTVZ(0),
144     fhDeltaPhiDeltaEtaNeutral(0),
145     fhPhiNeutral(0),                fhEtaNeutral(0), 
146     fhDeltaPhiNeutral(0),           fhDeltaEtaNeutral(0),
147     fhDeltaPhiNeutralPt(0),         fhDeltaPhiUeNeutralPt(0), 
148     fhXENeutral(0),                 fhXEUeNeutral(0),
149     fhPtHbpXENeutral(0),            fhPtHbpXEUeNeutral(0),
150     fhZTNeutral(0),                 fhZTUeNeutral(0),
151     fhPtHbpZTNeutral(0),            fhPtHbpZTUeNeutral(0),
152     fhDeltaPhiUeLeftNeutral(0),     fhXEUeLeftNeutral(0),
153     fhPtHbpXEUeLeftNeutral(0),      fhZTUeLeftNeutral(0),
154     fhPtHbpZTUeLeftNeutral(0),      fhPtPi0DecayRatio(0),
155     fhDeltaPhiPi0DecayCharged(0),   fhXEPi0DecayCharged(0),        fhZTPi0DecayCharged(0),
156     fhDeltaPhiPi0DecayNeutral(0),   fhXEPi0DecayNeutral(0),        fhZTPi0DecayNeutral(0),
157     fhDeltaPhiDecayCharged(),       fhXEDecayCharged(),            fhZTDecayCharged(),
158     fhDeltaPhiDecayChargedAssocPtBin(),
159     fhMCPtTrigger(),                fhMCPhiTrigger(),              fhMCEtaTrigger(),
160     fhMCPtTriggerNotLeading(),      fhMCPhiTriggerNotLeading(),    fhMCEtaTriggerNotLeading(),
161     fhMCEtaCharged(),               fhMCPhiCharged(),
162     fhMCDeltaEtaCharged(),          fhMCDeltaPhiCharged(),
163     fhMCDeltaPhiDeltaEtaCharged(),  fhMCDeltaPhiChargedPt(),
164     fhMCPtXECharged(),              fhMCPtXEUeCharged(),
165     fhMCPtXEUeLeftCharged(),
166     fhMCPtHbpXECharged(),           fhMCPtHbpXEUeCharged(),
167     fhMCPtHbpXEUeLeftCharged(),
168     fhMCUePart(),
169     fhMCPtZTCharged(),              fhMCPtZTUeCharged(),
170     fhMCPtZTUeLeftCharged(),
171     fhMCPtHbpZTCharged(),           fhMCPtHbpZTUeCharged(),
172     fhMCPtHbpZTUeLeftCharged(),
173     fhMCPtTrigPout(),               fhMCPtAssocDeltaPhi(),
174     //Mixing
175     fhNEventsTrigger(0),            fhNtracksMB(0),                 fhNclustersMB(0),
176     fhMixDeltaPhiCharged(0),        fhMixDeltaPhiDeltaEtaCharged(0),
177     fhMixXECharged(0),              fhMixXEUeCharged(0),            fhMixHbpXECharged(0),
178     fhMixDeltaPhiChargedAssocPtBin(),
179     fhMixDeltaPhiChargedAssocPtBinDEta08(),
180     fhMixDeltaPhiChargedAssocPtBinDEta0(),
181     fhMixDeltaPhiDeltaEtaChargedAssocPtBin(),
182     fhEventBin(0),                  fhEventMixBin(0),               fhEventMBBin(0)
183
184   //Default Ctor 
185     
186   //Initialize parameters
187   InitParameters();
188   
189   for(Int_t i = 0; i < fgkNmcTypes; i++)
190   { 
191     fhPtTriggerMC[i] = 0;
192     fhXEChargedMC[i] = 0;
193     fhDeltaPhiChargedMC[i] = 0;
194     for(Int_t ib = 0; ib < 4; ib++)  fhPtDecayTriggerMC[ib][i] = 0;
195   }
196
197   for(Int_t ib = 0; ib < 4; ib++)  fhPtDecayTrigger[ib] = 0;
198
199   for(Int_t i = 0; i < 7; i++)
200   {
201     fhPtTriggerPileUp             [i] = 0 ;
202     fhDeltaPhiChargedPileUp       [i] = 0 ; fhDeltaEtaChargedPileUp       [i] = 0 ;
203     fhXEChargedPileUp             [i] = 0 ; fhXEUeChargedPileUp           [i] = 0 ;
204     fhZTChargedPileUp             [i] = 0 ; fhZTUeChargedPileUp           [i] = 0 ;
205     fhPtTrigChargedPileUp         [i] = 0 ;
206     fhDeltaPhiChargedPtA3GeVPileUp[i] = 0 ; fhDeltaEtaChargedPtA3GeVPileUp[i] = 0 ;
207   }
208   
209 }
210
211 //_________________________________________________________________
212 AliAnaParticleHadronCorrelation::~AliAnaParticleHadronCorrelation()
213 {
214   // Remove event containers
215   
216   if(DoOwnMix())
217   {
218     if(fListMixTrackEvents)
219     {      
220       for(Int_t iz=0; iz < GetNZvertBin(); iz++)
221       {
222         for(Int_t ic=0; ic < GetNCentrBin(); ic++)
223         {
224           for(Int_t irp=0; irp<GetNRPBin(); irp++)
225           {
226             Int_t bin = GetEventMixBin(ic, iz, irp);
227             fListMixTrackEvents[bin]->Delete() ;
228             delete fListMixTrackEvents[bin] ;
229           }
230         }
231       }
232     }
233
234     delete[] fListMixTrackEvents;
235
236     if(fListMixCaloEvents)
237     {      
238       for(Int_t iz=0; iz < GetNZvertBin(); iz++)
239       {
240         for(Int_t ic=0; ic < GetNCentrBin(); ic++)
241         {
242           for(Int_t irp=0; irp<GetNRPBin(); irp++)
243           {
244             Int_t bin = GetEventMixBin(ic, iz, irp);
245             fListMixCaloEvents[bin]->Delete() ;
246             delete fListMixCaloEvents[bin] ;
247           }
248         }
249       }
250     }
251   
252     delete[] fListMixCaloEvents;
253     
254   }
255 }
256
257 //__________________________________________________________________________________________________________________________________________
258 void AliAnaParticleHadronCorrelation::FillChargedAngularCorrelationHistograms(Float_t ptAssoc,  Float_t ptTrig,      Int_t   bin,
259                                                                               Float_t phiAssoc, Float_t phiTrig,     Float_t deltaPhi,
260                                                                               Float_t etaAssoc, Float_t etaTrig,  
261                                                                               Int_t   decayTag, Float_t hmpidSignal, Int_t  outTOF,
262                                                                               Int_t   cen,      Int_t   mcTag)
263 {
264   // Fill angular correlation related histograms
265   
266   Float_t deltaEta    = etaTrig-etaAssoc;
267   Float_t deltaPhiOrg = phiTrig-phiAssoc;
268   
269   fhEtaCharged       ->Fill(ptAssoc,etaAssoc);
270   fhPhiCharged       ->Fill(ptAssoc,phiAssoc);
271   fhDeltaEtaCharged  ->Fill(ptTrig ,deltaEta);
272   fhDeltaPhiCharged  ->Fill(ptTrig ,deltaPhi);
273   fhDeltaPhiChargedPt->Fill(ptAssoc, deltaPhi);
274   fhDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta);
275   
276   if(ptAssoc > 3 )
277   {
278     fhDeltaEtaChargedPtA3GeV        ->Fill(ptTrig  ,deltaEta);
279     fhDeltaPhiChargedPtA3GeV        ->Fill(ptTrig  ,deltaPhi);
280     fhDeltaPhiDeltaEtaChargedPtA3GeV->Fill(deltaPhi, deltaEta);    
281   }   
282   
283   // Pile up studies
284   
285   if(IsPileUpAnalysisOn())
286   {
287     if     (outTOF==1)
288     {
289       fhDeltaPhiChargedOtherBC->Fill(ptTrig ,deltaPhi) ;
290       if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVOtherBC->Fill(ptTrig ,deltaPhi) ;
291     }
292     else if(outTOF==0)
293     {
294         fhDeltaPhiChargedBC0->Fill(ptTrig ,deltaPhi) ;
295         if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVBC0->Fill(ptTrig ,deltaPhi) ;
296     }
297     
298     Int_t vtxBC = GetReader()->GetVertexBC();
299     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
300     {
301       fhDeltaPhiChargedVtxBC0->Fill(ptTrig ,deltaPhi) ;
302       if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVVtxBC0->Fill(ptTrig ,deltaPhi) ;
303     }
304     
305     if(GetReader()->IsPileUpFromSPD())               { fhDeltaEtaChargedPileUp[0]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[0]->Fill(ptTrig ,deltaPhi) ; }
306     if(GetReader()->IsPileUpFromEMCal())             { fhDeltaEtaChargedPileUp[1]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[1]->Fill(ptTrig ,deltaPhi) ; }
307     if(GetReader()->IsPileUpFromSPDOrEMCal())        { fhDeltaEtaChargedPileUp[2]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[2]->Fill(ptTrig ,deltaPhi) ; }
308     if(GetReader()->IsPileUpFromSPDAndEMCal())       { fhDeltaEtaChargedPileUp[3]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[3]->Fill(ptTrig ,deltaPhi) ; }
309     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    { fhDeltaEtaChargedPileUp[4]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[4]->Fill(ptTrig ,deltaPhi) ; }
310     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    { fhDeltaEtaChargedPileUp[5]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[5]->Fill(ptTrig ,deltaPhi) ; }
311     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhDeltaEtaChargedPileUp[6]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[6]->Fill(ptTrig ,deltaPhi) ; }
312     
313     if(ptAssoc > 3 )
314     {
315       if(GetReader()->IsPileUpFromSPD())               { fhDeltaEtaChargedPtA3GeVPileUp[0]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[0]->Fill(ptTrig ,deltaPhi) ; }
316       if(GetReader()->IsPileUpFromEMCal())             { fhDeltaEtaChargedPtA3GeVPileUp[1]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[1]->Fill(ptTrig ,deltaPhi) ; }
317       if(GetReader()->IsPileUpFromSPDOrEMCal())        { fhDeltaEtaChargedPtA3GeVPileUp[2]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[2]->Fill(ptTrig ,deltaPhi) ; }
318       if(GetReader()->IsPileUpFromSPDAndEMCal())       { fhDeltaEtaChargedPtA3GeVPileUp[3]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[3]->Fill(ptTrig ,deltaPhi) ; }
319       if(GetReader()->IsPileUpFromSPDAndNotEMCal())    { fhDeltaEtaChargedPtA3GeVPileUp[4]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[4]->Fill(ptTrig ,deltaPhi) ; }
320       if(GetReader()->IsPileUpFromEMCalAndNotSPD())    { fhDeltaEtaChargedPtA3GeVPileUp[5]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[5]->Fill(ptTrig ,deltaPhi) ; }
321       if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhDeltaEtaChargedPtA3GeVPileUp[6]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[6]->Fill(ptTrig ,deltaPhi) ; }
322     }
323   }
324   
325   if(IsDataMC())
326   {
327     Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
328     fhDeltaPhiChargedMC[mcIndex]->Fill(ptTrig , deltaPhi);
329     if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) && mcIndex==2 )
330       fhDeltaPhiChargedMC[7]->Fill(ptTrig , deltaPhi);
331   }
332   
333   if(fDecayTrigger && decayTag > 0)
334   {
335     for(Int_t ibit = 0; ibit<fNDecayBits; ibit++)
336     {
337       if(GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[ibit])) fhDeltaPhiDecayCharged[ibit]->Fill(ptTrig,deltaPhi);
338     }
339   }
340   
341   Double_t  dphiBrad = -100;
342   if(fFillBradHisto)
343   {
344     dphiBrad = atan2(sin(deltaPhiOrg), cos(deltaPhiOrg))/TMath::Pi();//-1 to 1
345     if( TMath::Abs(dphiBrad) > 0.325 && TMath::Abs(dphiBrad) < 0.475 )  //Hardcoded values, BAD, FIXME
346     {
347       fhAssocPtBkg->Fill(ptTrig, ptAssoc);
348     }
349     
350     if( dphiBrad < -1./3 ) dphiBrad += 2;
351     fhDeltaPhiBrad->Fill(ptTrig, dphiBrad);
352   }
353   
354   // Fill histograms in bins of associated particle pT
355   if(bin>=0)
356   {
357       fhDeltaPhiDeltaEtaAssocPtBin    [bin]->Fill(deltaPhi,deltaEta);
358
359       fhDeltaPhiAssocPtBin            [bin]->Fill(ptTrig, deltaPhi);
360     
361     if(fFillEtaGapsHisto)
362     {
363       if(TMath::Abs(deltaEta)> 0.8)
364         fhDeltaPhiAssocPtBinDEta08      [bin]->Fill(ptTrig, deltaPhi);
365       
366       if(TMath::Abs(deltaEta)< 0.01)
367         fhDeltaPhiAssocPtBinDEta0       [bin]->Fill(ptTrig, deltaPhi);
368     }
369     
370     if (fFillBradHisto)
371       fhDeltaPhiBradAssocPtBin        [bin]->Fill(ptTrig, dphiBrad);
372     
373     if(fDecayTrigger && decayTag > 0 && fNDecayBits > 0 &&
374        GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[0]))
375       fhDeltaPhiDecayChargedAssocPtBin[bin]->Fill(ptTrig, deltaPhi);      
376     
377     if(fHMPIDCorrelation)
378     {
379       if( hmpidSignal > 0 )
380       {
381         //printf("Track pt %f with HMPID signal %f \n",pt,hmpidSignal);
382         fhDeltaPhiAssocPtBinHMPID[bin]->Fill(ptTrig, deltaPhi);        
383       }
384       
385       if(phiAssoc > 5*TMath::DegToRad() && phiAssoc < 20*TMath::DegToRad())
386       {
387         //printf("Track pt %f in HMPID acceptance phi %f \n ",pt,phi*TMath::RadToDeg() );
388         fhDeltaPhiAssocPtBinHMPIDAcc[bin]->Fill(ptTrig, deltaPhi);        
389       }
390     }
391   }
392   
393   //fill different multiplicity/centrality histogram
394   if(IsHighMultiplicityAnalysisOn() && cen >= 0 && cen < GetNCentrBin())
395   {
396     fhDeltaPhiChargedMult[cen]->Fill(ptTrig,deltaPhi);
397     fhDeltaEtaChargedMult[cen]->Fill(ptTrig,deltaEta);
398   }  
399 }
400
401 //___________________________________________________________________________________________________________________________________
402 Bool_t AliAnaParticleHadronCorrelation::FillChargedMCCorrelationHistograms(Float_t mcAssocPt, Float_t mcAssocPhi, Float_t mcAssocEta,
403                                                                            Float_t mcTrigPt,  Float_t mcTrigPhi,  Float_t mcTrigEta,
404                                                                            Int_t histoIndex, Bool_t lostDecayPair)
405 {
406   // Fill MC histograms independently of AOD or ESD
407   
408   Bool_t lead = kTRUE;
409   
410   // In case we requested the trigger to be a leading particle,
411   // check if this is true at the MC level.
412   // Not sure if it is correct to skip or not skip this.
413   // Absolute leading?
414   if( fMakeAbsoluteLeading && mcAssocPt > mcTrigPt )     lead = kFALSE; // skip event
415   
416   // Skip this event if near side associated particle pt larger than trigger
417   if( mcAssocPhi < 0 ) mcAssocPhi+=TMath::TwoPi();
418   
419   Float_t mcdeltaPhi= mcTrigPhi-mcAssocPhi;
420   if(mcdeltaPhi <= -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
421   if(mcdeltaPhi > 3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
422
423   if( fMakeNearSideLeading)
424   {
425     if( mcAssocPt > mcTrigPt && mcdeltaPhi < TMath::PiOver2() ) lead = kFALSE; // skip event
426   }
427   
428   //
429   // Select only hadrons in pt range
430   if ( mcAssocPt < fMinAssocPt || mcAssocPt > fMaxAssocPt ) return lead ; // exclude but continue
431   if ( mcAssocPt < GetReader()->GetCTSPtMin())              return lead ;
432
433   
434   //
435   // Remove trigger itself for correlation when use charged triggers
436   if(TMath::Abs(mcAssocPt -mcTrigPt ) < 1e-6 && 
437      mcdeltaPhi                       < 1e-6 &&
438      TMath::Abs(mcAssocEta-mcTrigEta) < 1e-6)            return lead ; // exclude but continue
439   
440   Float_t mcxE    =-mcAssocPt/mcTrigPt*TMath::Cos(mcdeltaPhi);// -(mcAssocPx*pxprim+mcAssocPy*pyprim)/(mcTrigPt*mcTrigPt);  
441   Float_t mchbpXE =-100 ;
442   if(mcxE > 0 ) mchbpXE = TMath::Log(1./mcxE);
443   
444   Float_t mczT    = mcAssocPt/mcTrigPt ;
445   Float_t mchbpZT =-100 ;
446   if(mczT > 0 ) mchbpZT = TMath::Log(1./mczT);
447   
448   Double_t mcpout = mcAssocPt*TMath::Sin(mcdeltaPhi) ; 
449   
450   AliDebug(1,Form("Charged hadron: track Pt %f, track Phi %f, phi trigger %f. Cuts:  delta phi  %2.2f < %2.2f < %2.2f",
451                   mcAssocPt,mcAssocPhi, mcTrigPhi,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut));
452   
453   
454   // Fill Histograms
455   fhMCEtaCharged     [histoIndex]->Fill(mcAssocPt, mcAssocEta);
456   fhMCPhiCharged     [histoIndex]->Fill(mcAssocPt, mcAssocPhi);
457   fhMCDeltaEtaCharged[histoIndex]->Fill(mcTrigPt , mcTrigEta-mcAssocEta);
458   fhMCDeltaPhiCharged[histoIndex]->Fill(mcTrigPt , mcdeltaPhi);
459   fhMCPtAssocDeltaPhi[histoIndex]->Fill(mcAssocPt, mcdeltaPhi);
460   
461   fhMCDeltaPhiDeltaEtaCharged[histoIndex]->Fill(mcdeltaPhi,mcTrigEta-mcAssocEta);
462   
463   //delta phi cut for correlation
464   if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) ) 
465   {
466     fhMCDeltaPhiChargedPt[histoIndex]->Fill(mcAssocPt,mcdeltaPhi);
467     fhMCPtXECharged      [histoIndex]->Fill(mcTrigPt, mcxE);
468     fhMCPtHbpXECharged   [histoIndex]->Fill(mcTrigPt, mchbpXE);
469     fhMCPtZTCharged      [histoIndex]->Fill(mcTrigPt, mczT);
470     fhMCPtHbpZTCharged   [histoIndex]->Fill(mcTrigPt, mchbpZT);
471     fhMCPtTrigPout       [histoIndex]->Fill(mcTrigPt, mcpout) ;
472   }
473
474   if(histoIndex==2 && lostDecayPair && 7 >= fMCGenTypeMin && 7 <= fMCGenTypeMax )
475   {
476     fhMCEtaCharged     [7]->Fill(mcAssocPt, mcAssocEta);
477     fhMCPhiCharged     [7]->Fill(mcAssocPt, mcAssocPhi);
478     fhMCDeltaEtaCharged[7]->Fill(mcTrigPt , mcTrigEta-mcAssocEta);
479     fhMCDeltaPhiCharged[7]->Fill(mcTrigPt , mcdeltaPhi);
480     fhMCPtAssocDeltaPhi[7]->Fill(mcAssocPt, mcdeltaPhi);
481     
482     fhMCDeltaPhiDeltaEtaCharged[7]->Fill(mcdeltaPhi,mcTrigEta-mcAssocEta);
483     
484     //delta phi cut for correlation
485     if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) )
486     {
487       fhMCDeltaPhiChargedPt[7]->Fill(mcAssocPt,mcdeltaPhi);
488       fhMCPtXECharged      [7]->Fill(mcTrigPt, mcxE);
489       fhMCPtHbpXECharged   [7]->Fill(mcTrigPt, mchbpXE);
490       fhMCPtZTCharged      [7]->Fill(mcTrigPt, mczT);
491       fhMCPtHbpZTCharged   [7]->Fill(mcTrigPt, mchbpZT);
492       fhMCPtTrigPout       [7]->Fill(mcTrigPt, mcpout) ;
493     }
494   }
495   
496   // Underlying event
497   
498   // Right
499   if ( (mcdeltaPhi > fUeDeltaPhiMinCut) && (mcdeltaPhi < fUeDeltaPhiMaxCut) )
500   {
501     //Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
502     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
503     Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
504     Double_t mcUezT =   mcAssocPt/mcTrigPt;
505     
506     if(mcUexE < 0.)
507       AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
508                       mcUexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
509
510     fhMCPtXEUeCharged[histoIndex]->Fill(mcTrigPt,mcUexE);
511     if(mcUexE > 0) fhMCPtHbpXEUeCharged[histoIndex]->Fill(mcTrigPt,TMath::Log(1/mcUexE));
512     
513     fhMCPtZTUeCharged[histoIndex]->Fill(mcTrigPt,mcUezT);
514     if(mcUezT > 0) fhMCPtHbpZTUeCharged[histoIndex]->Fill(mcTrigPt,TMath::Log(1/mcUezT));
515     
516     fhMCUePart[histoIndex]->Fill(mcTrigPt);
517     
518     if(histoIndex==2 && lostDecayPair && 7 >= fMCGenTypeMin && 7 <= fMCGenTypeMax )
519     {
520       fhMCPtXEUeCharged[7]->Fill(mcTrigPt,mcUexE);
521       if(mcUexE > 0) fhMCPtHbpXEUeCharged[histoIndex]->Fill(mcTrigPt,TMath::Log(1/mcUexE));
522       
523       fhMCPtZTUeCharged[7]->Fill(mcTrigPt,mcUezT);
524       if(mcUezT > 0) fhMCPtHbpZTUeCharged[histoIndex]->Fill(mcTrigPt,TMath::Log(1/mcUezT));
525       
526       fhMCUePart[7]->Fill(mcTrigPt);
527     }
528   }
529
530   if(fMakeSeveralUE)
531   {
532     // Left
533     if((mcdeltaPhi<-fUeDeltaPhiMinCut) || (mcdeltaPhi >2*fUeDeltaPhiMaxCut))
534     {
535       Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
536       Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
537       Double_t mcUezT =   mcAssocPt/mcTrigPt;
538       
539       if(mcUexE < 0.)
540         AliWarning(Form("Careful!!, negative xE %2.2f for left UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
541                         mcUexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
542       
543       fhMCPtXEUeLeftCharged[histoIndex]->Fill(mcTrigPt,mcUexE);
544       if(mcUexE > 0) fhMCPtHbpXEUeLeftCharged[histoIndex]->Fill(mcTrigPt,TMath::Log(1/mcUexE));
545       
546       fhMCPtZTUeLeftCharged[histoIndex]->Fill(mcTrigPt,mcUezT);
547       if(mcUezT > 0) fhMCPtHbpZTUeLeftCharged[histoIndex]->Fill(mcTrigPt,TMath::Log(1/mcUezT));
548       
549       if(histoIndex==2 && lostDecayPair && 7 >= fMCGenTypeMin && 7 <= fMCGenTypeMax )
550       {
551         fhMCPtXEUeLeftCharged[7]->Fill(mcTrigPt,mcUexE);
552         if(mcUexE > 0) fhMCPtHbpXEUeLeftCharged[7]->Fill(mcTrigPt,TMath::Log(1/mcUexE));
553         
554         fhMCPtZTUeLeftCharged[7]->Fill(mcTrigPt,mcUezT);
555         if(mcUezT > 0) fhMCPtHbpZTUeLeftCharged[7]->Fill(mcTrigPt,TMath::Log(1/mcUezT));
556       }
557     }
558   }
559   
560   return lead;
561 }
562
563 //___________________________________________________________________________________________________________________
564 void AliAnaParticleHadronCorrelation::FillChargedMomentumImbalanceHistograms(Float_t ptTrig,   Float_t ptAssoc, 
565                                                                              Float_t deltaPhi,
566                                                                              Int_t   cen,      Int_t   charge,
567                                                                              Int_t   bin,      Int_t   decayTag,
568                                                                              Int_t   outTOF,   Int_t   mcTag)
569
570 {
571   // Fill mostly momentum imbalance related histograms
572   
573   Float_t zT =   ptAssoc/ptTrig ;
574   Float_t xE   =-ptAssoc/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
575   Float_t pout = ptAssoc*TMath::Sin(deltaPhi) ;
576
577   if(xE < 0.)
578     AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
579                xE,deltaPhi*TMath::RadToDeg(),TMath::Cos(deltaPhi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
580
581   Float_t hbpXE = -100;
582   Float_t hbpZT = -100;
583
584   if(xE > 0 ) hbpXE = TMath::Log(1./xE);
585   if(zT > 0 ) hbpZT = TMath::Log(1./zT);
586   
587   fhXECharged         ->Fill(ptTrig , xE);
588   fhPtHbpXECharged    ->Fill(ptTrig , hbpXE);
589   fhZTCharged         ->Fill(ptTrig , zT); 
590   fhPtHbpZTCharged    ->Fill(ptTrig , hbpZT);
591   fhPtTrigPout        ->Fill(ptTrig , pout) ;
592   fhPtTrigCharged     ->Fill(ptTrig , ptAssoc) ;
593   if((deltaPhi > 5*TMath::Pi()/6.)   && (deltaPhi < 7*TMath::Pi()/6.))
594   {
595     fhXECharged_Cone2         ->Fill(ptTrig , xE);
596     fhPtHbpXECharged_Cone2    ->Fill(ptTrig , hbpXE);
597   }
598   
599   // MC
600   if(IsDataMC())
601   {
602     Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
603     fhXEChargedMC[mcIndex]->Fill(ptTrig , xE);
604     if(GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCDecayPairLost) && mcIndex==2 )
605       fhXEChargedMC[7]->Fill(ptTrig , xE);
606   }
607   
608   // Pile up studies
609   if(IsPileUpAnalysisOn())
610   {
611     if     (outTOF==1)
612     {
613       fhXEChargedOtherBC    ->Fill(ptTrig,xE);
614       fhZTChargedOtherBC    ->Fill(ptTrig,zT);
615       fhPtTrigChargedOtherBC->Fill(ptTrig,ptAssoc);
616     }
617     else if(outTOF==0)
618     {
619       fhXEChargedBC0    ->Fill(ptTrig,xE);
620       fhZTChargedBC0    ->Fill(ptTrig,zT);
621       fhPtTrigChargedBC0->Fill(ptTrig,ptAssoc);
622     }
623
624     Int_t vtxBC = GetReader()->GetVertexBC();
625     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
626     {
627       fhXEChargedVtxBC0    ->Fill(ptTrig,xE);
628       fhZTChargedVtxBC0    ->Fill(ptTrig,zT);
629       fhPtTrigChargedVtxBC0->Fill(ptTrig,ptAssoc);
630     }
631        
632     if(GetReader()->IsPileUpFromSPD())                { fhXEChargedPileUp[0]->Fill(ptTrig,xE); fhZTChargedPileUp[0]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[0]->Fill(ptTrig,ptAssoc); }
633     if(GetReader()->IsPileUpFromEMCal())              { fhXEChargedPileUp[1]->Fill(ptTrig,xE); fhZTChargedPileUp[1]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[1]->Fill(ptTrig,ptAssoc); }
634     if(GetReader()->IsPileUpFromSPDOrEMCal())         { fhXEChargedPileUp[2]->Fill(ptTrig,xE); fhZTChargedPileUp[2]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[2]->Fill(ptTrig,ptAssoc); }
635     if(GetReader()->IsPileUpFromSPDAndEMCal())        { fhXEChargedPileUp[3]->Fill(ptTrig,xE); fhZTChargedPileUp[3]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[3]->Fill(ptTrig,ptAssoc); }
636     if(GetReader()->IsPileUpFromSPDAndNotEMCal())     { fhXEChargedPileUp[4]->Fill(ptTrig,xE); fhZTChargedPileUp[4]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[4]->Fill(ptTrig,ptAssoc); }
637     if(GetReader()->IsPileUpFromEMCalAndNotSPD())     { fhXEChargedPileUp[5]->Fill(ptTrig,xE); fhZTChargedPileUp[5]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[5]->Fill(ptTrig,ptAssoc); }
638     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())  { fhXEChargedPileUp[6]->Fill(ptTrig,xE); fhZTChargedPileUp[6]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[6]->Fill(ptTrig,ptAssoc); }
639   }
640   
641   if(fDecayTrigger && decayTag > 0)
642   {
643     for(Int_t ibit = 0; ibit<fNDecayBits; ibit++)
644     {
645       if(GetNeutralMesonSelection()->CheckDecayBit(decayTag,fDecayBits[ibit]))
646       {
647         fhXEDecayCharged[ibit]->Fill(ptTrig,xE);
648         fhZTDecayCharged[ibit]->Fill(ptTrig,zT);
649       }
650     }
651   } // photon decay pi0/eta trigger
652   
653   if(bin >= 0 && fFillMomImbalancePtAssocBinsHisto)//away side
654   {
655     fhXEAssocPtBin[bin]->Fill(ptTrig, xE) ;
656     fhZTAssocPtBin[bin]->Fill(ptTrig, zT) ;
657   }        
658   
659   if(fCorrelVzBin)
660   {
661     Int_t vz = GetEventVzBin();
662     fhXEVZ[vz]->Fill(ptTrig, xE) ;
663     fhZTVZ[vz]->Fill(ptTrig, zT) ;
664   }
665   
666   if(charge > 0)
667   {
668     fhXEPosCharged->Fill(ptTrig,xE) ;
669     fhZTPosCharged->Fill(ptTrig,zT) ;
670   }
671   else  
672   {
673     fhXENegCharged->Fill(ptTrig,xE) ;
674     fhZTNegCharged->Fill(ptTrig,zT) ;
675   }
676   
677   //fill different multiplicity/centrality histogram
678   if(IsHighMultiplicityAnalysisOn() && cen >= 0 && cen < GetNCentrBin())
679   {
680     fhXEMult[cen]->Fill(ptTrig,xE);
681     fhZTMult[cen]->Fill(ptTrig,zT);
682   } //multiplicity/centrality events selection
683
684
685 //_______________________________________________________________________________________________________________________
686 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventHistograms(Float_t ptTrig,   Float_t ptAssoc,
687                                                                            Float_t deltaPhi, Int_t cen, Int_t outTOF)
688 {
689   // Fill underlying event histograms
690   
691   fhUePart->Fill(ptTrig);
692   
693   fhDeltaPhiUeChargedPt->Fill(ptAssoc,deltaPhi);
694   
695   Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
696   Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
697   Double_t uezT =   ptAssoc/ptTrig;
698   
699   if(uexE < 0.)
700     AliWarning(Form("Careful!!, negative xE %2.2f for right UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
701                uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
702   
703   fhXEUeCharged->Fill(ptTrig,uexE);
704   if(uexE > 0) fhPtHbpXEUeCharged->Fill(ptTrig,TMath::Log(1/uexE));
705   
706   fhZTUeCharged->Fill(ptTrig,uezT);
707   if(uezT > 0) fhPtHbpZTUeCharged->Fill(ptTrig,TMath::Log(1/uezT));
708   
709   // Pile up studies
710   
711   if(IsPileUpAnalysisOn())
712   {
713     if     (outTOF==1)
714     {
715       fhXEUeChargedOtherBC->Fill(ptTrig,uexE);
716       fhZTUeChargedOtherBC->Fill(ptTrig,uezT);
717     }
718     else if(outTOF==0)
719     {
720       fhXEUeChargedBC0->Fill(ptTrig,uexE);
721       fhZTUeChargedBC0->Fill(ptTrig,uezT);
722     }
723     
724     Int_t vtxBC = GetReader()->GetVertexBC();
725     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
726     {
727       fhXEUeChargedVtxBC0->Fill(ptTrig,uexE);
728       fhZTUeChargedVtxBC0->Fill(ptTrig,uezT);
729     }
730
731     if(GetReader()->IsPileUpFromSPD())               { fhXEUeChargedPileUp[0]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[0]->Fill(ptTrig,uezT);}
732     if(GetReader()->IsPileUpFromEMCal())             { fhXEUeChargedPileUp[1]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[1]->Fill(ptTrig,uezT);}
733     if(GetReader()->IsPileUpFromSPDOrEMCal())        { fhXEUeChargedPileUp[2]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[2]->Fill(ptTrig,uezT);}
734     if(GetReader()->IsPileUpFromSPDAndEMCal())       { fhXEUeChargedPileUp[3]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[3]->Fill(ptTrig,uezT);}
735     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    { fhXEUeChargedPileUp[4]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[4]->Fill(ptTrig,uezT);}
736     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    { fhXEUeChargedPileUp[5]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[5]->Fill(ptTrig,uezT);}
737     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhXEUeChargedPileUp[6]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[6]->Fill(ptTrig,uezT);}
738   }
739   
740   //fill different multiplicity/centrality histogram
741   if(IsHighMultiplicityAnalysisOn() && cen >= 0 && cen < GetNCentrBin())
742   {
743     fhXEUeMult[cen]->Fill(ptTrig,uexE);
744     fhZTUeMult[cen]->Fill(ptTrig,uezT);
745   }
746 }
747
748 //_____________________________________________________________________________________________________
749 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventSidesHistograms(Float_t ptTrig, 
750                                                                                 Float_t ptAssoc, 
751                                                                                 Float_t deltaPhi)
752 {
753   // Fill underlying event histograms to the left and right of trigger
754   // Right cone is the default UE.
755   
756   if((deltaPhi<-fUeDeltaPhiMinCut) || (deltaPhi >2*fUeDeltaPhiMaxCut))
757   {
758     fhDeltaPhiUeLeftCharged->Fill(ptAssoc,deltaPhi);
759     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
760     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
761     Double_t uezT =   ptAssoc/ptTrig;
762     
763     if(uexE < 0.)
764       AliWarning(Form("Careful!!, negative xE %2.2f for left UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
765                       uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
766     
767     fhXEUeLeftCharged->Fill(ptTrig,uexE);
768     if(uexE > 0) fhPtHbpXEUeLeftCharged->Fill(ptTrig,TMath::Log(1/uexE));
769     
770     fhZTUeLeftCharged->Fill(ptTrig,uezT);
771     if(uezT > 0) fhPtHbpZTUeLeftCharged->Fill(ptTrig,TMath::Log(1/uezT));
772     
773     fhDeltaPhiUeLeftCharged->Fill(ptAssoc, deltaPhi);
774   }
775   
776   if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-TMath::Pi()/2))
777   {
778     fhDeltaPhiUeLeftDownCharged->Fill(ptAssoc,deltaPhi);
779     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
780     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
781     
782     if(uexE < 0.)
783       AliWarning(Form("Careful!!, negative xE %2.2f for left-down UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
784                       uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
785     
786     fhXEUeLeftDownCharged->Fill(ptTrig,uexE);
787   }
788   
789   if((deltaPhi>2*fUeDeltaPhiMaxCut) && (deltaPhi <3*TMath::Pi()/2))
790   {
791     fhDeltaPhiUeLeftUpCharged->Fill(ptAssoc,deltaPhi);
792     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
793     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
794     
795     if(uexE < 0.)
796       AliWarning(Form("Careful!!, negative xE %2.2f for left-up UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
797                       uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
798     
799     fhXEUeLeftUpCharged->Fill(ptTrig,uexE);
800   }
801   
802   if((deltaPhi > TMath::Pi()/2) && (deltaPhi < fUeDeltaPhiMaxCut))
803   {
804     fhDeltaPhiUeRightUpCharged->Fill(ptAssoc,deltaPhi);
805     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
806     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
807     
808     if(uexE < 0.)
809       AliWarning(Form("Careful!!, negative xE %2.2f for right-up UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
810                       uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
811     
812     fhXEUeRightUpCharged->Fill(ptTrig,uexE);
813   }
814   
815   if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < TMath::Pi()/2))
816   {
817     fhDeltaPhiUeRightDownCharged->Fill(ptAssoc,deltaPhi);
818     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
819     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
820     
821     if(uexE < 0.)
822       AliWarning(Form("Careful!!, negative xE %2.2f for right-down UE cos(dPhi %2.2f) = %2.2f, check correlation dPhi limits %f to %f",
823              uexE,randomphi*TMath::RadToDeg(),TMath::Cos(randomphi),fDeltaPhiMinCut*TMath::RadToDeg(),fDeltaPhiMaxCut*TMath::RadToDeg()));
824     
825     fhXEUeRightDownCharged->Fill(ptTrig,uexE);
826   }
827
828
829 //_____________________________________________________________________________________________________________________________________
830 void AliAnaParticleHadronCorrelation::FillDecayPhotonCorrelationHistograms(Float_t ptAssoc, Float_t phiAssoc, Bool_t bChargedOrNeutral)
831 {
832   // Do correlation with decay photons of triggered pi0 or eta
833   
834   // Calculate the correlation parameters
835   Float_t ptDecay1 = fDecayMom1.Pt();
836   Float_t ptDecay2 = fDecayMom2.Pt();
837   
838   Float_t zTDecay1 = -100, zTDecay2 = -100;
839   if(ptDecay1 > 0) zTDecay1 = ptAssoc/ptDecay1 ;
840   if(ptDecay2 > 0) zTDecay2 = ptAssoc/ptDecay2 ;
841   
842   Float_t deltaPhiDecay1 = fDecayMom1.Phi()-phiAssoc;
843   if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
844   if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
845   
846   Float_t deltaPhiDecay2 = fDecayMom2.Phi()-phiAssoc;
847   if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
848   if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
849   
850   Float_t xEDecay1   =-zTDecay1*TMath::Cos(deltaPhiDecay1); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
851   Float_t xEDecay2   =-zTDecay2*TMath::Cos(deltaPhiDecay2); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
852   
853   if(bChargedOrNeutral) // correlate with charges
854   {
855     fhDeltaPhiPi0DecayCharged->Fill(ptDecay1, deltaPhiDecay1);
856     fhDeltaPhiPi0DecayCharged->Fill(ptDecay2, deltaPhiDecay2);
857     
858     AliDebug(2,Form("deltaPhoton1 = %f, deltaPhoton2 = %f", deltaPhiDecay1, deltaPhiDecay2));
859     
860     if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
861     {
862       fhZTPi0DecayCharged->Fill(ptDecay1,zTDecay1); 
863       fhXEPi0DecayCharged->Fill(ptDecay1,xEDecay1); 
864     }
865     if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
866     {
867       fhZTPi0DecayCharged->Fill(ptDecay2,zTDecay2);
868       fhXEPi0DecayCharged->Fill(ptDecay2,xEDecay2);
869     }
870   }
871   else // correlate with neutrals
872   {
873     fhDeltaPhiPi0DecayNeutral->Fill(ptDecay1, deltaPhiDecay1);
874     fhDeltaPhiPi0DecayNeutral->Fill(ptDecay2, deltaPhiDecay2);
875     
876     AliDebug(2,Form("deltaPhoton1 = %f, deltaPhoton2 = %f", deltaPhiDecay1, deltaPhiDecay2));
877     
878     if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
879     {
880       fhZTPi0DecayNeutral->Fill(ptDecay1,zTDecay1);
881       fhXEPi0DecayNeutral->Fill(ptDecay1,xEDecay1);
882     }
883     if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
884     {
885       fhZTPi0DecayNeutral->Fill(ptDecay2,zTDecay2);
886       fhXEPi0DecayNeutral->Fill(ptDecay2,xEDecay2);
887     }    
888   }
889 }
890
891 //_____________________________________________________________________________________________________________________________
892 void AliAnaParticleHadronCorrelation::FillNeutralUnderlyingEventSidesHistograms(Float_t ptTrig,   Float_t ptAssoc, 
893                                                                                 Float_t zT,       Float_t hbpZT,
894                                                                                 Float_t deltaPhi)
895 {
896   // Fill underlying event histograms to the left of trigger
897   // Right is the default case
898   
899   Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
900
901   Float_t xE  =-ptAssoc/ptTrig*TMath::Cos(randomphi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
902   Float_t hbpXE = -100;
903   if(xE > 0 ) hbpXE = TMath::Log(1./xE);
904   
905   if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut))
906   {
907     fhDeltaPhiUeLeftNeutral->Fill(ptAssoc, deltaPhi);
908     fhXEUeLeftNeutral      ->Fill(ptTrig , xE);
909     fhPtHbpXEUeLeftNeutral ->Fill(ptTrig , hbpXE);
910     fhZTUeLeftNeutral      ->Fill(ptTrig , zT);
911     fhPtHbpZTUeLeftNeutral ->Fill(ptTrig , hbpZT);
912   }
913
914
915 //______________________________________________________
916 void AliAnaParticleHadronCorrelation::FillEventMixPool()
917 {
918   // Fill the pool with tracks or clusters if requested
919     
920   if ( !DoOwnMix() ) return;
921   
922   FillChargedEventMixPool();
923   
924   // Do the cluster pool filling only if requested
925   // or in case of isolation cut using clusters in the cone.
926   Bool_t isoCase = OnlyIsolated() && (GetIsolationCut()->GetParticleTypeInCone() != AliIsolationCut::kOnlyCharged);
927   
928   if( !fFillNeutralEventMixPool && !isoCase) return;
929   
930   FillNeutralEventMixPool();
931 }
932
933 //_____________________________________________________________
934 void AliAnaParticleHadronCorrelation::FillChargedEventMixPool()
935 {
936   // Mixed event pool filling for tracks
937     
938   if(fUseMixStoredInReader && GetReader()->GetLastTracksMixedEvent() == GetEventNumber())
939   {
940     //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
941     return ; // pool filled previously for another trigger
942   }
943   
944   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
945   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
946   
947   if(!inputHandler) return ;
948     
949   // Do mixing only with MB event (or the chosen mask), if not skip
950   if( !(inputHandler->IsEventSelected( ) & GetReader()->GetMixEventTriggerMask()) ) return ;
951   
952   Int_t eventBin = GetEventMixBin();
953   
954   //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
955   if(eventBin < 0) return;
956   
957   fhEventMBBin->Fill(eventBin);
958   
959   TObjArray * mixEventTracks = new TObjArray;
960   
961   if(fUseMixStoredInReader) 
962   {
963     fListMixTrackEvents[eventBin] = GetReader()->GetListWithMixedEventsForTracks(eventBin);
964   }
965   
966   if(!fListMixTrackEvents[eventBin]) fListMixTrackEvents[eventBin] = new TList();
967   
968   //printf("%s ***** Pool Event bin : %d - nTracks %d\n",GetInputAODName().Data(),eventBin, GetCTSTracks()->GetEntriesFast());
969   
970   TList * pool = fListMixTrackEvents[eventBin];
971   
972   for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
973   {
974     AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
975     
976     fTrackVector.SetXYZ(track->Px(),track->Py(),track->Pz());
977     Float_t pt   = fTrackVector.Pt();
978     
979     //Select only hadrons in pt range
980     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
981     
982     AliAODPWG4Particle * mixedTrack = new AliAODPWG4Particle(track->Px(),track->Py(),track->Pz(),0);
983     mixedTrack->SetDetectorTag(kCTS);
984     mixedTrack->SetChargedBit(track->Charge()>0);
985     mixEventTracks->Add(mixedTrack);
986   }
987   
988   fhNtracksMB->Fill(mixEventTracks->GetEntriesFast(),eventBin);
989   
990   //Set the event number where the last event was added, to avoid double pool filling
991   GetReader()->SetLastTracksMixedEvent(GetEventNumber());
992   
993   //printf("Add event to pool with %d tracks \n ",mixEventTracks->GetEntries());
994   pool->AddFirst(mixEventTracks);
995   mixEventTracks = 0;
996   
997   //printf("Pool size %d, max %d\n",pool->GetSize(), GetNMaxEvMix());
998
999   if(pool->GetSize() > GetNMaxEvMix())
1000   {//Remove last event
1001     TClonesArray * tmp = static_cast<TClonesArray*>(pool->Last()) ;
1002     pool->RemoveLast() ;
1003     delete tmp ;
1004   }
1005 }
1006
1007 //_____________________________________________________________
1008 void AliAnaParticleHadronCorrelation::FillNeutralEventMixPool()
1009 {
1010   // Mixed event pool filling for neutral clusters
1011   // Right now only for EMCAL and in isolation case
1012   
1013   //printf("FillNeutralEventMixPool for %s\n",GetInputAODName().Data());
1014   
1015   if(fUseMixStoredInReader && GetReader()->GetLastCaloMixedEvent() == GetEventNumber())
1016   {
1017     //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
1018     return ; // pool filled previously for another trigger
1019   }
1020   
1021   TObjArray * pl = GetEMCALClusters();
1022   //if (GetAODObjArrayName.Contains("PHOS") )pl    = GetPHOSClusters();
1023   //else                                     pl    = GetEMCALClusters();
1024   
1025   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
1026   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
1027   
1028   if(!inputHandler) return ;
1029     
1030   // Do mixing only with MB event (or the chosen mask), if not skip
1031   if( !(inputHandler->IsEventSelected( ) & GetReader()->GetMixEventTriggerMask()) ) return ;
1032   
1033   Int_t eventBin = GetEventMixBin();
1034   
1035   //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
1036   if(eventBin < 0) return;
1037   
1038   TObjArray * mixEventCalo = new TObjArray;
1039   
1040   if(fUseMixStoredInReader) 
1041   {
1042     fListMixCaloEvents[eventBin] = GetReader()->GetListWithMixedEventsForCalo(eventBin);
1043   }
1044   
1045   if(!fListMixCaloEvents[eventBin]) fListMixCaloEvents[eventBin] = new TList();
1046   
1047   TList * poolCalo = fListMixCaloEvents[eventBin];
1048   
1049   for(Int_t ipr = 0;ipr <  pl->GetEntriesFast() ; ipr ++ )
1050   {
1051     AliVCluster * calo = (AliVCluster *) (pl->At(ipr)) ;
1052   
1053     // remove matched clusters
1054     if( IsTrackMatched( calo, GetReader()->GetInputEvent() ) ) continue ;
1055     
1056     //Cluster momentum calculation
1057     if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1058     {
1059       calo->GetMomentum(fMomentum,GetVertex(0)) ;
1060     }//Assume that come from vertex in straight line
1061     else
1062     {
1063       Double_t vertex[]={0,0,0};
1064       calo->GetMomentum(fMomentum,vertex) ;
1065     }
1066     
1067     Float_t pt = fMomentum.Pt();
1068     //Select only clusters in pt range
1069     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
1070     
1071     AliAODPWG4Particle * mixedCalo = new AliAODPWG4Particle(fMomentum);
1072     mixedCalo->SetDetectorTag(kEMCAL);
1073     mixEventCalo->Add(mixedCalo);
1074   }
1075   
1076   fhNclustersMB->Fill(mixEventCalo->GetEntriesFast(),eventBin);
1077   
1078   //Set the event number where the last event was added, to avoid double pool filling
1079   GetReader()->SetLastCaloMixedEvent(GetEventNumber());
1080   
1081   //printf("Add event to pool with %d clusters \n ",mixEventCalo->GetEntries());
1082   poolCalo->AddFirst(mixEventCalo);
1083   mixEventCalo = 0;
1084   
1085   //printf("Pool size %d, max %d\n",poolCalo->GetSize(), GetNMaxEvMix());
1086   
1087   if(poolCalo->GetSize() > GetNMaxEvMix())
1088   {//Remove last event
1089     TClonesArray * tmp = static_cast<TClonesArray*>(poolCalo->Last()) ;
1090     poolCalo->RemoveLast() ;
1091     delete tmp ;
1092   }  
1093 }
1094
1095 //_________________________________________________________________________________________________________________
1096 Bool_t AliAnaParticleHadronCorrelation::FindLeadingOppositeHadronInWindow(AliAODPWG4ParticleCorrelation * particle)
1097 {
1098   // Select events where the leading charged particle in the opposite hemisphere
1099   // to the trigger particle is in a window centered at 180 from the trigger
1100   
1101   Float_t etaTrig    = particle->Eta();
1102   Float_t ptTrig     = particle->Pt();
1103   Float_t phiTrig    = particle->Phi();
1104   if(phiTrig < 0 ) phiTrig+= TMath::TwoPi();
1105
1106   Float_t ptLeadHad  = 0 ;
1107   Float_t dphiLeadHad= -100 ;
1108   Float_t phiLeadHad = -100 ;
1109   Float_t etaLeadHad = -100 ;
1110   Int_t   nTrack     = 0;
1111
1112   for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
1113   {
1114     AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
1115     
1116     fTrackVector.SetXYZ(track->Px(),track->Py(),track->Pz());
1117     
1118     Float_t pt   = fTrackVector.Pt();
1119     Float_t phi  = fTrackVector.Phi() ;
1120     if(phi < 0 ) phi+= TMath::TwoPi();
1121     
1122     Float_t deltaPhi = phiTrig-phi;
1123     //
1124     // Calculate deltaPhi shift so that for the particles on the opposite side
1125     // it is defined between 90 and 270 degrees
1126     // Shift [-360,-90]  to [0, 270]
1127     // and [270,360] to [-90,0]
1128     if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
1129     if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
1130
1131     if(pt > ptLeadHad && deltaPhi > TMath::PiOver2()) // in opposite hemisphere
1132     {
1133       ptLeadHad  = pt ;
1134       phiLeadHad = phi;
1135       dphiLeadHad= deltaPhi;
1136       etaLeadHad = fTrackVector.Eta();
1137       nTrack++;
1138     }
1139   }// track loop
1140   
1141   if(fFillLeadHadOppositeHisto)
1142   {
1143     if(nTrack == 0)
1144     {
1145       fhPtNoLeadingOppositeHadron    ->Fill(ptTrig);
1146       fhEtaPhiNoLeadingOppositeHadron->Fill(etaTrig,phiTrig);
1147     }
1148     else
1149     {
1150       fhPtLeadingOppositeHadron       ->Fill(ptTrig,  ptLeadHad);
1151       fhPtDiffPhiLeadingOppositeHadron->Fill(ptTrig,dphiLeadHad);
1152       fhPtDiffEtaLeadingOppositeHadron->Fill(ptTrig, etaLeadHad-etaTrig);
1153     }
1154   }
1155   
1156   
1157     AliDebug(1,Form("pT %2.2f, phi %2.2f, eta %2.2f, nTracks away %d, total tracks %d",
1158                     ptLeadHad,phiLeadHad*TMath::RadToDeg(),etaLeadHad,nTrack, GetTrackMultiplicity()));
1159     AliDebug(1,Form("\t pT trig %2.2f, Dphi (trigger-hadron) %2.2f, Deta (trigger-hadron) %2.2f",
1160            ptTrig, dphiLeadHad*TMath::RadToDeg(), etaLeadHad-etaTrig));
1161     AliDebug(1,Form("\t cuts pT: min %2.2f, max %2.2f; DPhi: min %2.2f, max %2.2f",
1162                     fMinLeadHadPt,fMaxLeadHadPt,fMinLeadHadPhi*TMath::RadToDeg(),fMaxLeadHadPhi*TMath::RadToDeg()));
1163   
1164   
1165   // reject the trigger if the leading hadron is not in the requested pt or phi window and
1166   
1167   if( nTrack == 0 ) return kFALSE; // No track found in opposite hemisphere
1168   
1169   if( ptLeadHad < fMinLeadHadPt || ptLeadHad > fMaxLeadHadPt ) return kFALSE;
1170   
1171   //printf("Accept leading hadron pT \n");
1172   
1173   if( dphiLeadHad < fMinLeadHadPhi || dphiLeadHad > fMaxLeadHadPhi ) return kFALSE;
1174   
1175   //printf("Accept leading hadron phi \n");
1176   
1177   
1178   return kTRUE ;
1179 }
1180
1181 //____________________________________________________________
1182 TObjString* AliAnaParticleHadronCorrelation::GetAnalysisCuts()
1183 {
1184   //Save parameters used for analysis
1185   TString parList ; //this will be list of parameters used for this analysis.
1186   const Int_t buffersize = 560;
1187   char onePar[buffersize] ;
1188   
1189   snprintf(onePar,buffersize,"--- AliAnaPaticleHadronCorrelation ---:") ;
1190   parList+=onePar ;     
1191   snprintf(onePar,buffersize," %3.2f < Pt associated < %3.2f; ", fMinAssocPt,   fMaxAssocPt) ;
1192   parList+=onePar ;
1193   snprintf(onePar,buffersize," %3.2f < Phi trigger particle-Hadron < %3.2f; ",    fDeltaPhiMinCut,   fDeltaPhiMaxCut) ;
1194   parList+=onePar ;
1195   snprintf(onePar,buffersize," %3.2f < Phi trigger particle-UeHadron <  %3.2f; ", fUeDeltaPhiMinCut, fUeDeltaPhiMaxCut) ;
1196   parList+=onePar ;
1197   snprintf(onePar,buffersize,"Isolated Trigger?  %d;",    fSelectIsolated) ;
1198   parList+=onePar ;
1199   snprintf(onePar,buffersize,"Several UE?  %d;",          fMakeSeveralUE) ;
1200   parList+=onePar ;
1201   snprintf(onePar,buffersize,"Name of AOD Pi0 Branch %s;", fPi0AODBranchName.Data());
1202   parList+=onePar ;
1203   snprintf(onePar,buffersize,"Do Decay-hadron correlation ?  pi0 %d, decay %d;", fPi0Trigger, fDecayTrigger) ;
1204   parList+=onePar ;
1205   snprintf(onePar,buffersize,"Select absolute leading for cluster triggers ? %d or Near Side Leading %d;",
1206            fMakeAbsoluteLeading, fMakeNearSideLeading) ;
1207   parList+=onePar ;
1208   snprintf(onePar,buffersize,"Associated particle pt bins  %d: ", fNAssocPtBins) ;
1209   parList+=onePar ;
1210   for (Int_t ibin = 0; ibin<fNAssocPtBins; ibin++) {
1211     snprintf(onePar,buffersize,"bin %d = [%2.1f,%2.1f];", ibin, fAssocPtBinLimit[ibin], fAssocPtBinLimit[ibin+1]) ;
1212   }
1213   parList+=onePar ;
1214   
1215   //Get parameters set in base class.
1216   parList += GetBaseParametersList() ;
1217   
1218   //Get parameters set in FiducialCut class (not available yet)
1219   //parlist += GetFidCut()->GetFidCutParametersList() 
1220   
1221   return new TObjString(parList) ;  
1222   
1223
1224
1225 //________________________________________________________________
1226 TList *  AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
1227 {
1228   // Create histograms to be saved in output file and
1229   // store them in fOutputContainer
1230   
1231   TList * outputContainer = new TList() ;
1232   outputContainer->SetName("CorrelationHistos") ;
1233   
1234   Int_t   nptbins = GetHistogramRanges()->GetHistoPtBins(); Int_t  nphibins = GetHistogramRanges()->GetHistoPhiBins(); Int_t   netabins = GetHistogramRanges()->GetHistoEtaBins();
1235   Float_t ptmax   = GetHistogramRanges()->GetHistoPtMax();  Float_t phimax  = GetHistogramRanges()->GetHistoPhiMax();  Float_t etamax   = GetHistogramRanges()->GetHistoEtaMax();
1236   Float_t ptmin   = GetHistogramRanges()->GetHistoPtMin();  Float_t phimin  = GetHistogramRanges()->GetHistoPhiMin();  Float_t etamin   = GetHistogramRanges()->GetHistoEtaMin();
1237
1238   Int_t  ndeltaphibins = GetHistogramRanges()->GetHistoDeltaPhiBins(); Int_t   ndeltaetabins = GetHistogramRanges()->GetHistoDeltaEtaBins();
1239   Float_t deltaphimax  = GetHistogramRanges()->GetHistoDeltaPhiMax();  Float_t deltaetamax   = GetHistogramRanges()->GetHistoDeltaEtaMax();
1240   Float_t deltaphimin  = GetHistogramRanges()->GetHistoDeltaPhiMin();  Float_t deltaetamin   = GetHistogramRanges()->GetHistoDeltaEtaMin();
1241
1242   Int_t ntrbins = GetHistogramRanges()->GetHistoTrackMultiplicityBins(); Int_t nclbins = GetHistogramRanges()->GetHistoNClustersBins();
1243   Int_t trmax   = GetHistogramRanges()->GetHistoTrackMultiplicityMax();  Int_t clmax   = GetHistogramRanges()->GetHistoNClustersMax();
1244   Int_t trmin   = GetHistogramRanges()->GetHistoTrackMultiplicityMin();  Int_t clmin   = GetHistogramRanges()->GetHistoNClustersMin();
1245
1246   Int_t  nxeztbins = GetHistogramRanges()->GetHistoRatioBins();  Int_t  nhbpbins = GetHistogramRanges()->GetHistoHBPBins();
1247   Float_t xeztmax  = GetHistogramRanges()->GetHistoRatioMax();   Float_t hbpmax  = GetHistogramRanges()->GetHistoHBPMax();
1248   Float_t xeztmin  = GetHistogramRanges()->GetHistoRatioMin();   Float_t hbpmin  = GetHistogramRanges()->GetHistoHBPMin();
1249   
1250   Int_t nMixBins = GetNCentrBin()*GetNZvertBin()*GetNRPBin();
1251   
1252   TString nameMC[]     = {"Photon","Pi0","Pi0Decay","EtaDecay","OtherDecay","Electron","Hadron","Pi0DecayLostPair"};
1253   TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
1254   
1255   // For vz dependent histograms, if option ON
1256   Int_t   nz  = 1  ;
1257   if(fCorrelVzBin) nz = GetNZvertBin();
1258   TString sz  = "" ;
1259   TString tz  = "" ;
1260   
1261   // Fill histograms for neutral clusters in mixing?
1262   Bool_t isoCase = OnlyIsolated() && (GetIsolationCut()->GetParticleTypeInCone() != AliIsolationCut::kOnlyCharged);
1263   Bool_t neutralMix = fFillNeutralEventMixPool || isoCase ;
1264   
1265   fhPtTriggerInput  = new TH1F("hPtTriggerInput","Input trigger #it{p}_{T}", nptbins,ptmin,ptmax);
1266   fhPtTriggerInput->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1267   outputContainer->Add(fhPtTriggerInput);
1268   
1269   if( fM02MaxCut > 0 && fM02MinCut > 0 )
1270   {
1271     fhPtTriggerSSCut  = new TH1F("hPtTriggerSSCut","Trigger #it{p}_{T} after #lambda^{2}_{0} cut", nptbins,ptmin,ptmax);
1272     fhPtTriggerSSCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1273     outputContainer->Add(fhPtTriggerSSCut);
1274   }
1275   
1276   if( OnlyIsolated() )
1277   {
1278     fhPtTriggerIsoCut  = new TH1F("hPtTriggerIsoCut","Trigger #it{p}_{T} after isolation (and #lambda^{2}_{0} cut)", nptbins,ptmin,ptmax);
1279     fhPtTriggerIsoCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1280     outputContainer->Add(fhPtTriggerIsoCut);
1281   }
1282   
1283   fhPtTriggerFidCut  = new TH1F("hPtTriggerFidCut","Trigger #it{p}_{T} after fiducial (isolation and #lambda^{2}_{0}) cut", nptbins,ptmin,ptmax);
1284   fhPtTriggerFidCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1285   outputContainer->Add(fhPtTriggerFidCut);
1286   
1287   fhPtTrigger  = new TH1F("hPtTrigger","#it{p}_{T} distribution of trigger particles (after opposite hadron leading cut and rest)", nptbins,ptmin,ptmax);
1288   fhPtTrigger->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1289   outputContainer->Add(fhPtTrigger);
1290   
1291   if(IsDataMC())
1292   {
1293     for(Int_t i=0; i < fgkNmcTypes; i++)
1294     {
1295       fhPtTriggerMC[i]  = new TH1F(Form("hPtTrigger_MC%s",nameMC[i].Data()),
1296                                    Form("#it{p}_{T} distribution of trigger particles, trigger origin is %s",nameMC[i].Data()),
1297                                    nptbins,ptmin,ptmax);
1298       fhPtTriggerMC[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1299       outputContainer->Add(fhPtTriggerMC[i]);
1300     }
1301   }
1302   
1303   if(fDecayTrigger)
1304   {
1305     for(Int_t ibit = 0; ibit < fNDecayBits; ibit++)
1306     {
1307       fhPtDecayTrigger[ibit]  = new TH1F(Form("hPtDecayTrigger_bit%d",fDecayBits[ibit]),
1308                                          Form("#it{p}_{T} distribution of trigger particles, decay Bit %d",fDecayBits[ibit]),
1309                                          nptbins,ptmin,ptmax);
1310       fhPtDecayTrigger[ibit]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1311       outputContainer->Add(fhPtDecayTrigger[ibit]);
1312       
1313       if(IsDataMC())
1314       {
1315         for(Int_t i=0; i < fgkNmcTypes; i++)
1316         {
1317           fhPtDecayTriggerMC[ibit][i]  = new TH1F(Form("hPtDecayTrigger_bit%d_MC%s",fDecayBits[ibit], nameMC[i].Data()),
1318                                             Form("#it{p}_{T} distribution of trigger particles, decay Bit %d, trigger origin is %s",fDecayBits[ibit], nameMC[i].Data()),
1319                                             nptbins,ptmin,ptmax);
1320           fhPtDecayTriggerMC[ibit][i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1321           outputContainer->Add(fhPtDecayTriggerMC[ibit][i]);
1322         }
1323       }
1324     }
1325   }
1326   
1327   if(fCorrelVzBin)
1328   {
1329     fhPtTriggerVzBin  = new TH2F("hPtTriggerVzBin","#it{p}_{T} distribution of trigger particles vs vz bin", nptbins,ptmin,ptmax,GetNZvertBin(),0,GetNZvertBin());
1330     fhPtTriggerVzBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1331     fhPtTriggerVzBin->SetYTitle("#it{v}_{#it{z}} bin");
1332     outputContainer->Add(fhPtTriggerVzBin);
1333   }
1334   
1335   fhPtTriggerBin  = new TH2F ("hPtTriggerBin","#it{p}_{T} distribution of trigger particles", nptbins,ptmin,ptmax,nMixBins,0,nMixBins);
1336   fhPtTriggerBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1337   fhPtTriggerBin->SetYTitle("Bin");
1338   outputContainer->Add(fhPtTriggerBin);
1339   
1340   fhPhiTrigger  = new TH2F ("hPhiTrigger","#phi distribution of trigger Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
1341   fhPhiTrigger->SetYTitle("#phi (rad)");
1342   outputContainer->Add(fhPhiTrigger);
1343   
1344   fhEtaTrigger  = new TH2F ("hEtaTrigger","#eta distribution of trigger",nptbins,ptmin,ptmax, netabins,etamin,etamax);
1345   fhEtaTrigger->SetYTitle("#eta ");
1346   outputContainer->Add(fhEtaTrigger);
1347   
1348   if(IsHighMultiplicityAnalysisOn())
1349   {
1350     fhPtTriggerCentrality   = new TH2F("hPtTriggerCentrality","Trigger particle #it{p}_{T} vs centrality",nptbins,ptmin,ptmax,100,0.,100) ;
1351     fhPtTriggerCentrality->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1352     fhPtTriggerCentrality->SetYTitle("Centrality (%)");
1353     outputContainer->Add(fhPtTriggerCentrality) ;
1354     
1355     fhPtTriggerEventPlane  = new TH2F("hPtTriggerEventPlane","Trigger particle #it{p}_{T} vs event plane angle",nptbins,ptmin,ptmax, 100,0.,TMath::Pi()) ;
1356     fhPtTriggerEventPlane->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1357     fhPtTriggerEventPlane->SetXTitle("EP angle (rad)");
1358     outputContainer->Add(fhPtTriggerEventPlane) ;
1359     
1360     fhTriggerEventPlaneCentrality  = new TH2F("hTriggerEventPlaneCentrality","Trigger particle centrality vs event plane angle",100,0.,100,100,0.,TMath::Pi()) ;
1361     fhTriggerEventPlaneCentrality->SetXTitle("Centrality (%)");
1362     fhTriggerEventPlaneCentrality->SetYTitle("EP angle (rad)");
1363     outputContainer->Add(fhTriggerEventPlaneCentrality) ;
1364   }
1365
1366   // Leading hadron in oposite side
1367   if(fFillLeadHadOppositeHisto)
1368   {
1369     fhPtLeadingOppositeHadron  = new TH2F("hPtTriggerPtLeadingOppositeHadron","Leading hadron opposite to trigger vs trigger #it{p}_{T}",
1370                                           nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1371     fhPtLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1372     fhPtLeadingOppositeHadron->SetYTitle("#it{p}_{T}^{lead hadron} (GeV/#it{c})");
1373     outputContainer->Add(fhPtLeadingOppositeHadron);
1374
1375     fhPtNoLeadingOppositeHadron  = new TH1F("hPtTriggerNoLeadingOppositeHadron","No Leading hadron opposite to trigger #it{p}_{T}",
1376                                           nptbins,ptmin,ptmax);
1377     fhPtNoLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1378     outputContainer->Add(fhPtNoLeadingOppositeHadron);
1379
1380     fhEtaPhiNoLeadingOppositeHadron  = new TH2F("hEtaPhiTriggerNoLeadingOppositeHadron","No Leading hadron opposite to trigger #eta:#phi",
1381                                             netabins,etamin,etamax,nphibins,phimin,phimax);
1382     fhEtaPhiNoLeadingOppositeHadron->SetXTitle("#eta");
1383     fhEtaPhiNoLeadingOppositeHadron->SetYTitle("#phi");
1384     outputContainer->Add(fhEtaPhiNoLeadingOppositeHadron);
1385
1386     
1387     fhPtDiffPhiLeadingOppositeHadron  = new TH2F("hPtTriggerDiffPhiTriggerLeadingOppositeHadron","#phi_{trigger}-#phi_{leading opposite hadron} vs #it{p}_{T}^{trig}",
1388                                                  nptbins,ptmin,ptmax,ndeltaphibins,deltaphimin,deltaphimax);
1389     fhPtDiffPhiLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1390     fhPtDiffPhiLeadingOppositeHadron->SetYTitle("#phi_{trigger}-#phi_{leading opposite hadron} (rad)");
1391     outputContainer->Add(fhPtDiffPhiLeadingOppositeHadron);
1392     
1393     fhPtDiffEtaLeadingOppositeHadron  = new TH2F("hPtTriggerDiffEtaTriggerPhiLeadingOppositeHadron","#eta_{trigger}-#eta_{leading opposite hadron} vs #it{p}_{T}^{trig}",
1394                                                  nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1395     fhPtDiffEtaLeadingOppositeHadron->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1396     fhPtDiffEtaLeadingOppositeHadron->SetYTitle("#eta_{trigger}-#eta_{leading opposite hadron}");
1397     outputContainer->Add(fhPtDiffEtaLeadingOppositeHadron);
1398   }
1399   
1400   //Correlation with charged hadrons
1401   
1402   fhDeltaPhiDeltaEtaCharged  = new TH2F
1403   ("hDeltaPhiDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #phi_{trigger} - #phi_{h^{#pm}}",
1404    ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1405   fhDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi (rad)");
1406   fhDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
1407   
1408   fhDeltaPhiDeltaEtaChargedPtA3GeV  = new TH2F
1409   ("hDeltaPhiDeltaEtaChargedPtA3GeV","#eta_{trigger} - #eta_{h^{#pm}} vs #phi_{trigger} - #phi_{h^{#pm}, #it{p}_{TA}>3 GeV/#it{c}}",
1410    ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1411   fhDeltaPhiDeltaEtaChargedPtA3GeV->SetXTitle("#Delta #phi (rad)");
1412   fhDeltaPhiDeltaEtaChargedPtA3GeV->SetYTitle("#Delta #eta");
1413   
1414   fhPhiCharged  = new TH2F
1415   ("hPhiCharged","#phi_{h^{#pm}}  vs #it{p}_{T #pm}",
1416    nptbins,ptmin,ptmax,180,0,TMath::TwoPi());
1417   fhPhiCharged->SetYTitle("#phi_{h^{#pm}} (rad)");
1418   fhPhiCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
1419   
1420   fhEtaCharged  = new TH2F
1421   ("hEtaCharged","#eta_{h^{#pm}}  vs #it{p}_{T #pm}",
1422    nptbins,ptmin,ptmax,100,-1.,1.);
1423   fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
1424   fhEtaCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
1425   
1426   fhDeltaPhiCharged  = new TH2F
1427   ("hDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}",
1428    nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1429   fhDeltaPhiCharged->SetYTitle("#Delta #phi (rad)");
1430   fhDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1431   
1432   fhDeltaPhiChargedPtA3GeV  = new TH2F
1433   ("hDeltaPhiChargedPtA3GeV","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}",
1434    nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1435   fhDeltaPhiChargedPtA3GeV->SetYTitle("#Delta #phi (rad)");
1436   fhDeltaPhiChargedPtA3GeV->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1437   
1438   
1439   fhDeltaPhiChargedPt  = new TH2F
1440   ("hDeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs #it{p}_{T h^{#pm}}",
1441    nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1442   fhDeltaPhiChargedPt->SetYTitle("#Delta #phi (rad)");
1443   fhDeltaPhiChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1444   
1445   fhDeltaEtaCharged  = new TH2F
1446   ("hDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}",
1447    nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1448   fhDeltaEtaCharged->SetYTitle("#Delta #eta");
1449   fhDeltaEtaCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1450   
1451   fhDeltaEtaChargedPtA3GeV  = new TH2F
1452   ("hDeltaEtaChargedPtA3GeV","#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}",
1453    nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1454   fhDeltaEtaChargedPtA3GeV->SetYTitle("#Delta #eta");
1455   fhDeltaEtaChargedPtA3GeV->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1456   
1457   fhXECharged  =
1458   new TH2F("hXECharged","#it{x}_{#it{E}} for charged tracks",
1459            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1460   fhXECharged->SetYTitle("#it{x}_{#it{E}}");
1461   fhXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1462   
1463   fhXECharged_Cone2  =
1464   new TH2F("hXECharged_Cone2","#it{x}_{#it{E}} for charged tracks in cone 2 (5#pi/6-7#pi/6)",
1465            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1466   fhXECharged_Cone2->SetYTitle("#it{x}_{#it{E}}");
1467   fhXECharged_Cone2->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1468   
1469   fhXEPosCharged  =
1470   new TH2F("hXEPositiveCharged","#it{x}_{#it{E}} for positive charged tracks",
1471            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1472   fhXEPosCharged->SetYTitle("#it{x}_{#it{E}}");
1473   fhXEPosCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1474   
1475   fhXENegCharged  =
1476   new TH2F("hXENegativeCharged","#it{x}_{#it{E}} for negative charged tracks",
1477            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1478   fhXENegCharged->SetYTitle("#it{x}_{#it{E}}");
1479   fhXENegCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1480   
1481   fhPtHbpXECharged  =
1482   new TH2F("hHbpXECharged","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
1483            nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1484   fhPtHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1485   fhPtHbpXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1486   
1487   fhPtHbpXECharged_Cone2  =
1488   new TH2F("hHbpXECharged_Cone2","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons in cone 2 (5#pi/6-7#pi/6)",
1489            nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1490   fhPtHbpXECharged_Cone2->SetYTitle("ln(1/#it{x}_{#it{E}})");
1491   fhPtHbpXECharged_Cone2->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1492   
1493   fhZTCharged  =
1494   new TH2F("hZTCharged","#it{z}_{T} for charged tracks",
1495            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1496   fhZTCharged->SetYTitle("#it{z}_{T}");
1497   fhZTCharged->SetXTitle("#it{p}_{T trigger}");
1498   
1499   fhZTPosCharged  =
1500   new TH2F("hZTPositiveCharged","#it{z}_{T} for positive charged tracks",
1501            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1502   fhZTPosCharged->SetYTitle("#it{z}_{T}");
1503   fhZTPosCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1504   
1505   fhZTNegCharged  =
1506   new TH2F("hZTNegativeCharged","#it{z}_{T} for negative charged tracks",
1507            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1508   fhZTNegCharged->SetYTitle("#it{z}_{T}");
1509   fhZTNegCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1510   
1511   fhPtHbpZTCharged  =
1512   new TH2F("hHbpZTCharged","#xi = ln(1/#it{z}_{T}) with charged hadrons",
1513            nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1514   fhPtHbpZTCharged->SetYTitle("ln(1/#it{z}_{T})");
1515   fhPtHbpZTCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1516   
1517   fhPtTrigPout  =
1518   new TH2F("hPtTrigPout","Pout with triggers",
1519            nptbins,ptmin,ptmax,nptbins,-1.*ptmax/2.,ptmax/2.);
1520   fhPtTrigPout->SetYTitle("#it{p}_{out} (GeV/#it{c})");
1521   fhPtTrigPout->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1522   
1523   fhPtTrigCharged  =
1524   new TH2F("hPtTrigCharged","trigger and charged tracks pt distribution",
1525            nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1526   fhPtTrigCharged->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1527   fhPtTrigCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1528   
1529   outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
1530   outputContainer->Add(fhDeltaPhiDeltaEtaChargedPtA3GeV);
1531   outputContainer->Add(fhPhiCharged) ;
1532   outputContainer->Add(fhEtaCharged) ;
1533   outputContainer->Add(fhDeltaPhiCharged) ;
1534   outputContainer->Add(fhDeltaPhiChargedPtA3GeV) ;
1535   outputContainer->Add(fhDeltaEtaCharged) ;
1536   outputContainer->Add(fhDeltaEtaChargedPtA3GeV) ;
1537   outputContainer->Add(fhDeltaPhiChargedPt) ;
1538   
1539   outputContainer->Add(fhXECharged) ;
1540   outputContainer->Add(fhXECharged_Cone2) ;
1541   
1542   if(IsDataMC())
1543   {
1544     for(Int_t i=0; i < fgkNmcTypes; i++)
1545     {
1546       
1547       fhDeltaPhiChargedMC[i]  = new TH2F(Form("hDeltaPhiCharged_MC%s",nameMC[i].Data()),
1548                                          Form("#Delta #phi for charged tracks, trigger origin is %s",nameMC[i].Data()),
1549                                          nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
1550       fhDeltaPhiChargedMC[i]->SetYTitle("#it{x}_{#it{E}}");
1551       fhDeltaPhiChargedMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1552       outputContainer->Add(fhDeltaPhiChargedMC[i]) ;
1553       
1554       fhXEChargedMC[i]  = new TH2F(Form("hXECharged_MC%s",nameMC[i].Data()),
1555                                    Form("#it{x}_{#it{E}} for charged tracks, trigger origin is %s",nameMC[i].Data()),
1556                                    nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1557       fhXEChargedMC[i]->SetYTitle("#it{x}_{#it{E}}");
1558       fhXEChargedMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1559       outputContainer->Add(fhXEChargedMC[i]) ;
1560     }
1561   }
1562   
1563   outputContainer->Add(fhXEPosCharged) ;
1564   outputContainer->Add(fhXENegCharged) ;
1565   outputContainer->Add(fhPtHbpXECharged) ;
1566   outputContainer->Add(fhPtHbpXECharged_Cone2) ;
1567   
1568   outputContainer->Add(fhZTCharged) ;
1569   outputContainer->Add(fhZTPosCharged) ;
1570   outputContainer->Add(fhZTNegCharged) ;
1571   outputContainer->Add(fhPtHbpZTCharged) ;
1572   
1573   outputContainer->Add(fhPtTrigPout) ;
1574   outputContainer->Add(fhPtTrigCharged) ;
1575   
1576   TString right = "";
1577   if(fMakeSeveralUE) right = "Right";
1578   
1579   fhUePart  =  new TH1F("hUePart","UE particles distribution vs pt trig",
1580                         nptbins,ptmin,ptmax);
1581   fhUePart->SetYTitle("dNch");
1582   fhUePart->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1583   
1584   fhDeltaPhiUeChargedPt  = new TH2F
1585   (Form("hDeltaPhiUe%sChargedPt",right.Data()),"#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}}",
1586    nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1587   fhDeltaPhiUeChargedPt->SetYTitle("#Delta #phi (rad)");
1588   fhDeltaPhiUeChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1589   
1590   fhXEUeCharged  =
1591   new TH2F(Form("hXEUeCharged%s",right.Data()),"#it{x}_{#it{E}} for Underlying Event",
1592            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1593   fhXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
1594   fhXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1595   
1596   fhPtHbpXEUeCharged  =
1597   new TH2F(Form("hHbpXEUeCharged%s",right.Data()),"#xi = ln(1/#it{x}_{#it{E}}) for Underlying Event",
1598            nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1599   fhPtHbpXEUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1600   fhPtHbpXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1601   
1602   fhZTUeCharged  =
1603   new TH2F(Form("hZTUeCharged%s",right.Data()),"#it{z}_{T} for Underlying Event",
1604            nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1605   fhZTUeCharged->SetYTitle("#it{z}_{T}");
1606   fhZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1607   
1608   fhPtHbpZTUeCharged  =
1609   new TH2F(Form("hHbpZTUeCharged%s",right.Data()),"#xi = ln(1/#it{z}_{T}) for Underlying Event",
1610            nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1611   fhPtHbpZTUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1612   fhPtHbpZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1613   
1614   outputContainer->Add(fhUePart);
1615   outputContainer->Add(fhDeltaPhiUeChargedPt) ;
1616   outputContainer->Add(fhXEUeCharged) ;
1617   outputContainer->Add(fhPtHbpXEUeCharged) ;
1618   outputContainer->Add(fhZTUeCharged) ;
1619   outputContainer->Add(fhPtHbpZTUeCharged) ;
1620   
1621   if(fMakeSeveralUE)
1622   {
1623     fhDeltaPhiUeLeftCharged  = new TH2F
1624     ("hDeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left side range of trigger particles",
1625      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1626     fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #phi (rad)");
1627     fhDeltaPhiUeLeftCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1628     outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
1629     
1630     fhDeltaPhiUeLeftUpCharged  = new TH2F
1631     ("hDeltaPhiUeLeftUpChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left Up side range of trigger particles",
1632      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1633     fhDeltaPhiUeLeftUpCharged->SetYTitle("#Delta #phi (rad)");
1634     fhDeltaPhiUeLeftUpCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1635     outputContainer->Add(fhDeltaPhiUeLeftUpCharged) ;
1636     
1637     fhDeltaPhiUeRightUpCharged  = new TH2F
1638     ("hDeltaPhiUeRightUpChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right Up side range of trigger particles",
1639      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1640     fhDeltaPhiUeRightUpCharged->SetYTitle("#Delta #phi (rad)");
1641     fhDeltaPhiUeRightUpCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1642     outputContainer->Add(fhDeltaPhiUeRightUpCharged) ;
1643     
1644     fhDeltaPhiUeLeftDownCharged  = new TH2F
1645     ("hDeltaPhiUeLeftDownChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left Down side range of trigger particles",
1646      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1647     fhDeltaPhiUeLeftDownCharged->SetYTitle("#Delta #phi (rad)");
1648     fhDeltaPhiUeLeftDownCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1649     outputContainer->Add(fhDeltaPhiUeLeftDownCharged) ;
1650     
1651     fhDeltaPhiUeRightDownCharged  = new TH2F
1652     ("hDeltaPhiUeRightDownChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right Down side range of trigger particles",
1653      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1654     fhDeltaPhiUeRightDownCharged->SetYTitle("#Delta #phi (rad)");
1655     fhDeltaPhiUeRightDownCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1656     outputContainer->Add(fhDeltaPhiUeRightDownCharged) ;
1657     
1658     fhXEUeLeftCharged  =
1659     new TH2F("hXEUeChargedLeft","#it{x}_{#it{E}} with UE left side of trigger",
1660              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1661     fhXEUeLeftCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
1662     fhXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1663     outputContainer->Add(fhXEUeLeftCharged) ;
1664     
1665     fhXEUeLeftUpCharged  =
1666     new TH2F("hXEUeChargedLeftUp","#it{x}_{#it{E}} with UE left Up side of trigger",
1667              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1668     fhXEUeLeftUpCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
1669     fhXEUeLeftUpCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1670     outputContainer->Add(fhXEUeLeftUpCharged) ;
1671     
1672     fhXEUeRightUpCharged  =
1673     new TH2F("hXEUeChargedRightUp","#it{x}_{#it{E} h^{#pm}} with UE right Up side of trigger",
1674              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1675     fhXEUeRightUpCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
1676     fhXEUeRightUpCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1677     outputContainer->Add(fhXEUeRightUpCharged) ;
1678     
1679     fhXEUeLeftDownCharged  =
1680     new TH2F("hXEUeChargedLeftDown","#it{x}_{#it{E}} with UE left Down side of trigger",
1681              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1682     fhXEUeLeftDownCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
1683     fhXEUeLeftDownCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1684     outputContainer->Add(fhXEUeLeftDownCharged) ;
1685     
1686     fhXEUeRightDownCharged  =
1687     new TH2F("hXEUeChargedRightDown","#it{x}_{#it{E} h^{#pm}} with UE right Down side of trigger",
1688              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1689     fhXEUeRightDownCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
1690     fhXEUeRightDownCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1691     outputContainer->Add(fhXEUeRightDownCharged) ;
1692     
1693     fhPtHbpXEUeLeftCharged  =
1694     new TH2F("hHbpXEUeChargedLeft","#xi = ln(1/#it{x}_{#it{E}}) with charged UE left side of trigger",
1695              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1696     fhPtHbpXEUeLeftCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1697     fhPtHbpXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1698     outputContainer->Add(fhPtHbpXEUeLeftCharged) ;
1699     
1700     fhZTUeLeftCharged  =
1701     new TH2F("hZTUeChargedLeft","#it{z}_{trigger h^{#pm}} = #it{p}_{T Ueh^{#pm}} / #it{p}_{T trigger} with UE left side of trigger",
1702              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1703     fhZTUeLeftCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
1704     fhZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1705     outputContainer->Add(fhZTUeLeftCharged) ;
1706     
1707     fhPtHbpZTUeLeftCharged  =
1708     new TH2F("hHbpZTUeChargedLeft","#xi = ln(1/#it{z}_{T}) with charged UE left side of trigger",
1709              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
1710     fhPtHbpZTUeLeftCharged->SetYTitle("ln(1/#it{z}_{T})");
1711     fhPtHbpZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1712     outputContainer->Add(fhPtHbpZTUeLeftCharged) ;
1713   }
1714   
1715   if(IsPileUpAnalysisOn())
1716   {
1717     fhDeltaPhiChargedOtherBC  = new TH2F
1718     ("hDeltaPhiChargedOtherBC","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, track BC!=0",
1719      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1720     fhDeltaPhiChargedOtherBC->SetYTitle("#Delta #phi (rad)");
1721     fhDeltaPhiChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1722     
1723     fhDeltaPhiChargedPtA3GeVOtherBC  = new TH2F
1724     ("hDeltaPhiChargedPtA3GeVOtherBC","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC!=0",
1725      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1726     fhDeltaPhiChargedPtA3GeVOtherBC->SetYTitle("#Delta #phi (rad)");
1727     fhDeltaPhiChargedPtA3GeVOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1728     
1729     fhPtTrigChargedOtherBC  =
1730     new TH2F("hPtTrigChargedOtherBC","trigger and charged tracks pt distribution, track BC!=0",
1731              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1732     fhPtTrigChargedOtherBC->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1733     fhPtTrigChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1734     
1735     fhXEChargedOtherBC  =
1736     new TH2F("hXEChargedOtherBC","#it{x}_{#it{E}} for charged tracks, track BC!=0",
1737              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1738     fhXEChargedOtherBC->SetYTitle("#it{x}_{#it{E}}");
1739     fhXEChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1740     
1741     fhXEUeChargedOtherBC  =
1742     new TH2F("hXEUeChargedOtherBC","#it{x}_{#it{E}} for Underlying Event, track BC!=0",
1743              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1744     fhXEUeChargedOtherBC->SetYTitle("#it{x}_{#it{E}}");
1745     fhXEUeChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1746     
1747     fhZTChargedOtherBC  =
1748     new TH2F("hZTChargedOtherBC","#it{z}_{T} for charged tracks, track BC!=0",
1749              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1750     fhZTChargedOtherBC->SetYTitle("#it{z}_{T}");
1751     fhZTChargedOtherBC->SetXTitle("#it{p}_{T trigger}");
1752     
1753     fhZTUeChargedOtherBC  =
1754     new TH2F("hZTUeChargedOtherBC","#it{z}_{T} for Underlying Event, track BC!=0",
1755              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1756     fhZTUeChargedOtherBC->SetYTitle("#it{z}_{T}");
1757     fhZTUeChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1758     
1759     outputContainer->Add(fhDeltaPhiChargedOtherBC) ;
1760     outputContainer->Add(fhDeltaPhiChargedPtA3GeVOtherBC) ;
1761     outputContainer->Add(fhXEChargedOtherBC) ;
1762     outputContainer->Add(fhXEUeChargedOtherBC) ;
1763     outputContainer->Add(fhZTChargedOtherBC) ;
1764     outputContainer->Add(fhZTUeChargedOtherBC) ;
1765     outputContainer->Add(fhPtTrigChargedOtherBC) ;
1766     
1767     fhDeltaPhiChargedBC0  = new TH2F
1768     ("hDeltaPhiChargedBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, track BC==0",
1769      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1770     fhDeltaPhiChargedBC0->SetYTitle("#Delta #phi (rad)");
1771     fhDeltaPhiChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1772     
1773     fhDeltaPhiChargedPtA3GeVBC0  = new TH2F
1774     ("hDeltaPhiChargedPtA3GeVBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC==0",
1775      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1776     fhDeltaPhiChargedPtA3GeVBC0->SetYTitle("#Delta #phi (rad)");
1777     fhDeltaPhiChargedPtA3GeVBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1778     
1779     fhPtTrigChargedBC0  =
1780     new TH2F("hPtTrigChargedBC0","trigger and charged tracks pt distribution, track BC==0",
1781              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1782     fhPtTrigChargedBC0->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1783     fhPtTrigChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1784     
1785     fhXEChargedBC0  =
1786     new TH2F("hXEChargedBC0","#it{x}_{#it{E}} for charged tracks, track BC==0",
1787              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1788     fhXEChargedBC0->SetYTitle("#it{x}_{#it{E}}");
1789     fhXEChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1790     
1791     fhXEUeChargedBC0  =
1792     new TH2F("hXEUeChargedBC0","#it{x}_{#it{E}} for Underlying Event, track BC==0",
1793              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1794     fhXEUeChargedBC0->SetYTitle("#it{x}_{#it{E}}");
1795     fhXEUeChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1796     
1797     fhZTChargedBC0  =
1798     new TH2F("hZTChargedBC0","#it{z}_{T} for charged tracks, track BC==0",
1799              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1800     fhZTChargedBC0->SetYTitle("#it{z}_{T}");
1801     fhZTChargedBC0->SetXTitle("#it{p}_{T trigger}");
1802     
1803     fhZTUeChargedBC0  =
1804     new TH2F("hZTUeChargedBC0","#it{z}_{T} for Underlying Event, track BC==0",
1805              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1806     fhZTUeChargedBC0->SetYTitle("#it{z}_{T}");
1807     fhZTUeChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1808     
1809     outputContainer->Add(fhDeltaPhiChargedBC0) ;
1810     outputContainer->Add(fhDeltaPhiChargedPtA3GeVBC0) ;
1811     outputContainer->Add(fhXEChargedBC0) ;
1812     outputContainer->Add(fhXEUeChargedBC0) ;
1813     outputContainer->Add(fhZTChargedBC0) ;
1814     outputContainer->Add(fhZTUeChargedBC0) ;
1815     outputContainer->Add(fhPtTrigChargedBC0) ;
1816     
1817     fhPtTriggerVtxBC0  = new TH1F("hPtTriggerVtxBC0","#it{p}_{T} distribution of trigger particles", nptbins,ptmin,ptmax);
1818     fhPtTriggerVtxBC0->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1819     
1820     fhDeltaPhiChargedVtxBC0  = new TH2F
1821     ("hDeltaPhiChargedVtxBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, track BC==0",
1822      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1823     fhDeltaPhiChargedVtxBC0->SetYTitle("#Delta #phi (rad)");
1824     fhDeltaPhiChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1825     
1826     fhDeltaPhiChargedPtA3GeVVtxBC0  = new TH2F
1827     ("hDeltaPhiChargedPtA3GeVVtxBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC==0",
1828      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1829     fhDeltaPhiChargedPtA3GeVVtxBC0->SetYTitle("#Delta #phi (rad)");
1830     fhDeltaPhiChargedPtA3GeVVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1831     
1832     fhPtTrigChargedVtxBC0  =
1833     new TH2F("hPtTrigChargedVtxBC0","trigger and charged tracks pt distribution, track BC==0",
1834              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1835     fhPtTrigChargedVtxBC0->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1836     fhPtTrigChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1837     
1838     fhXEChargedVtxBC0  =
1839     new TH2F("hXEChargedVtxBC0","#it{x}_{#it{E}} for charged tracks, track BC==0",
1840              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1841     fhXEChargedVtxBC0->SetYTitle("#it{x}_{#it{E}}");
1842     fhXEChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1843     
1844     fhXEUeChargedVtxBC0  =
1845     new TH2F("hXEUeChargedVtxBC0","#it{x}_{#it{E}} for Underlying Event, track BC==0",
1846              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1847     fhXEUeChargedVtxBC0->SetYTitle("#it{x}_{#it{E}}");
1848     fhXEUeChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1849     
1850     fhZTChargedVtxBC0  =
1851     new TH2F("hZTChargedVtxBC0","#it{z}_{T} for charged tracks, track BC==0",
1852              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1853     fhZTChargedVtxBC0->SetYTitle("#it{z}_{T}");
1854     fhZTChargedVtxBC0->SetXTitle("#it{p}_{T trigger}");
1855     
1856     fhZTUeChargedVtxBC0  =
1857     new TH2F("hZTUeChargedVtxBC0","#it{z}_{T} for Underlying Event, track BC==0",
1858              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1859     fhZTUeChargedVtxBC0->SetYTitle("#it{z}_{T}");
1860     fhZTUeChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1861     
1862     outputContainer->Add(fhPtTriggerVtxBC0);
1863     outputContainer->Add(fhDeltaPhiChargedVtxBC0) ;
1864     outputContainer->Add(fhDeltaPhiChargedPtA3GeVVtxBC0) ;
1865     outputContainer->Add(fhXEChargedVtxBC0) ;
1866     outputContainer->Add(fhXEUeChargedVtxBC0) ;
1867     outputContainer->Add(fhZTChargedVtxBC0) ;
1868     outputContainer->Add(fhZTUeChargedVtxBC0) ;
1869     outputContainer->Add(fhPtTrigChargedVtxBC0) ;
1870     
1871     for(Int_t i = 0 ; i < 7 ; i++)
1872     {
1873       fhPtTriggerPileUp[i]  = new TH1F(Form("hPtTriggerPileUp%s",pileUpName[i].Data()),
1874                                        Form("#it{p}_{T} distribution of trigger particles, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
1875       fhPtTriggerPileUp[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1876       outputContainer->Add(fhPtTriggerPileUp[i]);
1877       
1878       fhDeltaPhiChargedPileUp[i]  = new TH2F(Form("hDeltaPhiChargedPileUp%s",pileUpName[i].Data()),
1879                                              Form("#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, %s Pile-Up event",pileUpName[i].Data()),
1880                                              nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1881       fhDeltaPhiChargedPileUp[i]->SetYTitle("#Delta #phi (rad)");
1882       fhDeltaPhiChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1883       outputContainer->Add(fhDeltaPhiChargedPileUp[i]) ;
1884       
1885       fhDeltaPhiChargedPtA3GeVPileUp[i]  = new TH2F(Form("hDeltaPhiChargedPtA3GeVPileUp%s",pileUpName[i].Data()),
1886                                                     Form("#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, %s Pile-Up event",pileUpName[i].Data()),
1887                                                     nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1888       fhDeltaPhiChargedPtA3GeVPileUp[i]->SetYTitle("#Delta #phi (rad)");
1889       fhDeltaPhiChargedPtA3GeVPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1890       outputContainer->Add(fhDeltaPhiChargedPtA3GeVPileUp[i]) ;
1891       
1892       fhDeltaEtaChargedPileUp[i]  = new TH2F(Form("hDeltaEtaChargedPileUp%s",pileUpName[i].Data()),
1893                                              Form("#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, %s Pile-Up event",pileUpName[i].Data()),
1894                                              nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1895       fhDeltaEtaChargedPileUp[i]->SetYTitle("#Delta #eta");
1896       fhDeltaEtaChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1897       outputContainer->Add(fhDeltaEtaChargedPileUp[i]) ;
1898       
1899       fhDeltaEtaChargedPtA3GeVPileUp[i]  = new TH2F(Form("hDeltaEtaChargedPtA3GeVPileUp%s",pileUpName[i].Data()),
1900                                                     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()),
1901                                                     nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1902       fhDeltaEtaChargedPtA3GeVPileUp[i]->SetYTitle("#Delta #eta");
1903       fhDeltaEtaChargedPtA3GeVPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1904       outputContainer->Add(fhDeltaEtaChargedPtA3GeVPileUp[i]) ;
1905       
1906       fhXEChargedPileUp[i]  = new TH2F(Form("hXEChargedPileUp%s",pileUpName[i].Data()),
1907                                        Form("#it{x}_{#it{E}} for charged tracks, %s Pile-Up event",pileUpName[i].Data()),
1908                                        nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1909       fhXEChargedPileUp[i]->SetYTitle("#it{x}_{#it{E}}");
1910       fhXEChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1911       outputContainer->Add(fhXEChargedPileUp[i]) ;
1912       
1913       fhXEUeChargedPileUp[i]  = new TH2F(Form("hXEUeChargedPileUp%s",pileUpName[i].Data()),
1914                                          Form("#it{x}_{#it{E}} for Underlying Event, %s Pile-Up event",pileUpName[i].Data()),
1915                                          nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1916       fhXEUeChargedPileUp[i]->SetYTitle("#it{x}_{#it{E}}");
1917       fhXEUeChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1918       outputContainer->Add(fhXEUeChargedPileUp[i]) ;
1919       
1920       fhZTChargedPileUp[i]  = new TH2F(Form("hZTChargedPileUp%s",pileUpName[i].Data()),
1921                                        Form("#it{z}_{T} for charged tracks, %s Pile-Up event",pileUpName[i].Data()),
1922                                        nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1923       fhZTChargedPileUp[i]->SetYTitle("#it{z}_{T}");
1924       fhZTChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1925       outputContainer->Add(fhZTChargedPileUp[i]) ;
1926       
1927       fhZTUeChargedPileUp[i]  = new TH2F(Form("hZTUeChargedPileUp%s",pileUpName[i].Data()),
1928                                          Form("#it{z}_{T} for Underlying Event, %s Pile-Up event",pileUpName[i].Data()),
1929                                          nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1930       fhZTUeChargedPileUp[i]->SetYTitle("#it{z}_{T}");
1931       fhZTUeChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1932       outputContainer->Add(fhZTUeChargedPileUp[i]) ;
1933       
1934       fhPtTrigChargedPileUp[i]  = new TH2F(Form("hPtTrigChargedPileUp%s",pileUpName[i].Data()),
1935                                            Form("trigger and charged tracks pt distribution, %s Pile-Up event",pileUpName[i].Data()),
1936                                            nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1937       fhPtTrigChargedPileUp[i]->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1938       fhPtTrigChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1939       outputContainer->Add(fhPtTrigChargedPileUp[i]) ;
1940       
1941     }
1942   }
1943
1944   if(IsHighMultiplicityAnalysisOn())
1945   {
1946     Int_t nMultiBins = GetNCentrBin();
1947     fhDeltaPhiChargedMult = new TH2F*[nMultiBins] ;
1948     fhDeltaEtaChargedMult = new TH2F*[nMultiBins] ;
1949     fhXEMult              = new TH2F*[nMultiBins] ;
1950     fhXEUeMult            = new TH2F*[nMultiBins] ;
1951     fhZTMult              = new TH2F*[nMultiBins] ;
1952     fhZTUeMult            = new TH2F*[nMultiBins] ;
1953     
1954     for(Int_t im=0; im<nMultiBins; im++)
1955     {
1956       fhDeltaPhiChargedMult[im]  = new TH2F
1957       (Form("hDeltaPhiCharged_Mult%d",im),Form("#Delta #phi charged Mult bin %d",im), nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1958       fhDeltaPhiChargedMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1959       fhDeltaPhiChargedMult[im]->SetYTitle("#Delta #phi (rad)");
1960       
1961       fhDeltaEtaChargedMult[im]  = new TH2F
1962       (Form("hDeltaEtaCharged_Mult%d",im),Form("#Delta #eta charged Mult bin %d",im), nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);
1963       fhDeltaEtaChargedMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1964       fhDeltaEtaChargedMult[im]->SetYTitle("#Delta #eta");
1965       
1966       fhXEMult[im]  = new TH2F
1967       (Form("hXECharged_Mult%d",im),Form("#it{x}_{E} charged Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1968       fhXEMult[im]->SetYTitle("#it{x}_{E}");
1969       fhXEMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1970       
1971       fhXEUeMult[im]  = new TH2F
1972       (Form("hXEUeCharged_Mult%d",im),Form("#it{x}_{E} UE charged Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1973       fhXEUeMult[im]->SetYTitle("#it{x}_{E}");
1974       fhXEUeMult[im]->SetXTitle("#it{p}_{T trigger}(GeV/#it{c})");
1975       
1976       fhZTMult[im]  = new TH2F
1977       (Form("hZTCharged_Mult%d",im),Form("#it{z}_{T} charged  Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1978       fhZTMult[im]->SetYTitle("#it{z}_{T}");
1979       fhZTMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1980       
1981       fhZTUeMult[im]  = new TH2F
1982       (Form("hZTUeCharged_Mult%d",im),Form("#it{z}_{T} UE charged  Mult bin %d",im), nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
1983       fhZTUeMult[im]->SetYTitle("#it{z}_{T}");
1984       fhZTUeMult[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1985       
1986       outputContainer->Add(fhDeltaPhiChargedMult[im]) ;
1987       outputContainer->Add(fhDeltaEtaChargedMult[im]) ;
1988       outputContainer->Add(fhXEMult  [im]);
1989       outputContainer->Add(fhXEUeMult[im]);
1990       outputContainer->Add(fhZTMult  [im]);
1991       outputContainer->Add(fhZTUeMult[im]);
1992     }
1993   }
1994   
1995   if(fFillBradHisto)
1996   {
1997     fhAssocPtBkg        = new TH2F("hAssocPtBkg", " Trigger #it{p}_{T} vs associated hadron #it{p}_{T} from background",
1998                                    nptbins, ptmin, ptmax,nptbins,ptmin,ptmax);
1999     fhAssocPtBkg->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2000     fhAssocPtBkg->SetYTitle("#it{p}_{T associated} (GeV/#it{c})");
2001     outputContainer->Add(fhAssocPtBkg) ;
2002     
2003     fhDeltaPhiBrad = new TH2F("hDeltaPhiBrad","atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi vs #it{p}_{T trigger} ",
2004                               nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
2005     fhDeltaPhiBrad->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2006     fhDeltaPhiBrad->SetYTitle("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi");
2007     outputContainer->Add(fhDeltaPhiBrad) ;
2008   }
2009   
2010   fhDeltaPhiDeltaEtaAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2011   fhDeltaPhiAssocPtBin       = new TH2F*[fNAssocPtBins*nz];
2012   fhDeltaPhiAssocPtBinDEta08 = new TH2F*[fNAssocPtBins*nz];
2013   fhDeltaPhiAssocPtBinDEta0  = new TH2F*[fNAssocPtBins*nz];
2014   if(fFillMomImbalancePtAssocBinsHisto)
2015   {
2016     fhXEAssocPtBin           = new TH2F*[fNAssocPtBins*nz];
2017     fhZTAssocPtBin           = new TH2F*[fNAssocPtBins*nz];
2018   }
2019
2020   if(fCorrelVzBin)
2021   {
2022     fhXEVZ = new TH2F*[nz];
2023     fhZTVZ = new TH2F*[nz];
2024   }
2025   
2026   if(fFillBradHisto)
2027     fhDeltaPhiBradAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2028   
2029
2030   fhDeltaPhiAssocPtBin       = new TH2F*[fNAssocPtBins*nz];
2031   if(fFillEtaGapsHisto)fhDeltaPhiAssocPtBinDEta08       = new TH2F*[fNAssocPtBins*nz];
2032   if(fDecayTrigger)    fhDeltaPhiDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2033   
2034   
2035   if(fHMPIDCorrelation)
2036   {
2037     fhDeltaPhiAssocPtBinHMPID   = new TH2F*[fNAssocPtBins*nz];
2038     fhDeltaPhiAssocPtBinHMPIDAcc= new TH2F*[fNAssocPtBins*nz];
2039   }
2040
2041   for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2042   {
2043     for(Int_t z = 0 ; z < nz ; z++)
2044     {
2045       Int_t bin = i*nz+z;
2046       
2047       if(fCorrelVzBin)
2048       {
2049         sz = Form("_vz%d",z);
2050         tz = Form(", #it{v}_{#it{z}} bin %d",z);
2051       }
2052       
2053       //printf("iAssoc %d, Vz %d, bin %d - sz %s, tz %s \n",i,z,bin,sz.Data(),tz.Data());
2054       
2055       fhDeltaPhiDeltaEtaAssocPtBin[bin]  = new TH2F(Form("hDeltaPhiDeltaEtaPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2056                                                     Form("#Delta #phi vs #Delta #eta for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2057                                                     ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
2058       fhDeltaPhiDeltaEtaAssocPtBin[bin]->SetXTitle("#Delta #phi (rad)");
2059       fhDeltaPhiDeltaEtaAssocPtBin[bin]->SetYTitle("#Delta #eta");
2060       
2061       fhDeltaPhiAssocPtBin[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2062                                            Form("#Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2063                                            nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2064       fhDeltaPhiAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2065       fhDeltaPhiAssocPtBin[bin]->SetYTitle("#Delta #phi (rad)");
2066       
2067       outputContainer->Add(fhDeltaPhiDeltaEtaAssocPtBin[bin]) ;
2068       outputContainer->Add(fhDeltaPhiAssocPtBin[bin]) ;
2069       
2070       if(fFillEtaGapsHisto)
2071       {
2072         fhDeltaPhiAssocPtBinDEta08[bin] = new TH2F(Form("hDeltaPhiDeltaEta0.8PtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2073                                                    Form("#Delta #phi 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()),
2074                                                    nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2075         fhDeltaPhiAssocPtBinDEta08[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2076         fhDeltaPhiAssocPtBinDEta08[bin]->SetYTitle("#Delta #phi (rad)");
2077         
2078         fhDeltaPhiAssocPtBinDEta0[bin] = new TH2F(Form("hDeltaPhiDeltaEta0PtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2079                                                   Form("#Delta #phi 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()),
2080                                                   nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2081         fhDeltaPhiAssocPtBinDEta0[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2082         fhDeltaPhiAssocPtBinDEta0[bin]->SetYTitle("#Delta #phi (rad)");
2083         
2084         outputContainer->Add(fhDeltaPhiAssocPtBinDEta08[bin]) ;
2085         outputContainer->Add(fhDeltaPhiAssocPtBinDEta0[bin]) ;
2086       }
2087       
2088       if(fDecayTrigger)
2089       {
2090         fhDeltaPhiDecayChargedAssocPtBin[bin] = new TH2F(Form("hDeltaPhiPtDecayChargedAssocPt%2.1f_%2.1f%s_bit%d", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data(),fDecayBits[0]),
2091                                                          Form("#Delta #phi 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]),
2092                                                          nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2093         fhDeltaPhiDecayChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2094         fhDeltaPhiDecayChargedAssocPtBin[bin]->SetYTitle("#Delta #phi (rad)");
2095         
2096         outputContainer->Add(fhDeltaPhiDecayChargedAssocPtBin[bin]) ;
2097       }
2098       
2099       if(fFillBradHisto)
2100       {
2101         fhDeltaPhiBradAssocPtBin[bin] = new TH2F(Form("hDeltaPhiBradPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2102                                                  Form("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2103                                                  nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
2104         fhDeltaPhiBradAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2105         fhDeltaPhiBradAssocPtBin[bin]->SetYTitle("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi");
2106         outputContainer->Add(fhDeltaPhiBradAssocPtBin[bin]) ;
2107       }
2108       
2109       if(fHMPIDCorrelation)
2110       {
2111         fhDeltaPhiAssocPtBinHMPID[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%sHMPID", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2112                                                   Form("#Delta #phi 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()),
2113                                                   nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2114         fhDeltaPhiAssocPtBinHMPID[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})" );
2115         fhDeltaPhiAssocPtBinHMPID[bin]->SetYTitle("#Delta #phi (rad)");
2116         
2117         fhDeltaPhiAssocPtBinHMPIDAcc[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%sHMPIDAcc", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2118                                                      Form("#Delta #phi 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()),
2119                                                      nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2120         fhDeltaPhiAssocPtBinHMPIDAcc[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2121         fhDeltaPhiAssocPtBinHMPIDAcc[bin]->SetYTitle("#Delta #phi (rad)");
2122         
2123         outputContainer->Add(fhDeltaPhiAssocPtBinHMPID   [bin]) ;
2124         outputContainer->Add(fhDeltaPhiAssocPtBinHMPIDAcc[bin]) ;
2125       }
2126     }
2127   }
2128   
2129   if(fFillMomImbalancePtAssocBinsHisto)
2130   {
2131     for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2132     {
2133       fhXEAssocPtBin[i]       = new TH2F(Form("hXEAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2134                                          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]),
2135                                          nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2136       fhXEAssocPtBin[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2137       fhXEAssocPtBin[i]->SetYTitle("#it{x}_{#it{E}}");
2138       
2139       fhZTAssocPtBin[i]       = new TH2F(Form("hZTAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2140                                          Form("#it{z}_{T} vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]),
2141                                          nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2142       fhZTAssocPtBin[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2143       fhZTAssocPtBin[i]->SetYTitle("#it{z}_{T}");
2144       
2145       
2146       outputContainer->Add(fhXEAssocPtBin[i]);
2147       outputContainer->Add(fhZTAssocPtBin[i]);
2148     }
2149   }
2150
2151   if(fCorrelVzBin)
2152   {
2153     for(Int_t z = 0 ; z < nz ; z++)
2154     {
2155       sz = Form("_vz%d",z);
2156       tz = Form(", #it{v}_{#it{z}} bin %d",z);
2157       
2158       fhXEVZ[z]       = new TH2F(Form("hXE%s", sz.Data()),
2159                                  Form("#it{x}_{#it{E}} vs #it{p}_{T trigger}%s", tz.Data()),
2160                                  nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2161       fhXEVZ[z]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2162       fhXEVZ[z]->SetYTitle("#it{x}_{#it{E}}");
2163       
2164       fhZTVZ[z]       = new TH2F(Form("hZT%s",sz.Data()),
2165                                  Form("#it{z}_{T} vs #it{p}_{T trigger}%s", tz.Data()),
2166                                  nptbins, ptmin, ptmax,nxeztbins,xeztmin,xeztmax);
2167       fhZTVZ[z]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2168       fhZTVZ[z]->SetYTitle("#it{z}_{T}");
2169       
2170       outputContainer->Add(fhXEVZ[z]);
2171       outputContainer->Add(fhZTVZ[z]);
2172     }
2173   }
2174
2175  
2176   if(fPi0Trigger)
2177   {
2178     fhPtPi0DecayRatio  = new TH2F
2179     ("hPtPi0DecayRatio","#it{p}_{T} of #pi^{0} and the ratio of pt for two decay",
2180      nptbins,ptmin,ptmax, 100,0.,2.);
2181     fhPtPi0DecayRatio->SetXTitle("#it{p}_{T}^{#pi^{0}} (GeV/#it{c})");
2182     fhPtPi0DecayRatio->SetYTitle("#it{p}_{T}^{Decay}/#it{p}_{T}^{#pi^{0}}");
2183     outputContainer->Add(fhPtPi0DecayRatio) ;
2184     
2185     fhDeltaPhiPi0DecayCharged  = new TH2F
2186     ("hDeltaPhiPi0DecayCharged","#phi_{Decay} - #phi_{h^{#pm}} vs #it{p}_{T Decay}",
2187      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2188     fhDeltaPhiPi0DecayCharged->SetYTitle("#Delta #phi (rad)");
2189     fhDeltaPhiPi0DecayCharged->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2190     
2191     fhXEPi0DecayCharged  =
2192     new TH2F("hXEPi0DecayCharged","#it{x}_{#it{E}}  Decay",
2193              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2194     fhXEPi0DecayCharged->SetYTitle("#it{x}_{#it{E}}");
2195     fhXEPi0DecayCharged->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2196     
2197     fhZTPi0DecayCharged  =
2198     new TH2F("hZTPi0DecayCharged","#it{z}_{trigger h^{#pm}} = #it{p}_{T h^{#pm}} / #it{p}_{T Decay}",
2199              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2200     fhZTPi0DecayCharged->SetYTitle("#it{z}_{decay h^{#pm}}");
2201     fhZTPi0DecayCharged->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2202     
2203     outputContainer->Add(fhDeltaPhiPi0DecayCharged) ;
2204     outputContainer->Add(fhXEPi0DecayCharged) ;
2205     outputContainer->Add(fhZTPi0DecayCharged) ;
2206   }
2207   
2208   if(fDecayTrigger)
2209   {
2210     for(Int_t ibit = 0; ibit< fNDecayBits; ibit++)
2211     {
2212       fhDeltaPhiDecayCharged[ibit]  = new TH2F
2213       (Form("hDeltaPhiDecayCharged_bit%d",fDecayBits[ibit]),
2214        Form("#phi_{Decay} - #phi_{h^{#pm}} vs #it{p}_{T Decay}, Bit %d",fDecayBits[ibit]),
2215        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2216       fhDeltaPhiDecayCharged[ibit]->SetYTitle("#Delta #phi (rad)");
2217       fhDeltaPhiDecayCharged[ibit]->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2218       
2219       fhXEDecayCharged[ibit]  =
2220       new TH2F(Form("hXEDecayCharged_bit%d",fDecayBits[ibit]),
2221                Form("#it{x}_{#it{E}}  Decay, Bit %d",fDecayBits[ibit]),
2222                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2223       fhXEDecayCharged[ibit]->SetYTitle("#it{x}_{#it{E}}");
2224       fhXEDecayCharged[ibit]->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2225       
2226       fhZTDecayCharged[ibit]  =
2227       new TH2F(Form("hZTDecayCharged_bit%d",fDecayBits[ibit]),
2228                Form("#it{z}_{trigger h^{#pm}} = #it{p}_{T h^{#pm}} / #it{p}_{T Decay}, Bit %d",fDecayBits[ibit]),
2229                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2230       fhZTDecayCharged[ibit]->SetYTitle("#it{z}_{decay h^{#pm}}");
2231       fhZTDecayCharged[ibit]->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
2232       
2233       outputContainer->Add(fhDeltaPhiDecayCharged[ibit]) ;
2234       outputContainer->Add(fhXEDecayCharged[ibit]) ;
2235       outputContainer->Add(fhZTDecayCharged[ibit]) ;
2236     }
2237   }
2238   
2239   //Correlation with neutral hadrons
2240   if(fNeutralCorr)
2241   {
2242     fhDeltaPhiDeltaEtaNeutral  = new TH2F
2243     ("hDeltaPhiDeltaEtaNeutral","#phi_{trigger} - #phi_{h^{0}} vs #eta_{trigger} - #eta_{h^{0}}",
2244      ndeltaphibins ,deltaphimin,deltaphimax, ndeltaetabins ,deltaetamin,deltaetamax);
2245     fhDeltaPhiDeltaEtaNeutral->SetXTitle("#Delta #phi (rad)");
2246     fhDeltaPhiDeltaEtaNeutral->SetYTitle("#Delta #eta");
2247           
2248     fhPhiNeutral  = new TH2F
2249     ("hPhiNeutral","#phi_{#pi^{0}}  vs #it{p}_{T #pi^{0}}",
2250      nptbins,ptmin,ptmax,180,0,TMath::TwoPi());
2251     fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
2252     fhPhiNeutral->SetXTitle("#it{p}_{T #pi^{0}} (GeV/#it{c})");
2253     
2254     fhEtaNeutral  = new TH2F
2255     ("hEtaNeutral","#eta_{#pi^{0}}  vs #it{p}_{T #pi^{0}}",
2256      nptbins,ptmin,ptmax,200,-1.,1.);
2257     fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
2258     fhEtaNeutral->SetXTitle("#it{p}_{T #pi^{0}} (GeV/#it{c})");
2259     
2260     fhDeltaPhiNeutral  = new TH2F
2261     ("hDeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs #it{p}_{T trigger}",
2262      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2263     fhDeltaPhiNeutral->SetYTitle("#Delta #phi (rad)");
2264     fhDeltaPhiNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2265     
2266     fhDeltaPhiNeutralPt  = new TH2F
2267     ("hDeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs #it{p}_{T #pi^{0}}}",
2268      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2269     fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi (rad)");
2270     fhDeltaPhiNeutralPt->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2271     
2272     fhDeltaEtaNeutral  = new TH2F
2273     ("hDeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs #it{p}_{T trigger}",
2274      nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);
2275     fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
2276     fhDeltaEtaNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2277     
2278     fhXENeutral  =
2279     new TH2F("hXENeutral","#it{x}_{#it{E}} for #pi^{0} associated",
2280              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2281     fhXENeutral->SetYTitle("#it{x}_{#it{E}}");
2282     fhXENeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2283     
2284     fhPtHbpXENeutral  =
2285     new TH2F("hHbpXENeutral","#xi = ln(1/#it{x}_{#it{E}})for #pi^{0} associated",
2286              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2287     fhPtHbpXENeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2288     fhPtHbpXENeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2289     
2290     fhZTNeutral  =
2291     new TH2F("hZTNeutral","#it{z}_{trigger #pi} = #it{p}_{T #pi^{0}} / #it{p}_{T trigger} for #pi^{0} associated",
2292              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2293     fhZTNeutral->SetYTitle("#it{z}_{trigger #pi^{0}}");
2294     fhZTNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2295     
2296     fhPtHbpZTNeutral  =
2297     new TH2F("hHbpZTNeutral","#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2298              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2299     fhPtHbpZTNeutral->SetYTitle("ln(1/#it{z}_{T})");
2300     fhPtHbpZTNeutral->SetXTitle("#it{p}_{T trigger}");
2301     
2302     outputContainer->Add(fhDeltaPhiDeltaEtaNeutral);
2303     outputContainer->Add(fhPhiNeutral) ;
2304     outputContainer->Add(fhEtaNeutral) ;
2305     outputContainer->Add(fhDeltaPhiNeutral) ;
2306     outputContainer->Add(fhDeltaPhiNeutralPt) ;
2307     outputContainer->Add(fhDeltaEtaNeutral) ;
2308     outputContainer->Add(fhXENeutral) ;
2309     outputContainer->Add(fhPtHbpXENeutral) ;
2310     outputContainer->Add(fhZTNeutral) ;
2311     outputContainer->Add(fhPtHbpZTNeutral) ;
2312     
2313     fhDeltaPhiUeNeutralPt  = new TH2F
2314     (Form("hDeltaPhiUe%sNeutralPt",right.Data()),"#phi_{trigger} - #phi_{#pi^{0}} vs #it{p}_{T #pi^{0}}}",
2315      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2316     fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi (rad)");
2317     fhDeltaPhiUeNeutralPt->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2318     
2319     fhXEUeNeutral  =
2320     new TH2F(Form("hXEUeNeutral%s",right.Data()),"#it{x}_{#it{E}} for #pi^{0} associated",
2321              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2322     fhXEUeNeutral->SetYTitle("#it{x}_{#it{E}}");
2323     fhXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2324     
2325     fhPtHbpXEUeNeutral  =
2326     new TH2F(Form("hHbpXEUeNeutral%s",right.Data()),"#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2327              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2328     fhPtHbpXEUeNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2329     fhPtHbpXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2330     
2331     fhZTUeNeutral  =
2332     new TH2F(Form("hZTUeNeutral%s",right.Data()),"#it{z}_{trigger #pi} = #it{p}_{T #pi^{0}} / #it{p}_{T trigger} for #pi^{0} associated",
2333              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2334     fhZTUeNeutral->SetYTitle("#it{z}_{trigger #pi^{0}}");
2335     fhZTUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2336     
2337     fhPtHbpZTUeNeutral  =
2338     new TH2F(Form("hHbpZTUeNeutral%s",right.Data()),"#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2339              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2340     fhPtHbpXEUeNeutral->SetYTitle("ln(1/#it{z}_{T})");
2341     fhPtHbpXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2342
2343     outputContainer->Add(fhDeltaPhiUeNeutralPt) ;
2344     outputContainer->Add(fhXEUeNeutral) ;
2345     outputContainer->Add(fhPtHbpXEUeNeutral) ;
2346     outputContainer->Add(fhZTUeNeutral) ;
2347     outputContainer->Add(fhPtHbpZTUeNeutral) ;
2348
2349     if(fMakeSeveralUE)
2350     {
2351       fhDeltaPhiUeLeftNeutral  = new TH2F
2352       ("hDeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs #it{p}_{T h^{0}} with neutral UE left side range of trigger particles",
2353        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2354       fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #phi (rad)");
2355       fhDeltaPhiUeLeftNeutral->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2356       outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
2357       
2358       fhXEUeLeftNeutral  =
2359       new TH2F("hXEUeNeutralLeft","#it{x}_{#it{E}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE left side of trigger",
2360                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2361       fhXEUeLeftNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2362       fhXEUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2363       outputContainer->Add(fhXEUeLeftNeutral) ;
2364       
2365       fhPtHbpXEUeLeftNeutral  =
2366       new TH2F("hHbpXEUeNeutralLeft","#xi = ln(1/#it{x}_{#it{E}}) with neutral UE left side of trigger",
2367                nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2368       fhPtHbpXEUeLeftNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2369       fhPtHbpXEUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2370       outputContainer->Add(fhPtHbpXEUeLeftNeutral) ;
2371       
2372       fhZTUeLeftNeutral  =
2373       new TH2F("hZTUeNeutralLeft","#it{z}_{trigger h^{0}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE left side of trigger",
2374                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2375       fhZTUeLeftNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2376       fhZTUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2377       outputContainer->Add(fhZTUeLeftNeutral) ;
2378       
2379       fhPtHbpZTUeLeftNeutral  =
2380       new TH2F("hHbpZTUeNeutralLeft","#xi = ln(1/#it{z}_{T}) with neutral UE left side of trigger",
2381                nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2382       fhPtHbpZTUeLeftNeutral->SetYTitle("ln(1/#it{z}_{T})");
2383       fhPtHbpZTUeLeftNeutral->SetXTitle("#it{p}_{T trigger}");
2384       outputContainer->Add(fhPtHbpZTUeLeftNeutral) ;
2385     }
2386     
2387     if(fPi0Trigger)
2388     {
2389       fhDeltaPhiPi0DecayNeutral  = new TH2F
2390       ("hDeltaPhiPi0DecayNeutral","#phi_{Decay} - #phi_{h^{0}} vs #it{p}_{T Decay}",
2391        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2392       fhDeltaPhiPi0DecayNeutral->SetYTitle("#Delta #phi (rad)");
2393       fhDeltaPhiPi0DecayNeutral->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2394       
2395       fhXEPi0DecayNeutral  =
2396       new TH2F("hXEPi0DecayNeutral","#it{x}_{#it{E}} for decay trigger",
2397                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2398       fhXEPi0DecayNeutral->SetYTitle("#it{x}_{#it{E}}");
2399       fhXEPi0DecayNeutral->SetXTitle("#it{p}_{T decay}");
2400       
2401       fhZTPi0DecayNeutral  =
2402       new TH2F("hZTPi0DecayNeutral","#it{z}_{trigger h^{0}} = #it{p}_{T h^{0}} / #it{p}_{T Decay}",
2403                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2404       fhZTPi0DecayNeutral->SetYTitle("#it{z}_{h^{0}}");
2405       fhZTPi0DecayNeutral->SetXTitle("#it{p}_{T decay}");
2406       
2407       outputContainer->Add(fhDeltaPhiPi0DecayNeutral) ;
2408       outputContainer->Add(fhXEPi0DecayNeutral) ;
2409       outputContainer->Add(fhZTPi0DecayNeutral) ;
2410     }
2411   }//Correlation with neutral hadrons
2412   
2413   // If data is MC, fill more histograms, depending on origin
2414   if(IsDataMC())
2415   {
2416     for(Int_t i= fMCGenTypeMin; i <= fMCGenTypeMax; i++)
2417     {
2418       fhMCPtTrigger[i]  = new TH1F (Form("hMCPtTrigger_%s",nameMC[i].Data()),
2419                                  Form("MC %s: generated trigger #it{p}_{T}",nameMC[i].Data()),
2420                                  nptbins,ptmin,ptmax);
2421       fhMCPtTrigger[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2422       
2423       fhMCPhiTrigger[i]  = new TH2F (Form("hMCPhiTrigger_%s",nameMC[i].Data()),
2424                                      Form("MC %s: generated trigger #phi",nameMC[i].Data()),
2425                                      nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2426       fhMCPhiTrigger[i]->SetYTitle("#phi (rad)");
2427       fhMCPhiTrigger[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2428       
2429       fhMCEtaTrigger[i]  = new TH2F (Form("hMCEtaTrigger_%s",nameMC[i].Data()),
2430                                      Form("MC %s: generated trigger #eta",nameMC[i].Data()),
2431                                      nptbins,ptmin,ptmax, netabins,etamin,etamax);
2432       fhMCEtaTrigger[i]->SetYTitle("#eta");
2433       fhMCEtaTrigger[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2434       
2435       if(fMakeAbsoluteLeading || fMakeNearSideLeading)
2436       {
2437         fhMCPtTriggerNotLeading[i]  = new TH1F (Form("hMCPtTriggerNotLeading_%s",nameMC[i].Data()),
2438                                                 Form("MC %s: generated trigger #it{p}_{T}, when not leading of primaries",nameMC[i].Data()),
2439                                                 nptbins,ptmin,ptmax);
2440         fhMCPtTriggerNotLeading[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2441         
2442         fhMCPhiTriggerNotLeading[i]  = new TH2F (Form("hMCPhiTriggerNotLeading_%s",nameMC[i].Data()),
2443                                                  Form("MC %s: generated trigger #phi, when not leading of primaries",nameMC[i].Data()),
2444                                                  nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2445         fhMCPhiTriggerNotLeading[i]->SetYTitle("#phi (rad)");
2446         fhMCPhiTriggerNotLeading[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2447         
2448         
2449         fhMCEtaTriggerNotLeading[i]  = new TH2F (Form("hMCEtaTriggerNotLeading_%s",nameMC[i].Data()),
2450                                                  Form("MC %s: generated triogger #eta, when not leading of primaries",nameMC[i].Data()),
2451                                                  nptbins,ptmin,ptmax, netabins,etamin,etamax);
2452         fhMCEtaTriggerNotLeading[i]->SetYTitle("#eta ");
2453         fhMCEtaTriggerNotLeading[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2454       }
2455       
2456       fhMCEtaCharged[i]  = new TH2F (Form("hMCEtaCharged_%s",nameMC[i].Data()),
2457                                      Form("MC %s: #eta_{h^{#pm}}  vs #it{p}_{T #pm}",nameMC[i].Data()),
2458                                      nptbins,ptmin,ptmax,100,-1.,1.);
2459       fhMCEtaCharged[i]->SetYTitle("#eta_{h^{#pm}} (rad)");
2460       fhMCEtaCharged[i]->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
2461       
2462       fhMCPhiCharged[i]  = new TH2F(Form("hMCPhiCharged_%s",nameMC[i].Data()),
2463                                     Form("MC %s: phi_{h^{#pm}}  vs #it{p}_{T #pm}",nameMC[i].Data()),
2464                                     nptbins,ptmin,ptmax,180,0,TMath::TwoPi());
2465       fhMCPhiCharged[i]->SetYTitle("MC #phi_{h^{#pm}} (rad)");
2466       fhMCPhiCharged[i]->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
2467       
2468       fhMCDeltaPhiDeltaEtaCharged[i]  = new TH2F (Form("hMCDeltaPhiDeltaEtaCharged_%s",nameMC[i].Data()),
2469                                                   Form("MC %s: phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",nameMC[i].Data()),
2470                                                   ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
2471       fhMCDeltaPhiDeltaEtaCharged[i]->SetXTitle("#Delta #phi (rad)");
2472       fhMCDeltaPhiDeltaEtaCharged[i]->SetYTitle("#Delta #eta");
2473       
2474       fhMCDeltaEtaCharged[i]  = new TH2F (Form("hMCDeltaEtaCharged_%s",nameMC[i].Data()),
2475                                           Form("MC %s: #eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger} and #it{p}_{T assoc}",nameMC[i].Data()),
2476                                           nptbins,ptmin,ptmax,ndeltaetabins ,deltaetamin,deltaetamax);
2477       fhMCDeltaEtaCharged[i]->SetYTitle("#Delta #eta");
2478       fhMCDeltaEtaCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2479       
2480       fhMCDeltaPhiCharged[i]  = new TH2F  (Form("hMCDeltaPhiCharged_%s",nameMC[i].Data()),
2481                                            Form("MC %s: #phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}",nameMC[i].Data()),
2482                                            nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2483       fhMCDeltaPhiCharged[i]->SetYTitle("#Delta #phi (rad)");
2484       fhMCDeltaPhiCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2485
2486       fhMCDeltaPhiChargedPt[i]  = new TH2F (Form("hMCDeltaPhiChargedPt_%s",nameMC[i].Data()),
2487                                             Form("MC %s: #phi_{trigger} - #phi_{#h^{#pm}} vs #it{p}_{T h^{#pm}}",nameMC[i].Data()),
2488                                             nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2489       fhMCDeltaPhiChargedPt[i]->SetYTitle("#Delta #phi (rad)");
2490       fhMCDeltaPhiChargedPt[i]->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2491       
2492       fhMCPtXECharged[i]  = new TH2F (Form("hMCPtXECharged_%s",nameMC[i].Data()),
2493                                       Form("MC %s: #it{x}_{#it{E}} with charged hadrons",nameMC[i].Data()),
2494                                       nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2495       fhMCPtXECharged[i]->SetYTitle("#it{x}_{#it{E}}");
2496       fhMCPtXECharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2497       
2498       fhMCPtHbpXECharged[i]  = new TH2F(Form("hMCHbpXECharged_%s",nameMC[i].Data()),
2499                                         Form("MC %s: #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",nameMC[i].Data()),
2500                                         nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2501       fhMCPtHbpXECharged[i]->SetYTitle("ln(1/#it{x}_{#it{E}})");
2502       fhMCPtHbpXECharged[i]->SetXTitle("#it{p}_{T trigger}");
2503       
2504       fhMCPtZTCharged[i]  = new TH2F(Form("hMCPtZTCharged_%s",nameMC[i].Data()),
2505                                      Form("MC %s: #it{z}_{T} with charged hadrons",nameMC[i].Data()),
2506                                      nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2507       fhMCPtZTCharged[i]->SetYTitle("#it{z}_{T}");
2508       fhMCPtZTCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2509       
2510       fhMCPtHbpZTCharged[i]  = new TH2F(Form("hMCHbpZTCharged_%s",nameMC[i].Data()),
2511                                         Form("MC %s: #xi = ln(1/#it{z}_{T}) with charged hadrons",nameMC[i].Data()),
2512                                              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2513       fhMCPtHbpZTCharged[i]->SetYTitle("ln(1/#it{z}_{T})");
2514       fhMCPtHbpZTCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2515       
2516       fhMCPtTrigPout[i]  = new TH2F(Form("hMCPtTrigPout_%s",nameMC[i].Data()),
2517                                     Form("MC %s: #it{p}_{out} with triggers",nameMC[i].Data()),
2518                                     nptbins,ptmin,ptmax,nptbins,-1.*ptmax/2.,ptmax/2.);
2519       fhMCPtTrigPout[i]->SetYTitle("#it{p}_{out} (GeV/#it{c})");
2520       fhMCPtTrigPout[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2521       
2522       fhMCPtAssocDeltaPhi[i]  = new TH2F(Form("hMCPtAssocDeltaPhi_%s",nameMC[i].Data()),
2523                                          Form("MC %s: #Delta #phi with associated charged hadrons",nameMC[i].Data()),
2524                                          nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2525       fhMCPtAssocDeltaPhi[i]->SetYTitle("#Delta #phi (rad)");
2526       fhMCPtAssocDeltaPhi[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2527       
2528       outputContainer->Add(fhMCPtTrigger[i]);
2529       outputContainer->Add(fhMCPhiTrigger[i]);
2530       outputContainer->Add(fhMCEtaTrigger[i]);
2531       
2532       if(fMakeAbsoluteLeading || fMakeNearSideLeading)
2533       {
2534         outputContainer->Add(fhMCPtTriggerNotLeading[i]);
2535         outputContainer->Add(fhMCPhiTriggerNotLeading[i]);
2536         outputContainer->Add(fhMCEtaTriggerNotLeading[i]);
2537       }
2538       
2539       outputContainer->Add(fhMCDeltaPhiDeltaEtaCharged[i]);
2540       outputContainer->Add(fhMCPhiCharged[i]) ;
2541       outputContainer->Add(fhMCEtaCharged[i]) ;
2542       outputContainer->Add(fhMCDeltaEtaCharged[i]) ;
2543       outputContainer->Add(fhMCDeltaPhiCharged[i]) ;
2544       
2545       outputContainer->Add(fhMCDeltaPhiChargedPt[i]) ;
2546       outputContainer->Add(fhMCPtXECharged[i]) ;
2547       outputContainer->Add(fhMCPtZTCharged[i]) ;
2548       outputContainer->Add(fhMCPtHbpXECharged[i]) ;
2549       outputContainer->Add(fhMCPtHbpZTCharged[i]) ;
2550       outputContainer->Add(fhMCPtTrigPout[i]) ;
2551       outputContainer->Add(fhMCPtAssocDeltaPhi[i]) ;
2552
2553       // Underlying event
2554       
2555       fhMCUePart[i]  =
2556       new TH1F(Form("hMCUePart_%s",nameMC[i].Data()),
2557                Form("MC %s: UE particles distribution vs #it{p}_{T trigger}",nameMC[i].Data()),
2558                nptbins,ptmin,ptmax);
2559       fhMCUePart[i]->SetYTitle("#it{dN}^{ch}");
2560       fhMCUePart[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2561       
2562       fhMCPtXEUeCharged[i]  =
2563       new TH2F(Form("hMCPtXEUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2564                Form("MC %s: #it{x}_{#it{E}} with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2565                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2566       fhMCPtXEUeCharged[i]->SetYTitle("#it{x}_{#it{E}}");
2567       fhMCPtXEUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2568       
2569       fhMCPtHbpXEUeCharged[i] =
2570       new TH2F(Form("hMCPtHbpXEUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2571                Form("MC %s: #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2572                nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2573       fhMCPtHbpXEUeCharged[i]->SetYTitle("ln(1/#it{x}_{#it{E}})");
2574       fhMCPtHbpXEUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2575       
2576       fhMCPtZTUeCharged[i] =
2577       new TH2F(Form("hMCPtZTUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2578                Form("MC %s: #it{z}_{T} with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2579                nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2580       fhMCPtZTUeCharged[i]->SetYTitle("#it{z}_{T}");
2581       fhMCPtZTUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2582       
2583       fhMCPtHbpZTUeCharged[i] =
2584       new TH2F(Form("hMCPtHbpZTUeCharged%s_%s",right.Data(),nameMC[i].Data()),
2585                Form("MC %s: #xi = ln(1/#it{z}_{T}) with charged hadrons, Underlying Event %s",nameMC[i].Data(),right.Data()),
2586                nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2587       fhMCPtHbpZTUeCharged[i]->SetYTitle("ln(1/#it{z}_{T})");
2588       fhMCPtHbpZTUeCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2589       
2590       outputContainer->Add(fhMCUePart[i]);
2591       outputContainer->Add(fhMCPtXEUeCharged[i]) ;
2592       outputContainer->Add(fhMCPtZTUeCharged[i]) ;
2593       outputContainer->Add(fhMCPtHbpZTUeCharged[i]);
2594       outputContainer->Add(fhMCPtHbpXEUeCharged[i]);
2595
2596       if(fMakeSeveralUE)
2597       {
2598         fhMCPtXEUeLeftCharged[i]  = new TH2F(Form("hMCPtXEUeChargedLeft_%s",nameMC[i].Data()),
2599                                              Form("MC %s: #it{x}_{#it{E}} with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
2600                                              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2601         fhMCPtXEUeLeftCharged[i]->SetYTitle("#it{x}_{#it{E}}");
2602         fhMCPtXEUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2603         
2604         fhMCPtHbpXEUeLeftCharged[i] = new TH2F(Form("hMCPtHbpXEUeChargedLeft_%s",nameMC[i].Data()),
2605                                                Form("MC %s: #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
2606                                                nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2607         fhMCPtHbpXEUeLeftCharged[i]->SetYTitle("ln(1/#it{x}_{#it{E}})");
2608         fhMCPtHbpXEUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2609         
2610         fhMCPtZTUeLeftCharged[i]  = new TH2F(Form("hMCPtZTUeChargedLeft_%s",nameMC[i].Data()),
2611                                              Form("MC %s: #it{z}_{T} with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
2612                                              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2613         fhMCPtZTUeLeftCharged[i]->SetYTitle("#it{z}_{T}");
2614         fhMCPtZTUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2615         
2616         fhMCPtHbpZTUeLeftCharged[i] = new TH2F(Form("hMCPtHbpZTUeChargedLeft_%s",nameMC[i].Data()),
2617                                                Form("MC %s: #xi = ln(1/#it{z}_{T}) with charged hadrons, with UE left side range of trigger particles",nameMC[i].Data()),
2618                                                nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2619         fhMCPtHbpZTUeLeftCharged[i]->SetYTitle("ln(1/#it{z}_{T})");
2620         fhMCPtHbpZTUeLeftCharged[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2621         
2622         outputContainer->Add(fhMCPtXEUeLeftCharged[i]) ;
2623         outputContainer->Add(fhMCPtZTUeLeftCharged[i]) ;
2624         outputContainer->Add(fhMCPtHbpXEUeLeftCharged[i]);
2625         outputContainer->Add(fhMCPtHbpZTUeLeftCharged[i]) ;
2626         
2627       }
2628     }
2629   } //for MC histogram
2630   
2631   if(DoOwnMix())
2632   {
2633     //create event containers
2634     
2635     if(!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForTracksExists()))
2636     {
2637       Int_t nvz = GetNZvertBin();
2638       Int_t nrp = GetNRPBin();
2639       Int_t nce = GetNCentrBin();
2640       
2641       fListMixTrackEvents= new TList*[nvz*nrp*nce] ;
2642       
2643       for( Int_t ice = 0 ; ice < nce ; ice++ )
2644       {
2645         for( Int_t ivz = 0 ; ivz < nvz ; ivz++ )
2646         {
2647           for( Int_t irp = 0 ; irp < nrp ; irp++ )
2648           {
2649             Int_t bin = GetEventMixBin(ice,ivz,irp); //ic*nvz*nrp+iz*nrp+irp;
2650             
2651             //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
2652             //       ic,iz, irp, bin);
2653             
2654             fListMixTrackEvents[bin] = new TList() ;
2655             fListMixTrackEvents[bin]->SetOwner(kFALSE);
2656           }
2657         }
2658       }
2659     }
2660     
2661     fhPtTriggerMixed  = new TH1F ("hPtTriggerMixed","#it{p}_{T} distribution of trigger particles, used for mixing", nptbins,ptmin,ptmax);
2662     fhPtTriggerMixed->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2663     
2664     if(fCorrelVzBin)
2665     {
2666       fhPtTriggerMixedVzBin  = new TH2F ("hPtTriggerMixedVzBin","#it{p}_{T} distribution of trigger particles, used for mixing", nptbins,ptmin,ptmax,GetNZvertBin(),0,GetNZvertBin());
2667       fhPtTriggerMixedVzBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2668       fhPtTriggerMixedVzBin->SetYTitle("#it{v}_{#it{z}} bin");
2669       outputContainer->Add(fhPtTriggerMixedVzBin);
2670     }
2671     
2672     fhPtTriggerMixedBin  = new TH2F ("hPtTriggerMixedBin","#it{p}_{T} distribution of trigger particles vs mixing bin", nptbins,ptmin,ptmax,nMixBins,0,nMixBins);
2673     fhPtTriggerMixedBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2674     fhPtTriggerMixedBin->SetYTitle("Bin");
2675     
2676     fhPhiTriggerMixed  = new TH2F ("hPhiTriggerMixed","#phi distribution of trigger Particles, used for mixing",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2677     fhPhiTriggerMixed->SetYTitle("#phi (rad)");
2678     
2679     fhEtaTriggerMixed  = new TH2F ("hEtaTriggerMixed","#eta distribution of trigger, used for mixing",nptbins,ptmin,ptmax, netabins,etamin,etamax);
2680     fhEtaTriggerMixed->SetYTitle("#eta ");
2681     
2682     outputContainer->Add(fhPtTriggerMixed);
2683     outputContainer->Add(fhPtTriggerMixedBin);
2684     outputContainer->Add(fhPhiTriggerMixed);
2685     outputContainer->Add(fhEtaTriggerMixed);
2686     
2687     // Fill the cluster pool only in isolation analysis or if requested
2688     if( neutralMix && (!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForCaloExists())))
2689     {
2690       Int_t nvz = GetNZvertBin();
2691       Int_t nrp = GetNRPBin();
2692       Int_t nce = GetNCentrBin();
2693       
2694       fListMixCaloEvents= new TList*[nvz*nrp*nce] ;
2695       
2696       for( Int_t ice = 0 ; ice < nce ; ice++ )
2697       {
2698         for( Int_t ivz = 0 ; ivz < nvz ; ivz++ )
2699         {
2700           for( Int_t irp = 0 ; irp < nrp ; irp++ )
2701           {
2702             Int_t bin = GetEventMixBin(ice,ivz,irp); //ic*nvz*nrp+iz*nrp+irp;
2703             
2704             //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
2705             //       ic,iz, irp, bin);
2706             
2707             fListMixCaloEvents[bin] = new TList() ;
2708             fListMixCaloEvents[bin]->SetOwner(kFALSE);
2709           }
2710         }
2711       }
2712     }
2713     
2714     //Init the list in the reader if not done previously
2715     if(fUseMixStoredInReader)
2716     {
2717       if( !GetReader()->ListWithMixedEventsForTracksExists() )
2718         GetReader()->SetListWithMixedEventsForTracks(fListMixTrackEvents);
2719       
2720       if( !GetReader()->ListWithMixedEventsForCaloExists()   )
2721         GetReader()->SetListWithMixedEventsForCalo  (fListMixCaloEvents );
2722     }
2723     
2724     fhEventBin=new TH1I("hEventBin","Number of triggers per bin(cen,vz,rp)",
2725                         GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2726                         GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2727     fhEventBin->SetXTitle("event bin");
2728     outputContainer->Add(fhEventBin) ;
2729     
2730     fhEventMixBin=new TH1I("hEventMixBin","Number of triggers mixed per event bin(cen,vz,rp)",
2731                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2732                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2733     fhEventMixBin->SetXTitle("event bin");
2734     outputContainer->Add(fhEventMixBin) ;
2735
2736     fhEventMBBin=new TH1I("hEventMBBin","Number of min bias events per bin(cen,vz,rp)",
2737                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2738                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2739     fhEventMBBin->SetXTitle("event bin");
2740     outputContainer->Add(fhEventMBBin) ;
2741
2742     fhNtracksMB=new TH2F("hNtracksMBEvent","Number of filtered tracks in MB event per event bin",ntrbins,trmin,trmax,
2743                          GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2744                          GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2745     fhNtracksMB->SetYTitle("event bin");
2746     fhNtracksMB->SetXTitle("#it{N}_{track}");
2747     outputContainer->Add(fhNtracksMB);
2748
2749     if( neutralMix )
2750     {
2751       fhNclustersMB=new TH2F("hNclustersMBEvent","Number of filtered clusters in MB events per event bin",nclbins,clmin,clmax,
2752                              GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2753                              GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2754       fhNclustersMB->SetYTitle("event bin");
2755       fhNclustersMB->SetXTitle("#it{N}_{cluster}");
2756       outputContainer->Add(fhNclustersMB);
2757     }
2758     
2759     fhMixDeltaPhiCharged  = new TH2F
2760     ("hMixDeltaPhiCharged","Mixed event : #phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}",
2761      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2762     fhMixDeltaPhiCharged->SetYTitle("#Delta #phi (rad)");
2763     fhMixDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2764     outputContainer->Add(fhMixDeltaPhiCharged);
2765     
2766     fhMixDeltaPhiDeltaEtaCharged  = new TH2F
2767     ("hMixDeltaPhiDeltaEtaCharged","Mixed event : #phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
2768      ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
2769     fhMixDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi (rad)");
2770     fhMixDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
2771     outputContainer->Add(fhMixDeltaPhiDeltaEtaCharged);
2772     
2773     fhMixXECharged  =
2774     new TH2F("hMixXECharged","Mixed event : #it{x}_{#it{E}} for charged tracks",
2775              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2776     fhMixXECharged->SetYTitle("#it{x}_{#it{E}}");
2777     fhMixXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2778     outputContainer->Add(fhMixXECharged);
2779     
2780     fhMixXEUeCharged  =
2781     new TH2F("hMixXEUeCharged","Mixed event : #it{x}_{#it{E}} for charged tracks in Ue region",
2782              nptbins,ptmin,ptmax,nxeztbins,xeztmin,xeztmax);
2783     fhMixXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
2784     fhMixXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2785     outputContainer->Add(fhMixXEUeCharged);
2786     
2787     fhMixHbpXECharged  =
2788     new TH2F("hMixHbpXECharged","mixed event : #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
2789              nptbins,ptmin,ptmax,nhbpbins,hbpmin,hbpmax);
2790     fhMixHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2791     fhMixHbpXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2792     outputContainer->Add(fhMixHbpXECharged);
2793     
2794     fhMixDeltaPhiChargedAssocPtBin         = new TH2F*[fNAssocPtBins*nz];
2795     fhMixDeltaPhiChargedAssocPtBinDEta08   = new TH2F*[fNAssocPtBins*nz];
2796     fhMixDeltaPhiChargedAssocPtBinDEta0    = new TH2F*[fNAssocPtBins*nz];
2797     fhMixDeltaPhiDeltaEtaChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2798     
2799     for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2800     {
2801       for(Int_t z = 0 ; z < nz ; z++)
2802       {
2803         Int_t bin = i*nz+z;
2804         
2805         if(fCorrelVzBin)
2806         {
2807           sz = Form("_vz%d",z);
2808           tz = Form(", #it{v}_{#it{z}} bin %d",z);
2809         }
2810         
2811         //printf("MIX : iAssoc %d, Vz %d, bin %d - sz %s, tz %s \n",i,z,bin,sz.Data(),tz.Data());
2812         
2813         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin] = new TH2F(Form("hMixDeltaPhiDeltaEtaChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2814                                                                Form("Mixed event #Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2815                                                                ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
2816         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->SetXTitle("#Delta #phi (rad)");
2817         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->SetYTitle("#Delta #eta");
2818         
2819         outputContainer->Add(fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]);
2820         
2821         fhMixDeltaPhiChargedAssocPtBin[bin] = new TH2F(Form("hMixDeltaPhiChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2822                                                        Form("Mixed event #Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2823                                                        nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
2824         fhMixDeltaPhiChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2825         fhMixDeltaPhiChargedAssocPtBin[bin]->SetYTitle("#Delta #phi (rad)");
2826         
2827         outputContainer->Add(fhMixDeltaPhiChargedAssocPtBin[bin]);
2828         
2829         if(fFillEtaGapsHisto)
2830         {
2831           fhMixDeltaPhiChargedAssocPtBinDEta08[bin] = new TH2F(Form("hMixDeltaPhiDeltaEta0.8ChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2832                                                                Form("Mixed event #Delta #phi 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()),
2833                                                                nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
2834           fhMixDeltaPhiChargedAssocPtBinDEta08[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2835           fhMixDeltaPhiChargedAssocPtBinDEta08[bin]->SetYTitle("#Delta #phi (rad)");
2836           
2837           fhMixDeltaPhiChargedAssocPtBinDEta0[bin] = new TH2F(Form("hMixDeltaPhiDeltaEta0ChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2838                                                               Form("Mixed event #Delta #phi 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()),
2839                                                               nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
2840           fhMixDeltaPhiChargedAssocPtBinDEta0[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2841           fhMixDeltaPhiChargedAssocPtBinDEta0[bin]->SetYTitle("#Delta #phi (rad)");
2842           
2843           outputContainer->Add(fhMixDeltaPhiChargedAssocPtBinDEta08[bin]);
2844           outputContainer->Add(fhMixDeltaPhiChargedAssocPtBinDEta0[bin]);
2845         }
2846       }
2847     }
2848   }
2849   
2850   return outputContainer;
2851   
2852 }
2853
2854 //_____________________________________________________________________________________________________________________
2855 Bool_t AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum(Int_t indexPhoton1, Int_t indexPhoton2, Int_t idetector)
2856 {
2857   // Get the momentum of the pi0/eta assigned decay photons
2858   // In case of pi0/eta trigger, we may want to check their decay correlation,
2859   // get their decay children
2860   
2861   if(indexPhoton1!=-1 || indexPhoton2!=-1) return kFALSE;
2862   
2863   AliDebug(1,Form("indexPhoton1 = %d, indexPhoton2 = %d", indexPhoton1, indexPhoton2));
2864   
2865   TObjArray * clusters  = 0x0 ;
2866   if(idetector==kEMCAL) clusters = GetEMCALClusters() ;
2867   else                  clusters = GetPHOSClusters()  ;
2868