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