7d7446af1e1e6e6568599f52d15f33d34faf66a2
[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     fListMixEvents(),               fUseMixStoredInReader(0),
74     //Histograms
75     fhPtLeading(0),                 fhPhiLeading(0),       
76     fhEtaLeading(0),                
77     fhPtLeadingCentrality(0),       fhPtLeadingEventPlane(0), 
78     fhLeadingEventPlaneCentrality(0),fhDeltaPhiDeltaEtaCharged(0),
79     fhPhiCharged(0),                fhEtaCharged(0), 
80     fhDeltaPhiCharged(0),           fhDeltaEtaCharged(0), 
81     fhDeltaPhiChargedPt(0),         fhDeltaPhiUeChargedPt(0), 
82     fhUePart(0),
83     fhXECharged(0),                 fhXEUeCharged(0),
84     fhXEPosCharged(0),              fhXENegCharged(0),
85     fhPtHbpXECharged(0),            fhPtHbpXEUeCharged(0),
86     fhZTCharged(0),                 fhZTUeCharged(0),
87     fhZTPosCharged(0),              fhZTNegCharged(0),
88     fhPtHbpZTCharged(0),            fhPtHbpZTUeCharged(0),
89     fhDeltaPhiUeLeftCharged(0),     fhDeltaPhiUeRightCharged(0),
90     fhXEUeLeftCharged(0),           fhXEUeRightCharged(0),
91     fhPtHbpXEUeLeftCharged(0),      fhPtHbpXEUeRightCharged(0), 
92     fhZTUeLeftCharged(0),           fhZTUeRightCharged(0),
93     fhPtHbpZTUeLeftCharged(0),      fhPtHbpZTUeRightCharged(0), 
94     fhPtTrigPout(0),                fhPtTrigCharged(0),
95     fhTrigDeltaPhiCharged(0x0),     fhTrigDeltaEtaCharged(0x0),
96     fhTrigXECorr(0x0),              fhTrigXEUeCorr(0x0),
97     fhTrigZTCorr(0x0),              fhTrigZTUeCorr(0x0),
98     fhAssocPtBkg(0),  
99     fhDeltaPhiAssocPtBin(0),        fhDeltaPhiAssocPtBinHMPID(0), fhDeltaPhiAssocPtBinHMPIDAcc(0),
100     fhDeltaPhiBradAssocPtBin(0),    fhDeltaPhiBrad(0),
101     fhXEAssocPtBin(0),              fhZTAssocPtBin(0),             
102     fhDeltaPhiDeltaEtaNeutral(0), 
103     fhPhiNeutral(0),                fhEtaNeutral(0), 
104     fhDeltaPhiNeutral(0),           fhDeltaEtaNeutral(0),
105     fhDeltaPhiNeutralPt(0),         fhDeltaPhiUeNeutralPt(0), 
106     fhXENeutral(0),                 fhXEUeNeutral(0),
107     fhPtHbpXENeutral(0),            fhPtHbpXEUeNeutral(0),
108     fhZTNeutral(0),                 fhZTUeNeutral(0),
109     fhPtHbpZTNeutral(0),            fhPtHbpZTUeNeutral(0),
110     fhDeltaPhiUeLeftNeutral(0),     fhDeltaPhiUeRightNeutral(0),
111     fhXEUeLeftNeutral(0),           fhXEUeRightNeutral(0),
112     fhPtHbpXEUeLeftNeutral(0),      fhPtHbpXEUeRightNeutral(0),
113     fhZTUeLeftNeutral(0),           fhZTUeRightNeutral(0),
114     fhPtHbpZTUeLeftNeutral(0),      fhPtHbpZTUeRightNeutral(0),
115     fhPtPi0DecayRatio(0),
116     fhDeltaPhiDecayCharged(0),      fhXEDecayCharged(0), fhZTDecayCharged(0), 
117     fhDeltaPhiDecayNeutral(0),      fhXEDecayNeutral(0), fhZTDecayNeutral(0),
118     fhDeltaPhiDecayChargedAssocPtBin(0), 
119     fhXEDecayChargedAssocPtBin(0),  fhZTDecayChargedAssocPtBin(0),                
120     fh2phiLeadingParticle(0x0),     fhMCPtLeading(0),
121     fhMCEtaCharged(0),              fhMCPhiCharged(0), 
122     fhMCDeltaEtaCharged(0),         fhMCDeltaPhiCharged(0x0),
123     fhMCDeltaPhiDeltaEtaCharged(0), fhMCDeltaPhiChargedPt(0),
124     fhMCPtXECharged(0),             fhMCPtXEUeCharged(0),
125     fhMCPtHbpXECharged(0),          fhMCPtHbpXEUeCharged(0),
126     fhMCUePart(0),
127     fhMCPtZTCharged(0),             fhMCPtHbpZTCharged(0),
128     fhMCPtTrigPout(0),              fhMCPtAssocDeltaPhi(0),
129     //Mixing
130     fhNEventsTrigger(0),            
131     fhNtracksAll(0),                fhNtracksTrigger(0),            
132     fhNtracksMB(0),               
133     fhMixDeltaPhiCharged(0),        fhMixDeltaPhiDeltaEtaCharged(0),
134     fhMixDeltaPhiChargedAssocPtBin(),fhMixDeltaPhiDeltaEtaChargedAssocPtBin(0),
135     fhEventBin(0),                  fhEventMixBin(0)
136
137   //Default Ctor
138     
139   //Initialize parameters
140   InitParameters();
141 }
142
143 //_________________________________________________________________
144 AliAnaParticleHadronCorrelation::~AliAnaParticleHadronCorrelation() 
145 {
146   // Remove event containers
147   
148   if(DoOwnMix() && fListMixEvents)
149   {      
150     for(Int_t iz=0; iz < GetNZvertBin(); iz++)
151     {
152       for(Int_t ic=0; ic < GetNCentrBin(); ic++)
153       {
154          for(Int_t irp=0; irp<GetNRPBin(); irp++)
155          {
156            Int_t bin = GetEventMixBin(ic, iz, irp);
157            fListMixEvents[bin]->Delete() ;
158            delete fListMixEvents[bin] ;
159          }
160       }
161     }
162     
163     delete[] fListMixEvents;
164     
165   }
166 }
167
168 //______________________________________________________________________________________________________________________________________________________
169 void AliAnaParticleHadronCorrelation::FillChargedAngularCorrelationHistograms(const Float_t ptAssoc,  const Float_t ptTrig,   const Int_t   assocBin,
170                                                                               const Float_t phiAssoc, const Float_t phiTrig,  Float_t &     deltaPhi,
171                                                                               const Float_t etaAssoc, const Float_t etaTrig,  
172                                                                               const Bool_t  decay,    const Float_t hmpidSignal, const Int_t nTracks)
173 {
174   // Fill angular correlation related histograms
175   
176   Float_t deltaEta    = etaTrig-etaAssoc;
177           deltaPhi    = phiTrig-phiAssoc;
178   Float_t deltaPhiOrg = deltaPhi;
179   
180   if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
181   if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
182   
183   fhEtaCharged     ->Fill(ptAssoc,etaAssoc);
184   fhPhiCharged     ->Fill(ptAssoc,phiAssoc);
185   fhDeltaEtaCharged->Fill(ptTrig ,deltaEta);
186   fhDeltaPhiCharged->Fill(ptTrig ,deltaPhi);
187   
188   if(ptAssoc > 2 )           fhDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta);
189   
190   if(fDecayTrigger && decay) fhDeltaPhiDecayCharged   ->Fill(ptTrig  , deltaPhi);      
191   
192   Double_t  dphiBrad = -100;
193   if(fFillBradHisto)
194   {
195     dphiBrad = atan2(sin(deltaPhiOrg), cos(deltaPhiOrg))/TMath::Pi();//-1 to 1
196     if(TMath::Abs(dphiBrad)>0.325 && TMath::Abs(dphiBrad)<0.475)  //Hardcoded values, BAD, FIXME
197     {
198       fhAssocPtBkg->Fill(ptTrig, ptAssoc);
199     }
200     
201     if(dphiBrad<-1./3) dphiBrad += 2;
202     fhDeltaPhiBrad->Fill(ptTrig, dphiBrad);
203   }
204   
205   // Fill histograms in bins of associated particle pT
206   if(assocBin>=0)
207   {
208     fhDeltaPhiAssocPtBin[assocBin]->Fill(ptTrig, deltaPhi);
209     if (fFillBradHisto)        fhDeltaPhiBradAssocPtBin        [assocBin]->Fill(ptTrig, dphiBrad);
210     if(fDecayTrigger && decay) fhDeltaPhiDecayChargedAssocPtBin[assocBin]->Fill(ptTrig, deltaPhi);      
211     
212     if(fHMPIDCorrelation)
213     {
214       if( hmpidSignal > 0 )
215       {
216         //printf("Track pt %f with HMPID signal %f \n",pt,hmpidSignal);
217         fhDeltaPhiAssocPtBinHMPID[assocBin]->Fill(ptTrig, deltaPhi);        
218       }
219       
220       if(phiAssoc > 5*TMath::DegToRad() && phiAssoc < 20*TMath::DegToRad())
221       {
222         //printf("Track pt %f in HMPID acceptance phi %f \n ",pt,phi*TMath::RadToDeg() );
223         fhDeltaPhiAssocPtBinHMPIDAcc[assocBin]->Fill(ptTrig, deltaPhi);        
224       }
225     }
226   }
227   
228   //fill different multiplicity histogram
229   if(DoEventSelect())
230   {
231     for(Int_t im = 0; im<GetMultiBin(); im++)
232     {
233       if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
234       {
235         fhTrigDeltaPhiCharged[im]->Fill(ptTrig,deltaPhi);
236         fhTrigDeltaEtaCharged[im]->Fill(ptTrig,deltaEta);
237       }
238     }
239   }  
240 }
241
242 //____________________________________________________________________________________________________________________________________________________
243 Bool_t AliAnaParticleHadronCorrelation::FillChargedMCCorrelationHistograms(const Float_t mcAssocPt,      Float_t mcAssocPhi, const Float_t mcAssocEta,
244                                                                            const Float_t mcTrigPt, const Float_t mcTrigPhi,  const Float_t mcTrigEta)
245 {
246   // Fill MC histograms independently of AOD or ESD
247   
248   //Select only hadrons in pt range
249   if(mcAssocPt < fMinAssocPt || mcAssocPt > fMaxAssocPt) return kTRUE ; // exclude but continue
250   
251   if(mcAssocPhi < 0) mcAssocPhi+=TMath::TwoPi();    
252   
253   //remove trigger itself for correlation when use charged triggers 
254   if(TMath::Abs(mcAssocPt -mcTrigPt ) < 1e-6 && 
255      TMath::Abs(mcAssocPhi-mcTrigPhi) < 1e-6 && 
256      TMath::Abs(mcAssocEta-mcTrigEta) < 1e-6)            return kTRUE ; // exclude but continue       
257   
258   // Absolute leading?
259   if( fMakeAbsoluteLeading && mcAssocPt > mcTrigPt )     return kFALSE; // jump event
260   
261   //jump out this event if near side associated partile pt larger than trigger
262   if( fMakeNearSideLeading && mcAssocPt > mcTrigPt && 
263      TMath::Abs(mcAssocPhi-mcTrigPhi)<TMath::PiOver2() ) return kFALSE; // jump event
264   
265   Float_t mcdeltaPhi= mcTrigPhi-mcAssocPhi; 
266   if(mcdeltaPhi <= -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
267   if(mcdeltaPhi > 3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();            
268   
269   Float_t mcxE    =-mcAssocPt/mcTrigPt*TMath::Cos(mcdeltaPhi);// -(mcAssocPx*pxprim+mcAssocPy*pyprim)/(mcTrigPt*mcTrigPt);  
270   Float_t mchbpXE =-100 ;
271   if(mcxE > 0 ) mchbpXE = TMath::Log(1./mcxE);
272   
273   Float_t mczT    = mcAssocPt/mcTrigPt ;
274   Float_t mchbpZT =-100 ;
275   if(mczT > 0 ) mchbpZT = TMath::Log(1./mczT);
276   
277   //Selection within angular range
278   if( mcdeltaPhi< -TMath::PiOver2())  mcdeltaPhi+=TMath::TwoPi();
279   if( mcdeltaPhi>3*TMath::PiOver2())  mcdeltaPhi-=TMath::TwoPi();              
280   
281   Double_t mcpout = mcAssocPt*TMath::Sin(mcdeltaPhi) ; 
282   
283   if(GetDebug() > 0 )  
284     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",
285            mcAssocPt,mcAssocPhi, mcTrigPhi,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut, GetMinPt());              
286   
287   // Fill Histograms
288   fhMCEtaCharged     ->Fill(mcAssocPt, mcAssocEta);
289   fhMCPhiCharged     ->Fill(mcAssocPt, mcAssocPhi);
290   fhMCDeltaEtaCharged->Fill(mcTrigPt,    mcTrigEta-mcAssocEta);
291   fhMCDeltaPhiCharged->Fill(mcTrigPt,    mcdeltaPhi);
292   fhMCPtAssocDeltaPhi->Fill(mcAssocPt, mcdeltaPhi);
293   
294   fhMCDeltaPhiDeltaEtaCharged->Fill(mcdeltaPhi,mcTrigEta-mcAssocEta);
295   
296   //delta phi cut for correlation
297   if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) ) 
298   {
299     fhMCDeltaPhiChargedPt->Fill(mcAssocPt,mcdeltaPhi);
300     fhMCPtXECharged      ->Fill(mcTrigPt,mcxE); 
301     fhMCPtHbpXECharged   ->Fill(mcTrigPt,mchbpXE);
302     fhMCPtZTCharged      ->Fill(mcTrigPt,mczT); 
303     fhMCPtHbpZTCharged   ->Fill(mcTrigPt,mchbpZT);
304     fhMCPtTrigPout       ->Fill(mcTrigPt, mcpout) ;
305   }
306
307   //underlying event
308   if ( (mcdeltaPhi > fUeDeltaPhiMinCut) && (mcdeltaPhi < fUeDeltaPhiMaxCut) )   
309     {
310        Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
311        Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
312   
313        if(mcUexE < 0.) mcUexE = -mcUexE;
314     
315        fhMCPtXEUeCharged->Fill(mcTrigPt,mcUexE);
316        if(mcUexE > 0) fhMCPtHbpXEUeCharged->Fill(mcTrigPt,TMath::Log(1/mcUexE));
317
318        fhMCUePart->Fill(mcTrigPt);
319     }
320   
321   return kTRUE;
322
323
324 //___________________________________________________________________________________________________________________
325 void AliAnaParticleHadronCorrelation::FillChargedMomentumImbalanceHistograms(const Float_t ptTrig,   const Float_t ptAssoc, 
326                                                                              const Float_t xE,       const Float_t hbpXE, 
327                                                                              const Float_t zT,       const Float_t hbpZT, 
328                                                                              const Float_t pout,     const Float_t deltaPhi,
329                                                                              const Int_t   nTracks,  const Int_t   charge,
330                                                                              const Int_t   assocBin, const Bool_t  decay    )
331
332 {
333   // Fill mostly momentum imbalance related histograms
334   
335   fhDeltaPhiChargedPt ->Fill(ptAssoc, deltaPhi);
336   fhXECharged         ->Fill(ptTrig , xE); 
337   fhPtHbpXECharged    ->Fill(ptTrig , hbpXE);
338   fhZTCharged         ->Fill(ptTrig , zT); 
339   fhPtHbpZTCharged    ->Fill(ptTrig , hbpZT);
340   fhPtTrigPout        ->Fill(ptTrig , pout) ;
341   fhPtTrigCharged     ->Fill(ptTrig , ptAssoc) ;
342   
343   if(fDecayTrigger && decay)
344   {          
345     fhXEDecayCharged->Fill(ptTrig,xE); 
346     fhZTDecayCharged->Fill(ptTrig,zT);
347   } // photon decay pi0/eta trigger        
348   
349   if(assocBin >= 0 )//away side 
350   {
351     fhXEAssocPtBin[assocBin]->Fill(ptTrig, xE) ;
352     fhZTAssocPtBin[assocBin]->Fill(ptTrig, zT) ;
353     
354     if(fDecayTrigger && decay)
355     {          
356       fhXEDecayChargedAssocPtBin[assocBin]->Fill(ptTrig, xE); 
357       fhZTDecayChargedAssocPtBin[assocBin]->Fill(ptTrig, zT);
358     }
359   }        
360   
361   if(charge > 0)
362   {
363     fhXEPosCharged->Fill(ptTrig,xE) ;
364     fhZTPosCharged->Fill(ptTrig,zT) ;
365   }
366   else  
367   {
368     fhXENegCharged->Fill(ptTrig,xE) ;
369     fhZTNegCharged->Fill(ptTrig,zT) ;
370   }
371   
372   //fill different multiplicity histogram
373   if(DoEventSelect())
374   {
375     for(Int_t im=0; im<GetMultiBin(); im++)
376     {
377       if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
378       {
379         fhTrigXECorr[im]->Fill(ptTrig,xE);
380         fhTrigZTCorr[im]->Fill(ptTrig,zT);
381       }
382     }
383   } //multiplicity events selection
384
385
386 //_______________________________________________________________________________________________________________________
387 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventHistograms(const Float_t ptTrig,   const Float_t ptAssoc,
388                                                                            const Float_t deltaPhi, const Int_t nTracks)
389 {
390   // Fill underlying event histograms
391   
392   fhDeltaPhiUeChargedPt->Fill(ptAssoc,deltaPhi);
393   
394   Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
395   Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
396   Double_t uezT =   ptAssoc/ptTrig;
397   
398   if(uexE < 0.) uexE = -uexE;
399     
400   fhXEUeCharged->Fill(ptTrig,uexE);
401   if(uexE > 0) fhPtHbpXEUeCharged->Fill(ptTrig,TMath::Log(1/uexE));
402   
403   fhZTUeCharged->Fill(ptTrig,uezT);
404   if(uexE > 0) fhPtHbpZTUeCharged->Fill(ptTrig,TMath::Log(1/uezT));
405   
406   if(DoEventSelect())
407   {
408     for(Int_t im=0; im<GetMultiBin(); im++)
409     {
410       if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
411       {
412         fhTrigXEUeCorr[im]->Fill(ptTrig,uexE); // xE? CHECK
413         fhTrigZTUeCorr[im]->Fill(ptTrig,uezT); // zT? CHECK
414       }
415     }
416   } //multiplicity events selection
417 }
418
419 //___________________________________________________________________________________________________________________________
420 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventSidesHistograms(const Float_t ptTrig,   const Float_t ptAssoc, 
421                                                                                 const Float_t xE,       const Float_t hbpXE, 
422                                                                                 const Float_t zT,       const Float_t hbpZT, 
423                                                                                 const Float_t deltaPhi)
424 {
425   // Fill underlying event histograms to the left and right of trigger
426   
427   if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut))
428   {  
429     fhDeltaPhiUeLeftCharged->Fill(ptAssoc, deltaPhi);
430     fhXEUeLeftCharged      ->Fill(ptTrig , xE);
431     fhPtHbpXEUeLeftCharged ->Fill(ptTrig , hbpXE);
432     fhZTUeLeftCharged      ->Fill(ptTrig , zT);
433     fhPtHbpZTUeLeftCharged ->Fill(ptTrig , hbpZT);
434   }
435   
436   if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut))
437   {  
438     fhDeltaPhiUeRightCharged->Fill(ptAssoc, deltaPhi);
439     fhXEUeRightCharged      ->Fill(ptTrig , xE);
440     fhPtHbpXEUeRightCharged ->Fill(ptTrig , hbpXE);
441     fhZTUeRightCharged      ->Fill(ptTrig , zT);
442     fhPtHbpZTUeRightCharged ->Fill(ptTrig , hbpZT);
443   }
444
445
446 //______________________________________________________________________________________________________________________________
447 void AliAnaParticleHadronCorrelation::FillDecayPhotonCorrelationHistograms(const Float_t ptAssoc,     const Float_t phiAssoc, 
448                                                                            const TLorentzVector mom1, const TLorentzVector mom2,
449                                                                            const Bool_t bChargedOrNeutral)
450 {
451   // Do correlation with decay photons of triggered pi0 or eta
452   
453   // Calculate the correlation parameters
454   Float_t ptDecay1 = mom1.Pt();
455   Float_t ptDecay2 = mom2.Pt();
456   
457   Float_t zTDecay1 = -100, zTDecay2 = -100;
458   if(ptDecay1) zTDecay1 = ptAssoc/ptDecay1 ;
459   if(ptDecay2) zTDecay2 = ptAssoc/ptDecay2 ;
460   
461   Float_t deltaPhiDecay1 = mom1.Phi()-phiAssoc;
462   if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
463   if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
464   
465   Float_t deltaPhiDecay2 = mom2.Phi()-phiAssoc;
466   if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
467   if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
468   
469   Float_t xEDecay1   =-zTDecay1*TMath::Cos(deltaPhiDecay1); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
470   Float_t xEDecay2   =-zTDecay2*TMath::Cos(deltaPhiDecay2); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
471   
472   if(bChargedOrNeutral) // correlate with charges
473   {
474     fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaPhiDecay1);
475     fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaPhiDecay2);
476     
477     if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::FillDecayPhotonHistograms( Charged corr) - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
478     
479     if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
480     {
481       fhZTDecayCharged->Fill(ptDecay1,zTDecay1); 
482       fhXEDecayCharged->Fill(ptDecay1,xEDecay1); 
483     }
484     if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
485     {
486       fhZTDecayCharged->Fill(ptDecay2,zTDecay2);
487       fhXEDecayCharged->Fill(ptDecay2,xEDecay2);
488     }
489   }
490   else // correlate with neutrals
491   {
492     fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaPhiDecay1);
493     fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaPhiDecay2);
494     
495     if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::FillDecayPhotonHistograms(Neutral corr) - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
496     
497     if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
498     {
499       fhZTDecayCharged->Fill(ptDecay1,zTDecay1); 
500       fhXEDecayCharged->Fill(ptDecay1,xEDecay1); 
501     }
502     if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
503     {
504       fhZTDecayCharged->Fill(ptDecay2,zTDecay2);
505       fhXEDecayCharged->Fill(ptDecay2,xEDecay2);
506     }    
507   }
508 }
509
510 //______________________________________________________________________________________________________________________________________________________
511 void AliAnaParticleHadronCorrelation::FillNeutralAngularCorrelationHistograms(const Float_t ptAssoc,  const Float_t ptTrig,  
512                                                                               const Float_t phiAssoc, const Float_t phiTrig,  Float_t &     deltaPhi,
513                                                                               const Float_t etaAssoc, const Float_t etaTrig)
514 {
515   // Fill angular correlation related histograms
516   
517   Float_t deltaEta    = etaTrig-etaAssoc;
518   deltaPhi    = phiTrig-phiAssoc;
519   
520   if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
521   if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
522   
523   fhEtaNeutral     ->Fill(ptAssoc,etaAssoc);
524   fhPhiNeutral     ->Fill(ptAssoc,phiAssoc);
525   fhDeltaEtaNeutral->Fill(ptTrig ,deltaEta);
526   fhDeltaPhiNeutral->Fill(ptTrig ,deltaPhi);
527   
528   if(ptAssoc > 2 ) fhDeltaPhiDeltaEtaNeutral->Fill(deltaPhi, deltaEta);
529   
530 }
531
532 //_____________________________________________________________________________________________________________________________
533 void AliAnaParticleHadronCorrelation::FillNeutralUnderlyingEventSidesHistograms(const Float_t ptTrig,   const Float_t ptAssoc, 
534                                                                                 const Float_t xE,       const Float_t hbpXE, 
535                                                                                 const Float_t zT,       const Float_t hbpZT, 
536                                                                                 const Float_t deltaPhi)
537 {
538   // Fill underlying event histograms to the left and right of trigger
539   
540   if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut))
541   {  
542     fhDeltaPhiUeLeftNeutral->Fill(ptAssoc, deltaPhi);
543     fhXEUeLeftNeutral      ->Fill(ptTrig , xE);
544     fhPtHbpXEUeLeftNeutral ->Fill(ptTrig , hbpXE);
545     fhZTUeLeftNeutral      ->Fill(ptTrig , zT);
546     fhPtHbpZTUeLeftNeutral ->Fill(ptTrig , hbpZT);
547   }
548   
549   if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut))
550   {  
551     fhDeltaPhiUeRightNeutral->Fill(ptAssoc, deltaPhi);
552     fhXEUeRightNeutral      ->Fill(ptTrig , xE);
553     fhPtHbpXEUeRightNeutral ->Fill(ptTrig , hbpXE);
554     fhZTUeRightNeutral      ->Fill(ptTrig , zT);
555     fhPtHbpZTUeRightNeutral ->Fill(ptTrig , hbpZT);
556   }
557
558
559 //_____________________________________________________________
560 void AliAnaParticleHadronCorrelation::FillChargedEventMixPool()
561 {
562   // Mixed event init
563     
564   //printf("FillChargedEventMixPool for %s\n",GetInputAODName().Data());
565    
566   Int_t nTracks = GetCTSTracks()->GetEntriesFast();
567   
568   fhNtracksAll->Fill(nTracks);
569   
570   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
571   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
572   
573   if(!inputHandler) return ;
574   
575   if(inputHandler->IsEventSelected( ) & GetReader()->GetEventTriggerMask())
576   {
577     fhNtracksTrigger->Fill(nTracks);
578   }
579   
580   if(inputHandler->IsEventSelected( ) & AliVEvent::kMB)
581   {
582     
583     fhNtracksMB->Fill(nTracks);
584         
585     if(fUseMixStoredInReader && GetReader()->GetLastTracksMixedEvent() == GetEventNumber())
586     {
587       //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
588       return ; // pool filled previously for another trigger
589     }
590     Int_t eventBin = GetEventMixBin();
591     
592     //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
593     if(eventBin < 0) return;
594     
595     TObjArray * mixEventTracks = new TObjArray;
596     
597     if(fUseMixStoredInReader) 
598     {
599       fListMixEvents[eventBin] = GetReader()->GetListWithMixedEventsForTracks(eventBin);
600     }
601     
602     if(!fListMixEvents[eventBin]) fListMixEvents[eventBin] = new TList();
603     
604     //printf("%s ***** Pool Event bin : %d - nTracks %d\n",GetInputAODName().Data(),eventBin, GetCTSTracks()->GetEntriesFast());
605     
606     TList * pool = fListMixEvents[eventBin];
607     
608     TVector3 p3;  
609     for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
610     {
611       AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
612       
613       Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
614       p3.SetXYZ(mom[0],mom[1],mom[2]);
615       Float_t pt   = p3.Pt();
616       
617       //Select only hadrons in pt range
618       if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
619       
620       AliAODPWG4Particle * mixedTrack = new AliAODPWG4Particle(mom[0],mom[1],mom[2],0);
621       mixedTrack->SetDetector("CTS");
622       mixedTrack->SetChargedBit(track->Charge()>0);
623       
624       mixEventTracks->Add(mixedTrack);
625     }
626     
627     //Set the event number where the last event was added, to avoid double pool filling
628     GetReader()->SetLastTracksMixedEvent(GetEventNumber());
629     
630     pool->AddFirst(mixEventTracks);
631     mixEventTracks = 0;
632     //printf("Pool size %d, max %d\n",pool->GetSize(), GetNMaxEvMix());
633     
634     if(pool->GetSize() > GetNMaxEvMix())
635     {//Remove last event
636       TClonesArray * tmp = static_cast<TClonesArray*>(pool->Last()) ;
637       pool->RemoveLast() ;
638       delete tmp ;
639     }
640         
641   } // MB event
642   
643 }
644
645 //____________________________________________________________
646 TObjString* AliAnaParticleHadronCorrelation::GetAnalysisCuts()
647 {
648   //Save parameters used for analysis
649   TString parList ; //this will be list of parameters used for this analysis.
650   const Int_t buffersize = 560;
651   char onePar[buffersize] ;
652   
653   snprintf(onePar,buffersize,"--- AliAnaPaticleHadronCorrelation ---\n") ;
654   parList+=onePar ;     
655   snprintf(onePar,buffersize," Pt Trigger > %3.2f ",    fMinTriggerPt) ; 
656   parList+=onePar ;
657   snprintf(onePar,buffersize," %3.2f < Pt associated < %3.2f ", fMinAssocPt,   fMaxAssocPt) ; 
658   parList+=onePar ;
659   snprintf(onePar,buffersize," %3.2f < Phi trigger particle-Hadron < %3.2f ",    fDeltaPhiMinCut,   fDeltaPhiMaxCut) ; 
660   parList+=onePar ;
661   snprintf(onePar,buffersize," %3.2f < Phi trigger particle-UeHadron <  %3.2f ", fUeDeltaPhiMinCut, fUeDeltaPhiMaxCut) ; 
662   parList+=onePar ;
663   snprintf(onePar,buffersize,"Isolated Trigger?  %d\n",    fSelectIsolated) ;
664   parList+=onePar ;
665   snprintf(onePar,buffersize,"Several UE?  %d\n",          fMakeSeveralUE) ;
666   parList+=onePar ;
667   snprintf(onePar,buffersize,"Name of AOD Pi0 Branch %s ", fPi0AODBranchName.Data());
668   parList+=onePar ;
669   snprintf(onePar,buffersize,"Do Decay-hadron correlation ?  pi0 %d, decay %d", fPi0Trigger, fDecayTrigger) ;
670   parList+=onePar ;
671   snprintf(onePar,buffersize,"Select absolute leading for cluster triggers ? %d or Near Side Leading %d \n", 
672            fMakeAbsoluteLeading, fMakeNearSideLeading) ;
673   parList+=onePar ;
674   snprintf(onePar,buffersize,"Associated particle pt bins  %d: ", fNAssocPtBins) ;
675   parList+=onePar ;
676   for (Int_t ibin = 0; ibin<fNAssocPtBins; ibin++) {
677     snprintf(onePar,buffersize,"bin %d = [%2.1f,%2.1f];", ibin, fAssocPtBinLimit[ibin], fAssocPtBinLimit[ibin+1]) ;
678   }
679   parList+=onePar ;
680   
681   //Get parameters set in base class.
682   parList += GetBaseParametersList() ;
683   
684   //Get parameters set in FiducialCut class (not available yet)
685   //parlist += GetFidCut()->GetFidCutParametersList() 
686   
687   return new TObjString(parList) ;  
688   
689
690
691 //________________________________________________________________
692 TList *  AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
693 {  
694   
695   // Create histograms to be saved in output file and 
696   // store them in fOutputContainer
697   
698   TList * outputContainer = new TList() ; 
699   outputContainer->SetName("CorrelationHistos") ; 
700   
701   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();
702   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();
703   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();     
704   
705   fhPtLeading  = new TH1F ("hPtLeading","p_T distribution of leading particles", nptbins,ptmin,ptmax); 
706   fhPtLeading->SetXTitle("p_{T}^{trig} (GeV/c)");
707   
708   fhPhiLeading  = new TH2F ("hPhiLeading","#phi distribution of leading Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax); 
709   fhPhiLeading->SetYTitle("#phi (rad)");
710   
711   fhEtaLeading  = new TH2F ("hEtaLeading","#eta distribution of leading",nptbins,ptmin,ptmax, netabins,etamin,etamax); 
712   fhEtaLeading->SetYTitle("#eta ");  
713   
714   outputContainer->Add(fhPtLeading);
715   outputContainer->Add(fhPhiLeading);
716   outputContainer->Add(fhEtaLeading);
717   
718   fhPtLeadingCentrality   = new TH2F("hPtLeadingCentrality","Leading particle p_{T} vs centrality",nptbins,ptmin,ptmax,100,0.,100) ;
719   fhPtLeadingCentrality->SetXTitle("p_{T}^{trig} (GeV/c)");
720   fhPtLeadingCentrality->SetYTitle("Centrality (%)");
721   outputContainer->Add(fhPtLeadingCentrality) ;  
722   
723   fhPtLeadingEventPlane  = new TH2F("hPtLeadingEventPlane","Leading particle p_{T} vs event plane angle",nptbins,ptmin,ptmax, 100,0.,TMath::Pi()) ;
724   fhPtLeadingEventPlane->SetXTitle("p_{T}^{trig} (GeV/c)");
725   fhPtLeadingEventPlane->SetXTitle("EP angle (rad)");
726   outputContainer->Add(fhPtLeadingEventPlane) ;
727   
728   fhLeadingEventPlaneCentrality  = new TH2F("hLeadingEventPlane","Leading particle centrality vs event plane angle",100,0.,100,100,0.,TMath::Pi()) ;
729   fhLeadingEventPlaneCentrality->SetXTitle("Centrality (%)");
730   fhLeadingEventPlaneCentrality->SetYTitle("EP angle (rad)");
731   outputContainer->Add(fhLeadingEventPlaneCentrality) ;
732   
733   //Correlation with charged hadrons
734   if(GetReader()->IsCTSSwitchedOn()) 
735   {
736     fhDeltaPhiDeltaEtaCharged  = new TH2F
737     ("hDeltaPhiDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #phi_{trigger} - #phi_{h^{#pm}}",
738     ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax); 
739     fhDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
740     fhDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");    
741     
742     fhPhiCharged  = new TH2F
743     ("hPhiCharged","#phi_{h^{#pm}}  vs p_{T #pm}",
744      nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
745     fhPhiCharged->SetYTitle("#phi_{h^{#pm}} (rad)");
746     fhPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
747     
748     fhEtaCharged  = new TH2F
749     ("hEtaCharged","#eta_{h^{#pm}}  vs p_{T #pm}",
750      nptbins,ptmin,ptmax,netabins,etamin,etamax); 
751     fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
752     fhEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
753     
754     fhDeltaPhiCharged  = new TH2F
755     ("hDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
756      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
757     fhDeltaPhiCharged->SetYTitle("#Delta #phi");
758     fhDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
759     
760     fhDeltaPhiChargedPt  = new TH2F
761     ("hDeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
762      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
763     fhDeltaPhiChargedPt->SetYTitle("#Delta #phi");
764     fhDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
765     
766     fhDeltaPhiUeChargedPt  = new TH2F
767     ("hDeltaPhiUeChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}}",
768      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
769     fhDeltaPhiUeChargedPt->SetYTitle("#Delta #phi");
770     fhDeltaPhiUeChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
771     
772     fhUePart  =  new TH1F("hUePart","UE particles distribution vs pt trig",
773              nptbins,ptmin,ptmax); 
774     fhUePart->SetYTitle("dNch");
775     fhUePart->SetXTitle("p_{T trigger}");
776     
777     
778     fhDeltaEtaCharged  = new TH2F
779     ("hDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger}",
780      nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);  
781     fhDeltaEtaCharged->SetYTitle("#Delta #eta");
782     fhDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
783     
784     fhXECharged  = 
785     new TH2F("hXECharged","x_{E} for charged tracks",
786              nptbins,ptmin,ptmax,200,0.,2.); 
787     fhXECharged->SetYTitle("x_{E}");
788     fhXECharged->SetXTitle("p_{T trigger}");
789     
790     fhXEUeCharged  = 
791     new TH2F("hXEUeCharged","x_{E} for Underlying Event",
792              nptbins,ptmin,ptmax,200,0.,2.); 
793     fhXEUeCharged->SetYTitle("x_{E}");
794     fhXEUeCharged->SetXTitle("p_{T trigger}");
795     
796     fhXEPosCharged  = 
797     new TH2F("hXEPositiveCharged","x_{E} for positive charged tracks",
798              nptbins,ptmin,ptmax,200,0.,2.); 
799     fhXEPosCharged->SetYTitle("x_{E}");
800     fhXEPosCharged->SetXTitle("p_{T trigger}");
801     
802     fhXENegCharged  = 
803     new TH2F("hXENegativeCharged","x_{E} for negative charged tracks",
804              nptbins,ptmin,ptmax,200,0.,2.); 
805     fhXENegCharged->SetYTitle("x_{E}");
806     fhXENegCharged->SetXTitle("p_{T trigger}");
807     
808     fhPtHbpXECharged  = 
809     new TH2F("hHbpXECharged","#xi = ln(1/x_{E}) with charged hadrons",
810              nptbins,ptmin,ptmax,200,0.,10.); 
811     fhPtHbpXECharged->SetYTitle("ln(1/x_{E})");
812     fhPtHbpXECharged->SetXTitle("p_{T trigger}");
813     
814     fhPtHbpXEUeCharged  = 
815     new TH2F("hHbpXEUeCharged","#xi = ln(1/x_{E}) with charged hadrons,Underlying Event",
816              nptbins,ptmin,ptmax,200,0.,10.); 
817     fhPtHbpXEUeCharged->SetYTitle("ln(1/x_{E})");
818     fhPtHbpXEUeCharged->SetXTitle("p_{T trigger}");
819     
820     fhZTCharged  = 
821     new TH2F("hZTCharged","z_{T} for charged tracks",
822              nptbins,ptmin,ptmax,200,0.,2.); 
823     fhZTCharged->SetYTitle("z_{T}");
824     fhZTCharged->SetXTitle("p_{T trigger}");
825     
826     fhZTUeCharged  = 
827     new TH2F("hZTUeCharged","z_{T} for Underlying Event",
828              nptbins,ptmin,ptmax,200,0.,2.); 
829     fhZTUeCharged->SetYTitle("z_{T}");
830     fhZTUeCharged->SetXTitle("p_{T trigger}");
831     
832     fhZTPosCharged  = 
833     new TH2F("hZTPositiveCharged","z_{T} for positive charged tracks",
834              nptbins,ptmin,ptmax,200,0.,2.); 
835     fhZTPosCharged->SetYTitle("z_{T}");
836     fhZTPosCharged->SetXTitle("p_{T trigger}");
837     
838     fhZTNegCharged  = 
839     new TH2F("hZTNegativeCharged","z_{T} for negative charged tracks",
840              nptbins,ptmin,ptmax,200,0.,2.); 
841     fhZTNegCharged->SetYTitle("z_{T}");
842     fhZTNegCharged->SetXTitle("p_{T trigger}");
843     
844     fhPtHbpZTCharged  = 
845     new TH2F("hHbpZTCharged","#xi = ln(1/z_{T}) with charged hadrons",
846              nptbins,ptmin,ptmax,200,0.,10.); 
847     fhPtHbpZTCharged->SetYTitle("ln(1/z_{T})");
848     fhPtHbpZTCharged->SetXTitle("p_{T trigger}");
849     
850     fhPtHbpZTUeCharged  = 
851     new TH2F("hHbpZTUeCharged","#xi = ln(1/z_{T}) with charged hadrons,Underlying Event",
852              nptbins,ptmin,ptmax,200,0.,10.); 
853     fhPtHbpZTUeCharged->SetYTitle("ln(1/x_{E})");
854     fhPtHbpZTUeCharged->SetXTitle("p_{T trigger}");
855     
856     fhPtTrigPout  = 
857     new TH2F("hPtTrigPout","Pout with triggers",
858              nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax); 
859     fhPtTrigPout->SetYTitle("p_{out} (GeV/c)");
860     fhPtTrigPout->SetXTitle("p_{T trigger} (GeV/c)"); 
861     
862     fhPtTrigCharged  = 
863     new TH2F("hPtTrigCharged","trgger and charged tracks pt distribution",
864              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax); 
865     fhPtTrigCharged->SetYTitle("p_{T h^{#pm}} (GeV/c)");
866     fhPtTrigCharged->SetXTitle("p_{T trigger} (GeV/c)");    
867           
868     outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
869     outputContainer->Add(fhPhiCharged) ;
870     outputContainer->Add(fhEtaCharged) ;
871     outputContainer->Add(fhDeltaPhiCharged) ; 
872     outputContainer->Add(fhDeltaEtaCharged) ;
873     outputContainer->Add(fhDeltaPhiChargedPt) ;
874     outputContainer->Add(fhDeltaPhiUeChargedPt) ;
875     outputContainer->Add(fhUePart);
876
877     outputContainer->Add(fhXECharged) ;
878     outputContainer->Add(fhXEPosCharged) ;
879     outputContainer->Add(fhXENegCharged) ;
880     outputContainer->Add(fhXEUeCharged) ;
881     outputContainer->Add(fhPtHbpXECharged) ;
882     outputContainer->Add(fhPtHbpXEUeCharged) ;
883
884     outputContainer->Add(fhZTCharged) ;
885     outputContainer->Add(fhZTPosCharged) ;
886     outputContainer->Add(fhZTNegCharged) ;
887     outputContainer->Add(fhZTUeCharged) ;
888     outputContainer->Add(fhPtHbpZTCharged) ;
889     outputContainer->Add(fhPtHbpZTUeCharged) ;
890     
891     outputContainer->Add(fhPtTrigPout) ;
892     outputContainer->Add(fhPtTrigCharged) ;
893     
894     if(DoEventSelect())
895     { 
896       Int_t nMultiBins = GetMultiBin();
897       fhTrigDeltaPhiCharged = new TH2F*[nMultiBins] ;
898       fhTrigDeltaEtaCharged = new TH2F*[nMultiBins] ;
899       fhTrigXECorr          = new TH2F*[nMultiBins] ;
900       fhTrigXEUeCorr        = new TH2F*[nMultiBins] ;
901       fhTrigZTCorr          = new TH2F*[nMultiBins] ;
902       fhTrigZTUeCorr        = new TH2F*[nMultiBins] ;
903       
904       for(Int_t im=0; im<nMultiBins; im++)
905       {
906         fhTrigDeltaPhiCharged[im]  = new TH2F 
907         (Form("hTrigDeltaPhiCharged_%d",im),Form("hTrigDeltaPhiCharged_%d",im), nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax); 
908         fhTrigDeltaPhiCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
909         fhTrigDeltaPhiCharged[im]->SetYTitle("#Delta #phi");
910         
911         fhTrigDeltaEtaCharged[im]  = new TH2F 
912         (Form("hTrigDeltaEtaCharged_%d",im),Form("hTrigDeltaEtaCharged_%d",im), nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax); 
913         fhTrigDeltaEtaCharged[im]->SetXTitle("p_{T trigger} (GeV/c)");
914         fhTrigDeltaEtaCharged[im]->SetYTitle("#Delta #eta");
915         
916         fhTrigXECorr[im]  = new TH2F
917         (Form("hTrigXEPtCorr_%d",im),Form("hTrigXEPtCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.); 
918         fhTrigXECorr[im]->SetYTitle("x_{E trigger h^{#pm}}");
919         fhTrigXECorr[im]->SetXTitle("p_{T trigger}");
920         
921         fhTrigXEUeCorr[im]  = new TH2F
922         (Form("hTrigXEPtUeCorr_%d",im),Form("hTrigXEPtUeCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.); 
923         fhTrigXEUeCorr[im]->SetYTitle("x_{E trigger h^{#pm}}");
924         fhTrigXEUeCorr[im]->SetXTitle("p_{T trigger}");       
925         
926         fhTrigZTCorr[im]  = new TH2F
927         (Form("hTrigZTPtCorr_%d",im),Form("hTrigZTPtCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.); 
928         fhTrigZTCorr[im]->SetYTitle("z_{trigger h^{#pm}}");
929         fhTrigZTCorr[im]->SetXTitle("p_{T trigger}");
930         
931         fhTrigZTUeCorr[im]  = new TH2F
932         (Form("hTrigZTPtUeCorr_%d",im),Form("hTrigZTPtUeCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.); 
933         fhTrigZTUeCorr[im]->SetYTitle("z_{trigger h^{#pm}}");
934         fhTrigZTUeCorr[im]->SetXTitle("p_{T trigger}");               
935         
936         outputContainer->Add(fhTrigDeltaPhiCharged[im]) ;
937         outputContainer->Add(fhTrigDeltaEtaCharged[im]) ;
938         outputContainer->Add(fhTrigXECorr[im]);
939         outputContainer->Add(fhTrigXEUeCorr[im]);
940         outputContainer->Add(fhTrigZTCorr[im]);
941         outputContainer->Add(fhTrigZTUeCorr[im]);
942       }
943     }
944     
945     if(fFillBradHisto)
946     {
947       fhAssocPtBkg        = new TH2F("hAssocPtBkg", " Trigger p_{T} vs associated hadron p_{T} from background",
948                                      nptbins, ptmin, ptmax,nptbins,ptmin,ptmax);
949       fhAssocPtBkg->SetXTitle("p_{T trigger}");
950       fhAssocPtBkg->SetYTitle("p_{T associated}");
951       outputContainer->Add(fhAssocPtBkg) ;
952       
953       fhDeltaPhiBrad = new TH2F("hDeltaPhiBrad","atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi vs p_{T trigger} ", 
954                                 nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
955       fhDeltaPhiBrad->SetXTitle("p_{T trigger}");
956       fhDeltaPhiBrad->SetYTitle("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi");
957       outputContainer->Add(fhDeltaPhiBrad) ;
958     }
959     
960     fhDeltaPhiAssocPtBin     = new TH2F*[fNAssocPtBins] ;
961     fhXEAssocPtBin           = new TH2F*[fNAssocPtBins] ;
962     fhZTAssocPtBin           = new TH2F*[fNAssocPtBins] ;
963     
964     if(fFillBradHisto)  
965       fhDeltaPhiBradAssocPtBin = new TH2F*[fNAssocPtBins] ;
966     
967     if(fPi0Trigger || fDecayTrigger)
968     {
969       fhDeltaPhiAssocPtBin       = new TH2F*[fNAssocPtBins] ;
970       fhXEAssocPtBin             = new TH2F*[fNAssocPtBins] ;
971       fhZTAssocPtBin             = new TH2F*[fNAssocPtBins] ;
972       fhXEDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins] ;
973       fhZTDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins] ;
974       fhDeltaPhiDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins] ;
975     }
976
977     if(fHMPIDCorrelation)
978     {
979       fhDeltaPhiAssocPtBinHMPID   = new TH2F*[fNAssocPtBins] ;
980       fhDeltaPhiAssocPtBinHMPIDAcc= new TH2F*[fNAssocPtBins] ;
981     }
982     
983     for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
984     {
985       fhDeltaPhiAssocPtBin[i] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
986                                          Form("#Delta #phi vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
987                                          nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
988       fhDeltaPhiAssocPtBin[i]->SetXTitle("p_{T trigger}");
989       fhDeltaPhiAssocPtBin[i]->SetYTitle("#Delta #phi");
990        
991       fhXEAssocPtBin[i]       = new TH2F(Form("hXEAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
992                                          Form("x_{E} vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
993                                          nptbins, ptmin, ptmax,200, 0.0, 2.0);
994       fhXEAssocPtBin[i]->SetXTitle("p_{T trigger}");
995       fhXEAssocPtBin[i]->SetYTitle("x_{E}");
996       
997       fhZTAssocPtBin[i]       = new TH2F(Form("hZTAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
998                                          Form("z_{T} vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
999                                          nptbins, ptmin, ptmax,200, 0.0, 2.0);
1000       fhZTAssocPtBin[i]->SetXTitle("p_{T trigger}");
1001       fhZTAssocPtBin[i]->SetYTitle("z_{T}");
1002       
1003       outputContainer->Add(fhDeltaPhiAssocPtBin[i]) ;
1004       outputContainer->Add(fhXEAssocPtBin[i]);
1005       outputContainer->Add(fhZTAssocPtBin[i]);
1006       
1007       if(fPi0Trigger || fDecayTrigger) 
1008       {
1009         fhDeltaPhiDecayChargedAssocPtBin[i] = new TH2F(Form("hDeltaPhiPtDecayChargedAssocPt%2.1f_%2.1f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1010                                            Form("#Delta #phi vs p_{T trigger} tagged as decay for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1011                                            nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1012         fhDeltaPhiDecayChargedAssocPtBin[i]->SetXTitle("p_{T trigger}");
1013         fhDeltaPhiDecayChargedAssocPtBin[i]->SetYTitle("#Delta #phi");
1014         
1015         fhXEDecayChargedAssocPtBin[i]       = new TH2F(Form("hXEDecayChargedAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1016                                            Form("x_{E} vs p_{T trigger} tagged as decay for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1017                                            nptbins, ptmin, ptmax,200, 0.0, 2.0);
1018         fhXEDecayChargedAssocPtBin[i]->SetXTitle("p_{T trigger}");
1019         fhXEDecayChargedAssocPtBin[i]->SetYTitle("x_{E}");
1020         
1021         fhZTDecayChargedAssocPtBin[i]       = new TH2F(Form("hZTDecayChargedAssocPtBin%1.f_%1.f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1022                                            Form("z_{T} vs p_{T trigger} tagged as decay for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1023                                            nptbins, ptmin, ptmax,200, 0.0, 2.0);
1024         fhZTDecayChargedAssocPtBin[i]->SetXTitle("p_{T trigger}");
1025         fhZTDecayChargedAssocPtBin[i]->SetYTitle("z_{T}");
1026         
1027         outputContainer->Add(fhDeltaPhiDecayChargedAssocPtBin[i]) ;
1028         outputContainer->Add(fhXEDecayChargedAssocPtBin[i]);
1029         outputContainer->Add(fhZTDecayChargedAssocPtBin[i]);
1030         
1031       }
1032       
1033       if(fFillBradHisto) 
1034       {
1035         fhDeltaPhiBradAssocPtBin[i] = new TH2F(Form("hDeltaPhiBradPtAssocPt%2.1f_%2.1f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1036                                                Form("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1037                                                nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
1038         fhDeltaPhiBradAssocPtBin[i]->SetXTitle("p_{T trigger}");
1039         fhDeltaPhiBradAssocPtBin[i]->SetYTitle("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi");
1040         outputContainer->Add(fhDeltaPhiBradAssocPtBin[i]) ;
1041       }       
1042       
1043       if(fHMPIDCorrelation)
1044       {
1045         fhDeltaPhiAssocPtBinHMPID[i] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1fHMPID", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1046                                                 Form("#Delta #phi vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f], with track having HMPID signal", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1047                                                 nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1048         fhDeltaPhiAssocPtBinHMPID[i]->SetXTitle("p_{T trigger}");
1049         fhDeltaPhiAssocPtBinHMPID[i]->SetYTitle("#Delta #phi");      
1050         
1051         fhDeltaPhiAssocPtBinHMPIDAcc[i] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1fHMPIDAcc", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1052                                                    Form("#Delta #phi vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f], with track within 5<phi<20 deg", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1053                                                    nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1054         fhDeltaPhiAssocPtBinHMPIDAcc[i]->SetXTitle("p_{T trigger}");
1055         fhDeltaPhiAssocPtBinHMPIDAcc[i]->SetYTitle("#Delta #phi"); 
1056         
1057         outputContainer->Add(fhDeltaPhiAssocPtBinHMPID[i]) ;
1058         outputContainer->Add(fhDeltaPhiAssocPtBinHMPIDAcc[i]) ;
1059         
1060       }      
1061     }
1062     
1063     if(fPi0Trigger || fDecayTrigger)
1064     {
1065       if(fPi0Trigger)
1066       {
1067         fhPtPi0DecayRatio  = new TH2F
1068         ("hPtPi0DecayRatio","p_T of #pi^{0} and the ratio of pt for two decay", 
1069          nptbins,ptmin,ptmax, 100,0.,2.); 
1070         fhPtPi0DecayRatio->SetXTitle("p_{T}^{#pi^{0}} (GeV/c)");
1071         fhPtPi0DecayRatio->SetYTitle("p_{T}^{Decay}/p_{T}^{#pi^{0}}");
1072         outputContainer->Add(fhPtPi0DecayRatio) ; 
1073       }
1074       
1075       fhDeltaPhiDecayCharged  = new TH2F
1076       ("hDeltaPhiDecayCharged","#phi_{Decay} - #phi_{h^{#pm}} vs p_{T Decay}",
1077        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax); 
1078       fhDeltaPhiDecayCharged->SetYTitle("#Delta #phi");
1079       fhDeltaPhiDecayCharged->SetXTitle("p_{T Decay} (GeV/c)");
1080       
1081       fhXEDecayCharged  = 
1082       new TH2F("hXEDecayCharged","x_{E}  Decay",
1083                nptbins,ptmin,ptmax,200,0.,2.); 
1084       fhXEDecayCharged->SetYTitle("x_{E}");
1085       fhXEDecayCharged->SetXTitle("p_{T decay}");
1086       
1087       fhZTDecayCharged  = 
1088       new TH2F("hZTDecayCharged","z_{trigger h^{#pm}} = p_{T h^{#pm}} / p_{T Decay}",
1089                nptbins,ptmin,ptmax,200,0.,2.); 
1090       fhZTDecayCharged->SetYTitle("z_{decay h^{#pm}}");
1091       fhZTDecayCharged->SetXTitle("p_{T decay}");      
1092       
1093       outputContainer->Add(fhDeltaPhiDecayCharged) ; 
1094       outputContainer->Add(fhXEDecayCharged) ;
1095       outputContainer->Add(fhZTDecayCharged) ;
1096     }    
1097     
1098     if(fMakeSeveralUE)
1099     { 
1100       fhDeltaPhiUeLeftCharged  = new TH2F
1101       ("hDeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE left side range of trigger particles",
1102        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1103       fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #phi");
1104       fhDeltaPhiUeLeftCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
1105       outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
1106       
1107       fhDeltaPhiUeRightCharged  = new TH2F
1108       ("hDeltaPhiUeRightChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs p_{T Ueh^{#pm}} with UE right side range of trigger particles",
1109        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1110       fhDeltaPhiUeRightCharged->SetYTitle("#Delta #phi");
1111       fhDeltaPhiUeRightCharged->SetXTitle("p_{T h^{#pm}} (GeV/c)");
1112       outputContainer->Add(fhDeltaPhiUeRightCharged) ;
1113       
1114       fhXEUeLeftCharged  = 
1115       new TH2F("hXEUeChargedLeft","x_{E} with UE left side of trigger",
1116                nptbins,ptmin,ptmax,200,0.,2.); 
1117       fhXEUeLeftCharged->SetYTitle("x_{E Ueh^{#pm}}");
1118       fhXEUeLeftCharged->SetXTitle("p_{T trigger}");
1119       outputContainer->Add(fhXEUeLeftCharged) ;
1120       
1121       fhXEUeRightCharged  = 
1122       new TH2F("hXEUeChargedRight","x_{E h^{#pm}} with UE right side of trigger",
1123                nptbins,ptmin,ptmax,200,0.,2.); 
1124       fhXEUeRightCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
1125       fhXEUeRightCharged->SetXTitle("p_{T trigger}");
1126       outputContainer->Add(fhXEUeRightCharged) ;
1127       
1128       fhPtHbpXEUeLeftCharged  = 
1129       new TH2F("hHbpXEUeChargedLeft","#xi = ln(1/x_{E}) with charged UE left side of trigger",
1130                nptbins,ptmin,ptmax,200,0.,10.); 
1131       fhPtHbpXEUeLeftCharged->SetYTitle("ln(1/x_{E})");
1132       fhPtHbpXEUeLeftCharged->SetXTitle("p_{T trigger}");
1133       outputContainer->Add(fhPtHbpXEUeLeftCharged) ;
1134       
1135       fhPtHbpXEUeRightCharged  = 
1136       new TH2F("hHbpXEUeChargedRight","#xi = ln(1/x_{E}) with charged UE right side of trigger",
1137                nptbins,ptmin,ptmax,200,0.,10.); 
1138       fhPtHbpXEUeRightCharged->SetYTitle("ln(1/x_{E})");
1139       fhPtHbpXEUeRightCharged->SetXTitle("p_{T trigger}");
1140       outputContainer->Add(fhPtHbpXEUeRightCharged) ;
1141       
1142       fhZTUeLeftCharged  = 
1143       new TH2F("hZTUeChargedLeft","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE left side of trigger",
1144                nptbins,ptmin,ptmax,200,0.,2.); 
1145       fhZTUeLeftCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
1146       fhZTUeLeftCharged->SetXTitle("p_{T trigger}");
1147       outputContainer->Add(fhZTUeLeftCharged) ;
1148       
1149       fhZTUeRightCharged  = 
1150       new TH2F("hZTUeChargedRight","z_{trigger h^{#pm}} = p_{T Ueh^{#pm}} / p_{T trigger} with UE right side of trigger",
1151                nptbins,ptmin,ptmax,200,0.,2.); 
1152       fhZTUeRightCharged->SetYTitle("z_{trigger Ueh^{#pm}}");
1153       fhZTUeRightCharged->SetXTitle("p_{T trigger}");
1154       outputContainer->Add(fhZTUeRightCharged) ;      
1155       
1156       fhPtHbpZTUeLeftCharged  = 
1157       new TH2F("hHbpZTUeChargedLeft","#xi = ln(1/z_{T}) with charged UE left side of trigger",
1158                nptbins,ptmin,ptmax,200,0.,10.); 
1159       fhPtHbpZTUeLeftCharged->SetYTitle("ln(1/z_{T})");
1160       fhPtHbpZTUeLeftCharged->SetXTitle("p_{T trigger}");
1161       outputContainer->Add(fhPtHbpZTUeLeftCharged) ;
1162       
1163       fhPtHbpZTUeRightCharged  = 
1164       new TH2F("hHbpZTUeChargedRight","#xi = ln(1/z_{T}) with charged UE right side of trigger",
1165                nptbins,ptmin,ptmax,200,0.,10.); 
1166       fhPtHbpZTUeRightCharged->SetYTitle("ln(1/z_{T})");
1167       fhPtHbpZTUeRightCharged->SetXTitle("p_{T trigger}");
1168       outputContainer->Add(fhPtHbpZTUeRightCharged) ;
1169       
1170     }  
1171   }  //Correlation with charged hadrons
1172
1173   //Correlation with neutral hadrons
1174   if(fNeutralCorr)
1175   {
1176     fhDeltaPhiDeltaEtaNeutral  = new TH2F
1177     ("hDeltaPhiDeltaEtaNeutral","#phi_{trigger} - #phi_{h^{0}} vs #eta_{trigger} - #eta_{h^{0}}",
1178      ndeltaphibins ,deltaphimin,deltaphimax, ndeltaetabins ,deltaetamin,deltaetamax); 
1179     fhDeltaPhiDeltaEtaNeutral->SetXTitle("#Delta #phi");
1180     fhDeltaPhiDeltaEtaNeutral->SetYTitle("#Delta #eta");   
1181           
1182     fhPhiNeutral  = new TH2F
1183     ("hPhiNeutral","#phi_{#pi^{0}}  vs p_{T #pi^{0}}",
1184      nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
1185     fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
1186     fhPhiNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
1187     
1188     fhEtaNeutral  = new TH2F
1189     ("hEtaNeutral","#eta_{#pi^{0}}  vs p_{T #pi^{0}}",
1190      nptbins,ptmin,ptmax,netabins,etamin,etamax); 
1191     fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
1192     fhEtaNeutral->SetXTitle("p_{T #pi^{0}} (GeV/c)");
1193     
1194     fhDeltaPhiNeutral  = new TH2F
1195     ("hDeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T trigger}",
1196      nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
1197     fhDeltaPhiNeutral->SetYTitle("#Delta #phi");
1198     fhDeltaPhiNeutral->SetXTitle("p_{T trigger} (GeV/c)");
1199     
1200     fhDeltaPhiNeutralPt  = new TH2F
1201     ("hDeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
1202      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax); 
1203     fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi");
1204     fhDeltaPhiNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
1205     
1206     fhDeltaPhiUeNeutralPt  = new TH2F
1207     ("hDeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs p_{T #pi^{0}}}",
1208      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax); 
1209     fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi");
1210     fhDeltaPhiUeNeutralPt->SetXTitle("p_{T h^{0}} (GeV/c)");
1211     
1212     fhDeltaEtaNeutral  = new TH2F
1213     ("hDeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs p_{T trigger}",
1214      nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);  
1215     fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
1216     fhDeltaEtaNeutral->SetXTitle("p_{T trigger} (GeV/c)");
1217     
1218     fhXENeutral  = 
1219     new TH2F("hXENeutral","x_{E} for #pi^{0} associated",
1220              nptbins,ptmin,ptmax,200,0.,2.); 
1221     fhXENeutral->SetYTitle("x_{E}");
1222     fhXENeutral->SetXTitle("p_{T trigger}");
1223     
1224     fhXEUeNeutral  = 
1225     new TH2F("hXEUeNeutral","x_{E} for #pi^{0} associated",
1226              nptbins,ptmin,ptmax,200,0.,2.); 
1227     fhXEUeNeutral->SetYTitle("x_{E}");
1228     fhXEUeNeutral->SetXTitle("p_{T trigger}");
1229     
1230     fhPtHbpXENeutral  = 
1231     new TH2F("hHbpXENeutral","#xi = ln(1/x_{E})for #pi^{0} associated",
1232              nptbins,ptmin,ptmax,200,0.,10.); 
1233     fhPtHbpXENeutral->SetYTitle("ln(1/x_{E})");
1234     fhPtHbpXENeutral->SetXTitle("p_{T trigger}");
1235     
1236     fhPtHbpXEUeNeutral  = 
1237     new TH2F("hHbpXEUeNeutral","#xi = ln(1/x_{E}) for #pi^{0} associated",
1238              nptbins,ptmin,ptmax,200,0.,10.); 
1239     fhPtHbpXEUeNeutral->SetYTitle("ln(1/x_{E})");
1240     fhPtHbpXEUeNeutral->SetXTitle("p_{T trigger}");
1241     
1242     fhZTNeutral  = 
1243     new TH2F("hZTNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger} for #pi^{0} associated",
1244              nptbins,ptmin,ptmax,200,0.,2.); 
1245     fhZTNeutral->SetYTitle("z_{trigger #pi^{0}}");
1246     fhZTNeutral->SetXTitle("p_{T trigger}");
1247     
1248     fhZTUeNeutral  = 
1249     new TH2F("hZTUeNeutral","z_{trigger #pi} = p_{T #pi^{0}} / p_{T trigger} for #pi^{0} associated",
1250              nptbins,ptmin,ptmax,200,0.,2.); 
1251     fhZTUeNeutral->SetYTitle("z_{trigger #pi^{0}}");
1252     fhZTUeNeutral->SetXTitle("p_{T trigger}");
1253     
1254     fhPtHbpZTNeutral  = 
1255     new TH2F("hHbpZTNeutral","#xi = ln(1/x_{E}) for #pi^{0} associated",
1256              nptbins,ptmin,ptmax,200,0.,10.); 
1257     fhPtHbpZTNeutral->SetYTitle("ln(1/z_{T})");
1258     fhPtHbpZTNeutral->SetXTitle("p_{T trigger}");
1259     
1260     fhPtHbpZTUeNeutral  = 
1261     new TH2F("hHbpZTUeNeutral","#xi = ln(1/x_{E}) for #pi^{0} associated",
1262              nptbins,ptmin,ptmax,200,0.,10.); 
1263     fhPtHbpXEUeNeutral->SetYTitle("ln(1/z_{T})");
1264     fhPtHbpXEUeNeutral->SetXTitle("p_{T trigger}");    
1265     
1266     outputContainer->Add(fhDeltaPhiDeltaEtaNeutral); 
1267     outputContainer->Add(fhPhiNeutral) ;  
1268     outputContainer->Add(fhEtaNeutral) ;   
1269     outputContainer->Add(fhDeltaPhiNeutral) ; 
1270     outputContainer->Add(fhDeltaPhiNeutralPt) ; 
1271     outputContainer->Add(fhDeltaPhiUeNeutralPt) ; 
1272     outputContainer->Add(fhDeltaEtaNeutral) ; 
1273     outputContainer->Add(fhXENeutral) ;
1274     outputContainer->Add(fhXEUeNeutral) ;  
1275     outputContainer->Add(fhPtHbpXENeutral) ;
1276     outputContainer->Add(fhPtHbpXEUeNeutral) ;    
1277     outputContainer->Add(fhZTNeutral) ;
1278     outputContainer->Add(fhZTUeNeutral) ;  
1279     outputContainer->Add(fhPtHbpZTNeutral) ;
1280     outputContainer->Add(fhPtHbpZTUeNeutral) ;    
1281     
1282     if(fPi0Trigger || fDecayTrigger)
1283     {
1284       fhDeltaPhiDecayNeutral  = new TH2F
1285       ("hDeltaPhiDecayNeutral","#phi_{Decay} - #phi_{h^{0}} vs p_{T Decay}",
1286        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);  
1287       fhDeltaPhiDecayNeutral->SetYTitle("#Delta #phi");
1288       fhDeltaPhiDecayNeutral->SetXTitle("p_{T Decay} (GeV/c)");
1289       
1290       fhXEDecayNeutral  = 
1291       new TH2F("hXEDecayNeutral","x_{E} for decay trigger",
1292                nptbins,ptmin,ptmax,200,0.,2.); 
1293       fhXEDecayNeutral->SetYTitle("x_{E}");
1294       fhXEDecayNeutral->SetXTitle("p_{T decay}");
1295       
1296       fhZTDecayNeutral  = 
1297       new TH2F("hZTDecayNeutral","z_{trigger h^{0}} = p_{T h^{0}} / p_{T Decay}",
1298                nptbins,ptmin,ptmax,200,0.,2.); 
1299       fhZTDecayNeutral->SetYTitle("z_{h^{0}}");
1300       fhZTDecayNeutral->SetXTitle("p_{T decay}");      
1301       
1302       outputContainer->Add(fhDeltaPhiDecayNeutral) ; 
1303       outputContainer->Add(fhXEDecayNeutral) ;      
1304       outputContainer->Add(fhZTDecayNeutral) ;
1305
1306     }
1307     
1308     if(fMakeSeveralUE)
1309     { 
1310       fhDeltaPhiUeLeftNeutral  = new TH2F
1311       ("hDeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T h^{0}} with neutral UE left side range of trigger particles",
1312        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax); 
1313       fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #phi");
1314       fhDeltaPhiUeLeftNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
1315       outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
1316       
1317       fhDeltaPhiUeRightNeutral  = new TH2F
1318       ("hDeltaPhiUeRightNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs p_{T Ueh^{0}} with neutral UE right side range of trigger particles",
1319        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax); 
1320       fhDeltaPhiUeRightNeutral->SetYTitle("#Delta #phi");
1321       fhDeltaPhiUeRightNeutral->SetXTitle("p_{T h^{0}} (GeV/c)");
1322       outputContainer->Add(fhDeltaPhiUeRightNeutral) ;
1323       
1324       fhXEUeLeftNeutral  = 
1325       new TH2F("hXEUeNeutralLeft","x_{E} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE left side of trigger",
1326                nptbins,ptmin,ptmax,140,0.,2.); 
1327       fhXEUeLeftNeutral->SetYTitle("z_{trigger Ueh^{0}}");
1328       fhXEUeLeftNeutral->SetXTitle("p_{T trigger}");
1329       outputContainer->Add(fhXEUeLeftNeutral) ;
1330       
1331       fhXEUeRightNeutral  = 
1332       new TH2F("hXEUeNeutralRight","x_{E} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE right side of trigger",
1333                nptbins,ptmin,ptmax,200,0.,2.); 
1334       fhXEUeRightNeutral->SetYTitle("z_{trigger Ueh^{0}}");
1335       fhXEUeRightNeutral->SetXTitle("p_{T trigger}");
1336       outputContainer->Add(fhXEUeRightNeutral) ;
1337       
1338       fhPtHbpXEUeLeftNeutral  = 
1339       new TH2F("hHbpXEUeNeutralLeft","#xi = ln(1/x_{E}) with neutral UE left side of trigger",
1340                nptbins,ptmin,ptmax,200,0.,10.); 
1341       fhPtHbpXEUeLeftNeutral->SetYTitle("ln(1/x_{E})");
1342       fhPtHbpXEUeLeftNeutral->SetXTitle("p_{T trigger}");
1343       outputContainer->Add(fhPtHbpXEUeLeftNeutral) ;
1344       
1345       fhPtHbpXEUeRightNeutral  = 
1346       new TH2F("hHbpXEUeNeutralRight","#xi = ln(1/x_{E}) with neutral UE right side of trigger",
1347                nptbins,ptmin,ptmax,200,0.,10.); 
1348       fhPtHbpXEUeRightNeutral->SetYTitle("ln(1/x_{E})");
1349       fhPtHbpXEUeRightNeutral->SetXTitle("p_{T trigger}");
1350       outputContainer->Add(fhPtHbpXEUeRightNeutral) ;
1351       
1352       fhZTUeLeftNeutral  = 
1353       new TH2F("hZTUeNeutralLeft","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE left side of trigger",
1354                nptbins,ptmin,ptmax,140,0.,2.); 
1355       fhZTUeLeftNeutral->SetYTitle("z_{trigger Ueh^{0}}");
1356       fhZTUeLeftNeutral->SetXTitle("p_{T trigger}");
1357       outputContainer->Add(fhZTUeLeftNeutral) ;
1358       
1359       fhZTUeRightNeutral  = 
1360       new TH2F("hZTUeNeutralRight","z_{trigger h^{0}} = p_{T Ueh^{0}} / p_{T trigger} with neutral UE right side of trigger",
1361                nptbins,ptmin,ptmax,200,0.,2.); 
1362       fhZTUeRightNeutral->SetYTitle("z_{trigger Ueh^{0}}");
1363       fhZTUeRightNeutral->SetXTitle("p_{T trigger}");
1364       outputContainer->Add(fhZTUeRightNeutral) ;
1365       
1366       fhPtHbpZTUeLeftNeutral  = 
1367       new TH2F("hHbpZTUeNeutralLeft","#xi = ln(1/z_{T}) with neutral UE left side of trigger",
1368                nptbins,ptmin,ptmax,200,0.,10.); 
1369       fhPtHbpZTUeLeftNeutral->SetYTitle("ln(1/z_{T})");
1370       fhPtHbpZTUeLeftNeutral->SetXTitle("p_{T trigger}");
1371       outputContainer->Add(fhPtHbpZTUeLeftNeutral) ;
1372       
1373       fhPtHbpZTUeRightNeutral  = 
1374       new TH2F("hHbpZTUeNeutralRight","#xi = ln(1/z_{T}) with neutral UE right side of trigger",
1375                nptbins,ptmin,ptmax,200,0.,10.); 
1376       fhPtHbpZTUeRightNeutral->SetYTitle("ln(1/z_{T})");
1377       fhPtHbpZTUeRightNeutral->SetXTitle("p_{T trigger}");
1378       outputContainer->Add(fhPtHbpZTUeRightNeutral) ;
1379       
1380     }  
1381         
1382   }//Correlation with neutral hadrons
1383   
1384   //if data is MC, fill more histograms
1385   if(IsDataMC())
1386   {
1387     fh2phiLeadingParticle=new TH2F("h2phiLeadingParticle","#phi resolustion for trigger particles",nptbins,ptmin,ptmax,100,-1,1);
1388     fh2phiLeadingParticle->GetXaxis()->SetTitle("p_{T gen Leading} (GeV/c)");
1389     fh2phiLeadingParticle->GetYaxis()->SetTitle("(#phi_{rec}-#phi_{gen})/#phi_{gen}");
1390
1391     fhMCPtLeading  = new TH1F ("hMCPtLeading","MC : p_T distribution of leading particles", nptbins,ptmin,ptmax); 
1392     fhMCPtLeading->SetXTitle("p_{T}^{trig} (GeV/c)");
1393     
1394     fhMCEtaCharged  = new TH2F
1395     ("hMCEtaCharged","MC #eta_{h^{#pm}}  vs p_{T #pm}",
1396      nptbins,ptmin,ptmax,netabins,etamin,etamax); 
1397     fhMCEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
1398     fhMCEtaCharged->SetXTitle("p_{T #pm} (GeV/c)");
1399     
1400     fhMCPhiCharged  = new TH2F
1401     ("hMCPhiCharged","#MC phi_{h^{#pm}}  vs p_{T #pm}",
1402      200,ptmin,ptmax,nphibins,phimin,phimax); 
1403     fhMCPhiCharged->SetYTitle("MC #phi_{h^{#pm}} (rad)");
1404     fhMCPhiCharged->SetXTitle("p_{T #pm} (GeV/c)");
1405     
1406     fhMCDeltaPhiDeltaEtaCharged  = new TH2F
1407     ("hMCDeltaPhiDeltaEtaCharged","#MC phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
1408      140,-2.,5.,200,-2,2); 
1409     fhMCDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
1410     fhMCDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");    
1411     
1412     fhMCDeltaEtaCharged  = new TH2F
1413     ("hMCDeltaEtaCharged","MC #eta_{trigger} - #eta_{h^{#pm}} vs p_{T trigger} and p_{T assoc}",
1414      nptbins,ptmin,ptmax,200,-2,2); 
1415     fhMCDeltaEtaCharged->SetYTitle("#Delta #eta");
1416     fhMCDeltaEtaCharged->SetXTitle("p_{T trigger} (GeV/c)");
1417     
1418     fhMCDeltaPhiCharged  = new TH2F
1419     ("hMCDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
1420      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax); 
1421     fhMCDeltaPhiCharged->SetYTitle("#Delta #phi");
1422     fhMCDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
1423     
1424     fhMCDeltaPhiChargedPt  = new TH2F
1425     ("hMCDeltaPhiChargedPt","MC #phi_{trigger} - #phi_{#h^{#pm}} vs p_{T h^{#pm}}",
1426      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax); 
1427     fhMCDeltaPhiChargedPt->SetYTitle("#Delta #phi");
1428     fhMCDeltaPhiChargedPt->SetXTitle("p_{T h^{#pm}} (GeV/c)");
1429     
1430     fhMCPtXECharged  = 
1431     new TH2F("hMCPtXECharged","x_{E}",
1432              nptbins,ptmin,ptmax,200,0.,2.); 
1433     fhMCPtXECharged->SetYTitle("x_{E}");
1434     fhMCPtXECharged->SetXTitle("p_{T trigger}");
1435
1436     fhMCPtXEUeCharged  = 
1437     new TH2F("hMCPtXEUeCharged","x_{E}",
1438              nptbins,ptmin,ptmax,200,0.,2.); 
1439     fhMCPtXEUeCharged->SetYTitle("x_{E}");
1440     fhMCPtXEUeCharged->SetXTitle("p_{T trigger}");
1441     
1442     fhMCPtHbpXECharged  = 
1443     new TH2F("hMCHbpXECharged","MC #xi = ln(1/x_{E}) with charged hadrons",
1444              nptbins,ptmin,ptmax,200,0.,10.); 
1445     fhMCPtHbpXECharged->SetYTitle("ln(1/x_{E})");
1446     fhMCPtHbpXECharged->SetXTitle("p_{T trigger}");
1447
1448     fhMCPtHbpXEUeCharged =
1449     new TH2F("hMCPtHbpXEUeCharged","#xi = ln(1/x_{E}) with charged hadrons,Underlying Event",
1450              nptbins,ptmin,ptmax,200,0.,10.); 
1451     fhMCPtHbpXEUeCharged->SetYTitle("ln(1/x_{E})");
1452     fhMCPtHbpXEUeCharged->SetXTitle("p_{T trigger}");
1453
1454     fhMCUePart  = 
1455     new TH1F("hMCUePart","MC UE particles distribution vs pt trig",
1456              nptbins,ptmin,ptmax); 
1457     fhMCUePart->SetYTitle("dNch");
1458     fhMCUePart->SetXTitle("p_{T trigger}");
1459     
1460     fhMCPtZTCharged  = 
1461     new TH2F("hMCPtZTCharged","z_{T}",
1462              nptbins,ptmin,ptmax,200,0.,2.); 
1463     fhMCPtZTCharged->SetYTitle("z_{T}");
1464     fhMCPtZTCharged->SetXTitle("p_{T trigger}"); 
1465     
1466     fhMCPtHbpZTCharged  = 
1467     new TH2F("hMCHbpZTCharged","MC #xi = ln(1/z_{T}) with charged hadrons",
1468              nptbins,ptmin,ptmax,200,0.,10.); 
1469     fhMCPtHbpZTCharged->SetYTitle("ln(1/z_{T})");
1470     fhMCPtHbpZTCharged->SetXTitle("p_{T trigger}");
1471     
1472     fhMCPtTrigPout  = 
1473     new TH2F("hMCPtTrigPout","AOD MC Pout with triggers",
1474              nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax); 
1475     fhMCPtTrigPout->SetYTitle("p_{out} (GeV/c)");
1476     fhMCPtTrigPout->SetXTitle("p_{T trigger} (GeV/c)"); 
1477     
1478     fhMCPtAssocDeltaPhi  = 
1479     new TH2F("hMCPtAssocDeltaPhi","AOD MC delta phi with associated charged hadrons",
1480              nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax); 
1481     fhMCPtAssocDeltaPhi->SetYTitle("#Delta #phi");
1482     fhMCPtAssocDeltaPhi->SetXTitle("p_{T trigger} (GeV/c)"); 
1483         
1484     outputContainer->Add(fh2phiLeadingParticle);
1485     outputContainer->Add(fhMCPtLeading);
1486     outputContainer->Add(fhMCDeltaPhiDeltaEtaCharged);
1487     outputContainer->Add(fhMCPhiCharged) ;
1488     outputContainer->Add(fhMCEtaCharged) ;
1489     outputContainer->Add(fhMCDeltaEtaCharged) ;
1490     outputContainer->Add(fhMCDeltaPhiCharged) ; 
1491     
1492     outputContainer->Add(fhMCDeltaPhiChargedPt) ;
1493     outputContainer->Add(fhMCPtXECharged) ;
1494     outputContainer->Add(fhMCPtXEUeCharged) ;
1495     outputContainer->Add(fhMCPtZTCharged) ;
1496     outputContainer->Add(fhMCPtHbpXECharged) ;
1497     outputContainer->Add(fhMCPtHbpXEUeCharged);
1498     outputContainer->Add(fhMCUePart);
1499     outputContainer->Add(fhMCPtHbpZTCharged) ;
1500     outputContainer->Add(fhMCPtTrigPout) ;
1501     outputContainer->Add(fhMCPtAssocDeltaPhi) ;      
1502   } //for MC histogram
1503   
1504   if(DoOwnMix())
1505   {
1506     //create event containers
1507     
1508     if(!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForTracksExists())) 
1509     {
1510       fListMixEvents= new TList*[GetNCentrBin()*GetNZvertBin()*GetNRPBin()] ;
1511       
1512       for(Int_t ic=0; ic<GetNCentrBin(); ic++)
1513       {
1514         for(Int_t iz=0; iz<GetNZvertBin(); iz++)
1515         {
1516           for(Int_t irp=0; irp<GetNRPBin(); irp++)
1517           {
1518             //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
1519             //       ic,iz, irp, ic*GetNZvertBin()*GetNRPBin()+iz*GetNRPBin()+irp,GetNZvertBin()*GetNRPBin()*GetNCentrBin());
1520             fListMixEvents[ic*GetNZvertBin()*GetNRPBin()+iz*GetNRPBin()+irp] = new TList() ;
1521             fListMixEvents[ic*GetNZvertBin()*GetNRPBin()+iz*GetNRPBin()+irp]->SetOwner(kFALSE);
1522           }
1523         }
1524       }    
1525     }
1526     
1527     //Init the list in the reader if not done previously
1528     if(!GetReader()->ListWithMixedEventsForTracksExists() && fUseMixStoredInReader) 
1529     {
1530       //printf("%s : Set the list of events \n",GetInputAODName().Data());
1531       GetReader()->SetListWithMixedEventsForTracks(fListMixEvents);
1532     }
1533     
1534     fhEventBin=new TH1I("hEventBin","Number of real events per bin(cen,vz,rp)",
1535                         GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0, 
1536                         GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
1537     fhEventBin->SetXTitle("bin");
1538     outputContainer->Add(fhEventBin) ;
1539     
1540     fhEventMixBin=new TH1I("hEventMixBin","Number of events  per bin(cen,vz,rp)",
1541                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
1542                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
1543     fhEventMixBin->SetXTitle("bin");
1544     outputContainer->Add(fhEventMixBin) ;
1545     
1546     fhNtracksAll=new TH1F("hNtracksAll","Number of tracks w/o event trigger",2000,0,2000);
1547     outputContainer->Add(fhNtracksAll);
1548     
1549     fhNtracksTrigger=new TH1F("hNtracksTriggerEvent","Number of tracks w/ event trigger",2000,0,2000);
1550     outputContainer->Add(fhNtracksTrigger);
1551     
1552     fhNtracksMB=new TH1F("hNtracksMBEvent","Number of tracks w/ event trigger kMB",2000,0,2000);
1553     outputContainer->Add(fhNtracksMB);
1554     
1555     fhMixDeltaPhiCharged  = new TH2F
1556     ("hMixDeltaPhiCharged","Mixed event : #phi_{trigger} - #phi_{h^{#pm}} vs p_{T trigger}",
1557      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax); 
1558     fhMixDeltaPhiCharged->SetYTitle("#Delta #phi");
1559     fhMixDeltaPhiCharged->SetXTitle("p_{T trigger} (GeV/c)");
1560     outputContainer->Add(fhMixDeltaPhiCharged);
1561     
1562     fhMixDeltaPhiDeltaEtaCharged  = new TH2F
1563     ("hMixDeltaPhiDeltaEtaCharged","Mixed event : #phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
1564      ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax); 
1565     fhMixDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi");
1566     fhMixDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
1567     outputContainer->Add(fhMixDeltaPhiDeltaEtaCharged);
1568     
1569     fhMixDeltaPhiChargedAssocPtBin         = new TH2F*[fNAssocPtBins] ;
1570     fhMixDeltaPhiDeltaEtaChargedAssocPtBin = new TH2F*[fNAssocPtBins] ;
1571     
1572     for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
1573     {
1574       fhMixDeltaPhiChargedAssocPtBin[i] = new TH2F(Form("hMixDeltaPhiChargedAssocPtBin%2.1f_%2.1f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1575                                                    Form("Mixed event #Delta #phi vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1576                                                    nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
1577       fhMixDeltaPhiChargedAssocPtBin[i]->SetXTitle("p_{T trigger}");
1578       fhMixDeltaPhiChargedAssocPtBin[i]->SetYTitle("#Delta #phi");
1579       
1580       fhMixDeltaPhiDeltaEtaChargedAssocPtBin[i] = new TH2F(Form("hMixDeltaPhiDeltaEtaChargedAssocPtBin%2.1f_%2.1f", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1581                                                            Form("Mixed event #Delta #phi vs p_{T trigger} for associated p_{T} bin [%2.1f,%2.1f]", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1]), 
1582                                                            ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax); 
1583       fhMixDeltaPhiDeltaEtaChargedAssocPtBin[i]->SetXTitle("#Delta #phi");
1584       fhMixDeltaPhiDeltaEtaChargedAssocPtBin[i]->SetYTitle("#Delta #eta");
1585       
1586       outputContainer->Add(fhMixDeltaPhiChargedAssocPtBin[i]);
1587       outputContainer->Add(fhMixDeltaPhiDeltaEtaChargedAssocPtBin[i]);
1588       
1589     }          
1590   }
1591
1592   return outputContainer;
1593   
1594 }
1595
1596 //_________________________________________________________________________________________________
1597 Bool_t AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum(const AliAODPWG4Particle* trigger, 
1598                                                              TLorentzVector & mom1, 
1599                                                              TLorentzVector & mom2)
1600 {
1601   // Get the momentum of the pi0/eta assigned decay photons
1602   // In case of pi0/eta trigger, we may want to check their decay correlation, 
1603   // get their decay children
1604   
1605   Int_t indexPhoton1 = trigger->GetCaloLabel(0);
1606   Int_t indexPhoton2 = trigger->GetCaloLabel(1);
1607   Float_t ptTrig     = trigger->Pt();
1608   
1609   if(indexPhoton1!=-1 || indexPhoton2!=-1) return kFALSE;
1610   
1611   if(GetDebug() > 1) 
1612     printf("AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum() - indexPhoton1 = %d, indexPhoton2 = %d \n", indexPhoton1, indexPhoton2);
1613   
1614   TObjArray * clusters  = 0x0 ;  
1615   if(trigger->GetDetector()=="EMCAL") clusters = GetEMCALClusters() ;
1616   else                                clusters = GetPHOSClusters()  ;
1617   
1618   for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++)
1619   {
1620     AliVCluster * photon =  (AliVCluster*) (clusters->At(iclus));       
1621     if(photon->GetID()==indexPhoton1) 
1622     {
1623       photon->GetMomentum(mom1,GetVertex(0)) ;
1624       if(ptTrig) fhPtPi0DecayRatio->Fill(ptTrig, mom1.Pt()/ptTrig);
1625     }
1626     if(photon->GetID()==indexPhoton2) 
1627     {
1628       photon->GetMomentum(mom1,GetVertex(0)) ;
1629       if(ptTrig > 0) fhPtPi0DecayRatio->Fill(ptTrig, mom2.Pt()/ptTrig);
1630     } 
1631     
1632     if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum() - Photon1 = %f, Photon2 = %f \n", mom1.Pt(), mom2.Pt());
1633     
1634   } //cluster loop        
1635   
1636   return kTRUE;
1637   
1638
1639
1640 //____________________________________________________
1641 void AliAnaParticleHadronCorrelation::InitParameters()
1642 {
1643   
1644   //Initialize the parameters of the analysis.
1645   SetInputAODName("PWG4Particle");
1646   SetAODObjArrayName("Hadrons");  
1647   AddToHistogramsName("AnaHadronCorr_");
1648   
1649   SetPtCutRange(0.,300);
1650   fDeltaPhiMinCut       = 1.5 ;
1651   fDeltaPhiMaxCut       = 4.5 ;
1652   fSelectIsolated       = kFALSE;
1653   fMakeSeveralUE        = kFALSE;
1654   fUeDeltaPhiMinCut     = 1. ;
1655   fUeDeltaPhiMaxCut     = 1.5 ;
1656   
1657   fNeutralCorr          = kFALSE ;
1658   fPi0Trigger           = kFALSE ;
1659   fDecayTrigger         = kFALSE ;
1660   fHMPIDCorrelation     = kFALSE ;
1661   
1662   fMakeAbsoluteLeading  = kTRUE;
1663   fMakeNearSideLeading  = kFALSE;
1664
1665   fNAssocPtBins         = 9   ;
1666   fAssocPtBinLimit[0]   = 0.2 ; 
1667   fAssocPtBinLimit[1]   = 2.0 ; 
1668   fAssocPtBinLimit[2]   = 4.0 ; 
1669   fAssocPtBinLimit[3]   = 6.0 ; 
1670   fAssocPtBinLimit[4]   = 8.0 ; 
1671   fAssocPtBinLimit[5]   = 10. ; 
1672   fAssocPtBinLimit[6]   = 12. ;
1673   fAssocPtBinLimit[7]   = 15. ;
1674   fAssocPtBinLimit[8]   = 25. ;
1675   fAssocPtBinLimit[9]   = 50. ;
1676   
1677   fUseMixStoredInReader = kTRUE;
1678   
1679 }
1680
1681 //__________________________________________________________
1682 void  AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD()  
1683 {  
1684   //Particle-Hadron Correlation Analysis, fill AODs
1685   
1686   if(!GetInputAODBranch())
1687   {
1688     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
1689     abort();
1690   }
1691         
1692   if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation"))
1693   {
1694     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());
1695     abort();
1696   }
1697         
1698   if(GetDebug() > 1)
1699   {
1700     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Begin hadron correlation analysis, fill AODs \n");
1701     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
1702     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n",   GetCTSTracks()    ->GetEntriesFast());
1703     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetEMCALClusters()->GetEntriesFast());
1704     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n",  GetPHOSClusters() ->GetEntriesFast());
1705   }
1706   
1707   //Get the vertex and check it is not too large in z
1708   Double_t v[3] = {0,0,0}; //vertex ;
1709   GetReader()->GetVertex(v);
1710   if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;   
1711   
1712   // Fill the pool with tracks if requested
1713   if(DoOwnMix()) FillChargedEventMixPool();
1714
1715   //Loop on stored AOD particles, find leading trigger
1716   Double_t ptTrig      = fMinTriggerPt ;
1717   fLeadingTriggerIndex = -1 ;
1718   Int_t    naod        = GetInputAODBranch()->GetEntriesFast() ;
1719   for(Int_t iaod = 0; iaod < naod ; iaod++)
1720   {
1721     AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1722     
1723     // Vertex cut in case of mixing
1724     Int_t check = CheckMixedEventVertex(particle->GetCaloLabel(0), particle->GetTrackLabel(0));
1725     if(check ==  0) continue;
1726     if(check == -1) return;
1727         
1728     // find the leading particles with highest momentum
1729     if (particle->Pt() > ptTrig) 
1730     {
1731       ptTrig               = particle->Pt() ;
1732       fLeadingTriggerIndex = iaod ;
1733     }
1734   }// finish search of leading trigger particle
1735         
1736   
1737   //Do correlation with leading particle
1738   if(fLeadingTriggerIndex >= 0)
1739   {
1740           
1741     AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
1742     
1743     //check if the particle is isolated or if we want to take the isolation into account
1744     if(OnlyIsolated() && !particle->IsIsolated()) return;
1745     
1746     //Make correlation with charged hadrons
1747     Bool_t okcharged = kTRUE;
1748     Bool_t okneutral = kTRUE;
1749     if(GetReader()->IsCTSSwitchedOn() )
1750       okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kFALSE);
1751     
1752     TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
1753     if(fNeutralCorr && pi0list && pi0list->GetEntriesFast() > 0)
1754       okneutral = MakeNeutralCorrelation(particle, pi0list,kFALSE);
1755     
1756   }//Correlate leading
1757   
1758   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
1759   
1760 }
1761
1762 //_________________________________________________________________
1763 void  AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms()  
1764 {  
1765   //Particle-Hadron Correlation Analysis, fill histograms
1766   
1767   if(!GetInputAODBranch())
1768   {
1769     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
1770     abort();
1771   }
1772   
1773   if(GetDebug() > 1)
1774   {
1775     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Begin hadron correlation analysis, fill histograms \n");
1776     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
1777   }
1778     
1779   //Get the vertex and check it is not too large in z
1780   Double_t v[3] = {0,0,0}; //vertex ;
1781   GetReader()->GetVertex(v);
1782   if(!GetMixedEvent() && TMath::Abs(v[2]) > GetZvertexCut()) return ;  
1783   
1784   //Loop on stored AOD particles, find leading
1785   Double_t ptTrig    = fMinTriggerPt;
1786   if(fLeadingTriggerIndex < 0)
1787   {
1788     //Search leading if not done before
1789     Int_t    naod      = GetInputAODBranch()->GetEntriesFast() ;
1790     for(Int_t iaod = 0; iaod < naod ; iaod++)
1791     {    //loop on input trigger AOD file 
1792       AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
1793
1794       // Vertex cut in case of mixing
1795       Int_t check = CheckMixedEventVertex(particle->GetCaloLabel(0), particle->GetTrackLabel(0));
1796       if(check ==  0) continue;
1797       if(check == -1) return;
1798             
1799       //check if the particle is isolated or if we want to take the isolation into account
1800       if(OnlyIsolated() && !particle->IsIsolated()) continue;
1801       
1802       //find the leading particles with highest momentum
1803       if (particle->Pt() > ptTrig) 
1804       {
1805         ptTrig               = particle->Pt() ;
1806         fLeadingTriggerIndex = iaod ;
1807       }
1808       
1809     }// Finish search of leading trigger particle
1810   }// Search leading if not done before
1811   
1812   if(fLeadingTriggerIndex >= 0 )
1813   { //using trigger particle to do correlations
1814     
1815     AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
1816
1817     // Check if trigger is in fiducial region
1818     if(IsFiducialCutOn())
1819     {
1820       Bool_t in = GetFiducialCut()->IsInFiducialCut(*particle->Momentum(),particle->GetDetector()) ;
1821       if(! in ) return ;
1822     }
1823     
1824     // Check if the particle is isolated or if we want to take the isolation into account
1825     if(OnlyIsolated() && !particle->IsIsolated()) return;
1826     
1827     // Make correlation with charged hadrons
1828     Bool_t okcharged = kTRUE;
1829     Bool_t okneutral = kTRUE;
1830     if(GetReader()->IsCTSSwitchedOn() )
1831     {
1832       okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kTRUE);
1833       if(IsDataMC())
1834       {      
1835         MakeMCChargedCorrelation(particle);
1836       }
1837     }  
1838     
1839     TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
1840     if(fNeutralCorr && pi0list)
1841     {
1842       if(pi0list->GetEntriesFast() > 0)
1843         okneutral = MakeNeutralCorrelation(particle, pi0list,kTRUE);
1844     }
1845     
1846     // Fill leading particle histogram if correlation went well and 
1847     // no problem was found, like not absolute leading, or bad vertex in mixing.
1848     if(okcharged && okneutral)
1849     {
1850       Float_t pt = particle->Pt();
1851       fhPtLeading->Fill(pt);
1852       
1853       Float_t phi = particle->Phi();
1854       if(phi<0)phi+=TMath::TwoPi();
1855       fhPhiLeading->Fill(pt, phi);
1856       
1857       fhEtaLeading->Fill(pt, particle->Eta());
1858       //printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Leading particle : pt %f, eta %f, phi %f\n",particle->Pt(),particle->Eta(),phi);
1859       
1860       Float_t cen = GetEventCentrality();
1861       Float_t ep  = GetEventPlaneAngle();
1862       
1863       fhPtLeadingCentrality        ->Fill(pt,cen);
1864       fhPtLeadingEventPlane        ->Fill(pt,ep);
1865       fhLeadingEventPlaneCentrality->Fill(cen,ep);
1866       
1867     }//ok charged && neutral
1868   }//Aod branch loop
1869   
1870   //Reinit for next event
1871   fLeadingTriggerIndex = -1;
1872   
1873   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
1874 }
1875
1876 //___________________________________________________________________________________________________________
1877 Bool_t  AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle, 
1878                                                                 const TObjArray* pl, const Bool_t bFillHisto)
1879 {  
1880   // Charged Hadron Correlation Analysis
1881   if(GetDebug() > 1) 
1882     printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
1883     
1884   Float_t phiTrig = aodParticle->Phi();
1885   Float_t etaTrig = aodParticle->Eta();
1886   Float_t ptTrig  = aodParticle->Pt();  
1887   Bool_t   decay  = aodParticle->IsTagged();
1888
1889   Float_t pt       = -100. ;
1890   Float_t zT       = -100. ; 
1891   Float_t xE       = -100. ; 
1892   Float_t hbpXE    = -100. ; 
1893   Float_t hbpZT    = -100. ; 
1894   Float_t phi      = -100. ;
1895   Float_t eta      = -100. ;
1896   Float_t pout     = -100. ;
1897   Float_t deltaPhi = -100. ;
1898   
1899   TVector3 p3;  
1900   TLorentzVector photonMom ;    
1901   TObjArray * reftracks = 0x0;
1902   Int_t nrefs           = 0;
1903   Int_t nTracks         = GetCTSTracks()->GetEntriesFast() ;
1904   
1905   // Mixed event settings
1906   Int_t evtIndex11   = -1 ; // cluster trigger or pi0 trigger 
1907   Int_t evtIndex12   = -1 ; // pi0 trigger
1908   Int_t evtIndex13   = -1 ; // charged trigger
1909   
1910   Double_t v[3]      = {0,0,0}; //vertex ;
1911   GetReader()->GetVertex(v);
1912   
1913   if (GetMixedEvent()) 
1914   {
1915     evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
1916     evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;    
1917     evtIndex13 = GetMixedEvent()->EventIndex(aodParticle->GetTrackLabel(0)) ;
1918   }
1919   
1920   // In case of pi0/eta trigger, we may want to check their decay correlation, 
1921   // get their decay children
1922   TLorentzVector decayMom1;
1923   TLorentzVector decayMom2;
1924   Bool_t decayFound = kFALSE;
1925   if(fPi0Trigger && bFillHisto) decayFound = GetDecayPhotonMomentum(aodParticle,decayMom1, decayMom2);
1926
1927   //-----------------------------------------------------------------------
1928   //Track loop, select tracks with good pt, phi and fill AODs or histograms
1929   //-----------------------------------------------------------------------
1930
1931   for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ )
1932   {
1933     AliVTrack * track = (AliVTrack *) (pl->At(ipr)) ;
1934     
1935     Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
1936     p3.SetXYZ(mom[0],mom[1],mom[2]);
1937     pt   = p3.Pt();
1938     eta  = p3.Eta();
1939     phi  = p3.Phi() ;
1940     if(phi < 0) phi+=TMath::TwoPi();
1941     
1942     //Select only hadrons in pt range
1943     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
1944     
1945     //remove trigger itself for correlation when use charged triggers    
1946     if( track->GetID() == aodParticle->GetTrackLabel(0) || track->GetID() == aodParticle->GetTrackLabel(1) ||
1947         track->GetID() == aodParticle->GetTrackLabel(2) || track->GetID() == aodParticle->GetTrackLabel(3)   ) 
1948       continue ;
1949     
1950     //jump out this event if near side associated particle pt larger than trigger
1951     if (fMakeNearSideLeading)
1952     {
1953       if(pt > ptTrig && TMath::Abs(phi-phiTrig) < TMath::PiOver2())  return kFALSE;
1954     }
1955     //jump out this event if there is any other particle with pt larger than trigger
1956     else if(fMakeAbsoluteLeading)
1957     {
1958       if(pt > ptTrig)  return kFALSE;
1959     }
1960     
1961     //Only for mixed event
1962     Int_t evtIndex2 = 0 ; 
1963     if (GetMixedEvent()) 
1964     {
1965       evtIndex2 = GetMixedEvent()->EventIndex(track->GetID()) ;
1966       if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 || evtIndex13 == evtIndex2 ) // photon and track from different events
1967         continue ; 
1968       //vertex cut
1969       if (TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut()) 
1970         return kFALSE;
1971     }    
1972     
1973     // Fill Histograms
1974     if(bFillHisto)
1975     {      
1976
1977       if(GetDebug() > 2 ) 
1978         printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
1979             
1980       // Set the pt associated bin for the defined bins
1981       Int_t assocBin   = -1; 
1982       for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
1983       {
1984         if(pt > fAssocPtBinLimit[i] && pt < fAssocPtBinLimit[i+1]) assocBin= i; 
1985       }      
1986       
1987       // Azimuthal Angle
1988       // calculate deltaPhi for later, shift when needed
1989       FillChargedAngularCorrelationHistograms(pt,  ptTrig,  assocBin, phi, phiTrig,  deltaPhi,
1990                                               eta, etaTrig, decay, track->GetHMPIDsignal(),nTracks);
1991       
1992       // Imbalance zT/xE/pOut
1993       zT = pt/ptTrig ;
1994       if(zT > 0 ) hbpZT = TMath::Log(1./zT);
1995       else        hbpZT =-100;
1996       
1997       xE   =-pt/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
1998       //if(xE <0.)xE =-xE;
1999       if(xE > 0 ) hbpXE = TMath::Log(1./xE); 
2000       else        hbpXE =-100;
2001     
2002       pout = pt*TMath::Sin(deltaPhi) ;
2003       
2004       //delta phi cut for momentum imbalance correlation
2005       if      ( (deltaPhi > fDeltaPhiMinCut)   && (deltaPhi < fDeltaPhiMaxCut)   ) 
2006       {
2007         
2008         FillChargedMomentumImbalanceHistograms(ptTrig, pt, xE, hbpXE, zT, hbpZT, pout, deltaPhi, 
2009                                                nTracks, track->Charge(), assocBin, decay);
2010         
2011       } 
2012       if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) ) 
2013       { //UE study
2014         
2015         FillChargedUnderlyingEventHistograms(ptTrig, pt, deltaPhi, nTracks);
2016
2017         fhUePart->Fill(ptTrig);
2018         
2019       }
2020       
2021       if(fPi0Trigger && decayFound) 
2022         FillDecayPhotonCorrelationHistograms(pt, phi, decayMom1,decayMom2, kTRUE) ;
2023       
2024       //several UE calculation 
2025       if(fMakeSeveralUE) FillChargedUnderlyingEventSidesHistograms(ptTrig,pt,xE,hbpXE,zT,hbpZT,deltaPhi);
2026       
2027     } //Fill histogram 
2028     else
2029     {
2030       nrefs++;
2031       if(nrefs==1)
2032       {
2033         reftracks = new TObjArray(0);
2034         TString trackname = Form("%s+Tracks", GetAODObjArrayName().Data());
2035         reftracks->SetName(trackname.Data());
2036         reftracks->SetOwner(kFALSE);        
2037       }
2038       
2039       reftracks->Add(track);
2040       
2041     }//aod particle loop
2042   }// track loop
2043   
2044   //Fill AOD with reference tracks, if not filling histograms
2045   if(!bFillHisto && reftracks) 
2046   {
2047     aodParticle->AddObjArray(reftracks);
2048   }
2049
2050   //Own mixed event, add event and remove previous or fill the mixed histograms
2051   if(DoOwnMix() && bFillHisto)
2052   {
2053       MakeChargedMixCorrelation(aodParticle);
2054   }
2055   
2056   return kTRUE;
2057   
2058 }  
2059
2060
2061 //______________________________________________________________________________________________________________________
2062 void AliAnaParticleHadronCorrelation::MakeChargedMixCorrelation(AliAODPWG4ParticleCorrelation *aodParticle)
2063 {  
2064   // Mix current trigger with tracks in another MB event
2065   
2066   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation() - Make trigger particle - charged hadron mixed event correlation \n");
2067   
2068   if(GetMixedEvent()) return;  // This is not the mixed event from general mixing frame
2069   
2070   // Get the event with similar caracteristics
2071   //printf("MakeChargedMixCorrelation for %s\n",GetInputAODName().Data());
2072
2073   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
2074   
2075   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
2076   
2077   if(!inputHandler) return;
2078   
2079   if(!(inputHandler->IsEventSelected( ) & GetReader()->GetEventTriggerMask())) return;
2080   
2081   // Get the pool, check if it exits
2082   Int_t eventBin = GetEventMixBin();
2083
2084   fhEventBin->Fill(eventBin);
2085   
2086   //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
2087   if(eventBin < 0) return;
2088   
2089   TList * pool = 0;
2090   if(fUseMixStoredInReader) pool = GetReader()->GetListWithMixedEventsForTracks(eventBin);
2091   else                      pool = fListMixEvents[eventBin];
2092
2093   if(!pool) return ;
2094     
2095   Double_t ptTrig  = aodParticle->Pt();
2096   Double_t etaTrig = aodParticle->Eta();
2097   Double_t phiTrig = aodParticle->Phi();
2098   if(phiTrig < 0.) phiTrig+=TMath::TwoPi();
2099   
2100   if(GetDebug() > 1) 
2101     printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation() - Pool bin %d size %d, leading trigger pt=%f, phi=%f, eta=%f\n",
2102            eventBin,pool->GetSize(), ptTrig,phiTrig,etaTrig);
2103   
2104   Double_t ptAssoc  = -999.;
2105   Double_t phiAssoc = -999.;
2106   Double_t etaAssoc = -999.;
2107   Double_t deltaPhi = -999.;
2108   Double_t deltaEta = -999.;
2109   
2110   //Start from first event in pool except if in this same event the pool was filled
2111   Int_t ev0 = 0;
2112   if(GetReader()->GetLastTracksMixedEvent() == GetEventNumber()) ev0 = 1;
2113   for(Int_t ev=ev0; ev < pool->GetSize(); ev++)
2114   {
2115     TObjArray* bgTracks = static_cast<TObjArray*>(pool->At(ev));
2116     
2117     fhEventMixBin->Fill(eventBin);
2118     
2119     Int_t nTracks=bgTracks->GetEntriesFast();
2120     //printf("\t Read Pool event %d, nTracks %d\n",ev,nTracks);
2121
2122     for(Int_t j1 = 0;j1 <nTracks; j1++ )
2123     {
2124       AliAODPWG4Particle *track = (AliAODPWG4Particle*) bgTracks->At(j1) ;
2125       
2126       if(!track) continue;
2127       
2128       ptAssoc  = track->Pt();
2129       etaAssoc = track->Eta();
2130       phiAssoc = track->Phi() ;
2131       if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
2132             
2133       if(IsFiducialCutOn())
2134       {
2135         Bool_t in = GetFiducialCut()->IsInFiducialCut(*aodParticle->Momentum(),"CTS") ;
2136         if(!in) continue ;
2137       }
2138       
2139       deltaPhi = phiTrig-phiAssoc;
2140       if(deltaPhi < -TMath::PiOver2())  deltaPhi+=TMath::TwoPi();
2141       if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
2142       deltaEta = etaTrig-etaAssoc;
2143       
2144       if(GetDebug()>0)
2145         printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation(): deltaPhi= %f, deltaEta=%f\n",deltaPhi, deltaEta);
2146       
2147       // Set the pt associated bin for the defined bins
2148       Int_t assocBin   = -1; 
2149       for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2150       {
2151         if(ptAssoc > fAssocPtBinLimit[i] && ptAssoc < fAssocPtBinLimit[i+1]) assocBin= i; 
2152       }      
2153
2154       fhMixDeltaPhiCharged        ->Fill(ptTrig,  deltaPhi);
2155       fhMixDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta);
2156
2157       if(assocBin < 0) continue ; // this pt bin was not considered
2158       
2159       fhMixDeltaPhiChargedAssocPtBin[assocBin]        ->Fill(ptTrig, deltaPhi);
2160       fhMixDeltaPhiDeltaEtaChargedAssocPtBin[assocBin]->Fill(deltaPhi, deltaEta);
2161       
2162     } // track loop
2163   } // mixed event loop
2164 }
2165   
2166
2167 //________________________________________________________________________________________________________________
2168 Bool_t  AliAnaParticleHadronCorrelation::MakeNeutralCorrelation(AliAODPWG4ParticleCorrelation * const aodParticle, 
2169                                                                 const TObjArray* pi0list, const Bool_t bFillHisto)  
2170 {  
2171   // Neutral Pion Correlation Analysis
2172   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Make trigger particle - pi0 correlation, %d pi0's \n",
2173                             pi0list->GetEntriesFast());
2174   
2175   Int_t evtIndex11 = 0 ; 
2176   Int_t evtIndex12 = 0 ; 
2177   if (GetMixedEvent()) 
2178   {
2179     evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
2180     evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;    
2181   }  
2182   
2183   Float_t pt   = -100. ;
2184   Float_t zT   = -100. ; 
2185   Float_t phi  = -100. ;
2186   Float_t eta  = -100. ;
2187   Float_t xE   = -100. ; 
2188   Float_t hbpXE= -100. ; 
2189   Float_t hbpZT= -100. ; 
2190
2191   Float_t ptTrig  = aodParticle->Pt();
2192   Float_t phiTrig = aodParticle->Phi();
2193   Float_t etaTrig = aodParticle->Eta();
2194   Float_t deltaPhi= -100. ;
2195
2196   TLorentzVector photonMom ;
2197         
2198   // In case of pi0/eta trigger, we may want to check their decay correlation, 
2199   // get their decay children
2200   TLorentzVector decayMom1;
2201   TLorentzVector decayMom2;
2202   Bool_t decayFound = kFALSE;
2203   if(fPi0Trigger && bFillHisto) decayFound = GetDecayPhotonMomentum(aodParticle,decayMom1, decayMom2);  
2204   
2205   TObjArray * refpi0 = 0x0;
2206   Int_t nrefs        = 0;
2207   
2208   //Loop on stored AOD pi0
2209   
2210   Int_t naod = pi0list->GetEntriesFast();
2211   if(GetDebug() > 0) 
2212     printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() -  aod branch entries %d\n", naod);
2213   
2214   for(Int_t iaod = 0; iaod < naod ; iaod++)
2215   {
2216     AliAODPWG4Particle* pi0 =  (AliAODPWG4Particle*) (pi0list->At(iaod));
2217     
2218     Int_t evtIndex2 = 0 ; 
2219     Int_t evtIndex3 = 0 ; 
2220     if (GetMixedEvent()) 
2221     {
2222       evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(0)) ;
2223       evtIndex3 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(1)) ;
2224       
2225       if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 || 
2226           evtIndex11 == evtIndex3 || evtIndex12 == evtIndex3) // trigger and pi0 are not from different events
2227         continue ; 
2228     }      
2229
2230     pt  = pi0->Pt();
2231      
2232     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
2233     
2234     //jump out this event if near side associated particle pt larger than trigger
2235     if (fMakeNearSideLeading)
2236     {
2237       if(pt > ptTrig && TMath::Abs(phi-phiTrig) < TMath::PiOver2())  return kFALSE;
2238     }
2239     //jump out this event if there is any other particle with pt larger than trigger
2240     else if(fMakeAbsoluteLeading)
2241     {
2242       if(pt > ptTrig)  return kFALSE;
2243     }
2244     
2245     if(bFillHisto)
2246     {
2247       phi = pi0->Phi() ;
2248       eta = pi0->Eta() ;
2249       
2250       FillNeutralAngularCorrelationHistograms(pt, ptTrig, phi, phiTrig, deltaPhi, eta, etaTrig);
2251       
2252       zT  = pt/ptTrig ;
2253       xE  =-pt/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
2254       
2255       //if(xE <0.)xE =-xE;
2256       
2257       hbpXE = -100;
2258       hbpZT = -100;
2259       
2260       if(xE > 0 ) hbpXE = TMath::Log(1./xE); 
2261       if(zT > 0 ) hbpZT = TMath::Log(1./zT); 
2262       
2263       if(fPi0Trigger && decayFound)
2264         FillDecayPhotonCorrelationHistograms(pt, phi, decayMom1,decayMom2,kFALSE) ;
2265       
2266       //delta phi cut for correlation
2267       if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) ) 
2268       {
2269         fhDeltaPhiNeutralPt->Fill(pt,deltaPhi);
2270         fhXENeutral        ->Fill(ptTrig,xE); 
2271         fhPtHbpXENeutral   ->Fill(ptTrig,hbpXE); 
2272       }
2273       else if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )      
2274       {
2275         fhDeltaPhiUeNeutralPt->Fill(pt,deltaPhi);
2276         fhXEUeNeutral        ->Fill(ptTrig,xE);
2277         fhPtHbpXEUeNeutral   ->Fill(ptTrig,hbpXE); 
2278       }
2279       
2280       //several UE calculation 
2281       if(fMakeSeveralUE) FillChargedUnderlyingEventSidesHistograms(ptTrig,pt,xE,hbpXE,zT,hbpZT,deltaPhi);
2282
2283           }
2284     else
2285     {
2286       nrefs++;
2287       if(nrefs==1)
2288       {
2289         refpi0 = new TObjArray(0);
2290         refpi0->SetName(GetAODObjArrayName()+"Pi0s");
2291         refpi0->SetOwner(kFALSE);
2292       }
2293       refpi0->Add(pi0);
2294     }//put references in trigger AOD 
2295     
2296     if(GetDebug() > 2 ) 
2297       printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
2298     
2299   }//loop
2300   
2301   return kTRUE;
2302 }
2303   
2304 //_________________________________________________________________________________________________________
2305 void  AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle)
2306 {  
2307   // Charged Hadron Correlation Analysis with MC information
2308   
2309   if(GetDebug()>1)
2310     printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation() - Make trigger particle - charged hadron correlation in AOD MC level\n");
2311   
2312   AliStack         * stack        = 0x0 ;
2313   TParticle        * primary      = 0x0 ;   
2314   TClonesArray     * mcparticles0 = 0x0 ;
2315   TClonesArray     * mcparticles  = 0x0 ;
2316   AliAODMCParticle * aodprimary   = 0x0 ; 
2317   
2318   Double_t eprim   = 0 ;
2319   Double_t ptprim  = 0 ;
2320   Double_t phiprim = 0 ;
2321   Double_t etaprim = 0 ;
2322   Int_t    nTracks = 0 ;  
2323   Int_t iParticle  = 0 ;
2324   Double_t charge  = 0.;
2325
2326
2327   if(GetReader()->ReadStack())
2328   {
2329     nTracks = GetMCStack()->GetNtrack() ;
2330   }
2331   else 
2332   {
2333     nTracks = GetReader()->GetAODMCParticles()->GetEntriesFast() ;
2334   }
2335   //Int_t trackIndex[nTracks];
2336   
2337   Int_t label= aodParticle->GetLabel();
2338   if(label < 0)
2339   {
2340     if(GetDebug() > 0) printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** bad label ***:  label %d \n", label);
2341     return;
2342   }  
2343
2344   
2345   if(GetReader()->ReadStack())
2346   {
2347     stack =  GetMCStack() ;
2348     if(!stack) 
2349     {
2350       printf(" AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation- Stack not available, is the MC handler called? STOP\n");
2351       abort();
2352     }
2353     
2354     nTracks=stack->GetNprimary();
2355     if(label >=  stack->GetNtrack()) 
2356     {
2357       if(GetDebug() > 2)  printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** large label ***:  label %d, n tracks %d \n", label, stack->GetNtrack());
2358       return ;
2359     }
2360     
2361     primary = stack->Particle(label);
2362     if(!primary)
2363     {
2364       printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** no primary ***:  label %d \n", label);   
2365       return;
2366     }
2367     
2368     if(primary)
2369     {
2370       eprim    = primary->Energy();
2371       ptprim   = primary->Pt();
2372       phiprim  = primary->Phi();
2373       etaprim  = primary->Eta();
2374       
2375       if(ptprim < 0.01 || eprim < 0.01) return ;
2376       
2377       for (iParticle = 0 ; iParticle <  nTracks ; iParticle++) 
2378       {
2379         TParticle * particle = stack->Particle(iParticle);
2380         TLorentzVector momentum;
2381         
2382         //keep only final state particles
2383         if(particle->GetStatusCode()!=1) continue ;
2384         
2385         Int_t pdg = particle->GetPdgCode();     
2386         
2387         charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
2388         
2389         particle->Momentum(momentum);
2390         
2391         //---------- Charged particles ----------------------
2392         if(charge != 0)
2393         {   
2394           //Particles in CTS acceptance
2395           Bool_t inCTS =  GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
2396           
2397           if(TMath::Abs(pdg) == 11 && stack->Particle(particle->GetFirstMother())->GetPdgCode()==22) continue ;
2398           
2399           if(inCTS)
2400           {            
2401             if( label!=iParticle) // avoid trigger particle
2402             {
2403               if(!FillChargedMCCorrelationHistograms(particle->Pt(),particle->Phi(),particle->Eta(),ptprim,phiprim,etaprim)) return;
2404             }
2405           }// in CTS acceptance
2406         }// charged
2407       } //track loop
2408     } //when the leading particles could trace back to MC
2409   } //ESD MC
2410   else if(GetReader()->ReadAODMCParticles())
2411   {
2412     //Get the list of MC particles
2413     mcparticles0 = GetReader()->GetAODMCParticles(0);
2414     if(!mcparticles0) return;
2415     if(label >=mcparticles0->GetEntriesFast())
2416     {
2417       if(GetDebug() > 2)  
2418         printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** large label ***:  label %d, n tracks %d \n", label,mcparticles0->GetEntriesFast());
2419       return;
2420     }
2421     //Get the particle
2422     aodprimary = (AliAODMCParticle*) mcparticles0->At(label);
2423     if(!aodprimary)  
2424     {
2425       printf("AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation *** no AOD primary ***:  label %d \n", label);   
2426       return;
2427     }
2428
2429    
2430     if(aodprimary)
2431     {
2432       ptprim  = aodprimary->Pt();
2433       phiprim = aodprimary->Phi();
2434       etaprim = aodprimary->Eta();
2435       eprim   = aodprimary->E();
2436
2437       Bool_t lead = kFALSE;
2438       
2439       if(ptprim < 0.01 || eprim < 0.01) return ;
2440  
2441       mcparticles= GetReader()->GetAODMCParticles();
2442       for (Int_t i = 0; i < nTracks; i++) 
2443       {
2444         AliAODMCParticle *part = (AliAODMCParticle*) mcparticles->At(i);
2445          
2446         if (!part->IsPhysicalPrimary()) continue;        
2447         Int_t pdg = part->GetPdgCode(); 
2448         charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
2449         TLorentzVector momentum(part->Px(),part->Py(),part->Pz(),part->E());        
2450         if(charge != 0)
2451         {
2452           if(part->Pt()> GetReader()->GetCTSPtMin())
2453           {
2454             //Particles in CTS acceptance
2455             Bool_t inCTS =  GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
2456             Int_t indexmother=part->GetMother();
2457             if(indexmother>-1)
2458             {
2459               Int_t mPdg = ((AliAODMCParticle*) mcparticles->At(indexmother)) ->GetPdgCode();
2460               if(TMath::Abs(pdg) == 11 && mPdg == 22) continue;
2461             }
2462             
2463             if(inCTS)
2464             {            
2465               if( label!=iParticle) // avoid trigger particle
2466               {
2467                 if(!FillChargedMCCorrelationHistograms(part->Pt(),part->Phi(),part->Eta(),ptprim,phiprim,etaprim)) return;
2468                 else lead = kTRUE;
2469               }
2470             } // in acceptance
2471           } // min pt cut
2472         } //only charged particles
2473       }  //MC particle loop    
2474       if (lead) fhMCPtLeading->Fill(ptprim);
2475     } //when the leading particles could trace back to MC
2476   }// AOD MC
2477 }
2478
2479 //_____________________________________________________________________
2480 void AliAnaParticleHadronCorrelation::Print(const Option_t * opt) const
2481 {
2482   
2483   //Print some relevant parameters set for the analysis
2484   if(! opt)
2485     return;
2486   
2487   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
2488   AliAnaCaloTrackCorrBaseClass::Print(" ");
2489   printf("Pt trigger           >  %3.2f\n", fMinTriggerPt) ;
2490   printf("Pt associated hadron <  %3.2f\n", fMaxAssocPt) ; 
2491   printf("Pt associated hadron >  %3.2f\n", fMinAssocPt) ;
2492   printf("Phi trigger particle-Hadron   <  %3.2f\n", fDeltaPhiMaxCut) ; 
2493   printf("Phi trigger particle-Hadron   >  %3.2f\n", fDeltaPhiMinCut) ;
2494   printf("Phi trigger particle-UeHadron <  %3.2f\n", fUeDeltaPhiMaxCut) ; 
2495   printf("Phi trigger particle-UeHadron >  %3.2f\n", fUeDeltaPhiMinCut) ;
2496   printf("Isolated Trigger?  %d\n"     , fSelectIsolated) ;
2497   printf("Several UE?  %d\n"           , fMakeSeveralUE) ;
2498   printf("Name of AOD Pi0 Branch %s \n", fPi0AODBranchName.Data());
2499   printf("Do Decay-hadron correlation ?  %d\n", fPi0Trigger) ;
2500   printf("Select absolute leading for cluster triggers ?  %d or Near Side %d\n", 
2501          fMakeAbsoluteLeading, fMakeNearSideLeading) ;
2502   printf("Trigger pt bins  %d\n", fNAssocPtBins) ;
2503   for (Int_t ibin = 0; ibin<fNAssocPtBins; ibin++) {
2504     printf("\t bin %d = [%2.1f,%2.1f]\n", ibin, fAssocPtBinLimit[ibin], fAssocPtBinLimit[ibin+1]) ;
2505   }
2506   
2507
2508
2509 //____________________________________________________________
2510 void AliAnaParticleHadronCorrelation::SetNAssocPtBins(Int_t n)
2511 {
2512   // Set number of bins
2513   
2514   fNAssocPtBins  = n ; 
2515   
2516   
2517   if(n < 10 && n > 0)
2518   {
2519     fNAssocPtBins  = n ; 
2520   }
2521   else 
2522   {
2523     printf("n = larger than 9 or too small, set to 9 \n");
2524     fNAssocPtBins = 9;
2525   }
2526 }
2527
2528 //______________________________________________________________________________
2529 void AliAnaParticleHadronCorrelation::SetAssocPtBinLimit(Int_t ibin, Float_t pt)
2530
2531   // Set the list of limits for the trigger pt bins
2532   
2533   if(ibin <= fNAssocPtBins || ibin >= 0) 
2534   {
2535     fAssocPtBinLimit[ibin] = pt  ;
2536   }
2537   else 
2538   {
2539     printf("AliAnaParticleHadronCorrelation::SetAssocPtBinLimit() - bin  number too large %d > %d or small, nothing done\n", ibin, fNAssocPtBins) ; 
2540   }
2541 }
2542