PWGGA/CaloTrackCorrelations/AliAnaParticleHadronCorrelation.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 //_________________________________________________________________________
  17 // Class for the analysis of particle - hadron correlations
  18 // Particle (for example direct gamma) must be found in a previous analysis
  19 //-- Author: Gustavo Conesa (LNF-INFN)
  20
  21 //  Modified by Yaxian Mao:
  22 // 1. add the UE subtraction for corrlation study
  23 // 2. change the correlation variable
  24 // 3. Only use leading particle(cluster/track) as trigger for correlation (2010/07/02)
  25 // 4. Make decay photon-hadron correlations where decay contribute pi0 mass (2010/09/09)
  26 // 5. fill the pout to extract kt at the end, also to study charge asymmetry(2010/10/06)
  27 // 6. Add the possibility for event selection analysis based on vertex and multiplicity bins (10/10/2010)
  28 // 7. change the way of delta phi cut for UE study due to memory issue (reduce histograms)
  29 // 8. Add the possibility to request the absolute leading particle at the near side or not, set trigger bins, general clean-up (08/2011)
  30 //////////////////////////////////////////////////////////////////////////////
  31
  32
  33 // --- ROOT system ---
  34 //#include "TClonesArray.h"
  35 #include <TClass.h>
  36 #include <TMath.h>
  37 #include <TH2F.h>
  38 #include <TDatabasePDG.h>
  39
  40 //---- ANALYSIS system ----
  41 #include "AliNeutralMesonSelection.h"
  42 #include "AliAnaParticleHadronCorrelation.h"
  43 #include "AliCaloTrackReader.h"
  44 #include "AliAODPWG4ParticleCorrelation.h"
  45 #include "AliFiducialCut.h"
  46 #include "AliVTrack.h"
  47 #include "AliVCluster.h"
  48 #include "AliMCAnalysisUtils.h"
  49 #include "TParticle.h"
  50 #include "AliStack.h"
  51 #include "AliAODMCParticle.h"
  52 #include "AliMixedEvent.h"
  53 #include "AliAnalysisManager.h"
  54 #include "AliInputEventHandler.h"
  55 #include "AliEventplane.h"
  56
  57 ClassImp(AliAnaParticleHadronCorrelation)
  58
  59
  60 //___________________________________________________________________
  61   AliAnaParticleHadronCorrelation::AliAnaParticleHadronCorrelation():
  62     AliAnaCaloTrackCorrBaseClass(),
  63     fMinTriggerPt(0.),
  64     fMaxAssocPt(1000.),             fMinAssocPt(0.),
  65     fDeltaPhiMaxCut(0.),            fDeltaPhiMinCut(0.),
  66     fSelectIsolated(0),             fMakeSeveralUE(0),
  67     fUeDeltaPhiMaxCut(0.),          fUeDeltaPhiMinCut(0.),
  68     fPi0AODBranchName(""),          fNeutralCorr(0),
  69     fPi0Trigger(0),                 fDecayTrigger(0),
  70     fMakeAbsoluteLeading(0),        fMakeNearSideLeading(0),
  71     fLeadingTriggerIndex(-1),       fHMPIDCorrelation(0),  fFillBradHisto(0),
  72     fNAssocPtBins(0),               fAssocPtBinLimit(),
  73     fCorrelVzBin(0),
  74     fListMixTrackEvents(),          fListMixCaloEvents(),
  75     fUseMixStoredInReader(0),       fFillNeutralEventMixPool(0),
  76     fM02MaxCut(0),                  fM02MinCut(0),
  77     fFillPileUpHistograms(0),
  78     fSelectLeadingHadronAngle(0),
  79     fMinLeadHadPhi(0),              fMaxLeadHadPhi(0),
  80     fMinLeadHadPt(0),               fMaxLeadHadPt(0),
  81
  82     //Histograms
  83     fhPtTriggerInput(0),            fhPtTriggerSSCut(0),
  84     fhPtTriggerIsoCut(0),           fhPtTriggerFidCut(0),
  85     fhPtLeading(0),                 fhPtLeadingVtxBC0(0),
  86     fhPtLeadingVzBin(0),            fhPtLeadingBin(0),
  87     fhPhiLeading(0),                fhEtaLeading(0),
  88     fhPtLeadingMC(),
  89     fhPtLeadingCentrality(0),       fhPtLeadingEventPlane(0),
  90     fhLeadingEventPlaneCentrality(0),
  91     fhPtLeadingMixed(0),            fhPtLeadingMixedVzBin(0), fhPtLeadingMixedBin(0),
  92     fhPhiLeadingMixed(0),           fhEtaLeadingMixed(0),
  93     fhDeltaPhiDeltaEtaCharged(0),
  94     fhPhiCharged(0),                fhEtaCharged(0),
  95     fhDeltaPhiCharged(0),           fhDeltaEtaCharged(0),
  96     fhDeltaPhiChargedPt(0),         fhDeltaPhiUeChargedPt(0),
  97     fhUePart(0),
  98     fhXECharged(0),                 fhXECharged_Cone2(0),      fhXEUeCharged(0),
  99     fhXEPosCharged(0),              fhXENegCharged(0),
 100     fhPtHbpXECharged(0),            fhPtHbpXECharged_Cone2(0), fhPtHbpXEUeCharged(0),
 101     fhZTCharged(0),                 fhZTUeCharged(0),
 102     fhZTPosCharged(0),              fhZTNegCharged(0),
 103     fhPtHbpZTCharged(0),            fhPtHbpZTUeCharged(0),
 104     fhXEChargedMC(),                fhDeltaPhiChargedMC(),
 105     fhDeltaPhiDeltaEtaChargedPtA3GeV(0),
 106     fhDeltaPhiChargedPtA3GeV(0),    fhDeltaEtaChargedPtA3GeV(0),
 107     //Pile-Up
 108     fhDeltaPhiChargedPileUp(),      fhDeltaEtaChargedPileUp(),
 109     fhDeltaPhiChargedPtA3GeVPileUp(), fhDeltaEtaChargedPtA3GeVPileUp(),
 110     fhXEChargedPileUp(),            fhXEUeChargedPileUp(),
 111     fhZTChargedPileUp(),            fhZTUeChargedPileUp(),
 112     fhPtTrigChargedPileUp(),
 113     fhDeltaPhiChargedOtherBC(),     fhDeltaPhiChargedPtA3GeVOtherBC(),
 114     fhXEChargedOtherBC(),           fhXEUeChargedOtherBC(),
 115     fhZTChargedOtherBC(),           fhZTUeChargedOtherBC(),
 116     fhPtTrigChargedOtherBC(),
 117     fhDeltaPhiChargedBC0(),         fhDeltaPhiChargedPtA3GeVBC0(),
 118     fhXEChargedBC0(),               fhXEUeChargedBC0(),
 119     fhZTChargedBC0(),               fhZTUeChargedBC0(),
 120     fhPtTrigChargedBC0(),
 121     fhDeltaPhiChargedVtxBC0(),      fhDeltaPhiChargedPtA3GeVVtxBC0(),
 122     fhXEChargedVtxBC0(),            fhXEUeChargedVtxBC0(),
 123     fhZTChargedVtxBC0(),            fhZTUeChargedVtxBC0(),
 124     fhPtTrigChargedVtxBC0(),
 125     fhDeltaPhiUeLeftCharged(0),     fhDeltaPhiUeRightCharged(0),
 126     fhDeltaPhiUeLeftUpCharged(0),   fhDeltaPhiUeRightUpCharged(0),
 127     fhDeltaPhiUeLeftDownCharged(0), fhDeltaPhiUeRightDownCharged(0),
 128     fhXEUeLeftCharged(0),           fhXEUeRightCharged(0),
 129     fhXEUeLeftUpCharged(0),         fhXEUeRightUpCharged(0),
 130     fhXEUeLeftDownCharged(0),       fhXEUeRightDownCharged(0),
 131     fhPtHbpXEUeLeftCharged(0),      fhPtHbpXEUeRightCharged(0),
 132     fhZTUeLeftCharged(0),           fhZTUeRightCharged(0),
 133     fhPtHbpZTUeLeftCharged(0),      fhPtHbpZTUeRightCharged(0),
 134     fhPtTrigPout(0),                fhPtTrigCharged(0),
 135     fhTrigDeltaPhiCharged(0x0),     fhTrigDeltaEtaCharged(0x0),
 136     fhTrigXECorr(0x0),              fhTrigXEUeCorr(0x0),
 137     fhTrigZTCorr(0x0),              fhTrigZTUeCorr(0x0),
 138     fhAssocPtBkg(0),                fhDeltaPhiDeltaEtaAssocPtBin(0),
 139     fhDeltaPhiAssocPtBin(0),
 140     fhDeltaPhiAssocPtBinDEta08(0),  fhDeltaPhiAssocPtBinDEta0(0),
 141     fhDeltaPhiAssocPtBinHMPID(0),   fhDeltaPhiAssocPtBinHMPIDAcc(0),
 142     fhDeltaPhiBradAssocPtBin(0),    fhDeltaPhiBrad(0),
 143     fhXEAssocPtBin(0),              fhZTAssocPtBin(0),
 144     fhDeltaPhiDeltaEtaNeutral(0),
 145     fhPhiNeutral(0),                fhEtaNeutral(0),
 146     fhDeltaPhiNeutral(0),           fhDeltaEtaNeutral(0),
 147     fhDeltaPhiNeutralPt(0),         fhDeltaPhiUeNeutralPt(0),
 148     fhXENeutral(0),                 fhXEUeNeutral(0),
 149     fhPtHbpXENeutral(0),            fhPtHbpXEUeNeutral(0),
 150     fhZTNeutral(0),                 fhZTUeNeutral(0),
 151     fhPtHbpZTNeutral(0),            fhPtHbpZTUeNeutral(0),
 152     fhDeltaPhiUeLeftNeutral(0),     fhDeltaPhiUeRightNeutral(0),
 153     fhXEUeLeftNeutral(0),           fhXEUeRightNeutral(0),
 154     fhPtHbpXEUeLeftNeutral(0),      fhPtHbpXEUeRightNeutral(0),
 155     fhZTUeLeftNeutral(0),           fhZTUeRightNeutral(0),
 156     fhPtHbpZTUeLeftNeutral(0),      fhPtHbpZTUeRightNeutral(0),
 157     fhPtPi0DecayRatio(0),
 158     fhDeltaPhiDecayCharged(0),      fhXEDecayCharged(0), fhZTDecayCharged(0),
 159     fhDeltaPhiDecayNeutral(0),      fhXEDecayNeutral(0), fhZTDecayNeutral(0),
 160     fhDeltaPhiDecayChargedAssocPtBin(0),
 161     fhXEDecayChargedAssocPtBin(0),  fhZTDecayChargedAssocPtBin(0),
 162     fh2phiLeadingParticle(0x0),     fhMCPtLeading(0),
 163     fhMCPhiLeading(0),              fhMCEtaLeading(0),
 164     fhMCEtaCharged(0),              fhMCPhiCharged(0),
 165     fhMCDeltaEtaCharged(0),         fhMCDeltaPhiCharged(0x0),
 166     fhMCDeltaPhiDeltaEtaCharged(0), fhMCDeltaPhiChargedPt(0),
 167     fhMCPtXECharged(0),             fhMCPtXEUeCharged(0),
 168     fhMCPtXEUeLeftCharged(0),       fhMCPtXEUeRightCharged(0),
 169     fhMCPtHbpXECharged(0),          fhMCPtHbpXEUeCharged(0),
 170     fhMCPtHbpXEUeLeftCharged(0),    fhMCPtHbpXEUeRightCharged(0),
 171     fhMCUePart(0),
 172     fhMCPtZTCharged(0),             fhMCPtZTUeCharged(0),
 173     fhMCPtZTUeLeftCharged(0),       fhMCPtZTUeRightCharged(0),
 174     fhMCPtHbpZTCharged(0),          fhMCPtHbpZTUeCharged(0),
 175     fhMCPtHbpZTUeLeftCharged(0),    fhMCPtHbpZTUeRightCharged(0),
 176     fhMCPtTrigPout(0),              fhMCPtAssocDeltaPhi(0),
 177     //Mixing
 178     fhNEventsTrigger(0),            fhNtracksMB(0),                 fhNclustersMB(0),
 179     fhMixDeltaPhiCharged(0),        fhMixDeltaPhiDeltaEtaCharged(0),
 180     fhMixXECharged(0),              fhMixXEUeCharged(0),            fhMixHbpXECharged(0),
 181     fhMixDeltaPhiChargedAssocPtBin(),
 182     fhMixDeltaPhiChargedAssocPtBinDEta08(),
 183     fhMixDeltaPhiChargedAssocPtBinDEta0(),
 184     fhMixDeltaPhiDeltaEtaChargedAssocPtBin(),
 185     fhEventBin(0),                  fhEventMixBin(0)
 186 {
 187   //Default Ctor
 188
 189   //Initialize parameters
 190   InitParameters();
 191
 192   for(Int_t i = 0; i < 7; i++)
 193   {
 194     fhPtLeadingMC[i] = 0;
 195     fhXEChargedMC[i] = 0;
 196     fhDeltaPhiChargedMC[i] = 0;
 197   }
 198
 199   for(Int_t i = 0; i < 7; i++)
 200   {
 201     fhPtLeadingPileUp             [i] = 0 ;
 202     fhDeltaPhiChargedPileUp       [i] = 0 ; fhDeltaEtaChargedPileUp       [i] = 0 ;
 203     fhXEChargedPileUp             [i] = 0 ; fhXEUeChargedPileUp           [i] = 0 ;
 204     fhZTChargedPileUp             [i] = 0 ; fhZTUeChargedPileUp           [i] = 0 ;
 205     fhPtTrigChargedPileUp         [i] = 0 ;
 206     fhDeltaPhiChargedPtA3GeVPileUp[i] = 0 ; fhDeltaEtaChargedPtA3GeVPileUp[i] = 0 ;
 207   }
 208
 209 }
 210
 211 //_________________________________________________________________
 212 AliAnaParticleHadronCorrelation::~AliAnaParticleHadronCorrelation()
 213 {
 214   // Remove event containers
 215
 216   if(DoOwnMix())
 217   {
 218     if(fListMixTrackEvents)
 219     {
 220       for(Int_t iz=0; iz < GetNZvertBin(); iz++)
 221       {
 222         for(Int_t ic=0; ic < GetNCentrBin(); ic++)
 223         {
 224           for(Int_t irp=0; irp<GetNRPBin(); irp++)
 225           {
 226             Int_t bin = GetEventMixBin(ic, iz, irp);
 227             fListMixTrackEvents[bin]->Delete() ;
 228             delete fListMixTrackEvents[bin] ;
 229           }
 230         }
 231       }
 232     }
 233
 234     delete[] fListMixTrackEvents;
 235
 236     if(fListMixCaloEvents)
 237     {
 238       for(Int_t iz=0; iz < GetNZvertBin(); iz++)
 239       {
 240         for(Int_t ic=0; ic < GetNCentrBin(); ic++)
 241         {
 242           for(Int_t irp=0; irp<GetNRPBin(); irp++)
 243           {
 244             Int_t bin = GetEventMixBin(ic, iz, irp);
 245             fListMixCaloEvents[bin]->Delete() ;
 246             delete fListMixCaloEvents[bin] ;
 247           }
 248         }
 249       }
 250     }
 251
 252     delete[] fListMixCaloEvents;
 253
 254   }
 255 }
 256
 257 //______________________________________________________________________________________________________________________________________________________
 258 void AliAnaParticleHadronCorrelation::FillChargedAngularCorrelationHistograms(Float_t ptAssoc,  Float_t ptTrig,      Int_t   bin,
 259                                                                               Float_t phiAssoc, Float_t phiTrig,     Float_t &     deltaPhi,
 260                                                                               Float_t etaAssoc, Float_t etaTrig,
 261                                                                               Bool_t  decay,    Float_t hmpidSignal, Int_t  outTOF,
 262                                                                               Int_t nTracks,       Int_t   mcTag)
 263 {
 264   // Fill angular correlation related histograms
 265
 266   Float_t deltaEta    = etaTrig-etaAssoc;
 267           deltaPhi    = phiTrig-phiAssoc;
 268   Float_t deltaPhiOrg = deltaPhi;
 269
 270   if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
 271   if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
 272
 273   fhEtaCharged       ->Fill(ptAssoc,etaAssoc);
 274   fhPhiCharged       ->Fill(ptAssoc,phiAssoc);
 275   fhDeltaEtaCharged  ->Fill(ptTrig ,deltaEta);
 276   fhDeltaPhiCharged  ->Fill(ptTrig ,deltaPhi);
 277   fhDeltaPhiChargedPt->Fill(ptAssoc, deltaPhi);
 278   fhDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta);
 279
 280   if(ptAssoc > 3 )
 281   {
 282     fhDeltaEtaChargedPtA3GeV        ->Fill(ptTrig  ,deltaEta);
 283     fhDeltaPhiChargedPtA3GeV        ->Fill(ptTrig  ,deltaPhi);
 284     fhDeltaPhiDeltaEtaChargedPtA3GeV->Fill(deltaPhi, deltaEta);
 285   }
 286
 287   // Pile up studies
 288
 289   if(fFillPileUpHistograms)
 290   {
 291     if     (outTOF==1)
 292     {
 293       fhDeltaPhiChargedOtherBC->Fill(ptTrig ,deltaPhi) ;
 294       if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVOtherBC->Fill(ptTrig ,deltaPhi) ;
 295     }
 296     else if(outTOF==0)
 297     {
 298         fhDeltaPhiChargedBC0->Fill(ptTrig ,deltaPhi) ;
 299         if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVBC0->Fill(ptTrig ,deltaPhi) ;
 300     }
 301
 302     Int_t vtxBC = GetReader()->GetVertexBC();
 303     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
 304     {
 305       fhDeltaPhiChargedVtxBC0->Fill(ptTrig ,deltaPhi) ;
 306       if(ptAssoc > 3 ) fhDeltaPhiChargedPtA3GeVVtxBC0->Fill(ptTrig ,deltaPhi) ;
 307     }
 308
 309     if(GetReader()->IsPileUpFromSPD())               { fhDeltaEtaChargedPileUp[0]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[0]->Fill(ptTrig ,deltaPhi) ; }
 310     if(GetReader()->IsPileUpFromEMCal())             { fhDeltaEtaChargedPileUp[1]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[1]->Fill(ptTrig ,deltaPhi) ; }
 311     if(GetReader()->IsPileUpFromSPDOrEMCal())        { fhDeltaEtaChargedPileUp[2]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[2]->Fill(ptTrig ,deltaPhi) ; }
 312     if(GetReader()->IsPileUpFromSPDAndEMCal())       { fhDeltaEtaChargedPileUp[3]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[3]->Fill(ptTrig ,deltaPhi) ; }
 313     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    { fhDeltaEtaChargedPileUp[4]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[4]->Fill(ptTrig ,deltaPhi) ; }
 314     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    { fhDeltaEtaChargedPileUp[5]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[5]->Fill(ptTrig ,deltaPhi) ; }
 315     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhDeltaEtaChargedPileUp[6]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPileUp[6]->Fill(ptTrig ,deltaPhi) ; }
 316
 317     if(ptAssoc > 3 )
 318     {
 319       if(GetReader()->IsPileUpFromSPD())               { fhDeltaEtaChargedPtA3GeVPileUp[0]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[0]->Fill(ptTrig ,deltaPhi) ; }
 320       if(GetReader()->IsPileUpFromEMCal())             { fhDeltaEtaChargedPtA3GeVPileUp[1]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[1]->Fill(ptTrig ,deltaPhi) ; }
 321       if(GetReader()->IsPileUpFromSPDOrEMCal())        { fhDeltaEtaChargedPtA3GeVPileUp[2]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[2]->Fill(ptTrig ,deltaPhi) ; }
 322       if(GetReader()->IsPileUpFromSPDAndEMCal())       { fhDeltaEtaChargedPtA3GeVPileUp[3]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[3]->Fill(ptTrig ,deltaPhi) ; }
 323       if(GetReader()->IsPileUpFromSPDAndNotEMCal())    { fhDeltaEtaChargedPtA3GeVPileUp[4]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[4]->Fill(ptTrig ,deltaPhi) ; }
 324       if(GetReader()->IsPileUpFromEMCalAndNotSPD())    { fhDeltaEtaChargedPtA3GeVPileUp[5]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[5]->Fill(ptTrig ,deltaPhi) ; }
 325       if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhDeltaEtaChargedPtA3GeVPileUp[6]->Fill(ptTrig ,deltaEta) ; fhDeltaPhiChargedPtA3GeVPileUp[6]->Fill(ptTrig ,deltaPhi) ; }
 326     }
 327   }
 328
 329   if(IsDataMC())
 330   {
 331     Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
 332     fhDeltaPhiChargedMC[mcIndex]->Fill(ptTrig , deltaPhi);
 333   }
 334
 335   if(fDecayTrigger && decay) fhDeltaPhiDecayCharged   ->Fill(ptTrig  , deltaPhi);
 336
 337   Double_t  dphiBrad = -100;
 338   if(fFillBradHisto)
 339   {
 340     dphiBrad = atan2(sin(deltaPhiOrg), cos(deltaPhiOrg))/TMath::Pi();//-1 to 1
 341     if(TMath::Abs(dphiBrad)>0.325 && TMath::Abs(dphiBrad)<0.475)  //Hardcoded values, BAD, FIXME
 342     {
 343       fhAssocPtBkg->Fill(ptTrig, ptAssoc);
 344     }
 345
 346     if(dphiBrad<-1./3) dphiBrad += 2;
 347     fhDeltaPhiBrad->Fill(ptTrig, dphiBrad);
 348   }
 349
 350   // Fill histograms in bins of associated particle pT
 351   if(bin>=0)
 352   {
 353       fhDeltaPhiDeltaEtaAssocPtBin    [bin]->Fill(deltaPhi,deltaEta);
 354
 355       fhDeltaPhiAssocPtBin            [bin]->Fill(ptTrig, deltaPhi);
 356
 357     if(TMath::Abs(deltaEta)> 0.8)
 358       fhDeltaPhiAssocPtBinDEta08      [bin]->Fill(ptTrig, deltaPhi);
 359
 360     if(TMath::Abs(deltaEta)< 0.01)
 361       fhDeltaPhiAssocPtBinDEta0       [bin]->Fill(ptTrig, deltaPhi);
 362
 363     if (fFillBradHisto)
 364       fhDeltaPhiBradAssocPtBin        [bin]->Fill(ptTrig, dphiBrad);
 365
 366     if(fDecayTrigger && decay)
 367       fhDeltaPhiDecayChargedAssocPtBin[bin]->Fill(ptTrig, deltaPhi);
 368
 369     if(fHMPIDCorrelation)
 370     {
 371       if( hmpidSignal > 0 )
 372       {
 373         //printf("Track pt %f with HMPID signal %f \n",pt,hmpidSignal);
 374         fhDeltaPhiAssocPtBinHMPID[bin]->Fill(ptTrig, deltaPhi);
 375       }
 376
 377       if(phiAssoc > 5*TMath::DegToRad() && phiAssoc < 20*TMath::DegToRad())
 378       {
 379         //printf("Track pt %f in HMPID acceptance phi %f \n ",pt,phi*TMath::RadToDeg() );
 380         fhDeltaPhiAssocPtBinHMPIDAcc[bin]->Fill(ptTrig, deltaPhi);
 381       }
 382     }
 383   }
 384
 385   //fill different multiplicity histogram
 386   if(DoEventSelect())
 387   {
 388     for(Int_t im = 0; im<GetMultiBin(); im++)
 389     {
 390       if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
 391       {
 392         fhTrigDeltaPhiCharged[im]->Fill(ptTrig,deltaPhi);
 393         fhTrigDeltaEtaCharged[im]->Fill(ptTrig,deltaEta);
 394       }
 395     }
 396   }
 397 }
 398
 399 //____________________________________________________________________________________________________________________________________________________
 400 Bool_t AliAnaParticleHadronCorrelation::FillChargedMCCorrelationHistograms(Float_t mcAssocPt,      Float_t mcAssocPhi, Float_t mcAssocEta,
 401                                                                            Float_t mcTrigPt, Float_t mcTrigPhi,  Float_t mcTrigEta)
 402 {
 403   // Fill MC histograms independently of AOD or ESD
 404
 405   //Select only hadrons in pt range
 406   if( mcAssocPt < fMinAssocPt || mcAssocPt > fMaxAssocPt ) return kTRUE ; // exclude but continue
 407
 408   if( mcAssocPhi < 0 ) mcAssocPhi+=TMath::TwoPi();
 409
 410   //remove trigger itself for correlation when use charged triggers
 411   if(TMath::Abs(mcAssocPt -mcTrigPt ) < 1e-6 &&
 412      TMath::Abs(mcAssocPhi-mcTrigPhi) < 1e-6 &&
 413      TMath::Abs(mcAssocEta-mcTrigEta) < 1e-6)            return kTRUE ; // exclude but continue
 414
 415   // Absolute leading?
 416   if( fMakeAbsoluteLeading && mcAssocPt > mcTrigPt )     return kFALSE; // skip event
 417
 418   // Skip this event if near side associated particle pt larger than trigger
 419   if( fMakeNearSideLeading && mcAssocPt > mcTrigPt &&
 420      TMath::Abs(mcAssocPhi-mcTrigPhi)<TMath::PiOver2() ) return kFALSE; // skip event
 421
 422   Float_t mcdeltaPhi= mcTrigPhi-mcAssocPhi;
 423   if(mcdeltaPhi <= -TMath::PiOver2()) mcdeltaPhi+=TMath::TwoPi();
 424   if(mcdeltaPhi > 3*TMath::PiOver2()) mcdeltaPhi-=TMath::TwoPi();
 425
 426   Float_t mcxE    =-mcAssocPt/mcTrigPt*TMath::Cos(mcdeltaPhi);// -(mcAssocPx*pxprim+mcAssocPy*pyprim)/(mcTrigPt*mcTrigPt);
 427   Float_t mchbpXE =-100 ;
 428   if(mcxE > 0 ) mchbpXE = TMath::Log(1./mcxE);
 429
 430   Float_t mczT    = mcAssocPt/mcTrigPt ;
 431   Float_t mchbpZT =-100 ;
 432   if(mczT > 0 ) mchbpZT = TMath::Log(1./mczT);
 433
 434   //Selection within angular range
 435   if( mcdeltaPhi< -TMath::PiOver2())  mcdeltaPhi+=TMath::TwoPi();
 436   if( mcdeltaPhi>3*TMath::PiOver2())  mcdeltaPhi-=TMath::TwoPi();
 437
 438   Double_t mcpout = mcAssocPt*TMath::Sin(mcdeltaPhi) ;
 439
 440   if(GetDebug() > 0 )
 441   {
 442     AliInfo(Form("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",
 443                  mcAssocPt,mcAssocPhi, mcTrigPhi,fDeltaPhiMinCut, mcdeltaPhi, fDeltaPhiMaxCut, GetMinPt()));
 444   }
 445
 446   // Fill Histograms
 447   fhMCEtaCharged     ->Fill(mcAssocPt, mcAssocEta);
 448   fhMCPhiCharged     ->Fill(mcAssocPt, mcAssocPhi);
 449   fhMCDeltaEtaCharged->Fill(mcTrigPt , mcTrigEta-mcAssocEta);
 450   fhMCDeltaPhiCharged->Fill(mcTrigPt , mcdeltaPhi);
 451   fhMCPtAssocDeltaPhi->Fill(mcAssocPt, mcdeltaPhi);
 452
 453   fhMCDeltaPhiDeltaEtaCharged->Fill(mcdeltaPhi,mcTrigEta-mcAssocEta);
 454
 455   //delta phi cut for correlation
 456   if( (mcdeltaPhi > fDeltaPhiMinCut) && (mcdeltaPhi < fDeltaPhiMaxCut) )
 457   {
 458     fhMCDeltaPhiChargedPt->Fill(mcAssocPt,mcdeltaPhi);
 459     fhMCPtXECharged      ->Fill(mcTrigPt, mcxE);
 460     fhMCPtHbpXECharged   ->Fill(mcTrigPt, mchbpXE);
 461     fhMCPtZTCharged      ->Fill(mcTrigPt, mczT);
 462     fhMCPtHbpZTCharged   ->Fill(mcTrigPt, mchbpZT);
 463     fhMCPtTrigPout       ->Fill(mcTrigPt, mcpout) ;
 464   }
 465
 466   //underlying event
 467
 468   if ( (mcdeltaPhi > fUeDeltaPhiMinCut) && (mcdeltaPhi < fUeDeltaPhiMaxCut) )
 469   {
 470     //Double_t randomphi = gRandom->Uniform(TMath::Pi()/2,3*TMath::Pi()/2);
 471     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 472     Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
 473     Double_t mcUezT =   mcAssocPt/mcTrigPt;
 474
 475     if(mcUexE < 0.) mcUexE = -mcUexE;
 476
 477     fhMCPtXEUeCharged->Fill(mcTrigPt,mcUexE);
 478     if(mcUexE > 0) fhMCPtHbpXEUeCharged->Fill(mcTrigPt,TMath::Log(1/mcUexE));
 479
 480     fhMCPtZTUeCharged->Fill(mcTrigPt,mcUezT);
 481     if(mcUezT > 0) fhMCPtHbpZTUeCharged->Fill(mcTrigPt,TMath::Log(1/mcUezT));
 482
 483     fhMCUePart->Fill(mcTrigPt);
 484   }
 485
 486   //left
 487   if((mcdeltaPhi<-fUeDeltaPhiMinCut) || (mcdeltaPhi >2*fUeDeltaPhiMaxCut))
 488   {
 489     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 490     Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
 491     Double_t mcUezT =   mcAssocPt/mcTrigPt;
 492
 493     if(mcUexE < 0.) mcUexE = -mcUexE;
 494
 495     fhMCPtXEUeLeftCharged->Fill(mcTrigPt,mcUexE);
 496     if(mcUexE > 0) fhMCPtHbpXEUeLeftCharged->Fill(mcTrigPt,TMath::Log(1/mcUexE));
 497
 498     fhMCPtZTUeLeftCharged->Fill(mcTrigPt,mcUezT);
 499     if(mcUexE > 0) fhMCPtHbpZTUeLeftCharged->Fill(mcTrigPt,TMath::Log(1/mcUezT));
 500
 501   }
 502
 503   //right
 504   if((mcdeltaPhi > fUeDeltaPhiMinCut) && (mcdeltaPhi < fUeDeltaPhiMaxCut))
 505   {
 506     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 507     Double_t mcUexE = -(mcAssocPt/mcTrigPt)*TMath::Cos(randomphi);
 508     Double_t mcUezT =   mcAssocPt/mcTrigPt;
 509
 510     if(mcUexE < 0.) mcUexE = -mcUexE;
 511
 512     fhMCPtXEUeRightCharged->Fill(mcTrigPt,mcUexE);
 513     if(mcUexE > 0) fhMCPtHbpXEUeRightCharged->Fill(mcTrigPt,TMath::Log(1/mcUexE));
 514
 515     fhMCPtZTUeRightCharged->Fill(mcTrigPt,mcUezT);
 516     if(mcUexE > 0) fhMCPtHbpZTUeRightCharged->Fill(mcTrigPt,TMath::Log(1/mcUezT));
 517   }
 518
 519   return kTRUE;
 520 }
 521
 522 //___________________________________________________________________________________________________________________
 523 void AliAnaParticleHadronCorrelation::FillChargedMomentumImbalanceHistograms(Float_t ptTrig,   Float_t ptAssoc,
 524                                                                              Float_t xE,       Float_t hbpXE,
 525                                                                              Float_t zT,       Float_t hbpZT,
 526                                                                              Float_t pout,     Float_t deltaPhi,
 527                                                                              Int_t   nTracks,  Int_t   charge,
 528                                                                              Int_t   bin,      Bool_t  decay,
 529                                                                              Int_t   outTOF,   Int_t   mcTag)
 530
 531 {
 532   // Fill mostly momentum imbalance related histograms
 533
 534   fhXECharged         ->Fill(ptTrig , xE);
 535   fhPtHbpXECharged    ->Fill(ptTrig , hbpXE);
 536   fhZTCharged         ->Fill(ptTrig , zT);
 537   fhPtHbpZTCharged    ->Fill(ptTrig , hbpZT);
 538   fhPtTrigPout        ->Fill(ptTrig , pout) ;
 539   fhPtTrigCharged     ->Fill(ptTrig , ptAssoc) ;
 540   if((deltaPhi > 5*TMath::Pi()/6.)   && (deltaPhi < 7*TMath::Pi()/6.)) {
 541     fhXECharged_Cone2         ->Fill(ptTrig , xE);
 542     fhPtHbpXECharged_Cone2    ->Fill(ptTrig , hbpXE);
 543   }
 544
 545   // Pile up studies
 546   if(fFillPileUpHistograms)
 547   {
 548     if     (outTOF==1)
 549     {
 550       fhXEChargedOtherBC    ->Fill(ptTrig,xE);
 551       fhZTChargedOtherBC    ->Fill(ptTrig,zT);
 552       fhPtTrigChargedOtherBC->Fill(ptTrig,ptAssoc);
 553     }
 554     else if(outTOF==0)
 555     {
 556       fhXEChargedBC0    ->Fill(ptTrig,xE);
 557       fhZTChargedBC0    ->Fill(ptTrig,zT);
 558       fhPtTrigChargedBC0->Fill(ptTrig,ptAssoc);
 559     }
 560
 561     Int_t vtxBC = GetReader()->GetVertexBC();
 562     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
 563     {
 564       fhXEChargedVtxBC0    ->Fill(ptTrig,xE);
 565       fhZTChargedVtxBC0    ->Fill(ptTrig,zT);
 566       fhPtTrigChargedVtxBC0->Fill(ptTrig,ptAssoc);
 567     }
 568
 569     if(GetReader()->IsPileUpFromSPD())                { fhXEChargedPileUp[0]->Fill(ptTrig,xE); fhZTChargedPileUp[0]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[0]->Fill(ptTrig,ptAssoc); }
 570     if(GetReader()->IsPileUpFromEMCal())              { fhXEChargedPileUp[1]->Fill(ptTrig,xE); fhZTChargedPileUp[1]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[1]->Fill(ptTrig,ptAssoc); }
 571     if(GetReader()->IsPileUpFromSPDOrEMCal())         { fhXEChargedPileUp[2]->Fill(ptTrig,xE); fhZTChargedPileUp[2]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[2]->Fill(ptTrig,ptAssoc); }
 572     if(GetReader()->IsPileUpFromSPDAndEMCal())        { fhXEChargedPileUp[3]->Fill(ptTrig,xE); fhZTChargedPileUp[3]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[3]->Fill(ptTrig,ptAssoc); }
 573     if(GetReader()->IsPileUpFromSPDAndNotEMCal())     { fhXEChargedPileUp[4]->Fill(ptTrig,xE); fhZTChargedPileUp[4]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[4]->Fill(ptTrig,ptAssoc); }
 574     if(GetReader()->IsPileUpFromEMCalAndNotSPD())     { fhXEChargedPileUp[5]->Fill(ptTrig,xE); fhZTChargedPileUp[5]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[5]->Fill(ptTrig,ptAssoc); }
 575     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())  { fhXEChargedPileUp[6]->Fill(ptTrig,xE); fhZTChargedPileUp[6]->Fill(ptTrig,zT); fhPtTrigChargedPileUp[6]->Fill(ptTrig,ptAssoc); }
 576   }
 577
 578   if(IsDataMC())
 579   {
 580     Int_t mcIndex = GetMCTagHistogramIndex(mcTag);
 581     fhXEChargedMC      [mcIndex]->Fill(ptTrig , xE      );
 582   }
 583
 584   if(fDecayTrigger && decay)
 585   {
 586     fhXEDecayCharged->Fill(ptTrig,xE);
 587     fhZTDecayCharged->Fill(ptTrig,zT);
 588   } // photon decay pi0/eta trigger
 589
 590   if(bin >= 0 )//away side
 591   {
 592     fhXEAssocPtBin[bin]->Fill(ptTrig, xE) ;
 593     fhZTAssocPtBin[bin]->Fill(ptTrig, zT) ;
 594
 595     if(fDecayTrigger && decay)
 596     {
 597       fhXEDecayChargedAssocPtBin[bin]->Fill(ptTrig, xE);
 598       fhZTDecayChargedAssocPtBin[bin]->Fill(ptTrig, zT);
 599     }
 600   }
 601
 602   if(charge > 0)
 603   {
 604     fhXEPosCharged->Fill(ptTrig,xE) ;
 605     fhZTPosCharged->Fill(ptTrig,zT) ;
 606   }
 607   else
 608   {
 609     fhXENegCharged->Fill(ptTrig,xE) ;
 610     fhZTNegCharged->Fill(ptTrig,zT) ;
 611   }
 612
 613   //fill different multiplicity histogram
 614   if(DoEventSelect())
 615   {
 616     for(Int_t im=0; im<GetMultiBin(); im++)
 617     {
 618       if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
 619       {
 620         fhTrigXECorr[im]->Fill(ptTrig,xE);
 621         fhTrigZTCorr[im]->Fill(ptTrig,zT);
 622       }
 623     }
 624   } //multiplicity events selection
 625 }
 626
 627 //_______________________________________________________________________________________________________________________
 628 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventHistograms(Float_t ptTrig,   Float_t ptAssoc,
 629                                                                            Float_t deltaPhi, Int_t nTracks, Int_t outTOF)
 630 {
 631   // Fill underlying event histograms
 632
 633   fhDeltaPhiUeChargedPt->Fill(ptAssoc,deltaPhi);
 634
 635   Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 636   Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 637   Double_t uezT =   ptAssoc/ptTrig;
 638
 639   if(uexE < 0.) uexE = -uexE;
 640
 641   fhXEUeCharged->Fill(ptTrig,uexE);
 642   if(uexE > 0) fhPtHbpXEUeCharged->Fill(ptTrig,TMath::Log(1/uexE));
 643
 644   fhZTUeCharged->Fill(ptTrig,uezT);
 645   if(uezT > 0) fhPtHbpZTUeCharged->Fill(ptTrig,TMath::Log(1/uezT));
 646
 647   // Pile up studies
 648
 649   if(fFillPileUpHistograms)
 650   {
 651     if     (outTOF==1)
 652     {
 653       fhXEUeChargedOtherBC->Fill(ptTrig,uexE);
 654       fhZTUeChargedOtherBC->Fill(ptTrig,uezT);
 655     }
 656     else if(outTOF==0)
 657     {
 658       fhXEUeChargedBC0->Fill(ptTrig,uexE);
 659       fhZTUeChargedBC0->Fill(ptTrig,uezT);
 660     }
 661
 662     Int_t vtxBC = GetReader()->GetVertexBC();
 663     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)
 664     {
 665       fhXEUeChargedVtxBC0->Fill(ptTrig,uexE);
 666       fhZTUeChargedVtxBC0->Fill(ptTrig,uezT);
 667     }
 668
 669     if(GetReader()->IsPileUpFromSPD())               { fhXEUeChargedPileUp[0]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[0]->Fill(ptTrig,uezT);}
 670     if(GetReader()->IsPileUpFromEMCal())             { fhXEUeChargedPileUp[1]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[1]->Fill(ptTrig,uezT);}
 671     if(GetReader()->IsPileUpFromSPDOrEMCal())        { fhXEUeChargedPileUp[2]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[2]->Fill(ptTrig,uezT);}
 672     if(GetReader()->IsPileUpFromSPDAndEMCal())       { fhXEUeChargedPileUp[3]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[3]->Fill(ptTrig,uezT);}
 673     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    { fhXEUeChargedPileUp[4]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[4]->Fill(ptTrig,uezT);}
 674     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    { fhXEUeChargedPileUp[5]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[5]->Fill(ptTrig,uezT);}
 675     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) { fhXEUeChargedPileUp[6]->Fill(ptTrig,uexE); fhZTUeChargedPileUp[6]->Fill(ptTrig,uezT);}
 676   }
 677
 678   if(DoEventSelect())
 679   {
 680     for(Int_t im=0; im<GetMultiBin(); im++)
 681     {
 682       if(nTracks < ( GetMaxMulti() - GetMinMulti() )/GetMultiBin()*(im+1))
 683       {
 684         fhTrigXEUeCorr[im]->Fill(ptTrig,uexE); // xE? CHECK
 685         fhTrigZTUeCorr[im]->Fill(ptTrig,uezT); // zT? CHECK
 686       }
 687     }
 688   } //multiplicity events selection
 689 }
 690
 691 //_____________________________________________________________________________________________________
 692 void AliAnaParticleHadronCorrelation::FillChargedUnderlyingEventSidesHistograms(Float_t ptTrig,
 693                                                                                 Float_t ptAssoc,
 694                                                                                 Float_t deltaPhi)
 695 {
 696  // Fill underlying event histograms to the left and right of trigger
 697   if((deltaPhi<-fUeDeltaPhiMinCut) || (deltaPhi >2*fUeDeltaPhiMaxCut))
 698   {
 699     fhDeltaPhiUeLeftCharged->Fill(ptAssoc,deltaPhi);
 700     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 701     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 702     Double_t uezT =   ptAssoc/ptTrig;
 703
 704     if(uexE < 0.) uexE = -uexE;
 705
 706     fhXEUeLeftCharged->Fill(ptTrig,uexE);
 707     if(uexE > 0) fhPtHbpXEUeLeftCharged->Fill(ptTrig,TMath::Log(1/uexE));
 708
 709     fhZTUeLeftCharged->Fill(ptTrig,uezT);
 710     if(uexE > 0) fhPtHbpZTUeLeftCharged->Fill(ptTrig,TMath::Log(1/uezT));
 711     fhDeltaPhiUeLeftCharged->Fill(ptAssoc, deltaPhi);
 712   }
 713
 714   if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut))
 715   {
 716     fhDeltaPhiUeRightCharged->Fill(ptAssoc,deltaPhi);
 717     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 718     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 719     Double_t uezT =   ptAssoc/ptTrig;
 720
 721     if(uexE < 0.) uexE = -uexE;
 722
 723     fhXEUeRightCharged->Fill(ptTrig,uexE);
 724     if(uexE > 0) fhPtHbpXEUeRightCharged->Fill(ptTrig,TMath::Log(1/uexE));
 725
 726     fhZTUeRightCharged->Fill(ptTrig,uezT);
 727     if(uexE > 0) fhPtHbpZTUeRightCharged->Fill(ptTrig,TMath::Log(1/uezT));
 728     fhDeltaPhiUeRightCharged->Fill(ptAssoc, deltaPhi);
 729   }
 730
 731   if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-TMath::Pi()/2))
 732   {
 733     fhDeltaPhiUeLeftDownCharged->Fill(ptAssoc,deltaPhi);
 734     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 735     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 736
 737     if(uexE < 0.) uexE = -uexE;
 738
 739     fhXEUeLeftDownCharged->Fill(ptTrig,uexE);
 740   }
 741
 742   if((deltaPhi>2*fUeDeltaPhiMaxCut) && (deltaPhi <3*TMath::Pi()/2))
 743   {
 744     fhDeltaPhiUeLeftUpCharged->Fill(ptAssoc,deltaPhi);
 745     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 746     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 747
 748     if(uexE < 0.) uexE = -uexE;
 749
 750     fhXEUeLeftUpCharged->Fill(ptTrig,uexE);
 751   }
 752
 753   if((deltaPhi > TMath::Pi()/2) && (deltaPhi < fUeDeltaPhiMaxCut))
 754   {
 755     fhDeltaPhiUeRightUpCharged->Fill(ptAssoc,deltaPhi);
 756     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 757     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 758
 759     if(uexE < 0.) uexE = -uexE;
 760
 761     fhXEUeRightUpCharged->Fill(ptTrig,uexE);
 762   }
 763
 764   if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < TMath::Pi()/2))
 765   {
 766     fhDeltaPhiUeRightDownCharged->Fill(ptAssoc,deltaPhi);
 767     Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
 768     Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
 769
 770     if(uexE < 0.) uexE = -uexE;
 771
 772     fhXEUeRightDownCharged->Fill(ptTrig,uexE);
 773   }
 774 }
 775
 776 //______________________________________________________________________________________________________________________________
 777 void AliAnaParticleHadronCorrelation::FillDecayPhotonCorrelationHistograms(Float_t ptAssoc,     Float_t phiAssoc,
 778                                                                            TLorentzVector mom1, TLorentzVector mom2,
 779                                                                            Bool_t bChargedOrNeutral)
 780 {
 781   // Do correlation with decay photons of triggered pi0 or eta
 782
 783   // Calculate the correlation parameters
 784   Float_t ptDecay1 = mom1.Pt();
 785   Float_t ptDecay2 = mom2.Pt();
 786
 787   Float_t zTDecay1 = -100, zTDecay2 = -100;
 788   if(ptDecay1) zTDecay1 = ptAssoc/ptDecay1 ;
 789   if(ptDecay2) zTDecay2 = ptAssoc/ptDecay2 ;
 790
 791   Float_t deltaPhiDecay1 = mom1.Phi()-phiAssoc;
 792   if(deltaPhiDecay1< -TMath::PiOver2()) deltaPhiDecay1+=TMath::TwoPi();
 793   if(deltaPhiDecay1>3*TMath::PiOver2()) deltaPhiDecay1-=TMath::TwoPi();
 794
 795   Float_t deltaPhiDecay2 = mom2.Phi()-phiAssoc;
 796   if(deltaPhiDecay2< -TMath::PiOver2()) deltaPhiDecay2+=TMath::TwoPi();
 797   if(deltaPhiDecay2>3*TMath::PiOver2()) deltaPhiDecay2-=TMath::TwoPi();
 798
 799   Float_t xEDecay1   =-zTDecay1*TMath::Cos(deltaPhiDecay1); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
 800   Float_t xEDecay2   =-zTDecay2*TMath::Cos(deltaPhiDecay2); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
 801
 802   if(bChargedOrNeutral) // correlate with charges
 803   {
 804     fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaPhiDecay1);
 805     fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaPhiDecay2);
 806
 807     if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::FillDecayPhotonHistograms( Charged corr) - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
 808
 809     if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
 810     {
 811       fhZTDecayCharged->Fill(ptDecay1,zTDecay1);
 812       fhXEDecayCharged->Fill(ptDecay1,xEDecay1);
 813     }
 814     if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
 815     {
 816       fhZTDecayCharged->Fill(ptDecay2,zTDecay2);
 817       fhXEDecayCharged->Fill(ptDecay2,xEDecay2);
 818     }
 819   }
 820   else // correlate with neutrals
 821   {
 822     fhDeltaPhiDecayCharged->Fill(ptDecay1, deltaPhiDecay1);
 823     fhDeltaPhiDecayCharged->Fill(ptDecay2, deltaPhiDecay2);
 824
 825     if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::FillDecayPhotonHistograms(Neutral corr) - deltaPhoton1 = %f, deltaPhoton2 = %f \n", deltaPhiDecay1, deltaPhiDecay2);
 826
 827     if( (deltaPhiDecay1 > fDeltaPhiMinCut) && ( deltaPhiDecay1 < fDeltaPhiMaxCut) )
 828     {
 829       fhZTDecayCharged->Fill(ptDecay1,zTDecay1);
 830       fhXEDecayCharged->Fill(ptDecay1,xEDecay1);
 831     }
 832     if( (deltaPhiDecay2 > fDeltaPhiMinCut) && ( deltaPhiDecay2 < fDeltaPhiMaxCut) )
 833     {
 834       fhZTDecayCharged->Fill(ptDecay2,zTDecay2);
 835       fhXEDecayCharged->Fill(ptDecay2,xEDecay2);
 836     }
 837   }
 838 }
 839
 840 //______________________________________________________________________________________________________________________________________________________
 841 void AliAnaParticleHadronCorrelation::FillNeutralAngularCorrelationHistograms(Float_t ptAssoc,  Float_t ptTrig,
 842                                                                               Float_t phiAssoc, Float_t phiTrig,  Float_t &     deltaPhi,
 843                                                                               Float_t etaAssoc, Float_t etaTrig)
 844 {
 845   // Fill angular correlation related histograms
 846
 847   Float_t deltaEta    = etaTrig-etaAssoc;
 848   deltaPhi    = phiTrig-phiAssoc;
 849
 850   if(deltaPhi <= -TMath::PiOver2()) deltaPhi+=TMath::TwoPi();
 851   if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
 852
 853   fhEtaNeutral     ->Fill(ptAssoc,etaAssoc);
 854   fhPhiNeutral     ->Fill(ptAssoc,phiAssoc);
 855   fhDeltaEtaNeutral->Fill(ptTrig ,deltaEta);
 856   fhDeltaPhiNeutral->Fill(ptTrig ,deltaPhi);
 857
 858   if(ptAssoc > 2 ) fhDeltaPhiDeltaEtaNeutral->Fill(deltaPhi, deltaEta);
 859
 860 }
 861
 862 //_____________________________________________________________________________________________________________________________
 863 void AliAnaParticleHadronCorrelation::FillNeutralUnderlyingEventSidesHistograms(Float_t ptTrig,   Float_t ptAssoc,
 864                                                                                 Float_t xE,       Float_t hbpXE,
 865                                                                                 Float_t zT,       Float_t hbpZT,
 866                                                                                 Float_t deltaPhi)
 867 {
 868   // Fill underlying event histograms to the left and right of trigger
 869
 870   if((deltaPhi<-fUeDeltaPhiMinCut) && (deltaPhi >-fUeDeltaPhiMaxCut))
 871   {
 872     fhDeltaPhiUeLeftNeutral->Fill(ptAssoc, deltaPhi);
 873     fhXEUeLeftNeutral      ->Fill(ptTrig , xE);
 874     fhPtHbpXEUeLeftNeutral ->Fill(ptTrig , hbpXE);
 875     fhZTUeLeftNeutral      ->Fill(ptTrig , zT);
 876     fhPtHbpZTUeLeftNeutral ->Fill(ptTrig , hbpZT);
 877   }
 878
 879   if((deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut))
 880   {
 881     fhDeltaPhiUeRightNeutral->Fill(ptAssoc, deltaPhi);
 882     fhXEUeRightNeutral      ->Fill(ptTrig , xE);
 883     fhPtHbpXEUeRightNeutral ->Fill(ptTrig , hbpXE);
 884     fhZTUeRightNeutral      ->Fill(ptTrig , zT);
 885     fhPtHbpZTUeRightNeutral ->Fill(ptTrig , hbpZT);
 886   }
 887 }
 888
 889 //______________________________________________________
 890 void AliAnaParticleHadronCorrelation::FillEventMixPool()
 891 {
 892   // Fill the pool with tracks if requested
 893
 894   if(DoOwnMix())
 895   {
 896     FillChargedEventMixPool();
 897
 898     if(OnlyIsolated() || fFillNeutralEventMixPool)
 899       FillNeutralEventMixPool();
 900   }
 901
 902 }
 903
 904 //_____________________________________________________________
 905 void AliAnaParticleHadronCorrelation::FillChargedEventMixPool()
 906 {
 907   // Mixed event pool filling for tracks
 908
 909   if(fUseMixStoredInReader && GetReader()->GetLastTracksMixedEvent() == GetEventNumber())
 910   {
 911     //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
 912     return ; // pool filled previously for another trigger
 913   }
 914
 915   Int_t nTracks = GetCTSTracks()->GetEntriesFast();
 916
 917   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
 918   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
 919
 920   if(!inputHandler) return ;
 921
 922   // Do mixing only with MB event (or the chosen mask), if not skip
 923   if( !(inputHandler->IsEventSelected( ) & GetReader()->GetMixEventTriggerMask()) ) return ;
 924
 925   fhNtracksMB->Fill(nTracks);
 926
 927   Int_t eventBin = GetEventMixBin();
 928
 929   //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
 930   if(eventBin < 0) return;
 931
 932   TObjArray * mixEventTracks = new TObjArray;
 933
 934   if(fUseMixStoredInReader)
 935   {
 936     fListMixTrackEvents[eventBin] = GetReader()->GetListWithMixedEventsForTracks(eventBin);
 937   }
 938
 939   if(!fListMixTrackEvents[eventBin]) fListMixTrackEvents[eventBin] = new TList();
 940
 941   //printf("%s ***** Pool Event bin : %d - nTracks %d\n",GetInputAODName().Data(),eventBin, GetCTSTracks()->GetEntriesFast());
 942
 943   TList * pool = fListMixTrackEvents[eventBin];
 944
 945   TVector3 p3;
 946   for(Int_t ipr = 0;ipr < GetCTSTracks()->GetEntriesFast() ; ipr ++ )
 947   {
 948     AliVTrack * track = (AliVTrack *) (GetCTSTracks()->At(ipr)) ;
 949
 950     Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
 951     p3.SetXYZ(mom[0],mom[1],mom[2]);
 952     Float_t pt   = p3.Pt();
 953
 954     //Select only hadrons in pt range
 955     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
 956
 957     AliAODPWG4Particle * mixedTrack = new AliAODPWG4Particle(mom[0],mom[1],mom[2],0);
 958     mixedTrack->SetDetector("CTS");
 959     mixedTrack->SetChargedBit(track->Charge()>0);
 960     mixEventTracks->Add(mixedTrack);
 961   }
 962
 963   //Set the event number where the last event was added, to avoid double pool filling
 964   GetReader()->SetLastTracksMixedEvent(GetEventNumber());
 965
 966   //printf("Add event to pool with %d tracks \n ",mixEventTracks->GetEntries());
 967   pool->AddFirst(mixEventTracks);
 968   mixEventTracks = 0;
 969
 970   //printf("Pool size %d, max %d\n",pool->GetSize(), GetNMaxEvMix());
 971
 972   if(pool->GetSize() > GetNMaxEvMix())
 973   {//Remove last event
 974     TClonesArray * tmp = static_cast<TClonesArray*>(pool->Last()) ;
 975     pool->RemoveLast() ;
 976     delete tmp ;
 977   }
 978 }
 979
 980 //_____________________________________________________________
 981 void AliAnaParticleHadronCorrelation::FillNeutralEventMixPool()
 982 {
 983   // Mixed event pool filling for neutral clusters
 984   // Right now only for EMCAL and in isolation case
 985
 986   //printf("FillNeutralEventMixPool for %s\n",GetInputAODName().Data());
 987
 988   TObjArray * pl = GetEMCALClusters();
 989   //if (GetAODObjArrayName.Contains("PHOS") )pl    = GetPHOSClusters();
 990   //else                                     pl    = GetEMCALClusters();
 991
 992   Int_t nClusters   = pl->GetEntriesFast();
 993
 994   if(fUseMixStoredInReader && GetReader()->GetLastCaloMixedEvent() == GetEventNumber())
 995   {
 996     //printf("%s : Pool already filled for this event !!!\n",GetInputAODName().Data());
 997     return ; // pool filled previously for another trigger
 998   }
 999
1000   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
1001   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
1002
1003   if(!inputHandler) return ;
1004
1005   // Do mixing only with MB event (or the chosen mask), if not skip
1006   if( !(inputHandler->IsEventSelected( ) & GetReader()->GetMixEventTriggerMask()) ) return ;
1007
1008   fhNclustersMB->Fill(nClusters);
1009
1010   Int_t eventBin = GetEventMixBin();
1011
1012   //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
1013   if(eventBin < 0) return;
1014
1015   TObjArray * mixEventCalo = new TObjArray;
1016
1017   if(fUseMixStoredInReader)
1018   {
1019     fListMixCaloEvents[eventBin] = GetReader()->GetListWithMixedEventsForCalo(eventBin);
1020   }
1021
1022   if(!fListMixCaloEvents[eventBin]) fListMixCaloEvents[eventBin] = new TList();
1023
1024   TList * poolCalo = fListMixCaloEvents[eventBin];
1025
1026   TLorentzVector mom;
1027
1028   for(Int_t ipr = 0;ipr <  nClusters ; ipr ++ )
1029   {
1030     AliVCluster * calo = (AliVCluster *) (pl->At(ipr)) ;
1031
1032     // remove matched clusters
1033     if( IsTrackMatched( calo, GetReader()->GetInputEvent() ) ) continue ;
1034
1035     //Cluster momentum calculation
1036     if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
1037     {
1038       calo->GetMomentum(mom,GetVertex(0)) ;
1039     }//Assume that come from vertex in straight line
1040     else
1041     {
1042       Double_t vertex[]={0,0,0};
1043       calo->GetMomentum(mom,vertex) ;
1044     }
1045
1046     Float_t pt = mom.Pt();
1047     //Select only clusters in pt range
1048     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
1049
1050     AliAODPWG4Particle * mixedCalo = new AliAODPWG4Particle(mom);
1051     mixedCalo->SetDetector("EMCAL");
1052     mixEventCalo->Add(mixedCalo);
1053   }
1054
1055   //Set the event number where the last event was added, to avoid double pool filling
1056   GetReader()->SetLastCaloMixedEvent(GetEventNumber());
1057
1058   //printf("Add event to pool with %d clusters \n ",mixEventCalo->GetEntries());
1059   poolCalo->AddFirst(mixEventCalo);
1060   mixEventCalo = 0;
1061
1062   //printf("Pool size %d, max %d\n",poolCalo->GetSize(), GetNMaxEvMix());
1063
1064   if(poolCalo->GetSize() > GetNMaxEvMix())
1065   {//Remove last event
1066     TClonesArray * tmp = static_cast<TClonesArray*>(poolCalo->Last()) ;
1067     poolCalo->RemoveLast() ;
1068     delete tmp ;
1069   }
1070 }
1071
1072 //____________________________________________________________
1073 TObjString* AliAnaParticleHadronCorrelation::GetAnalysisCuts()
1074 {
1075   //Save parameters used for analysis
1076   TString parList ; //this will be list of parameters used for this analysis.
1077   const Int_t buffersize = 560;
1078   char onePar[buffersize] ;
1079
1080   snprintf(onePar,buffersize,"--- AliAnaPaticleHadronCorrelation ---\n") ;
1081   parList+=onePar ;
1082   snprintf(onePar,buffersize," Pt Trigger > %3.2f ",    fMinTriggerPt) ;
1083   parList+=onePar ;
1084   snprintf(onePar,buffersize," %3.2f < Pt associated < %3.2f ", fMinAssocPt,   fMaxAssocPt) ;
1085   parList+=onePar ;
1086   snprintf(onePar,buffersize," %3.2f < Phi trigger particle-Hadron < %3.2f ",    fDeltaPhiMinCut,   fDeltaPhiMaxCut) ;
1087   parList+=onePar ;
1088   snprintf(onePar,buffersize," %3.2f < Phi trigger particle-UeHadron <  %3.2f ", fUeDeltaPhiMinCut, fUeDeltaPhiMaxCut) ;
1089   parList+=onePar ;
1090   snprintf(onePar,buffersize,"Isolated Trigger?  %d\n",    fSelectIsolated) ;
1091   parList+=onePar ;
1092   snprintf(onePar,buffersize,"Several UE?  %d\n",          fMakeSeveralUE) ;
1093   parList+=onePar ;
1094   snprintf(onePar,buffersize,"Name of AOD Pi0 Branch %s ", fPi0AODBranchName.Data());
1095   parList+=onePar ;
1096   snprintf(onePar,buffersize,"Do Decay-hadron correlation ?  pi0 %d, decay %d", fPi0Trigger, fDecayTrigger) ;
1097   parList+=onePar ;
1098   snprintf(onePar,buffersize,"Select absolute leading for cluster triggers ? %d or Near Side Leading %d \n",
1099            fMakeAbsoluteLeading, fMakeNearSideLeading) ;
1100   parList+=onePar ;
1101   snprintf(onePar,buffersize,"Associated particle pt bins  %d: ", fNAssocPtBins) ;
1102   parList+=onePar ;
1103   for (Int_t ibin = 0; ibin<fNAssocPtBins; ibin++) {
1104     snprintf(onePar,buffersize,"bin %d = [%2.1f,%2.1f];", ibin, fAssocPtBinLimit[ibin], fAssocPtBinLimit[ibin+1]) ;
1105   }
1106   parList+=onePar ;
1107
1108   //Get parameters set in base class.
1109   parList += GetBaseParametersList() ;
1110
1111   //Get parameters set in FiducialCut class (not available yet)
1112   //parlist += GetFidCut()->GetFidCutParametersList()
1113
1114   return new TObjString(parList) ;
1115
1116 }
1117
1118 //________________________________________________________________
1119 TList *  AliAnaParticleHadronCorrelation::GetCreateOutputObjects()
1120 {
1121
1122   // Create histograms to be saved in output file and
1123   // store them in fOutputContainer
1124
1125   TList * outputContainer = new TList() ;
1126   outputContainer->SetName("CorrelationHistos") ;
1127
1128   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();
1129   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();
1130   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();
1131
1132   Int_t nMixBins = GetNCentrBin()*GetNZvertBin()*GetNRPBin();
1133
1134   TString nameMC[]     = {"Photon","Pi0","Pi0Decay","EtaDecay","OtherDecay","Electron","Hadron"};
1135   TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
1136
1137   // For vz dependent histograms, if option ON
1138   Int_t   nz  = 1  ;
1139   if(fCorrelVzBin) nz = GetNZvertBin();
1140   TString sz  = "" ;
1141   TString tz  = "" ;
1142
1143   fhPtTriggerInput  = new TH1F("hPtInput","Input trigger #it{p}_{T}", nptbins,ptmin,ptmax);
1144   fhPtTriggerInput->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1145   outputContainer->Add(fhPtTriggerInput);
1146
1147   if( fM02MaxCut > 0 && fM02MinCut > 0 )
1148   {
1149     fhPtTriggerSSCut  = new TH1F("hPtTriggerSSCut","Trigger #it{p}_{T} after #lambda^{2}_{0} cut", nptbins,ptmin,ptmax);
1150     fhPtTriggerSSCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1151     outputContainer->Add(fhPtTriggerSSCut);
1152   }
1153
1154   if( OnlyIsolated() )
1155   {
1156     fhPtTriggerIsoCut  = new TH1F("hPtTriggerIsoCut","Trigger #it{p}_{T} after isolation (and #lambda^{2}_{0}) cut", nptbins,ptmin,ptmax);
1157     fhPtTriggerIsoCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1158     outputContainer->Add(fhPtTriggerIsoCut);
1159   }
1160
1161   fhPtTriggerFidCut  = new TH1F("hPtTriggerFidCut","Trigger #it{p}_{T} after fiducial (isolation and #lambda^{2}_{0}) cut", nptbins,ptmin,ptmax);
1162   fhPtTriggerFidCut->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1163   outputContainer->Add(fhPtTriggerFidCut);
1164
1165   fhPtLeading  = new TH1F("hPtLeading","#it{p}_{T} distribution of leading particles", nptbins,ptmin,ptmax);
1166   fhPtLeading->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1167   outputContainer->Add(fhPtLeading);
1168
1169   if(IsDataMC())
1170   {
1171     for(Int_t i=0; i < 7; i++)
1172     {
1173       fhPtLeadingMC[i]  = new TH1F(Form("hPtLeading_MC%s",nameMC[i].Data()),
1174                                    Form("#it{p}_{T} distribution of leading particles, trigger origin is %s",nameMC[i].Data()),
1175                                    nptbins,ptmin,ptmax);
1176       fhPtLeadingMC[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1177       outputContainer->Add(fhPtLeadingMC[i]);
1178     }
1179   }
1180
1181   if(fCorrelVzBin)
1182   {
1183     fhPtLeadingVzBin  = new TH2F("hPtLeadingVzBin","#it{p}_{T} distribution of leading particles vs vz bin", nptbins,ptmin,ptmax,GetNZvertBin(),0,GetNZvertBin());
1184     fhPtLeadingVzBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1185     fhPtLeadingVzBin->SetYTitle("#it{v}_{#it{z}} bin");
1186     outputContainer->Add(fhPtLeadingVzBin);
1187   }
1188
1189   fhPtLeadingBin  = new TH2F ("hPtLeadingBin","#it{p}_{T} distribution of leading particles", nptbins,ptmin,ptmax,nMixBins,0,nMixBins);
1190   fhPtLeadingBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1191   fhPtLeadingBin->SetYTitle("Bin");
1192   outputContainer->Add(fhPtLeadingBin);
1193
1194   fhPhiLeading  = new TH2F ("hPhiLeading","#phi distribution of leading Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
1195   fhPhiLeading->SetYTitle("#phi (rad)");
1196   outputContainer->Add(fhPhiLeading);
1197
1198   fhEtaLeading  = new TH2F ("hEtaLeading","#eta distribution of leading",nptbins,ptmin,ptmax, netabins,etamin,etamax);
1199   fhEtaLeading->SetYTitle("#eta ");
1200   outputContainer->Add(fhEtaLeading);
1201
1202   fhPtLeadingCentrality   = new TH2F("hPtLeadingCentrality","Leading particle #it{p}_{T} vs centrality",nptbins,ptmin,ptmax,100,0.,100) ;
1203   fhPtLeadingCentrality->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1204   fhPtLeadingCentrality->SetYTitle("Centrality (%)");
1205   outputContainer->Add(fhPtLeadingCentrality) ;
1206
1207   fhPtLeadingEventPlane  = new TH2F("hPtLeadingEventPlane","Leading particle #it{p}_{T} vs event plane angle",nptbins,ptmin,ptmax, 100,0.,TMath::Pi()) ;
1208   fhPtLeadingEventPlane->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1209   fhPtLeadingEventPlane->SetXTitle("EP angle (rad)");
1210   outputContainer->Add(fhPtLeadingEventPlane) ;
1211
1212   fhLeadingEventPlaneCentrality  = new TH2F("hLeadingEventPlane","Leading particle centrality vs event plane angle",100,0.,100,100,0.,TMath::Pi()) ;
1213   fhLeadingEventPlaneCentrality->SetXTitle("Centrality (%)");
1214   fhLeadingEventPlaneCentrality->SetYTitle("EP angle (rad)");
1215   outputContainer->Add(fhLeadingEventPlaneCentrality) ;
1216
1217   //Correlation with charged hadrons
1218   if(GetReader()->IsCTSSwitchedOn())
1219   {
1220     fhDeltaPhiDeltaEtaCharged  = new TH2F
1221     ("hDeltaPhiDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #phi_{trigger} - #phi_{h^{#pm}}",
1222     ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1223     fhDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi (rad)");
1224     fhDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
1225
1226     fhDeltaPhiDeltaEtaChargedPtA3GeV  = new TH2F
1227     ("hDeltaPhiDeltaEtaChargedPtA3GeV","#eta_{trigger} - #eta_{h^{#pm}} vs #phi_{trigger} - #phi_{h^{#pm}, #it{p}_{TA}>3 GeV/#it{c}}",
1228      ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1229     fhDeltaPhiDeltaEtaChargedPtA3GeV->SetXTitle("#Delta #phi (rad)");
1230     fhDeltaPhiDeltaEtaChargedPtA3GeV->SetYTitle("#Delta #eta");
1231
1232     fhPhiCharged  = new TH2F
1233     ("hPhiCharged","#phi_{h^{#pm}}  vs #it{p}_{T #pm}",
1234      nptbins,ptmin,ptmax,nphibins,phimin,phimax);
1235     fhPhiCharged->SetYTitle("#phi_{h^{#pm}} (rad)");
1236     fhPhiCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
1237
1238     fhEtaCharged  = new TH2F
1239     ("hEtaCharged","#eta_{h^{#pm}}  vs #it{p}_{T #pm}",
1240      nptbins,ptmin,ptmax,netabins,etamin,etamax);
1241     fhEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
1242     fhEtaCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
1243
1244     fhDeltaPhiCharged  = new TH2F
1245     ("hDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}",
1246      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1247     fhDeltaPhiCharged->SetYTitle("#Delta #phi (rad)");
1248     fhDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1249
1250     fhDeltaPhiChargedPtA3GeV  = new TH2F
1251     ("hDeltaPhiChargedPtA3GeV","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}",
1252      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1253     fhDeltaPhiChargedPtA3GeV->SetYTitle("#Delta #phi (rad)");
1254     fhDeltaPhiChargedPtA3GeV->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1255
1256
1257     fhDeltaPhiChargedPt  = new TH2F
1258     ("hDeltaPhiChargedPt","#phi_{trigger} - #phi_{#h^{#pm}} vs #it{p}_{T h^{#pm}}",
1259      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1260     fhDeltaPhiChargedPt->SetYTitle("#Delta #phi (rad)");
1261     fhDeltaPhiChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1262
1263     fhDeltaPhiUeChargedPt  = new TH2F
1264     ("hDeltaPhiUeChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}}",
1265      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1266     fhDeltaPhiUeChargedPt->SetYTitle("#Delta #phi (rad)");
1267     fhDeltaPhiUeChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1268
1269     fhUePart  =  new TH1F("hUePart","UE particles distribution vs pt trig",
1270              nptbins,ptmin,ptmax);
1271     fhUePart->SetYTitle("dNch");
1272     fhUePart->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1273
1274
1275     fhDeltaEtaCharged  = new TH2F
1276     ("hDeltaEtaCharged","#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}",
1277      nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1278     fhDeltaEtaCharged->SetYTitle("#Delta #eta");
1279     fhDeltaEtaCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1280
1281     fhDeltaEtaChargedPtA3GeV  = new TH2F
1282     ("hDeltaEtaChargedPtA3GeV","#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}",
1283      nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1284     fhDeltaEtaChargedPtA3GeV->SetYTitle("#Delta #eta");
1285     fhDeltaEtaChargedPtA3GeV->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1286
1287     fhXECharged  =
1288     new TH2F("hXECharged","#it{x}_{#it{E}} for charged tracks",
1289              nptbins,ptmin,ptmax,200,0.,2.);
1290     fhXECharged->SetYTitle("#it{x}_{#it{E}}");
1291     fhXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1292
1293     fhXECharged_Cone2  =
1294     new TH2F("hXECharged_Cone2","#it{x}_{#it{E}} for charged tracks in cone 2 (5#pi/6-7#pi/6)",
1295              nptbins,ptmin,ptmax,200,0.,2.);
1296     fhXECharged_Cone2->SetYTitle("#it{x}_{#it{E}}");
1297     fhXECharged_Cone2->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1298
1299     fhXEUeCharged  =
1300     new TH2F("hXEUeCharged","#it{x}_{#it{E}} for Underlying Event",
1301              nptbins,ptmin,ptmax,200,0.,2.);
1302     fhXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
1303     fhXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1304
1305     fhXEPosCharged  =
1306     new TH2F("hXEPositiveCharged","#it{x}_{#it{E}} for positive charged tracks",
1307              nptbins,ptmin,ptmax,200,0.,2.);
1308     fhXEPosCharged->SetYTitle("#it{x}_{#it{E}}");
1309     fhXEPosCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1310
1311     fhXENegCharged  =
1312     new TH2F("hXENegativeCharged","#it{x}_{#it{E}} for negative charged tracks",
1313              nptbins,ptmin,ptmax,200,0.,2.);
1314     fhXENegCharged->SetYTitle("#it{x}_{#it{E}}");
1315     fhXENegCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1316
1317     fhPtHbpXECharged  =
1318     new TH2F("hHbpXECharged","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
1319              nptbins,ptmin,ptmax,200,0.,10.);
1320     fhPtHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1321     fhPtHbpXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1322
1323     fhPtHbpXECharged_Cone2  =
1324     new TH2F("hHbpXECharged_Cone2","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons in cone 2 (5#pi/6-7#pi/6)",
1325              nptbins,ptmin,ptmax,200,0.,10.);
1326     fhPtHbpXECharged_Cone2->SetYTitle("ln(1/#it{x}_{#it{E}})");
1327     fhPtHbpXECharged_Cone2->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1328
1329     fhPtHbpXEUeCharged  =
1330     new TH2F("hHbpXEUeCharged","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons,Underlying Event",
1331              nptbins,ptmin,ptmax,200,0.,10.);
1332     fhPtHbpXEUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1333     fhPtHbpXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1334
1335     fhZTCharged  =
1336     new TH2F("hZTCharged","#it{z}_{T} for charged tracks",
1337              nptbins,ptmin,ptmax,200,0.,2.);
1338     fhZTCharged->SetYTitle("#it{z}_{T}");
1339     fhZTCharged->SetXTitle("#it{p}_{T trigger}");
1340
1341     fhZTUeCharged  =
1342     new TH2F("hZTUeCharged","#it{z}_{T} for Underlying Event",
1343              nptbins,ptmin,ptmax,200,0.,2.);
1344     fhZTUeCharged->SetYTitle("#it{z}_{T}");
1345     fhZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1346
1347     fhZTPosCharged  =
1348     new TH2F("hZTPositiveCharged","#it{z}_{T} for positive charged tracks",
1349              nptbins,ptmin,ptmax,200,0.,2.);
1350     fhZTPosCharged->SetYTitle("#it{z}_{T}");
1351     fhZTPosCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1352
1353     fhZTNegCharged  =
1354     new TH2F("hZTNegativeCharged","#it{z}_{T} for negative charged tracks",
1355              nptbins,ptmin,ptmax,200,0.,2.);
1356     fhZTNegCharged->SetYTitle("#it{z}_{T}");
1357     fhZTNegCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1358
1359     fhPtHbpZTCharged  =
1360     new TH2F("hHbpZTCharged","#xi = ln(1/#it{z}_{T}) with charged hadrons",
1361              nptbins,ptmin,ptmax,200,0.,10.);
1362     fhPtHbpZTCharged->SetYTitle("ln(1/#it{z}_{T})");
1363     fhPtHbpZTCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1364
1365     fhPtHbpZTUeCharged  =
1366     new TH2F("hHbpZTUeCharged","#xi = ln(1/#it{z}_{T}) with charged hadrons,Underlying Event",
1367              nptbins,ptmin,ptmax,200,0.,10.);
1368     fhPtHbpZTUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1369     fhPtHbpZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1370
1371     fhPtTrigPout  =
1372     new TH2F("hPtTrigPout","Pout with triggers",
1373              nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax);
1374     fhPtTrigPout->SetYTitle("#it{p}_{out} (GeV/#it{c})");
1375     fhPtTrigPout->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1376
1377     fhPtTrigCharged  =
1378     new TH2F("hPtTrigCharged","trigger and charged tracks pt distribution",
1379              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1380     fhPtTrigCharged->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1381     fhPtTrigCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1382
1383     outputContainer->Add(fhDeltaPhiDeltaEtaCharged);
1384     outputContainer->Add(fhDeltaPhiDeltaEtaChargedPtA3GeV);
1385     outputContainer->Add(fhPhiCharged) ;
1386     outputContainer->Add(fhEtaCharged) ;
1387     outputContainer->Add(fhDeltaPhiCharged) ;
1388     outputContainer->Add(fhDeltaPhiChargedPtA3GeV) ;
1389     outputContainer->Add(fhDeltaEtaCharged) ;
1390     outputContainer->Add(fhDeltaEtaChargedPtA3GeV) ;
1391     outputContainer->Add(fhDeltaPhiChargedPt) ;
1392     outputContainer->Add(fhDeltaPhiUeChargedPt) ;
1393     outputContainer->Add(fhUePart);
1394
1395     outputContainer->Add(fhXECharged) ;
1396     outputContainer->Add(fhXECharged_Cone2) ;
1397     if(IsDataMC())
1398     {
1399       for(Int_t i=0; i < 7; i++)
1400       {
1401
1402         fhDeltaPhiChargedMC[i]  = new TH2F(Form("hDeltaPhiCharged_MC%s",nameMC[i].Data()),
1403                                      Form("#Delta #phi for charged tracks, trigger origin is %s",nameMC[i].Data()),
1404                                      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
1405         fhDeltaPhiChargedMC[i]->SetYTitle("#it{x}_{#it{E}}");
1406         fhDeltaPhiChargedMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1407         outputContainer->Add(fhDeltaPhiChargedMC[i]) ;
1408
1409         fhXEChargedMC[i]  = new TH2F(Form("hXECharged_MC%s",nameMC[i].Data()),
1410                                      Form("#it{x}_{#it{E}} for charged tracks, trigger origin is %s",nameMC[i].Data()),
1411          nptbins,ptmin,ptmax,200,0.,2.);
1412         fhXEChargedMC[i]->SetYTitle("#it{x}_{#it{E}}");
1413         fhXEChargedMC[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1414         outputContainer->Add(fhXEChargedMC[i]) ;
1415       }
1416     }
1417
1418     outputContainer->Add(fhXEPosCharged) ;
1419     outputContainer->Add(fhXENegCharged) ;
1420     outputContainer->Add(fhXEUeCharged) ;
1421     outputContainer->Add(fhPtHbpXECharged) ;
1422     outputContainer->Add(fhPtHbpXECharged_Cone2) ;
1423     outputContainer->Add(fhPtHbpXEUeCharged) ;
1424
1425     outputContainer->Add(fhZTCharged) ;
1426     outputContainer->Add(fhZTPosCharged) ;
1427     outputContainer->Add(fhZTNegCharged) ;
1428     outputContainer->Add(fhZTUeCharged) ;
1429     outputContainer->Add(fhPtHbpZTCharged) ;
1430     outputContainer->Add(fhPtHbpZTUeCharged) ;
1431
1432     outputContainer->Add(fhPtTrigPout) ;
1433     outputContainer->Add(fhPtTrigCharged) ;
1434
1435     if(fFillPileUpHistograms)
1436     {
1437       fhDeltaPhiChargedOtherBC  = new TH2F
1438       ("hDeltaPhiChargedOtherBC","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, track BC!=0",
1439        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1440       fhDeltaPhiChargedOtherBC->SetYTitle("#Delta #phi (rad)");
1441       fhDeltaPhiChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1442
1443       fhDeltaPhiChargedPtA3GeVOtherBC  = new TH2F
1444       ("hDeltaPhiChargedPtA3GeVOtherBC","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC!=0",
1445        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1446       fhDeltaPhiChargedPtA3GeVOtherBC->SetYTitle("#Delta #phi (rad)");
1447       fhDeltaPhiChargedPtA3GeVOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1448
1449       fhPtTrigChargedOtherBC  =
1450       new TH2F("hPtTrigChargedOtherBC","trigger and charged tracks pt distribution, track BC!=0",
1451                nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1452       fhPtTrigChargedOtherBC->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1453       fhPtTrigChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1454
1455       fhXEChargedOtherBC  =
1456       new TH2F("hXEChargedOtherBC","#it{x}_{#it{E}} for charged tracks, track BC!=0",
1457                nptbins,ptmin,ptmax,200,0.,2.);
1458       fhXEChargedOtherBC->SetYTitle("#it{x}_{#it{E}}");
1459       fhXEChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1460
1461       fhXEUeChargedOtherBC  =
1462       new TH2F("hXEUeChargedOtherBC","#it{x}_{#it{E}} for Underlying Event, track BC!=0",
1463                nptbins,ptmin,ptmax,200,0.,2.);
1464       fhXEUeChargedOtherBC->SetYTitle("#it{x}_{#it{E}}");
1465       fhXEUeChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1466
1467       fhZTChargedOtherBC  =
1468       new TH2F("hZTChargedOtherBC","#it{z}_{T} for charged tracks, track BC!=0",
1469                nptbins,ptmin,ptmax,200,0.,2.);
1470       fhZTChargedOtherBC->SetYTitle("#it{z}_{T}");
1471       fhZTChargedOtherBC->SetXTitle("#it{p}_{T trigger}");
1472
1473       fhZTUeChargedOtherBC  =
1474       new TH2F("hZTUeChargedOtherBC","#it{z}_{T} for Underlying Event, track BC!=0",
1475                nptbins,ptmin,ptmax,200,0.,2.);
1476       fhZTUeChargedOtherBC->SetYTitle("#it{z}_{T}");
1477       fhZTUeChargedOtherBC->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1478
1479       outputContainer->Add(fhDeltaPhiChargedOtherBC) ;
1480       outputContainer->Add(fhDeltaPhiChargedPtA3GeVOtherBC) ;
1481       outputContainer->Add(fhXEChargedOtherBC) ;
1482       outputContainer->Add(fhXEUeChargedOtherBC) ;
1483       outputContainer->Add(fhZTChargedOtherBC) ;
1484       outputContainer->Add(fhZTUeChargedOtherBC) ;
1485       outputContainer->Add(fhPtTrigChargedOtherBC) ;
1486
1487       fhDeltaPhiChargedBC0  = new TH2F
1488       ("hDeltaPhiChargedBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, track BC==0",
1489        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1490       fhDeltaPhiChargedBC0->SetYTitle("#Delta #phi (rad)");
1491       fhDeltaPhiChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1492
1493       fhDeltaPhiChargedPtA3GeVBC0  = new TH2F
1494       ("hDeltaPhiChargedPtA3GeVBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC==0",
1495        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1496       fhDeltaPhiChargedPtA3GeVBC0->SetYTitle("#Delta #phi (rad)");
1497       fhDeltaPhiChargedPtA3GeVBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1498
1499       fhPtTrigChargedBC0  =
1500       new TH2F("hPtTrigChargedBC0","trigger and charged tracks pt distribution, track BC==0",
1501                nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1502       fhPtTrigChargedBC0->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1503       fhPtTrigChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1504
1505       fhXEChargedBC0  =
1506       new TH2F("hXEChargedBC0","#it{x}_{#it{E}} for charged tracks, track BC==0",
1507                nptbins,ptmin,ptmax,200,0.,2.);
1508       fhXEChargedBC0->SetYTitle("#it{x}_{#it{E}}");
1509       fhXEChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1510
1511       fhXEUeChargedBC0  =
1512       new TH2F("hXEUeChargedBC0","#it{x}_{#it{E}} for Underlying Event, track BC==0",
1513                nptbins,ptmin,ptmax,200,0.,2.);
1514       fhXEUeChargedBC0->SetYTitle("#it{x}_{#it{E}}");
1515       fhXEUeChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1516
1517       fhZTChargedBC0  =
1518       new TH2F("hZTChargedBC0","#it{z}_{T} for charged tracks, track BC==0",
1519                nptbins,ptmin,ptmax,200,0.,2.);
1520       fhZTChargedBC0->SetYTitle("#it{z}_{T}");
1521       fhZTChargedBC0->SetXTitle("#it{p}_{T trigger}");
1522
1523       fhZTUeChargedBC0  =
1524       new TH2F("hZTUeChargedBC0","#it{z}_{T} for Underlying Event, track BC==0",
1525                nptbins,ptmin,ptmax,200,0.,2.);
1526       fhZTUeChargedBC0->SetYTitle("#it{z}_{T}");
1527       fhZTUeChargedBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1528
1529       outputContainer->Add(fhDeltaPhiChargedBC0) ;
1530       outputContainer->Add(fhDeltaPhiChargedPtA3GeVBC0) ;
1531       outputContainer->Add(fhXEChargedBC0) ;
1532       outputContainer->Add(fhXEUeChargedBC0) ;
1533       outputContainer->Add(fhZTChargedBC0) ;
1534       outputContainer->Add(fhZTUeChargedBC0) ;
1535       outputContainer->Add(fhPtTrigChargedBC0) ;
1536
1537       fhPtLeadingVtxBC0  = new TH1F("hPtLeadingVtxBC0","#it{p}_{T} distribution of leading particles", nptbins,ptmin,ptmax);
1538       fhPtLeadingVtxBC0->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1539
1540       fhDeltaPhiChargedVtxBC0  = new TH2F
1541       ("hDeltaPhiChargedVtxBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, track BC==0",
1542        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1543       fhDeltaPhiChargedVtxBC0->SetYTitle("#Delta #phi (rad)");
1544       fhDeltaPhiChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1545
1546       fhDeltaPhiChargedPtA3GeVVtxBC0  = new TH2F
1547       ("hDeltaPhiChargedPtA3GeVVtxBC0","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, track BC==0",
1548        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1549       fhDeltaPhiChargedPtA3GeVVtxBC0->SetYTitle("#Delta #phi (rad)");
1550       fhDeltaPhiChargedPtA3GeVVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1551
1552       fhPtTrigChargedVtxBC0  =
1553       new TH2F("hPtTrigChargedVtxBC0","trigger and charged tracks pt distribution, track BC==0",
1554                nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1555       fhPtTrigChargedVtxBC0->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1556       fhPtTrigChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1557
1558       fhXEChargedVtxBC0  =
1559       new TH2F("hXEChargedVtxBC0","#it{x}_{#it{E}} for charged tracks, track BC==0",
1560                nptbins,ptmin,ptmax,200,0.,2.);
1561       fhXEChargedVtxBC0->SetYTitle("#it{x}_{#it{E}}");
1562       fhXEChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1563
1564       fhXEUeChargedVtxBC0  =
1565       new TH2F("hXEUeChargedVtxBC0","#it{x}_{#it{E}} for Underlying Event, track BC==0",
1566                nptbins,ptmin,ptmax,200,0.,2.);
1567       fhXEUeChargedVtxBC0->SetYTitle("#it{x}_{#it{E}}");
1568       fhXEUeChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1569
1570       fhZTChargedVtxBC0  =
1571       new TH2F("hZTChargedVtxBC0","#it{z}_{T} for charged tracks, track BC==0",
1572                nptbins,ptmin,ptmax,200,0.,2.);
1573       fhZTChargedVtxBC0->SetYTitle("#it{z}_{T}");
1574       fhZTChargedVtxBC0->SetXTitle("#it{p}_{T trigger}");
1575
1576       fhZTUeChargedVtxBC0  =
1577       new TH2F("hZTUeChargedVtxBC0","#it{z}_{T} for Underlying Event, track BC==0",
1578                nptbins,ptmin,ptmax,200,0.,2.);
1579       fhZTUeChargedVtxBC0->SetYTitle("#it{z}_{T}");
1580       fhZTUeChargedVtxBC0->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1581
1582       outputContainer->Add(fhPtLeadingVtxBC0);
1583       outputContainer->Add(fhDeltaPhiChargedVtxBC0) ;
1584       outputContainer->Add(fhDeltaPhiChargedPtA3GeVVtxBC0) ;
1585       outputContainer->Add(fhXEChargedVtxBC0) ;
1586       outputContainer->Add(fhXEUeChargedVtxBC0) ;
1587       outputContainer->Add(fhZTChargedVtxBC0) ;
1588       outputContainer->Add(fhZTUeChargedVtxBC0) ;
1589       outputContainer->Add(fhPtTrigChargedVtxBC0) ;
1590
1591       for(Int_t i = 0 ; i < 7 ; i++)
1592       {
1593         fhPtLeadingPileUp[i]  = new TH1F(Form("hPtLeadingPileUp%s",pileUpName[i].Data()),
1594                                          Form("#it{p}_{T} distribution of leading particles, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
1595         fhPtLeadingPileUp[i]->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
1596         outputContainer->Add(fhPtLeadingPileUp[i]);
1597
1598         fhDeltaPhiChargedPileUp[i]  = new TH2F(Form("hDeltaPhiChargedPileUp%s",pileUpName[i].Data()),
1599                                                     Form("#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, %s Pile-Up event",pileUpName[i].Data()),
1600          nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1601         fhDeltaPhiChargedPileUp[i]->SetYTitle("#Delta #phi (rad)");
1602         fhDeltaPhiChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1603         outputContainer->Add(fhDeltaPhiChargedPileUp[i]) ;
1604
1605         fhDeltaPhiChargedPtA3GeVPileUp[i]  = new TH2F(Form("hDeltaPhiChargedPtA3GeVPileUp%s",pileUpName[i].Data()),
1606                                                            Form("#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, %s Pile-Up event",pileUpName[i].Data()),
1607          nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1608         fhDeltaPhiChargedPtA3GeVPileUp[i]->SetYTitle("#Delta #phi (rad)");
1609         fhDeltaPhiChargedPtA3GeVPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1610         outputContainer->Add(fhDeltaPhiChargedPtA3GeVPileUp[i]) ;
1611
1612         fhDeltaEtaChargedPileUp[i]  = new TH2F(Form("hDeltaEtaChargedPileUp%s",pileUpName[i].Data()),
1613                                                     Form("#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, %s Pile-Up event",pileUpName[i].Data()),
1614          nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1615         fhDeltaEtaChargedPileUp[i]->SetYTitle("#Delta #eta");
1616         fhDeltaEtaChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1617         outputContainer->Add(fhDeltaEtaChargedPileUp[i]) ;
1618
1619         fhDeltaEtaChargedPtA3GeVPileUp[i]  = new TH2F(Form("hDeltaEtaChargedPtA3GeVPileUp%s",pileUpName[i].Data()),
1620                                                            Form("#eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger}, #it{p}_{TA}>3 GeV/#it{c}, %s Pile-Up event",pileUpName[i].Data()),
1621          nptbins,ptmin,ptmax,ndeltaetabins,deltaetamin,deltaetamax);
1622         fhDeltaEtaChargedPtA3GeVPileUp[i]->SetYTitle("#Delta #eta");
1623         fhDeltaEtaChargedPtA3GeVPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1624         outputContainer->Add(fhDeltaEtaChargedPtA3GeVPileUp[i]) ;
1625
1626         fhXEChargedPileUp[i]  = new TH2F(Form("hXEChargedPileUp%s",pileUpName[i].Data()),
1627                                               Form("#it{x}_{#it{E}} for charged tracks, %s Pile-Up event",pileUpName[i].Data()),
1628                  nptbins,ptmin,ptmax,200,0.,2.);
1629         fhXEChargedPileUp[i]->SetYTitle("#it{x}_{#it{E}}");
1630         fhXEChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1631         outputContainer->Add(fhXEChargedPileUp[i]) ;
1632
1633         fhXEUeChargedPileUp[i]  = new TH2F(Form("hXEUeChargedPileUp%s",pileUpName[i].Data()),
1634                                                 Form("#it{x}_{#it{E}} for Underlying Event, %s Pile-Up event",pileUpName[i].Data()),
1635                  nptbins,ptmin,ptmax,200,0.,2.);
1636         fhXEUeChargedPileUp[i]->SetYTitle("#it{x}_{#it{E}}");
1637         fhXEUeChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1638         outputContainer->Add(fhXEUeChargedPileUp[i]) ;
1639
1640         fhZTChargedPileUp[i]  = new TH2F(Form("hZTChargedPileUp%s",pileUpName[i].Data()),
1641                                               Form("#it{z}_{T} for charged tracks, %s Pile-Up event",pileUpName[i].Data()),
1642                  nptbins,ptmin,ptmax,200,0.,2.);
1643         fhZTChargedPileUp[i]->SetYTitle("#it{z}_{T}");
1644         fhZTChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1645         outputContainer->Add(fhZTChargedPileUp[i]) ;
1646
1647         fhZTUeChargedPileUp[i]  = new TH2F(Form("hZTUeChargedPileUp%s",pileUpName[i].Data()),
1648                                                 Form("#it{z}_{T} for Underlying Event, %s Pile-Up event",pileUpName[i].Data()),
1649                  nptbins,ptmin,ptmax,200,0.,2.);
1650         fhZTUeChargedPileUp[i]->SetYTitle("#it{z}_{T}");
1651         fhZTUeChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1652         outputContainer->Add(fhZTUeChargedPileUp[i]) ;
1653
1654         fhPtTrigChargedPileUp[i]  = new TH2F(Form("hPtTrigChargedPileUp%s",pileUpName[i].Data()),
1655                                                   Form("trigger and charged tracks pt distribution, %s Pile-Up event",pileUpName[i].Data()),
1656                  nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
1657         fhPtTrigChargedPileUp[i]->SetYTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1658         fhPtTrigChargedPileUp[i]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1659         outputContainer->Add(fhPtTrigChargedPileUp[i]) ;
1660
1661       }
1662     }
1663
1664     if(DoEventSelect())
1665     {
1666       Int_t nMultiBins = GetMultiBin();
1667       fhTrigDeltaPhiCharged = new TH2F*[nMultiBins] ;
1668       fhTrigDeltaEtaCharged = new TH2F*[nMultiBins] ;
1669       fhTrigXECorr          = new TH2F*[nMultiBins] ;
1670       fhTrigXEUeCorr        = new TH2F*[nMultiBins] ;
1671       fhTrigZTCorr          = new TH2F*[nMultiBins] ;
1672       fhTrigZTUeCorr        = new TH2F*[nMultiBins] ;
1673
1674       for(Int_t im=0; im<nMultiBins; im++)
1675       {
1676         fhTrigDeltaPhiCharged[im]  = new TH2F
1677         (Form("hTrigDeltaPhiCharged_%d",im),Form("hTrigDeltaPhiCharged_%d",im), nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1678         fhTrigDeltaPhiCharged[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1679         fhTrigDeltaPhiCharged[im]->SetYTitle("#Delta #phi (rad)");
1680
1681         fhTrigDeltaEtaCharged[im]  = new TH2F
1682         (Form("hTrigDeltaEtaCharged_%d",im),Form("hTrigDeltaEtaCharged_%d",im), nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);
1683         fhTrigDeltaEtaCharged[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1684         fhTrigDeltaEtaCharged[im]->SetYTitle("#Delta #eta");
1685
1686         fhTrigXECorr[im]  = new TH2F
1687         (Form("hTrigXEPtCorr_%d",im),Form("hTrigXEPtCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
1688         fhTrigXECorr[im]->SetYTitle("#it{x}_{#it{E} trigger h^{#pm}}");
1689         fhTrigXECorr[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1690
1691         fhTrigXEUeCorr[im]  = new TH2F
1692         (Form("hTrigXEPtUeCorr_%d",im),Form("hTrigXEPtUeCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
1693         fhTrigXEUeCorr[im]->SetYTitle("#it{x}_{#it{E} trigger h^{#pm}}");
1694         fhTrigXEUeCorr[im]->SetXTitle("#it{p}_{T trigger}(GeV/#it{c})");
1695
1696         fhTrigZTCorr[im]  = new TH2F
1697         (Form("hTrigZTPtCorr_%d",im),Form("hTrigZTPtCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
1698         fhTrigZTCorr[im]->SetYTitle("#it{z}_{trigger h^{#pm}}");
1699         fhTrigZTCorr[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1700
1701         fhTrigZTUeCorr[im]  = new TH2F
1702         (Form("hTrigZTPtUeCorr_%d",im),Form("hTrigZTPtUeCorr_%d",im), nptbins,ptmin,ptmax,200,0.,2.);
1703         fhTrigZTUeCorr[im]->SetYTitle("#it{z}_{trigger h^{#pm}}");
1704         fhTrigZTUeCorr[im]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1705
1706         outputContainer->Add(fhTrigDeltaPhiCharged[im]) ;
1707         outputContainer->Add(fhTrigDeltaEtaCharged[im]) ;
1708         outputContainer->Add(fhTrigXECorr[im]);
1709         outputContainer->Add(fhTrigXEUeCorr[im]);
1710         outputContainer->Add(fhTrigZTCorr[im]);
1711         outputContainer->Add(fhTrigZTUeCorr[im]);
1712       }
1713     }
1714
1715     if(fFillBradHisto)
1716     {
1717       fhAssocPtBkg        = new TH2F("hAssocPtBkg", " Trigger #it{p}_{T} vs associated hadron #it{p}_{T} from background",
1718                                      nptbins, ptmin, ptmax,nptbins,ptmin,ptmax);
1719       fhAssocPtBkg->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1720       fhAssocPtBkg->SetYTitle("#it{p}_{T associated} (GeV/#it{c})");
1721       outputContainer->Add(fhAssocPtBkg) ;
1722
1723       fhDeltaPhiBrad = new TH2F("hDeltaPhiBrad","atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi vs #it{p}_{T trigger} ",
1724                                 nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
1725       fhDeltaPhiBrad->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1726       fhDeltaPhiBrad->SetYTitle("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi");
1727       outputContainer->Add(fhDeltaPhiBrad) ;
1728     }
1729
1730     fhDeltaPhiDeltaEtaAssocPtBin = new TH2F*[fNAssocPtBins*nz];
1731     fhDeltaPhiAssocPtBin       = new TH2F*[fNAssocPtBins*nz];
1732     fhDeltaPhiAssocPtBinDEta08 = new TH2F*[fNAssocPtBins*nz];
1733     fhDeltaPhiAssocPtBinDEta0  = new TH2F*[fNAssocPtBins*nz];
1734     fhXEAssocPtBin             = new TH2F*[fNAssocPtBins*nz];
1735     fhZTAssocPtBin             = new TH2F*[fNAssocPtBins*nz];
1736
1737     if(fFillBradHisto)
1738       fhDeltaPhiBradAssocPtBin = new TH2F*[fNAssocPtBins*nz];
1739
1740     if(fPi0Trigger || fDecayTrigger)
1741     {
1742       fhDeltaPhiAssocPtBin       = new TH2F*[fNAssocPtBins*nz];
1743       fhDeltaPhiAssocPtBinDEta08 = new TH2F*[fNAssocPtBins*nz];
1744       fhXEAssocPtBin             = new TH2F*[fNAssocPtBins*nz];
1745       fhZTAssocPtBin             = new TH2F*[fNAssocPtBins*nz];
1746       fhXEDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
1747       fhZTDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
1748       fhDeltaPhiDecayChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
1749     }
1750
1751     if(fHMPIDCorrelation)
1752     {
1753       fhDeltaPhiAssocPtBinHMPID   = new TH2F*[fNAssocPtBins*nz];
1754       fhDeltaPhiAssocPtBinHMPIDAcc= new TH2F*[fNAssocPtBins*nz];
1755     }
1756
1757     for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
1758     {
1759       for(Int_t z = 0 ; z < nz ; z++)
1760       {
1761         Int_t bin = i*nz+z;
1762
1763         if(fCorrelVzBin)
1764         {
1765           sz = Form("_vz%d",z);
1766           tz = Form(", #it{v}_{#it{z}} bin %d",z);
1767         }
1768
1769         //printf("iAssoc %d, Vz %d, bin %d - sz %s, tz %s       \n",i,z,bin,sz.Data(),tz.Data());
1770
1771         fhDeltaPhiDeltaEtaAssocPtBin[bin]  = new TH2F(Form("hDeltaPhiDeltaEtaPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1772                                                  Form("#Delta #phi vs #Delta #eta for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1773                                                  ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins,deltaetamin,deltaetamax);
1774         fhDeltaPhiDeltaEtaAssocPtBin[bin]->SetXTitle("#Delta #phi (rad)");
1775         fhDeltaPhiDeltaEtaAssocPtBin[bin]->SetYTitle("#Delta #eta");
1776
1777         fhDeltaPhiAssocPtBin[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1778                                            Form("#Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1779                                            nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1780         fhDeltaPhiAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1781         fhDeltaPhiAssocPtBin[bin]->SetYTitle("#Delta #phi (rad)");
1782
1783         fhDeltaPhiAssocPtBinDEta08[bin] = new TH2F(Form("hDeltaPhiDeltaEta0.8PtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1784                                                  Form("#Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta > 0.8", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1785                                                  nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1786         fhDeltaPhiAssocPtBinDEta08[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1787         fhDeltaPhiAssocPtBinDEta08[bin]->SetYTitle("#Delta #phi (rad)");
1788
1789         fhDeltaPhiAssocPtBinDEta0[bin] = new TH2F(Form("hDeltaPhiDeltaEta0PtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1790                                                    Form("#Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta = 0.", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1791                                                    nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1792         fhDeltaPhiAssocPtBinDEta0[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1793         fhDeltaPhiAssocPtBinDEta0[bin]->SetYTitle("#Delta #phi (rad)");
1794
1795         fhXEAssocPtBin[bin]       = new TH2F(Form("hXEAssocPtBin%1.f_%1.f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1796                                            Form("#it{x}_{#it{E}} vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1797                                            nptbins, ptmin, ptmax,200, 0.0, 2.0);
1798         fhXEAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1799         fhXEAssocPtBin[bin]->SetYTitle("#it{x}_{#it{E}}");
1800
1801         fhZTAssocPtBin[bin]       = new TH2F(Form("hZTAssocPtBin%1.f_%1.f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1802                                            Form("#it{z}_{T} vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1803                                            nptbins, ptmin, ptmax,200, 0.0, 2.0);
1804         fhZTAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1805         fhZTAssocPtBin[bin]->SetYTitle("#it{z}_{T}");
1806
1807         outputContainer->Add(fhDeltaPhiDeltaEtaAssocPtBin[bin]) ;
1808         outputContainer->Add(fhDeltaPhiAssocPtBin[bin]) ;
1809         outputContainer->Add(fhDeltaPhiAssocPtBinDEta08[bin]) ;
1810         outputContainer->Add(fhDeltaPhiAssocPtBinDEta0[bin]) ;
1811         outputContainer->Add(fhXEAssocPtBin[bin]);
1812         outputContainer->Add(fhZTAssocPtBin[bin]);
1813
1814         if(fPi0Trigger || fDecayTrigger)
1815         {
1816           fhDeltaPhiDecayChargedAssocPtBin[bin] = new TH2F(Form("hDeltaPhiPtDecayChargedAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1817                                                          Form("#Delta #phi vs #it{p}_{T trigger} tagged as decay for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1818                                                          nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1819           fhDeltaPhiDecayChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1820           fhDeltaPhiDecayChargedAssocPtBin[bin]->SetYTitle("#Delta #phi (rad)");
1821
1822           fhXEDecayChargedAssocPtBin[bin]       = new TH2F(Form("hXEDecayChargedAssocPtBin%1.f_%1.f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1823                                                          Form("#it{x}_{#it{E}} vs #it{p}_{T trigger} tagged as decay for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1824                                                          nptbins, ptmin, ptmax,200, 0.0, 2.0);
1825           fhXEDecayChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1826           fhXEDecayChargedAssocPtBin[bin]->SetYTitle("#it{x}_{#it{E}}");
1827
1828           fhZTDecayChargedAssocPtBin[bin]       = new TH2F(Form("hZTDecayChargedAssocPtBin%1.f_%1.f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1829                                                          Form("#it{z}_{T} vs #it{p}_{T trigger} tagged as decay for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1830                                                          nptbins, ptmin, ptmax,200, 0.0, 2.0);
1831           fhZTDecayChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1832           fhZTDecayChargedAssocPtBin[bin]->SetYTitle("#it{z}_{T}");
1833
1834           outputContainer->Add(fhDeltaPhiDecayChargedAssocPtBin[bin]) ;
1835           outputContainer->Add(fhXEDecayChargedAssocPtBin[bin]);
1836           outputContainer->Add(fhZTDecayChargedAssocPtBin[bin]);
1837
1838         }
1839
1840         if(fFillBradHisto)
1841         {
1842           fhDeltaPhiBradAssocPtBin[bin] = new TH2F(Form("hDeltaPhiBradPtAssocPt%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1843                                                  Form("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1844                                                  nptbins, ptmin, ptmax,288, -1.0/3.0, 5.0/3.0);
1845           fhDeltaPhiBradAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1846           fhDeltaPhiBradAssocPtBin[bin]->SetYTitle("atan2(sin(#Delta #phi), cos(#Delta #phi))/#pi");
1847           outputContainer->Add(fhDeltaPhiBradAssocPtBin[bin]) ;
1848         }
1849
1850         if(fHMPIDCorrelation)
1851         {
1852           fhDeltaPhiAssocPtBinHMPID[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%sHMPID", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1853                                                   Form("#Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, with track having HMPID signal", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1854                                                   nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1855           fhDeltaPhiAssocPtBinHMPID[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})" );
1856           fhDeltaPhiAssocPtBinHMPID[bin]->SetYTitle("#Delta #phi (rad)");
1857
1858           fhDeltaPhiAssocPtBinHMPIDAcc[bin] = new TH2F(Form("hDeltaPhiPtAssocPt%2.1f_%2.1f%sHMPIDAcc", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
1859                                                      Form("#Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, with track within 5<phi<20 deg", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
1860                                                      nptbins, ptmin, ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1861           fhDeltaPhiAssocPtBinHMPIDAcc[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1862           fhDeltaPhiAssocPtBinHMPIDAcc[bin]->SetYTitle("#Delta #phi (rad)");
1863
1864           outputContainer->Add(fhDeltaPhiAssocPtBinHMPID[bin]) ;
1865           outputContainer->Add(fhDeltaPhiAssocPtBinHMPIDAcc[bin]) ;
1866
1867         }
1868       }
1869     }
1870
1871     if(fPi0Trigger || fDecayTrigger)
1872     {
1873       if(fPi0Trigger)
1874       {
1875         fhPtPi0DecayRatio  = new TH2F
1876         ("hPtPi0DecayRatio","#it{p}_{T} of #pi^{0} and the ratio of pt for two decay",
1877          nptbins,ptmin,ptmax, 100,0.,2.);
1878         fhPtPi0DecayRatio->SetXTitle("#it{p}_{T}^{#pi^{0}} (GeV/#it{c})");
1879         fhPtPi0DecayRatio->SetYTitle("#it{p}_{T}^{Decay}/#it{p}_{T}^{#pi^{0}}");
1880         outputContainer->Add(fhPtPi0DecayRatio) ;
1881       }
1882
1883       fhDeltaPhiDecayCharged  = new TH2F
1884       ("hDeltaPhiDecayCharged","#phi_{Decay} - #phi_{h^{#pm}} vs #it{p}_{T Decay}",
1885        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1886       fhDeltaPhiDecayCharged->SetYTitle("#Delta #phi (rad)");
1887       fhDeltaPhiDecayCharged->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
1888
1889       fhXEDecayCharged  =
1890       new TH2F("hXEDecayCharged","#it{x}_{#it{E}}  Decay",
1891                nptbins,ptmin,ptmax,200,0.,2.);
1892       fhXEDecayCharged->SetYTitle("#it{x}_{#it{E}}");
1893       fhXEDecayCharged->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
1894
1895       fhZTDecayCharged  =
1896       new TH2F("hZTDecayCharged","#it{z}_{trigger h^{#pm}} = #it{p}_{T h^{#pm}} / #it{p}_{T Decay}",
1897                nptbins,ptmin,ptmax,200,0.,2.);
1898       fhZTDecayCharged->SetYTitle("#it{z}_{decay h^{#pm}}");
1899       fhZTDecayCharged->SetXTitle("#it{p}_{T decay} (GeV/#it{c})");
1900
1901       outputContainer->Add(fhDeltaPhiDecayCharged) ;
1902       outputContainer->Add(fhXEDecayCharged) ;
1903       outputContainer->Add(fhZTDecayCharged) ;
1904     }
1905
1906     if(fMakeSeveralUE)
1907     {
1908       fhDeltaPhiUeLeftCharged  = new TH2F
1909       ("hDeltaPhiUeLeftChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left side range of trigger particles",
1910        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1911       fhDeltaPhiUeLeftCharged->SetYTitle("#Delta #phi (rad)");
1912       fhDeltaPhiUeLeftCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1913       outputContainer->Add(fhDeltaPhiUeLeftCharged) ;
1914
1915       fhDeltaPhiUeRightCharged  = new TH2F
1916       ("hDeltaPhiUeRightChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right side range of trigger particles",
1917        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1918       fhDeltaPhiUeRightCharged->SetYTitle("#Delta #phi (rad)");
1919       fhDeltaPhiUeRightCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1920       outputContainer->Add(fhDeltaPhiUeRightCharged) ;
1921
1922       fhDeltaPhiUeLeftUpCharged  = new TH2F
1923       ("hDeltaPhiUeLeftUpChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left Up side range of trigger particles",
1924        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1925       fhDeltaPhiUeLeftUpCharged->SetYTitle("#Delta #phi (rad)");
1926       fhDeltaPhiUeLeftUpCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1927       outputContainer->Add(fhDeltaPhiUeLeftUpCharged) ;
1928
1929       fhDeltaPhiUeRightUpCharged  = new TH2F
1930       ("hDeltaPhiUeRightUpChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right Up side range of trigger particles",
1931        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1932       fhDeltaPhiUeRightUpCharged->SetYTitle("#Delta #phi (rad)");
1933       fhDeltaPhiUeRightUpCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1934       outputContainer->Add(fhDeltaPhiUeRightUpCharged) ;
1935
1936       fhDeltaPhiUeLeftDownCharged  = new TH2F
1937       ("hDeltaPhiUeLeftDownChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE left Down side range of trigger particles",
1938        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1939       fhDeltaPhiUeLeftDownCharged->SetYTitle("#Delta #phi (rad)");
1940       fhDeltaPhiUeLeftDownCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1941       outputContainer->Add(fhDeltaPhiUeLeftDownCharged) ;
1942
1943       fhDeltaPhiUeRightDownCharged  = new TH2F
1944       ("hDeltaPhiUeRightDownChargedPt","#phi_{trigger} - #phi_{#Ueh^{#pm}} vs #it{p}_{T Ueh^{#pm}} with UE right Down side range of trigger particles",
1945        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
1946       fhDeltaPhiUeRightDownCharged->SetYTitle("#Delta #phi (rad)");
1947       fhDeltaPhiUeRightDownCharged->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
1948       outputContainer->Add(fhDeltaPhiUeRightDownCharged) ;
1949
1950       fhXEUeLeftCharged  =
1951       new TH2F("hXEUeChargedLeft","#it{x}_{#it{E}} with UE left side of trigger",
1952                nptbins,ptmin,ptmax,200,0.,2.);
1953       fhXEUeLeftCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
1954       fhXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1955       outputContainer->Add(fhXEUeLeftCharged) ;
1956
1957       fhXEUeRightCharged  =
1958       new TH2F("hXEUeChargedRight","#it{x}_{#it{E} h^{#pm}} with UE right side of trigger",
1959                nptbins,ptmin,ptmax,200,0.,2.);
1960       fhXEUeRightCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
1961       fhXEUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1962       outputContainer->Add(fhXEUeRightCharged) ;
1963
1964       fhXEUeLeftUpCharged  =
1965       new TH2F("hXEUeChargedLeftUp","#it{x}_{#it{E}} with UE left Up side of trigger",
1966                nptbins,ptmin,ptmax,200,0.,2.);
1967       fhXEUeLeftUpCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
1968       fhXEUeLeftUpCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1969       outputContainer->Add(fhXEUeLeftUpCharged) ;
1970
1971       fhXEUeRightUpCharged  =
1972       new TH2F("hXEUeChargedRightUp","#it{x}_{#it{E} h^{#pm}} with UE right Up side of trigger",
1973                nptbins,ptmin,ptmax,200,0.,2.);
1974       fhXEUeRightUpCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
1975       fhXEUeRightUpCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1976       outputContainer->Add(fhXEUeRightUpCharged) ;
1977
1978       fhXEUeLeftDownCharged  =
1979       new TH2F("hXEUeChargedLeftDown","#it{x}_{#it{E}} with UE left Down side of trigger",
1980                nptbins,ptmin,ptmax,200,0.,2.);
1981       fhXEUeLeftDownCharged->SetYTitle("#it{x}_{#it{E} Ueh^{#pm}}");
1982       fhXEUeLeftDownCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1983       outputContainer->Add(fhXEUeLeftDownCharged) ;
1984
1985       fhXEUeRightDownCharged  =
1986       new TH2F("hXEUeChargedRightDown","#it{x}_{#it{E} h^{#pm}} with UE right Down side of trigger",
1987                nptbins,ptmin,ptmax,200,0.,2.);
1988       fhXEUeRightDownCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
1989       fhXEUeRightDownCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1990       outputContainer->Add(fhXEUeRightDownCharged) ;
1991
1992       fhPtHbpXEUeLeftCharged  =
1993       new TH2F("hHbpXEUeChargedLeft","#xi = ln(1/#it{x}_{#it{E}}) with charged UE left side of trigger",
1994                nptbins,ptmin,ptmax,200,0.,10.);
1995       fhPtHbpXEUeLeftCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
1996       fhPtHbpXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
1997       outputContainer->Add(fhPtHbpXEUeLeftCharged) ;
1998
1999       fhPtHbpXEUeRightCharged  =
2000       new TH2F("hHbpXEUeChargedRight","#xi = ln(1/#it{x}_{#it{E}}) with charged UE right side of trigger",
2001                nptbins,ptmin,ptmax,200,0.,10.);
2002       fhPtHbpXEUeRightCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2003       fhPtHbpXEUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2004       outputContainer->Add(fhPtHbpXEUeRightCharged) ;
2005
2006       fhZTUeLeftCharged  =
2007       new TH2F("hZTUeChargedLeft","#it{z}_{trigger h^{#pm}} = #it{p}_{T Ueh^{#pm}} / #it{p}_{T trigger} with UE left side of trigger",
2008                nptbins,ptmin,ptmax,200,0.,2.);
2009       fhZTUeLeftCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
2010       fhZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2011       outputContainer->Add(fhZTUeLeftCharged) ;
2012
2013       fhZTUeRightCharged  =
2014       new TH2F("hZTUeChargedRight","#it{z}_{trigger h^{#pm}} = #it{p}_{T Ueh^{#pm}} / #it{p}_{T trigger} with UE right side of trigger",
2015                nptbins,ptmin,ptmax,200,0.,2.);
2016       fhZTUeRightCharged->SetYTitle("#it{z}_{trigger Ueh^{#pm}}");
2017       fhZTUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2018       outputContainer->Add(fhZTUeRightCharged) ;
2019
2020       fhPtHbpZTUeLeftCharged  =
2021       new TH2F("hHbpZTUeChargedLeft","#xi = ln(1/#it{z}_{T}) with charged UE left side of trigger",
2022                nptbins,ptmin,ptmax,200,0.,10.);
2023       fhPtHbpZTUeLeftCharged->SetYTitle("ln(1/#it{z}_{T})");
2024       fhPtHbpZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2025       outputContainer->Add(fhPtHbpZTUeLeftCharged) ;
2026
2027       fhPtHbpZTUeRightCharged  =
2028       new TH2F("hHbpZTUeChargedRight","#xi = ln(1/#it{z}_{T}) with charged UE right side of trigger",
2029                nptbins,ptmin,ptmax,200,0.,10.);
2030       fhPtHbpZTUeRightCharged->SetYTitle("ln(1/#it{z}_{T})");
2031       fhPtHbpZTUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2032       outputContainer->Add(fhPtHbpZTUeRightCharged) ;
2033
2034     }
2035   }  //Correlation with charged hadrons
2036
2037   //Correlation with neutral hadrons
2038   if(fNeutralCorr)
2039   {
2040     fhDeltaPhiDeltaEtaNeutral  = new TH2F
2041     ("hDeltaPhiDeltaEtaNeutral","#phi_{trigger} - #phi_{h^{0}} vs #eta_{trigger} - #eta_{h^{0}}",
2042      ndeltaphibins ,deltaphimin,deltaphimax, ndeltaetabins ,deltaetamin,deltaetamax);
2043     fhDeltaPhiDeltaEtaNeutral->SetXTitle("#Delta #phi (rad)");
2044     fhDeltaPhiDeltaEtaNeutral->SetYTitle("#Delta #eta");
2045
2046     fhPhiNeutral  = new TH2F
2047     ("hPhiNeutral","#phi_{#pi^{0}}  vs #it{p}_{T #pi^{0}}",
2048      nptbins,ptmin,ptmax,nphibins,phimin,phimax);
2049     fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
2050     fhPhiNeutral->SetXTitle("#it{p}_{T #pi^{0}} (GeV/#it{c})");
2051
2052     fhEtaNeutral  = new TH2F
2053     ("hEtaNeutral","#eta_{#pi^{0}}  vs #it{p}_{T #pi^{0}}",
2054      nptbins,ptmin,ptmax,netabins,etamin,etamax);
2055     fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
2056     fhEtaNeutral->SetXTitle("#it{p}_{T #pi^{0}} (GeV/#it{c})");
2057
2058     fhDeltaPhiNeutral  = new TH2F
2059     ("hDeltaPhiNeutral","#phi_{trigger} - #phi_{#pi^{0}} vs #it{p}_{T trigger}",
2060      nptbins,ptmin,ptmax,nphibins,phimin,phimax);
2061     fhDeltaPhiNeutral->SetYTitle("#Delta #phi (rad)");
2062     fhDeltaPhiNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2063
2064     fhDeltaPhiNeutralPt  = new TH2F
2065     ("hDeltaPhiNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs #it{p}_{T #pi^{0}}}",
2066      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2067     fhDeltaPhiNeutralPt->SetYTitle("#Delta #phi (rad)");
2068     fhDeltaPhiNeutralPt->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2069
2070     fhDeltaPhiUeNeutralPt  = new TH2F
2071     ("hDeltaPhiUeNeutralPt","#phi_{trigger} - #phi_{#pi^{0}} vs #it{p}_{T #pi^{0}}}",
2072      nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2073     fhDeltaPhiUeNeutralPt->SetYTitle("#Delta #phi (rad)");
2074     fhDeltaPhiUeNeutralPt->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2075
2076     fhDeltaEtaNeutral  = new TH2F
2077     ("hDeltaEtaNeutral","#eta_{trigger} - #eta_{#pi^{0}} vs #it{p}_{T trigger}",
2078      nptbins,ptmin,ptmax, ndeltaetabins ,deltaetamin,deltaetamax);
2079     fhDeltaEtaNeutral->SetYTitle("#Delta #eta");
2080     fhDeltaEtaNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2081
2082     fhXENeutral  =
2083     new TH2F("hXENeutral","#it{x}_{#it{E}} for #pi^{0} associated",
2084              nptbins,ptmin,ptmax,200,0.,2.);
2085     fhXENeutral->SetYTitle("#it{x}_{#it{E}}");
2086     fhXENeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2087
2088     fhXEUeNeutral  =
2089     new TH2F("hXEUeNeutral","#it{x}_{#it{E}} for #pi^{0} associated",
2090              nptbins,ptmin,ptmax,200,0.,2.);
2091     fhXEUeNeutral->SetYTitle("#it{x}_{#it{E}}");
2092     fhXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2093
2094     fhPtHbpXENeutral  =
2095     new TH2F("hHbpXENeutral","#xi = ln(1/#it{x}_{#it{E}})for #pi^{0} associated",
2096              nptbins,ptmin,ptmax,200,0.,10.);
2097     fhPtHbpXENeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2098     fhPtHbpXENeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2099
2100     fhPtHbpXEUeNeutral  =
2101     new TH2F("hHbpXEUeNeutral","#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2102              nptbins,ptmin,ptmax,200,0.,10.);
2103     fhPtHbpXEUeNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2104     fhPtHbpXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2105
2106     fhZTNeutral  =
2107     new TH2F("hZTNeutral","#it{z}_{trigger #pi} = #it{p}_{T #pi^{0}} / #it{p}_{T trigger} for #pi^{0} associated",
2108              nptbins,ptmin,ptmax,200,0.,2.);
2109     fhZTNeutral->SetYTitle("#it{z}_{trigger #pi^{0}}");
2110     fhZTNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2111
2112     fhZTUeNeutral  =
2113     new TH2F("hZTUeNeutral","#it{z}_{trigger #pi} = #it{p}_{T #pi^{0}} / #it{p}_{T trigger} for #pi^{0} associated",
2114              nptbins,ptmin,ptmax,200,0.,2.);
2115     fhZTUeNeutral->SetYTitle("#it{z}_{trigger #pi^{0}}");
2116     fhZTUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2117
2118     fhPtHbpZTNeutral  =
2119     new TH2F("hHbpZTNeutral","#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2120              nptbins,ptmin,ptmax,200,0.,10.);
2121     fhPtHbpZTNeutral->SetYTitle("ln(1/#it{z}_{T})");
2122     fhPtHbpZTNeutral->SetXTitle("#it{p}_{T trigger}");
2123
2124     fhPtHbpZTUeNeutral  =
2125     new TH2F("hHbpZTUeNeutral","#xi = ln(1/#it{x}_{#it{E}}) for #pi^{0} associated",
2126              nptbins,ptmin,ptmax,200,0.,10.);
2127     fhPtHbpXEUeNeutral->SetYTitle("ln(1/#it{z}_{T})");
2128     fhPtHbpXEUeNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2129
2130     outputContainer->Add(fhDeltaPhiDeltaEtaNeutral);
2131     outputContainer->Add(fhPhiNeutral) ;
2132     outputContainer->Add(fhEtaNeutral) ;
2133     outputContainer->Add(fhDeltaPhiNeutral) ;
2134     outputContainer->Add(fhDeltaPhiNeutralPt) ;
2135     outputContainer->Add(fhDeltaPhiUeNeutralPt) ;
2136     outputContainer->Add(fhDeltaEtaNeutral) ;
2137     outputContainer->Add(fhXENeutral) ;
2138     outputContainer->Add(fhXEUeNeutral) ;
2139     outputContainer->Add(fhPtHbpXENeutral) ;
2140     outputContainer->Add(fhPtHbpXEUeNeutral) ;
2141     outputContainer->Add(fhZTNeutral) ;
2142     outputContainer->Add(fhZTUeNeutral) ;
2143     outputContainer->Add(fhPtHbpZTNeutral) ;
2144     outputContainer->Add(fhPtHbpZTUeNeutral) ;
2145
2146     if(fPi0Trigger || fDecayTrigger)
2147     {
2148       fhDeltaPhiDecayNeutral  = new TH2F
2149       ("hDeltaPhiDecayNeutral","#phi_{Decay} - #phi_{h^{0}} vs #it{p}_{T Decay}",
2150        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2151       fhDeltaPhiDecayNeutral->SetYTitle("#Delta #phi (rad)");
2152       fhDeltaPhiDecayNeutral->SetXTitle("#it{p}_{T Decay} (GeV/#it{c})");
2153
2154       fhXEDecayNeutral  =
2155       new TH2F("hXEDecayNeutral","#it{x}_{#it{E}} for decay trigger",
2156                nptbins,ptmin,ptmax,200,0.,2.);
2157       fhXEDecayNeutral->SetYTitle("#it{x}_{#it{E}}");
2158       fhXEDecayNeutral->SetXTitle("#it{p}_{T decay}");
2159
2160       fhZTDecayNeutral  =
2161       new TH2F("hZTDecayNeutral","#it{z}_{trigger h^{0}} = #it{p}_{T h^{0}} / #it{p}_{T Decay}",
2162                nptbins,ptmin,ptmax,200,0.,2.);
2163       fhZTDecayNeutral->SetYTitle("#it{z}_{h^{0}}");
2164       fhZTDecayNeutral->SetXTitle("#it{p}_{T decay}");
2165
2166       outputContainer->Add(fhDeltaPhiDecayNeutral) ;
2167       outputContainer->Add(fhXEDecayNeutral) ;
2168       outputContainer->Add(fhZTDecayNeutral) ;
2169
2170     }
2171
2172     if(fMakeSeveralUE)
2173     {
2174       fhDeltaPhiUeLeftNeutral  = new TH2F
2175       ("hDeltaPhiUeLeftNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs #it{p}_{T h^{0}} with neutral UE left side range of trigger particles",
2176        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2177       fhDeltaPhiUeLeftNeutral->SetYTitle("#Delta #phi (rad)");
2178       fhDeltaPhiUeLeftNeutral->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2179       outputContainer->Add(fhDeltaPhiUeLeftNeutral) ;
2180
2181       fhDeltaPhiUeRightNeutral  = new TH2F
2182       ("hDeltaPhiUeRightNeutralPt","#phi_{trigger} - #phi_{#Ueh^{0}} vs #it{p}_{T Ueh^{0}} with neutral UE right side range of trigger particles",
2183        nptbins,ptmin,ptmax, ndeltaphibins ,deltaphimin,deltaphimax);
2184       fhDeltaPhiUeRightNeutral->SetYTitle("#Delta #phi (rad)");
2185       fhDeltaPhiUeRightNeutral->SetXTitle("#it{p}_{T h^{0}} (GeV/#it{c})");
2186       outputContainer->Add(fhDeltaPhiUeRightNeutral) ;
2187
2188       fhXEUeLeftNeutral  =
2189       new TH2F("hXEUeNeutralLeft","#it{x}_{#it{E}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE left side of trigger",
2190                nptbins,ptmin,ptmax,140,0.,2.);
2191       fhXEUeLeftNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2192       fhXEUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2193       outputContainer->Add(fhXEUeLeftNeutral) ;
2194
2195       fhXEUeRightNeutral  =
2196       new TH2F("hXEUeNeutralRight","#it{x}_{#it{E}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE right side of trigger",
2197                nptbins,ptmin,ptmax,200,0.,2.);
2198       fhXEUeRightNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2199       fhXEUeRightNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2200       outputContainer->Add(fhXEUeRightNeutral) ;
2201
2202       fhPtHbpXEUeLeftNeutral  =
2203       new TH2F("hHbpXEUeNeutralLeft","#xi = ln(1/#it{x}_{#it{E}}) with neutral UE left side of trigger",
2204                nptbins,ptmin,ptmax,200,0.,10.);
2205       fhPtHbpXEUeLeftNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2206       fhPtHbpXEUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2207       outputContainer->Add(fhPtHbpXEUeLeftNeutral) ;
2208
2209       fhPtHbpXEUeRightNeutral  =
2210       new TH2F("hHbpXEUeNeutralRight","#xi = ln(1/#it{x}_{#it{E}}) with neutral UE right side of trigger",
2211                nptbins,ptmin,ptmax,200,0.,10.);
2212       fhPtHbpXEUeRightNeutral->SetYTitle("ln(1/#it{x}_{#it{E}})");
2213       fhPtHbpXEUeRightNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2214       outputContainer->Add(fhPtHbpXEUeRightNeutral) ;
2215
2216       fhZTUeLeftNeutral  =
2217       new TH2F("hZTUeNeutralLeft","#it{z}_{trigger h^{0}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE left side of trigger",
2218                nptbins,ptmin,ptmax,140,0.,2.);
2219       fhZTUeLeftNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2220       fhZTUeLeftNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2221       outputContainer->Add(fhZTUeLeftNeutral) ;
2222
2223       fhZTUeRightNeutral  =
2224       new TH2F("hZTUeNeutralRight","#it{z}_{trigger h^{0}} = #it{p}_{T Ueh^{0}} / #it{p}_{T trigger} with neutral UE right side of trigger",
2225                nptbins,ptmin,ptmax,200,0.,2.);
2226       fhZTUeRightNeutral->SetYTitle("#it{z}_{trigger Ueh^{0}}");
2227       fhZTUeRightNeutral->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2228       outputContainer->Add(fhZTUeRightNeutral) ;
2229
2230       fhPtHbpZTUeLeftNeutral  =
2231       new TH2F("hHbpZTUeNeutralLeft","#xi = ln(1/#it{z}_{T}) with neutral UE left side of trigger",
2232                nptbins,ptmin,ptmax,200,0.,10.);
2233       fhPtHbpZTUeLeftNeutral->SetYTitle("ln(1/#it{z}_{T})");
2234       fhPtHbpZTUeLeftNeutral->SetXTitle("#it{p}_{T trigger}");
2235       outputContainer->Add(fhPtHbpZTUeLeftNeutral) ;
2236
2237       fhPtHbpZTUeRightNeutral  =
2238       new TH2F("hHbpZTUeNeutralRight","#xi = ln(1/#it{z}_{T}) with neutral UE right side of trigger",
2239                nptbins,ptmin,ptmax,200,0.,10.);
2240       fhPtHbpZTUeRightNeutral->SetYTitle("ln(1/#it{z}_{T})");
2241       fhPtHbpZTUeRightNeutral->SetXTitle("#it{p}_{T trigger}");
2242       outputContainer->Add(fhPtHbpZTUeRightNeutral) ;
2243
2244     }
2245
2246   }//Correlation with neutral hadrons
2247
2248   //if data is MC, fill more histograms
2249   if(IsDataMC())
2250   {
2251     fh2phiLeadingParticle=new TH2F("h2phiLeadingParticle","#phi resolustion for trigger particles",nptbins,ptmin,ptmax,100,-1,1);
2252     fh2phiLeadingParticle->GetXaxis()->SetTitle("#it{p}_{T gen Leading} (GeV/#it{c})");
2253     fh2phiLeadingParticle->GetYaxis()->SetTitle("(#phi_{rec}-#phi_{gen})/#phi_{gen}");
2254
2255     fhMCPtLeading  = new TH1F ("hMCPtLeading","MC : #it{p}_{T} distribution of leading particles", nptbins,ptmin,ptmax);
2256     fhMCPtLeading->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2257
2258     fhMCPhiLeading  = new TH2F ("hMCPhiLeading","MC : #phi distribution of leading Particles",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2259     fhMCPhiLeading->SetYTitle("#phi (rad)");
2260
2261     fhMCEtaLeading  = new TH2F ("hMCEtaLeading","MC : #eta distribution of leading",nptbins,ptmin,ptmax, netabins,etamin,etamax);
2262     fhMCEtaLeading->SetYTitle("#eta ");
2263
2264
2265     fhMCEtaCharged  = new TH2F
2266     ("hMCEtaCharged","MC #eta_{h^{#pm}}  vs #it{p}_{T #pm}",
2267      nptbins,ptmin,ptmax,netabins,etamin,etamax);
2268     fhMCEtaCharged->SetYTitle("#eta_{h^{#pm}} (rad)");
2269     fhMCEtaCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
2270
2271     fhMCPhiCharged  = new TH2F
2272     ("hMCPhiCharged","#MC phi_{h^{#pm}}  vs #it{p}_{T #pm}",
2273      200,ptmin,ptmax,nphibins,phimin,phimax);
2274     fhMCPhiCharged->SetYTitle("MC #phi_{h^{#pm}} (rad)");
2275     fhMCPhiCharged->SetXTitle("#it{p}_{T #pm} (GeV/#it{c})");
2276
2277     fhMCDeltaPhiDeltaEtaCharged  = new TH2F
2278     ("hMCDeltaPhiDeltaEtaCharged","#MC phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
2279      140,-2.,5.,200,-2,2);
2280     fhMCDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi (rad)");
2281     fhMCDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
2282
2283     fhMCDeltaEtaCharged  = new TH2F
2284     ("hMCDeltaEtaCharged","MC #eta_{trigger} - #eta_{h^{#pm}} vs #it{p}_{T trigger} and #it{p}_{T assoc}",
2285      nptbins,ptmin,ptmax,200,-2,2);
2286     fhMCDeltaEtaCharged->SetYTitle("#Delta #eta");
2287     fhMCDeltaEtaCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2288
2289     fhMCDeltaPhiCharged  = new TH2F
2290     ("hMCDeltaPhiCharged","#phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}",
2291      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2292     fhMCDeltaPhiCharged->SetYTitle("#Delta #phi (rad)");
2293     fhMCDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2294
2295     fhMCDeltaPhiChargedPt  = new TH2F
2296     ("hMCDeltaPhiChargedPt","MC #phi_{trigger} - #phi_{#h^{#pm}} vs #it{p}_{T h^{#pm}}",
2297      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2298     fhMCDeltaPhiChargedPt->SetYTitle("#Delta #phi (rad)");
2299     fhMCDeltaPhiChargedPt->SetXTitle("#it{p}_{T h^{#pm}} (GeV/#it{c})");
2300
2301     fhMCPtXECharged  =
2302     new TH2F("hMCPtXECharged","#it{x}_{#it{E}} with charged hadrons",
2303              nptbins,ptmin,ptmax,200,0.,2.);
2304     fhMCPtXECharged->SetYTitle("#it{x}_{#it{E}}");
2305     fhMCPtXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2306
2307     fhMCPtXEUeCharged  =
2308     new TH2F("hMCPtXEUeCharged","#it{x}_{#it{E}} with charged hadrons, Underlying Event",
2309              nptbins,ptmin,ptmax,200,0.,2.);
2310     fhMCPtXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
2311     fhMCPtXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2312
2313     fhMCPtXEUeLeftCharged  =
2314     new TH2F("hMCPtXEUeChargedLeft","#it{x}_{#it{E}} with charged hadrons, with UE left side range of trigger particles",
2315              nptbins,ptmin,ptmax,200,0.,2.);
2316     fhMCPtXEUeLeftCharged->SetYTitle("#it{x}_{#it{E}}");
2317     fhMCPtXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2318
2319     fhMCPtXEUeRightCharged  =
2320     new TH2F("hMCPtXEUeChargedRight","#it{x}_{#it{E}} with charged hadrons, with UE left side range of trigger particles",
2321              nptbins,ptmin,ptmax,200,0.,2.);
2322     fhMCPtXEUeRightCharged->SetYTitle("#it{x}_{#it{E}}");
2323     fhMCPtXEUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2324
2325
2326     fhMCPtHbpXECharged  =
2327     new TH2F("hMCHbpXECharged","MC #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
2328              nptbins,ptmin,ptmax,200,0.,10.);
2329     fhMCPtHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2330     fhMCPtHbpXECharged->SetXTitle("#it{p}_{T trigger}");
2331
2332     fhMCPtHbpXEUeCharged =
2333     new TH2F("hMCPtHbpXEUeCharged","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, Underlying Event",
2334              nptbins,ptmin,ptmax,200,0.,10.);
2335     fhMCPtHbpXEUeCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2336     fhMCPtHbpXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2337
2338     fhMCPtHbpXEUeLeftCharged =
2339     new TH2F("hMCPtHbpXEUeChargedLeft","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, with UE left side range of trigger particles",
2340              nptbins,ptmin,ptmax,200,0.,10.);
2341     fhMCPtHbpXEUeLeftCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2342     fhMCPtHbpXEUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2343
2344     fhMCPtHbpXEUeRightCharged =
2345     new TH2F("hMCPtHbpXEUeChargedRight","#xi = ln(1/#it{x}_{#it{E}}) with charged hadrons, with UE right side range of trigger particles",
2346              nptbins,ptmin,ptmax,200,0.,10.);
2347     fhMCPtHbpXEUeRightCharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2348     fhMCPtHbpXEUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2349
2350
2351     fhMCUePart  =
2352     new TH1F("hMCUePart","MC UE particles distribution vs pt trig",
2353              nptbins,ptmin,ptmax);
2354     fhMCUePart->SetYTitle("#it{dN}^{ch}");
2355     fhMCUePart->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2356
2357     fhMCPtZTCharged  =
2358     new TH2F("hMCPtZTCharged","#it{z}_{T} with charged hadrons",
2359              nptbins,ptmin,ptmax,200,0.,2.);
2360     fhMCPtZTCharged->SetYTitle("#it{z}_{T}");
2361     fhMCPtZTCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2362
2363     fhMCPtZTUeCharged  =
2364     new TH2F("hMCPtZTUeCharged","#it{z}_{T} with charged hadrons, Underlying Event",
2365              nptbins,ptmin,ptmax,200,0.,2.);
2366     fhMCPtZTUeCharged->SetYTitle("#it{z}_{T}");
2367     fhMCPtZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2368
2369     fhMCPtZTUeLeftCharged  =
2370     new TH2F("hMCPtZTUeChargedLeft","#it{z}_{T} with charged hadrons, with UE left side range of trigger particles",
2371              nptbins,ptmin,ptmax,200,0.,2.);
2372     fhMCPtZTUeLeftCharged->SetYTitle("#it{z}_{T}");
2373     fhMCPtZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2374
2375     fhMCPtZTUeRightCharged  =
2376     new TH2F("hMCPtZTUeChargedRight","#it{z}_{T} with charged hadrons, with UE right side range of trigger particles",
2377              nptbins,ptmin,ptmax,200,0.,2.);
2378     fhMCPtZTUeRightCharged->SetYTitle("#it{z}_{T}");
2379     fhMCPtZTUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2380
2381     fhMCPtHbpZTCharged  =
2382     new TH2F("hMCHbpZTCharged","MC #xi = ln(1/#it{z}_{T}) with charged hadrons",
2383              nptbins,ptmin,ptmax,200,0.,10.);
2384     fhMCPtHbpZTCharged->SetYTitle("ln(1/#it{z}_{T})");
2385     fhMCPtHbpZTCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2386
2387     fhMCPtHbpZTUeCharged =
2388     new TH2F("hMCPtHbpZTUeCharged","#xi = ln(1/#it{z}_{T}) with charged hadrons, Underlying Event",
2389              nptbins,ptmin,ptmax,200,0.,10.);
2390     fhMCPtHbpZTUeCharged->SetYTitle("ln(1/#it{z}_{T})");
2391     fhMCPtHbpZTUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2392
2393     fhMCPtHbpZTUeLeftCharged =
2394     new TH2F("hMCPtHbpZTUeChargedLeft","#xi = ln(1/#it{z}_{T}) with charged hadrons, with UE left side range of trigger particles",
2395              nptbins,ptmin,ptmax,200,0.,10.);
2396     fhMCPtHbpZTUeLeftCharged->SetYTitle("ln(1/#it{z}_{T})");
2397     fhMCPtHbpZTUeLeftCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2398
2399     fhMCPtHbpZTUeRightCharged =
2400     new TH2F("hMCPtHbpZTUeChargedRight","#xi = ln(1/#it{z}_{T}) with charged hadrons, with UE right side range of trigger particles",
2401              nptbins,ptmin,ptmax,200,0.,10.);
2402     fhMCPtHbpZTUeRightCharged->SetYTitle("ln(1/#it{z}_{T})");
2403     fhMCPtHbpZTUeRightCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2404
2405     fhMCPtTrigPout  =
2406     new TH2F("hMCPtTrigPout","AOD MC Pout with triggers",
2407              nptbins,ptmin,ptmax,2*nptbins,-ptmax,ptmax);
2408     fhMCPtTrigPout->SetYTitle("#it{p}_{out} (GeV/#it{c})");
2409     fhMCPtTrigPout->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2410
2411     fhMCPtAssocDeltaPhi  =
2412     new TH2F("hMCPtAssocDeltaPhi","AOD MC delta phi with associated charged hadrons",
2413              nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2414     fhMCPtAssocDeltaPhi->SetYTitle("#Delta #phi (rad)");
2415     fhMCPtAssocDeltaPhi->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2416
2417     outputContainer->Add(fh2phiLeadingParticle);
2418     outputContainer->Add(fhMCPtLeading);
2419     outputContainer->Add(fhMCPhiLeading);
2420     outputContainer->Add(fhMCEtaLeading);
2421     outputContainer->Add(fhMCDeltaPhiDeltaEtaCharged);
2422     outputContainer->Add(fhMCPhiCharged) ;
2423     outputContainer->Add(fhMCEtaCharged) ;
2424     outputContainer->Add(fhMCDeltaEtaCharged) ;
2425     outputContainer->Add(fhMCDeltaPhiCharged) ;
2426
2427     outputContainer->Add(fhMCDeltaPhiChargedPt) ;
2428     outputContainer->Add(fhMCPtXECharged) ;
2429     outputContainer->Add(fhMCPtXEUeCharged) ;
2430     outputContainer->Add(fhMCPtXEUeLeftCharged) ;
2431     outputContainer->Add(fhMCPtXEUeRightCharged) ;
2432     outputContainer->Add(fhMCPtZTCharged) ;
2433     outputContainer->Add(fhMCPtZTUeCharged) ;
2434     outputContainer->Add(fhMCPtZTUeLeftCharged) ;
2435     outputContainer->Add(fhMCPtZTUeRightCharged) ;
2436     outputContainer->Add(fhMCPtHbpXECharged) ;
2437     outputContainer->Add(fhMCPtHbpXEUeCharged);
2438     outputContainer->Add(fhMCPtHbpXEUeLeftCharged);
2439     outputContainer->Add(fhMCPtHbpXEUeRightCharged);
2440     outputContainer->Add(fhMCUePart);
2441     outputContainer->Add(fhMCPtHbpZTCharged) ;
2442     outputContainer->Add(fhMCPtHbpZTUeCharged) ;
2443     outputContainer->Add(fhMCPtHbpZTUeLeftCharged) ;
2444     outputContainer->Add(fhMCPtHbpZTUeRightCharged) ;
2445     outputContainer->Add(fhMCPtTrigPout) ;
2446     outputContainer->Add(fhMCPtAssocDeltaPhi) ;
2447   } //for MC histogram
2448
2449   if(DoOwnMix())
2450   {
2451     //create event containers
2452
2453     if(!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForTracksExists()))
2454     {
2455       Int_t nvz = GetNZvertBin();
2456       Int_t nrp = GetNRPBin();
2457       Int_t nce = GetNCentrBin();
2458
2459       fListMixTrackEvents= new TList*[nvz*nrp*nce] ;
2460
2461       for( Int_t ice = 0 ; ice < nce ; ice++ )
2462       {
2463         for( Int_t ivz = 0 ; ivz < nvz ; ivz++ )
2464         {
2465           for( Int_t irp = 0 ; irp < nrp ; irp++ )
2466           {
2467             Int_t bin = GetEventMixBin(ice,ivz,irp); //ic*nvz*nrp+iz*nrp+irp;
2468
2469             //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
2470             //       ic,iz, irp, bin);
2471
2472             fListMixTrackEvents[bin] = new TList() ;
2473             fListMixTrackEvents[bin]->SetOwner(kFALSE);
2474           }
2475         }
2476       }
2477     }
2478
2479     fhPtLeadingMixed  = new TH1F ("hPtLeadingMixed","#it{p}_{T} distribution of leading particles, used for mixing", nptbins,ptmin,ptmax);
2480     fhPtLeadingMixed->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2481
2482     if(fCorrelVzBin)
2483     {
2484       fhPtLeadingMixedVzBin  = new TH2F ("hPtLeadingMixedVzBin","#it{p}_{T} distribution of leading particles, used for mixing", nptbins,ptmin,ptmax,GetNZvertBin(),0,GetNZvertBin());
2485       fhPtLeadingMixedVzBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2486       fhPtLeadingMixedVzBin->SetYTitle("#it{v}_{#it{z}} bin");
2487       outputContainer->Add(fhPtLeadingMixedVzBin);
2488     }
2489
2490     fhPtLeadingMixedBin  = new TH2F ("hPtLeadingMixedBin","#it{p}_{T} distribution of leading particles vs mixing bin", nptbins,ptmin,ptmax,nMixBins,0,nMixBins);
2491     fhPtLeadingMixedBin->SetXTitle("#it{p}_{T}^{trig} (GeV/#it{c})");
2492     fhPtLeadingMixedBin->SetYTitle("Bin");
2493
2494     fhPhiLeadingMixed  = new TH2F ("hPhiLeadingMixed","#phi distribution of leading Particles, used for mixing",nptbins,ptmin,ptmax, nphibins,phimin,phimax);
2495     fhPhiLeadingMixed->SetYTitle("#phi (rad)");
2496
2497     fhEtaLeadingMixed  = new TH2F ("hEtaLeadingMixed","#eta distribution of leading, used for mixing",nptbins,ptmin,ptmax, netabins,etamin,etamax);
2498     fhEtaLeadingMixed->SetYTitle("#eta ");
2499
2500     outputContainer->Add(fhPtLeadingMixed);
2501     outputContainer->Add(fhPtLeadingMixedBin);
2502     outputContainer->Add(fhPhiLeadingMixed);
2503     outputContainer->Add(fhEtaLeadingMixed);
2504
2505     // Fill the cluster pool only in isolation analysis or if requested
2506     if( ( OnlyIsolated()        ||  fFillNeutralEventMixPool) &&
2507         (!fUseMixStoredInReader || (fUseMixStoredInReader && !GetReader()->ListWithMixedEventsForCaloExists())))
2508     {
2509       Int_t nvz = GetNZvertBin();
2510       Int_t nrp = GetNRPBin();
2511       Int_t nce = GetNCentrBin();
2512
2513       fListMixCaloEvents= new TList*[nvz*nrp*nce] ;
2514
2515       for( Int_t ice = 0 ; ice < nce ; ice++ )
2516       {
2517         for( Int_t ivz = 0 ; ivz < nvz ; ivz++ )
2518         {
2519           for( Int_t irp = 0 ; irp < nrp ; irp++ )
2520           {
2521             Int_t bin = GetEventMixBin(ice,ivz,irp); //ic*nvz*nrp+iz*nrp+irp;
2522
2523             //printf("GetCreateOutputObjects - Bins : cent %d, vz %d, RP %d, event %d/%d\n",
2524             //       ic,iz, irp, bin);
2525
2526             fListMixCaloEvents[bin] = new TList() ;
2527             fListMixCaloEvents[bin]->SetOwner(kFALSE);
2528           }
2529         }
2530       }
2531     }
2532
2533     //Init the list in the reader if not done previously
2534     if(fUseMixStoredInReader)
2535     {
2536       if( !GetReader()->ListWithMixedEventsForTracksExists() )
2537         GetReader()->SetListWithMixedEventsForTracks(fListMixTrackEvents);
2538
2539       if( !GetReader()->ListWithMixedEventsForCaloExists()   )
2540         GetReader()->SetListWithMixedEventsForCalo  (fListMixCaloEvents );
2541     }
2542
2543     fhEventBin=new TH1I("hEventBin","Number of real events per bin(cen,vz,rp)",
2544                         GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2545                         GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2546     fhEventBin->SetXTitle("bin");
2547     outputContainer->Add(fhEventBin) ;
2548
2549     fhEventMixBin=new TH1I("hEventMixBin","Number of events  per bin(cen,vz,rp)",
2550                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1,0,
2551                            GetNCentrBin()*GetNZvertBin()*GetNRPBin()+1) ;
2552     fhEventMixBin->SetXTitle("bin");
2553     outputContainer->Add(fhEventMixBin) ;
2554
2555     fhNtracksMB=new TH1F("hNtracksMBEvent","Number of tracks w/ event trigger kMB",2000,0,2000);
2556     outputContainer->Add(fhNtracksMB);
2557
2558     if(fFillNeutralEventMixPool || OnlyIsolated())
2559     {
2560       fhNclustersMB=new TH1F("hNclustersMBEvent","Number of clusters w/ event trigger kMB",2000,0,2000);
2561       outputContainer->Add(fhNclustersMB);
2562     }
2563
2564     fhMixDeltaPhiCharged  = new TH2F
2565     ("hMixDeltaPhiCharged","Mixed event : #phi_{trigger} - #phi_{h^{#pm}} vs #it{p}_{T trigger}",
2566      nptbins,ptmin,ptmax,ndeltaphibins ,deltaphimin,deltaphimax);
2567     fhMixDeltaPhiCharged->SetYTitle("#Delta #phi (rad)");
2568     fhMixDeltaPhiCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2569     outputContainer->Add(fhMixDeltaPhiCharged);
2570
2571     fhMixDeltaPhiDeltaEtaCharged  = new TH2F
2572     ("hMixDeltaPhiDeltaEtaCharged","Mixed event : #phi_{trigger} - #phi_{h^{#pm}} vs #eta_{trigger} - #eta_{h^{#pm}}",
2573      ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
2574     fhMixDeltaPhiDeltaEtaCharged->SetXTitle("#Delta #phi (rad)");
2575     fhMixDeltaPhiDeltaEtaCharged->SetYTitle("#Delta #eta");
2576     outputContainer->Add(fhMixDeltaPhiDeltaEtaCharged);
2577
2578     fhMixXECharged  =
2579     new TH2F("hMixXECharged","Mixed event : #it{x}_{#it{E}} for charged tracks",
2580              nptbins,ptmin,ptmax,200,0.,2.);
2581     fhMixXECharged->SetYTitle("#it{x}_{#it{E}}");
2582     fhMixXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2583     outputContainer->Add(fhMixXECharged);
2584
2585     fhMixXEUeCharged  =
2586     new TH2F("hMixXEUeCharged","Mixed event : #it{x}_{#it{E}} for charged tracks in Ue region",
2587              nptbins,ptmin,ptmax,200,0.,2.);
2588     fhMixXEUeCharged->SetYTitle("#it{x}_{#it{E}}");
2589     fhMixXEUeCharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2590     outputContainer->Add(fhMixXEUeCharged);
2591
2592     fhMixHbpXECharged  =
2593     new TH2F("hMixHbpXECharged","mixed event : #xi = ln(1/#it{x}_{#it{E}}) with charged hadrons",
2594              nptbins,ptmin,ptmax,200,0.,10.);
2595     fhMixHbpXECharged->SetYTitle("ln(1/#it{x}_{#it{E}})");
2596     fhMixHbpXECharged->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2597     outputContainer->Add(fhMixHbpXECharged);
2598
2599     fhMixDeltaPhiChargedAssocPtBin         = new TH2F*[fNAssocPtBins*nz];
2600     fhMixDeltaPhiChargedAssocPtBinDEta08   = new TH2F*[fNAssocPtBins*nz];
2601     fhMixDeltaPhiChargedAssocPtBinDEta0    = new TH2F*[fNAssocPtBins*nz];
2602     fhMixDeltaPhiDeltaEtaChargedAssocPtBin = new TH2F*[fNAssocPtBins*nz];
2603
2604     for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
2605     {
2606       for(Int_t z = 0 ; z < nz ; z++)
2607       {
2608         Int_t bin = i*nz+z;
2609
2610         if(fCorrelVzBin)
2611         {
2612           sz = Form("_vz%d",z);
2613           tz = Form(", #it{v}_{#it{z}} bin %d",z);
2614         }
2615
2616         //printf("MIX : iAssoc %d, Vz %d, bin %d - sz %s, tz %s \n",i,z,bin,sz.Data(),tz.Data());
2617
2618         fhMixDeltaPhiChargedAssocPtBin[bin] = new TH2F(Form("hMixDeltaPhiChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2619                                                      Form("Mixed event #Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2620                                                      nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
2621         fhMixDeltaPhiChargedAssocPtBin[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2622         fhMixDeltaPhiChargedAssocPtBin[bin]->SetYTitle("#Delta #phi (rad)");
2623
2624         fhMixDeltaPhiChargedAssocPtBinDEta08[bin] = new TH2F(Form("hMixDeltaPhiDeltaEta0.8ChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2625                                                            Form("Mixed event #Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta > 0.8", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2626                                                            nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
2627         fhMixDeltaPhiChargedAssocPtBinDEta08[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2628         fhMixDeltaPhiChargedAssocPtBinDEta08[bin]->SetYTitle("#Delta #phi (rad)");
2629
2630         fhMixDeltaPhiChargedAssocPtBinDEta0[bin] = new TH2F(Form("hMixDeltaPhiDeltaEta0ChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2631                                                              Form("Mixed event #Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s, for #Delta #eta = 0", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2632                                                              nptbins, ptmin, ptmax,  ndeltaphibins ,deltaphimin,deltaphimax);
2633         fhMixDeltaPhiChargedAssocPtBinDEta0[bin]->SetXTitle("#it{p}_{T trigger} (GeV/#it{c})");
2634         fhMixDeltaPhiChargedAssocPtBinDEta0[bin]->SetYTitle("#Delta #phi (rad)");
2635
2636         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin] = new TH2F(Form("hMixDeltaPhiDeltaEtaChargedAssocPtBin%2.1f_%2.1f%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],sz.Data()),
2637                                                              Form("Mixed event #Delta #phi vs #it{p}_{T trigger} for associated #it{p}_{T} bin [%2.1f,%2.1f]%s", fAssocPtBinLimit[i], fAssocPtBinLimit[i+1],tz.Data()),
2638                                                              ndeltaphibins ,deltaphimin,deltaphimax,ndeltaetabins ,deltaetamin,deltaetamax);
2639         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->SetXTitle("#Delta #phi (rad)");
2640         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->SetYTitle("#Delta #eta");
2641
2642         outputContainer->Add(fhMixDeltaPhiChargedAssocPtBin[bin]);
2643         outputContainer->Add(fhMixDeltaPhiChargedAssocPtBinDEta08[bin]);
2644         outputContainer->Add(fhMixDeltaPhiChargedAssocPtBinDEta0[bin]);
2645         outputContainer->Add(fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]);
2646       }
2647     }
2648   }
2649
2650   return outputContainer;
2651
2652 }
2653
2654 //_________________________________________________________________________________________________
2655 Bool_t AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum(const AliAODPWG4Particle* trigger,
2656                                                              TLorentzVector & mom1,
2657                                                              TLorentzVector & mom2)
2658 {
2659   // Get the momentum of the pi0/eta assigned decay photons
2660   // In case of pi0/eta trigger, we may want to check their decay correlation,
2661   // get their decay children
2662
2663   Int_t indexPhoton1 = trigger->GetCaloLabel(0);
2664   Int_t indexPhoton2 = trigger->GetCaloLabel(1);
2665   Float_t ptTrig     = trigger->Pt();
2666
2667   if(indexPhoton1!=-1 || indexPhoton2!=-1) return kFALSE;
2668
2669   if(GetDebug() > 1)
2670     printf("AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum() - indexPhoton1 = %d, indexPhoton2 = %d \n", indexPhoton1, indexPhoton2);
2671
2672   TObjArray * clusters  = 0x0 ;
2673   if(trigger->GetDetector()=="EMCAL") clusters = GetEMCALClusters() ;
2674   else                                clusters = GetPHOSClusters()  ;
2675
2676   for(Int_t iclus = 0; iclus < clusters->GetEntriesFast(); iclus++)
2677   {
2678     AliVCluster * photon =  (AliVCluster*) (clusters->At(iclus));
2679     if(photon->GetID()==indexPhoton1)
2680     {
2681       photon->GetMomentum(mom1,GetVertex(0)) ;
2682       if(ptTrig) fhPtPi0DecayRatio->Fill(ptTrig, mom1.Pt()/ptTrig);
2683     }
2684     if(photon->GetID()==indexPhoton2)
2685     {
2686       photon->GetMomentum(mom1,GetVertex(0)) ;
2687       if(ptTrig > 0) fhPtPi0DecayRatio->Fill(ptTrig, mom2.Pt()/ptTrig);
2688     }
2689
2690     if(GetDebug() > 1)printf("AliAnaParticleHadronCorrelation::GetDecayPhotonMomentum() - Photon1 = %f, Photon2 = %f \n", mom1.Pt(), mom2.Pt());
2691
2692   } //cluster loop
2693
2694   return kTRUE;
2695
2696 }
2697
2698 //_____________________________________________________________
2699 Int_t AliAnaParticleHadronCorrelation::GetMCTagHistogramIndex(Int_t mcTag)
2700 {
2701   // Index of MC histograms depending on MC origin
2702
2703   if     ( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPrompt) ||
2704            GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCFragmentation)) return 0;
2705   else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0))           return 1;
2706   else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCPi0Decay))      return 2;
2707   else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCEtaDecay))      return 3;
2708   else if( GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCOtherDecay))    return 4;
2709   else if(!GetMCAnalysisUtils()->CheckTagBit(mcTag,AliMCAnalysisUtils::kMCElectron))      return 5;
2710   else                                                                                    return 6;
2711
2712 }
2713
2714 //____________________________________________________
2715 void AliAnaParticleHadronCorrelation::InitParameters()
2716 {
2717
2718   //Initialize the parameters of the analysis.
2719   SetInputAODName("Particle");
2720   SetAODObjArrayName("Hadrons");
2721   AddToHistogramsName("AnaHadronCorr_");
2722
2723   SetPtCutRange(0.,300);
2724   fDeltaPhiMinCut       = 1.5 ;
2725   fDeltaPhiMaxCut       = 4.5 ;
2726   fSelectIsolated       = kFALSE;
2727   fMakeSeveralUE        = kFALSE;
2728   fUeDeltaPhiMinCut     = 1. ;
2729   fUeDeltaPhiMaxCut     = 1.5 ;
2730
2731   fNeutralCorr          = kFALSE ;
2732   fPi0Trigger           = kFALSE ;
2733   fDecayTrigger         = kFALSE ;
2734   fHMPIDCorrelation     = kFALSE ;
2735
2736   fMakeAbsoluteLeading  = kTRUE;
2737   fMakeNearSideLeading  = kFALSE;
2738
2739   fNAssocPtBins         = 9   ;
2740   fAssocPtBinLimit[0]   = 0.2 ;
2741   fAssocPtBinLimit[1]   = 0.5 ;
2742   fAssocPtBinLimit[2]   = 1.0 ;
2743   fAssocPtBinLimit[3]   = 2.0 ;
2744   fAssocPtBinLimit[4]   = 3.0 ;
2745   fAssocPtBinLimit[5]   = 4.0 ;
2746   fAssocPtBinLimit[6]   = 5.0 ;
2747   fAssocPtBinLimit[7]   = 6.0 ;
2748   fAssocPtBinLimit[8]   = 7.0 ;
2749   fAssocPtBinLimit[9]   = 8.0 ;
2750   fAssocPtBinLimit[10]  = 9.0 ;
2751   fAssocPtBinLimit[11]  = 10.0 ;
2752   fAssocPtBinLimit[12]  = 12.0 ;
2753   fAssocPtBinLimit[13]  = 14.0 ;
2754   fAssocPtBinLimit[14]  = 16.0 ;
2755   fAssocPtBinLimit[15]  = 20.0 ;
2756   fAssocPtBinLimit[16]  = 30.0 ;
2757   fAssocPtBinLimit[17]  = 40.0 ;
2758   fAssocPtBinLimit[18]  = 50.0 ;
2759   fAssocPtBinLimit[19]  = 200.0 ;
2760
2761   fUseMixStoredInReader = kTRUE;
2762
2763   fM02MinCut   = -1 ;
2764   fM02MaxCut   = -1 ;
2765
2766   fSelectLeadingHadronAngle = kFALSE;
2767   fMinLeadHadPhi = 150*TMath::DegToRad();
2768   fMaxLeadHadPhi = 210*TMath::DegToRad();
2769
2770   fMinLeadHadPt  = 1;
2771   fMaxLeadHadPt  = 100;
2772
2773 }
2774
2775 //__________________________________________________________
2776 void  AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD()
2777 {
2778   //Particle-Hadron Correlation Analysis, fill AODs
2779
2780   if(!GetInputAODBranch())
2781   {
2782     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data());
2783     abort();
2784   }
2785
2786   if(strcmp(GetInputAODBranch()->GetClass()->GetName(), "AliAODPWG4ParticleCorrelation"))
2787   {
2788     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());
2789     abort();
2790   }
2791
2792   if(GetDebug() > 1)
2793   {
2794     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - Begin hadron correlation analysis, fill AODs \n");
2795     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
2796     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In CTS aod entries %d\n",   GetCTSTracks()    ->GetEntriesFast());
2797     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In EMCAL aod entries %d\n", GetEMCALClusters()->GetEntriesFast());
2798     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - In PHOS aod entries %d\n",  GetPHOSClusters() ->GetEntriesFast());
2799   }
2800
2801   //Loop on stored AOD particles, find leading trigger
2802   Double_t ptTrig      = fMinTriggerPt ;
2803   fLeadingTriggerIndex = -1 ;
2804   Int_t    naod        = GetInputAODBranch()->GetEntriesFast() ;
2805   for(Int_t iaod = 0; iaod < naod ; iaod++)
2806   {
2807     AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
2808
2809     // Vertex cut in case of mixing
2810     Int_t check = CheckMixedEventVertex(particle->GetCaloLabel(0), particle->GetTrackLabel(0));
2811     if(check ==  0) continue;
2812     if(check == -1) return;
2813
2814     // find the leading particles with highest momentum
2815     if (particle->Pt() > ptTrig)
2816     {
2817       ptTrig               = particle->Pt() ;
2818       fLeadingTriggerIndex = iaod ;
2819     }
2820   }// finish search of leading trigger particle
2821
2822
2823   //Do correlation with leading particle
2824   if(fLeadingTriggerIndex < 0) return ;
2825
2826   AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
2827
2828   //check if the particle is isolated or if we want to take the isolation into account
2829   if(OnlyIsolated() && !particle->IsIsolated()) return;
2830
2831   //Make correlation with charged hadrons
2832   Bool_t okcharged = kTRUE;
2833   Bool_t okneutral = kTRUE;
2834   if(GetReader()->IsCTSSwitchedOn() )
2835     okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kFALSE);
2836
2837   TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
2838   if(fNeutralCorr && pi0list && pi0list->GetEntriesFast() > 0)
2839     okneutral = MakeNeutralCorrelation(particle, pi0list,kFALSE);
2840
2841   if( GetDebug() > 1 ) printf("Charged correlation OK %d; Neutral correlation OK %d",okcharged,okneutral);
2842
2843   if( GetDebug() > 1 ) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillAOD() - End fill AODs \n");
2844
2845 }
2846
2847 //_________________________________________________________________
2848 void  AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms()
2849 {
2850   //Particle-Hadron Correlation Analysis, fill histograms
2851
2852   if(!GetInputAODBranch())
2853   {
2854     AliFatal(Form("No input particles in AOD with name branch < %s >, STOP \n",GetInputAODName().Data()));
2855     return ; // coverity
2856   }
2857
2858   if(GetDebug() > 1)
2859   {
2860     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Begin hadron correlation analysis, fill histograms \n");
2861     printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - In particle branch aod entries %d\n", GetInputAODBranch()->GetEntriesFast());
2862   }
2863
2864   //---------------------------------------------------
2865   // Loop on stored AOD particles, find leading trigger
2866
2867   Double_t ptTrig    = fMinTriggerPt;
2868   if(fLeadingTriggerIndex < 0)
2869   {
2870     //Search leading if not done before
2871     Int_t    naod      = GetInputAODBranch()->GetEntriesFast() ;
2872     for(Int_t iaod = 0; iaod < naod ; iaod++)
2873     {    //loop on input trigger AOD file
2874       AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(iaod));
2875
2876       // Vertex cut in case of mixing
2877       Int_t check = CheckMixedEventVertex(particle->GetCaloLabel(0), particle->GetTrackLabel(0));
2878       if(check ==  0) continue;
2879       if(check == -1) return;
2880
2881       //check if the particle is isolated or if we want to take the isolation into account
2882       if(OnlyIsolated() && !particle->IsIsolated()) continue;
2883
2884       //find the leading particles with highest momentum
2885       if (particle->Pt() > ptTrig)
2886       {
2887         ptTrig               = particle->Pt() ;
2888         fLeadingTriggerIndex = iaod ;
2889       }
2890
2891     }// Finish search of leading trigger particle
2892   }// Search leading if not done before
2893
2894   //-----------------------------------------
2895   // Fill Leading trigger related histograms
2896   //-----------------------------------------
2897   if(fLeadingTriggerIndex < 0 ) return ;
2898
2899   AliAODPWG4ParticleCorrelation* particle =  (AliAODPWG4ParticleCorrelation*) (GetInputAODBranch()->At(fLeadingTriggerIndex));
2900
2901   Float_t pt = particle->Pt();
2902   fhPtTriggerInput->Fill(pt);
2903
2904   // check if it was a calorimeter cluster and if the SS cut was requested, if so, apply it
2905   Int_t clID1  = particle->GetCaloLabel(0) ;
2906   Int_t clID2  = particle->GetCaloLabel(1) ; // for photon clusters should not be set.
2907   //printf("Leading for for %s: id1 %d, id2 %d, min %f, max %f, det %s\n",
2908   //       GetInputAODName().Data(),clID1,clID2,fM02MinCut,fM02MaxCut,(particle->GetDetector()).Data());
2909
2910
2911   // if requested, do the rejection of clusters out of a shower shape window
2912   // not needed if already done at the particle identification level, but for isolation
2913   // studies, it is preferred not to remove so we do it here
2914   if(clID1 > 0 && clID2 < 0 && fM02MaxCut > 0 && fM02MinCut > 0)
2915   {
2916     Int_t iclus = -1;
2917     TObjArray* clusters = 0x0;
2918     if     (particle->GetDetector() == "EMCAL") clusters = GetEMCALClusters();
2919     else if(particle->GetDetector() == "PHOS" ) clusters = GetPHOSClusters();
2920
2921     if(clusters)
2922     {
2923       AliVCluster *cluster = FindCluster(clusters,clID1,iclus);
2924       Float_t m02 = cluster->GetM02();
2925       if(m02 > fM02MaxCut || m02 < fM02MinCut) return ;
2926     }
2927
2928     fhPtTriggerSSCut->Fill(pt);
2929   }
2930
2931   // Check if the particle is isolated or if we want to take the isolation into account
2932   if(OnlyIsolated())
2933   {
2934     if( !particle->IsIsolated() ) return;
2935     fhPtTriggerIsoCut->Fill(pt);
2936   }
2937
2938   // Check if trigger is in fiducial region
2939   if(IsFiducialCutOn())
2940   {
2941     Bool_t in = GetFiducialCut()->IsInFiducialCut(*particle->Momentum(),particle->GetDetector()) ;
2942     if(! in ) return ;
2943   }
2944
2945   fhPtTriggerFidCut->Fill(pt);
2946
2947   // Make correlation with charged hadrons
2948   Bool_t okcharged = kTRUE;
2949   Bool_t okneutral = kTRUE;
2950   if(GetReader()->IsCTSSwitchedOn() )
2951   {
2952     okcharged = MakeChargedCorrelation(particle, GetCTSTracks(),kTRUE);
2953     if(IsDataMC())
2954     {
2955       MakeMCChargedCorrelation(particle);
2956     }
2957   }
2958
2959   TObjArray * pi0list = (TObjArray*) GetAODBranch(fPi0AODBranchName); //For the future, foresee more possible pi0 lists
2960   if(fNeutralCorr && pi0list)
2961   {
2962     if(pi0list->GetEntriesFast() > 0)
2963       okneutral = MakeNeutralCorrelation(particle, pi0list,kTRUE);
2964   }
2965
2966   // If the correlation or the finding of the leading did not succeed.
2967   if(!okcharged || !okneutral) return ;
2968
2969   // Fill leading particle histogram if correlation went well and
2970   // no problem was found, like not absolute leading, or bad vertex in mixing.
2971
2972   // pT of the leading, vs leading origin if MC
2973   fhPtLeading->Fill(pt);
2974
2975   if(IsDataMC())
2976   {
2977     Int_t mcIndex = GetMCTagHistogramIndex(particle->GetTag());
2978     fhPtLeadingMC[mcIndex]->Fill(pt);
2979   }
2980
2981   // Acceptance of the leading
2982   Float_t phi = particle->Phi();
2983   if( phi<0 ) phi+=TMath::TwoPi();
2984   fhPhiLeading->Fill(pt, phi);
2985
2986   fhEtaLeading->Fill(pt, particle->Eta());
2987   //printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - Leading particle : pt %f, eta %f, phi %f\n",particle->Pt(),particle->Eta(),phi);
2988
2989   //----------------------------------
2990   // Leading particle pT vs event bins
2991
2992   fhPtLeadingBin->Fill(pt,GetEventMixBin());
2993   if(fCorrelVzBin) fhPtLeadingVzBin->Fill(pt,GetEventVzBin());
2994
2995   Float_t cen = GetEventCentrality();
2996   Float_t ep  = GetEventPlaneAngle();
2997
2998   fhPtLeadingCentrality        ->Fill(pt,cen);
2999   fhPtLeadingEventPlane        ->Fill(pt,ep);
3000   fhLeadingEventPlaneCentrality->Fill(cen,ep);
3001
3002   //----------------------------------
3003   // Leading particle pT vs pile-up
3004
3005   if(fFillPileUpHistograms)
3006   {
3007     Int_t vtxBC = GetReader()->GetVertexBC();
3008     if(vtxBC == 0 || vtxBC==AliVTrack::kTOFBCNA)     fhPtLeadingVtxBC0->Fill(pt);
3009
3010     if(GetReader()->IsPileUpFromSPD())               fhPtLeadingPileUp[0]->Fill(pt);
3011     if(GetReader()->IsPileUpFromEMCal())             fhPtLeadingPileUp[1]->Fill(pt);
3012     if(GetReader()->IsPileUpFromSPDOrEMCal())        fhPtLeadingPileUp[2]->Fill(pt);
3013     if(GetReader()->IsPileUpFromSPDAndEMCal())       fhPtLeadingPileUp[3]->Fill(pt);
3014     if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhPtLeadingPileUp[4]->Fill(pt);
3015     if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhPtLeadingPileUp[5]->Fill(pt);
3016     if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhPtLeadingPileUp[6]->Fill(pt);
3017   }
3018
3019   //Reinit for next event
3020   fLeadingTriggerIndex = -1;
3021
3022   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeAnalysisFillHistograms() - End fill histograms \n");
3023 }
3024
3025 //___________________________________________________________________________________________________________
3026 Bool_t  AliAnaParticleHadronCorrelation::MakeChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle,
3027                                                                 const TObjArray* pl, Bool_t bFillHisto)
3028 {
3029   // Charged Hadron Correlation Analysis
3030   if(GetDebug() > 1)
3031     printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Make trigger particle - charged hadron correlation \n");
3032
3033   Float_t phiTrig = aodParticle->Phi();
3034   Float_t etaTrig = aodParticle->Eta();
3035   Float_t ptTrig  = aodParticle->Pt();
3036   Bool_t   decay  = aodParticle->IsTagged();
3037   Int_t    mcTag  = aodParticle->GetTag();
3038   Double_t bz     = GetReader()->GetInputEvent()->GetMagneticField();
3039
3040   Float_t pt       = -100. ;
3041   Float_t zT       = -100. ;
3042   Float_t xE       = -100. ;
3043   Float_t hbpXE    = -100. ;
3044   Float_t hbpZT    = -100. ;
3045   Float_t phi      = -100. ;
3046   Float_t eta      = -100. ;
3047   Float_t pout     = -100. ;
3048   Float_t deltaPhi = -100. ;
3049   Float_t ptLeadHad  = -100 ;
3050   Float_t phiLeadHad = -100 ;
3051
3052   TVector3 p3;
3053   TLorentzVector photonMom ;
3054   TObjArray * reftracks = 0x0;
3055   Int_t nrefs           = 0;
3056   Int_t nTracks         = GetCTSTracks()->GetEntriesFast() ;
3057
3058   // Mixed event settings
3059   Int_t evtIndex11   = -1 ; // cluster trigger or pi0 trigger
3060   Int_t evtIndex12   = -1 ; // pi0 trigger
3061   Int_t evtIndex13   = -1 ; // charged trigger
3062
3063   Double_t v[3]      = {0,0,0}; //vertex ;
3064   GetReader()->GetVertex(v);
3065
3066   if (GetMixedEvent())
3067   {
3068     evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
3069     evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
3070     evtIndex13 = GetMixedEvent()->EventIndex(aodParticle->GetTrackLabel(0)) ;
3071   }
3072
3073   // In case of pi0/eta trigger, we may want to check their decay correlation,
3074   // get their decay children
3075   TLorentzVector decayMom1;
3076   TLorentzVector decayMom2;
3077   Bool_t decayFound = kFALSE;
3078   if(fPi0Trigger && bFillHisto) decayFound = GetDecayPhotonMomentum(aodParticle,decayMom1, decayMom2);
3079
3080   //-----------------------------------------------------------------------
3081   // Track loop, select tracks with good pt, phi and fill AODs or histograms
3082   //-----------------------------------------------------------------------
3083
3084   // select events where the leading particle in the opposite hemisphere to the trigger particle is
3085   // is in a window centered at 180 from the trigger
3086   // Find the leading hadron
3087   if(fSelectLeadingHadronAngle)
3088   {
3089     for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ )
3090     {
3091       AliVTrack * track = (AliVTrack *) (pl->At(ipr)) ;
3092       Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
3093       p3.SetXYZ(mom[0],mom[1],mom[2]);
3094       pt   = p3.Pt();
3095       phi  = p3.Phi() ;
3096
3097       if(pt > ptLeadHad && TMath::Abs(phi-phiTrig) > TMath::PiOver2())
3098       {
3099         ptLeadHad = pt ;
3100         phiLeadHad = phi;
3101       }
3102
3103     }// track loop
3104
3105     if( ptLeadHad < fMinLeadHadPt ||
3106         ptLeadHad > fMaxLeadHadPt ) return kFALSE;
3107
3108     if( TMath::Abs(phiLeadHad-phiTrig) < fMinLeadHadPhi ||
3109         TMath::Abs(phiLeadHad-phiTrig) > fMaxLeadHadPhi ) return kFALSE;
3110
3111   }// select leading hadron
3112
3113   for(Int_t ipr = 0;ipr < pl->GetEntriesFast() ; ipr ++ )
3114   {
3115     AliVTrack * track = (AliVTrack *) (pl->At(ipr)) ;
3116
3117     Double_t mom[3] = {track->Px(),track->Py(),track->Pz()};
3118     p3.SetXYZ(mom[0],mom[1],mom[2]);
3119     pt   = p3.Pt();
3120     eta  = p3.Eta();
3121     phi  = p3.Phi() ;
3122     if(phi < 0) phi+=TMath::TwoPi();
3123
3124     //Select only hadrons in pt range
3125     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
3126
3127     //remove trigger itself for correlation when use charged triggers
3128     if( track->GetID() == aodParticle->GetTrackLabel(0) || track->GetID() == aodParticle->GetTrackLabel(1) ||
3129         track->GetID() == aodParticle->GetTrackLabel(2) || track->GetID() == aodParticle->GetTrackLabel(3)   )
3130       continue ;
3131
3132     //jump out this event if near side associated particle pt larger than trigger
3133     if (fMakeNearSideLeading)
3134     {
3135       if(pt > ptTrig && TMath::Abs(phi-phiTrig) < TMath::PiOver2())  return kFALSE;
3136     }
3137     //jump out this event if there is any other particle with pt larger than trigger
3138     else if(fMakeAbsoluteLeading)
3139     {
3140       if(pt > ptTrig)  return kFALSE;
3141     }
3142
3143     //Only for mixed event
3144     Int_t evtIndex2 = 0 ;
3145     if (GetMixedEvent())
3146     {
3147       evtIndex2 = GetMixedEvent()->EventIndex(track->GetID()) ;
3148       if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 || evtIndex13 == evtIndex2 ) // photon and track from different events
3149         continue ;
3150       //vertex cut
3151       if (TMath::Abs(GetVertex(evtIndex2)[2]) > GetZvertexCut())
3152         return kFALSE;
3153     }
3154
3155     // Fill Histograms
3156     if(bFillHisto)
3157     {
3158
3159       if(GetDebug() > 2 )
3160         printf("AliAnaParticleHadronCorrelation::MakeChargedCorrelation() - Selected charge for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
3161
3162       // Set the pt associated bin for the defined bins
3163       Int_t assocBin   = -1;
3164       for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
3165       {
3166         if(pt > fAssocPtBinLimit[i] && pt < fAssocPtBinLimit[i+1]) assocBin= i;
3167       }
3168
3169       // Assign to the histogram array a bin corresponding to a combination of pTa and vz bins
3170       Int_t nz = 1;
3171       Int_t vz = 0;
3172
3173       if(fCorrelVzBin)
3174       {
3175         nz = GetNZvertBin();
3176         vz = GetEventVzBin();
3177       }
3178
3179       Int_t bin = assocBin*nz+vz;
3180
3181       //printf("assoc Bin = %d, vZ bin  = %d, bin = %d \n", assocBin,GetEventVzBin(),bin);
3182
3183       ULong_t status = track->GetStatus();
3184       Bool_t okTOF = ( (status & AliVTrack::kTOFout) == AliVTrack::kTOFout ) ;
3185       //Double32_t tof = track->GetTOFsignal()*1e-3;
3186       Int_t trackBC = track->GetTOFBunchCrossing(bz);
3187
3188       Int_t outTOF = -1;
3189       if     (okTOF && trackBC!=0) outTOF = 1;
3190       else if(okTOF && trackBC==0) outTOF = 0;
3191
3192       // Azimuthal Angle
3193       // calculate deltaPhi for later, shift when needed
3194       FillChargedAngularCorrelationHistograms(pt,  ptTrig,  bin, phi, phiTrig,  deltaPhi,
3195                                               eta, etaTrig, decay, track->GetHMPIDsignal(),outTOF,nTracks,mcTag);
3196
3197       // Imbalance zT/xE/pOut
3198       zT = pt/ptTrig ;
3199       if(zT > 0 ) hbpZT = TMath::Log(1./zT);
3200       else        hbpZT =-100;
3201
3202       xE   =-pt/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
3203       //if(xE <0.)xE =-xE;
3204       if(xE > 0 ) hbpXE = TMath::Log(1./xE);
3205       else        hbpXE =-100;
3206
3207       pout = pt*TMath::Sin(deltaPhi) ;
3208
3209       //delta phi cut for momentum imbalance correlation
3210       if  ( (deltaPhi > fDeltaPhiMinCut)   && (deltaPhi < fDeltaPhiMaxCut)   )
3211       {
3212
3213         FillChargedMomentumImbalanceHistograms(ptTrig, pt, xE, hbpXE, zT, hbpZT, pout, deltaPhi,
3214                                                nTracks, track->Charge(), bin, decay,outTOF,mcTag);
3215
3216       }
3217
3218       if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )
3219       { //UE study
3220
3221         FillChargedUnderlyingEventHistograms(ptTrig, pt, deltaPhi, nTracks,outTOF);
3222
3223         fhUePart->Fill(ptTrig);
3224
3225       }
3226
3227       if(fPi0Trigger && decayFound)
3228         FillDecayPhotonCorrelationHistograms(pt, phi, decayMom1,decayMom2, kTRUE) ;
3229
3230       //several UE calculation
3231       if(fMakeSeveralUE) FillChargedUnderlyingEventSidesHistograms(ptTrig,pt,deltaPhi);
3232
3233     } //Fill histogram
3234     else
3235     {
3236       nrefs++;
3237       if(nrefs==1)
3238       {
3239         reftracks = new TObjArray(0);
3240         TString trackname = Form("%s+Tracks", GetAODObjArrayName().Data());
3241         reftracks->SetName(trackname.Data());
3242         reftracks->SetOwner(kFALSE);
3243       }
3244
3245       reftracks->Add(track);
3246
3247     }//aod particle loop
3248   }// track loop
3249
3250   //Fill AOD with reference tracks, if not filling histograms
3251   if(!bFillHisto && reftracks)
3252   {
3253     aodParticle->AddObjArray(reftracks);
3254   }
3255
3256   //Own mixed event, add event and remove previous or fill the mixed histograms
3257   if(DoOwnMix() && bFillHisto)
3258   {
3259       MakeChargedMixCorrelation(aodParticle);
3260   }
3261
3262   return kTRUE;
3263
3264 }
3265
3266
3267 //_________________________________________________________________________________________________________
3268 void AliAnaParticleHadronCorrelation::MakeChargedMixCorrelation(AliAODPWG4ParticleCorrelation *aodParticle)
3269 {
3270   // Mix current trigger with tracks in another MB event
3271
3272   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation() - Make trigger particle - charged hadron mixed event correlation \n");
3273
3274   if(GetMixedEvent()) return;  // This is not the mixed event from general mixing frame
3275
3276   // Get the event with similar caracteristics
3277   //printf("MakeChargedMixCorrelation for %s\n",GetInputAODName().Data());
3278
3279   AliAnalysisManager   * manager      = AliAnalysisManager::GetAnalysisManager();
3280
3281   AliInputEventHandler * inputHandler = dynamic_cast<AliInputEventHandler*>(manager->GetInputEventHandler());
3282
3283   if(!inputHandler) return;
3284
3285   if(!(inputHandler->IsEventSelected( ) & GetReader()->GetEventTriggerMask())) return;
3286
3287   // Get the pool, check if it exits
3288   Int_t eventBin = GetEventMixBin();
3289
3290   fhEventBin->Fill(eventBin);
3291
3292   //Check that the bin exists, if not (bad determination of RP, centrality or vz bin) do nothing
3293   if(eventBin < 0) return;
3294
3295   TList * pool     = 0;
3296   TList * poolCalo = 0;
3297   if(fUseMixStoredInReader)
3298   {
3299     pool     = GetReader()->GetListWithMixedEventsForTracks(eventBin);
3300     if(OnlyIsolated() || fFillNeutralEventMixPool) poolCalo = GetReader()->GetListWithMixedEventsForCalo  (eventBin);
3301   }
3302   else
3303   {
3304     pool     = fListMixTrackEvents[eventBin];
3305     if(OnlyIsolated()  || fFillNeutralEventMixPool) poolCalo = fListMixCaloEvents [eventBin];
3306   }
3307
3308   if(!pool) return ;
3309
3310   if((OnlyIsolated()  || fFillNeutralEventMixPool ) && !poolCalo &&
3311      (GetIsolationCut()->GetParticleTypeInCone()!=AliIsolationCut::AliIsolationCut::kOnlyCharged))
3312     printf("AliAnaParticleHadronCorrelation::MakeChargedMixCorrelation() - Careful, cluster pool not available\n");
3313
3314   Double_t ptTrig  = aodParticle->Pt();
3315   Double_t etaTrig = aodParticle->Eta();
3316   Double_t phiTrig = aodParticle->Phi();
3317   if(phiTrig < 0.) phiTrig+=TMath::TwoPi();
3318
3319   if(GetDebug() > 1)
3320     printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation() - Pool bin %d size %d, leading trigger pt=%f, phi=%f, eta=%f\n",
3321            eventBin,pool->GetSize(), ptTrig,phiTrig,etaTrig);
3322
3323   Double_t ptAssoc  = -999.;
3324   Double_t phiAssoc = -999.;
3325   Double_t etaAssoc = -999.;
3326   Double_t deltaPhi = -999.;
3327   Double_t deltaEta = -999.;
3328   Double_t xE = -999.;
3329   Double_t hbpXE = -999.;
3330
3331   //Start from first event in pool except if in this same event the pool was filled
3332   Int_t ev0 = 0;
3333   if(GetReader()->GetLastTracksMixedEvent() == GetEventNumber()) ev0 = 1;
3334
3335   for(Int_t ev=ev0; ev < pool->GetSize(); ev++)
3336   {
3337     TObjArray* bgTracks = static_cast<TObjArray*>(pool->At(ev));
3338     TObjArray* bgCalo   = 0;
3339
3340     // Check if the particle is isolated in the mixed event, it not, do not fill the histograms
3341     if((OnlyIsolated() || fFillNeutralEventMixPool) && poolCalo)
3342     {
3343       if(pool->GetSize()!=poolCalo->GetSize())
3344         printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation() - Different size of calo and track pools\n");
3345
3346       bgCalo = static_cast<TObjArray*>(poolCalo->At(ev));
3347
3348       if(!bgCalo)
3349         printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation() - Event %d in calo pool not available?\n",ev);
3350
3351       if(OnlyIsolated())
3352       {
3353         Int_t n=0; Int_t nfrac = 0; Bool_t isolated = kFALSE; Float_t coneptsum = 0;
3354         GetIsolationCut()->MakeIsolationCut(bgTracks,bgCalo,
3355                                             GetReader(), GetCaloPID(),
3356                                             kFALSE, aodParticle, "",
3357                                             n,nfrac,coneptsum, isolated);
3358
3359         //printf("AliAnaParticleHadronCorrelation::MakeChargedMixCorrelation() - Isolated? %d - cone %f, ptthres %f",
3360         //       isolated,GetIsolationCut()->GetConeSize(),GetIsolationCut()->GetPtThreshold());
3361         //if(bgTracks)printf(" - n track %d", bgTracks->GetEntriesFast());
3362         //printf("\n");
3363
3364         if(!isolated) continue ;
3365       }
3366     }
3367
3368     fhEventMixBin->Fill(eventBin);
3369
3370     Int_t nTracks=bgTracks->GetEntriesFast();
3371     //printf("\t Read Pool event %d, nTracks %d\n",ev,nTracks);
3372
3373     //Check if it is leading if mixed event
3374     if(fMakeNearSideLeading || fMakeAbsoluteLeading)
3375     {
3376       Bool_t leading = kTRUE;
3377       for(Int_t jlead = 0;jlead <nTracks; jlead++ )
3378       {
3379         AliAODPWG4Particle *track = (AliAODPWG4Particle*) bgTracks->At(jlead) ;
3380
3381         ptAssoc  = track->Pt();
3382         phiAssoc = track->Phi() ;
3383
3384         if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
3385         if (fMakeNearSideLeading)
3386         {
3387           if(ptAssoc > ptTrig && TMath::Abs(phiAssoc-phiTrig) < TMath::PiOver2())
3388           {
3389             leading = kFALSE;
3390             break;
3391           }
3392         }
3393         //jump out this event if there is any other particle with pt larger than trigger
3394         else if(fMakeAbsoluteLeading)
3395         {
3396           if(ptAssoc > ptTrig)
3397           {
3398             leading = kFALSE;
3399             break;
3400           }
3401         }
3402       }
3403
3404       if(fFillNeutralEventMixPool && bgCalo)
3405       {
3406         Int_t nClusters=bgCalo->GetEntriesFast();
3407         TLorentzVector mom ;
3408         for(Int_t jlead = 0;jlead <nClusters; jlead++ )
3409         {
3410           AliAODPWG4Particle *cluster= (AliAODPWG4Particle*) bgCalo->At(jlead) ;
3411
3412           ptAssoc  = cluster->Pt();
3413           phiAssoc = cluster->Phi() ;
3414
3415           if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
3416           if (fMakeNearSideLeading)
3417           {
3418             if(ptAssoc > ptTrig && TMath::Abs(phiAssoc-phiTrig) < TMath::PiOver2())
3419             {
3420               leading = kFALSE;
3421               break;
3422             }
3423           }
3424           //jump out this event if there is any other particle with pt larger than trigger
3425           else if(fMakeAbsoluteLeading)
3426           {
3427             if(ptAssoc > ptTrig)
3428             {
3429               leading = kFALSE;
3430               break;
3431             }
3432           }
3433         }
3434       }
3435
3436       if(!leading) continue; // not leading, check the next event in pool
3437
3438     }
3439
3440     fhPtLeadingMixed   ->Fill(ptTrig);
3441     fhPhiLeadingMixed  ->Fill(ptTrig, phiTrig);
3442     fhEtaLeadingMixed  ->Fill(ptTrig, etaTrig);
3443     fhPtLeadingMixedBin->Fill(ptTrig,eventBin);
3444     if(fCorrelVzBin)fhPtLeadingMixedVzBin->Fill(ptTrig, GetEventVzBin());
3445
3446     for(Int_t j1 = 0;j1 <nTracks; j1++ )
3447     {
3448       AliAODPWG4Particle *track = (AliAODPWG4Particle*) bgTracks->At(j1) ;
3449
3450       if(!track) continue;
3451
3452       ptAssoc  = track->Pt();
3453       etaAssoc = track->Eta();
3454       phiAssoc = track->Phi() ;
3455       if(phiAssoc < 0) phiAssoc+=TMath::TwoPi();
3456
3457       if(IsFiducialCutOn())
3458       {
3459         Bool_t in = GetFiducialCut()->IsInFiducialCut(*aodParticle->Momentum(),"CTS") ;
3460         if(!in) continue ;
3461       }
3462
3463       deltaPhi = phiTrig-phiAssoc;
3464       if(deltaPhi < -TMath::PiOver2())  deltaPhi+=TMath::TwoPi();
3465       if(deltaPhi > 3*TMath::PiOver2()) deltaPhi-=TMath::TwoPi();
3466       deltaEta = etaTrig-etaAssoc;
3467
3468       if(GetDebug()>0)
3469         printf("AliAnaParticleHadronCorrelationNew::MakeChargedMixCorrelation(): deltaPhi= %f, deltaEta=%f\n",deltaPhi, deltaEta);
3470
3471       // Set the pt associated bin for the defined bins
3472       Int_t assocBin   = -1;
3473       for(Int_t i = 0 ; i < fNAssocPtBins ; i++)
3474       {
3475         if(ptAssoc > fAssocPtBinLimit[i] && ptAssoc < fAssocPtBinLimit[i+1]) assocBin= i;
3476       }
3477
3478       // Assign to the histogram array a bin corresponding to a combination of pTa and vz bins
3479       Int_t nz = 1;
3480       Int_t vz = 0;
3481
3482       if(fCorrelVzBin)
3483       {
3484         nz = GetNZvertBin();
3485         vz = GetEventVzBin();
3486       }
3487
3488       Int_t bin = assocBin*nz+vz;
3489
3490       fhMixDeltaPhiCharged        ->Fill(ptTrig,  deltaPhi);
3491       fhMixDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta);
3492
3493       fhMixDeltaPhiCharged        ->Fill(ptTrig,  deltaPhi);
3494       fhMixDeltaPhiDeltaEtaCharged->Fill(deltaPhi, deltaEta);
3495
3496       xE   =-ptAssoc/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
3497       //if(xE <0.)xE =-xE;
3498       if(xE > 0 ) hbpXE = TMath::Log(1./xE);
3499       else        hbpXE =-100;
3500
3501       if ( (deltaPhi > fDeltaPhiMinCut)   && (deltaPhi < fDeltaPhiMaxCut)   )
3502       {
3503         fhMixXECharged->Fill(ptTrig,xE);
3504         fhMixHbpXECharged->Fill(ptTrig,hbpXE);
3505       }
3506
3507       if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )
3508       {
3509         //Underlying event region
3510         Double_t randomphi = gRandom->Uniform(fDeltaPhiMinCut,fDeltaPhiMaxCut);
3511         Double_t uexE = -(ptAssoc/ptTrig)*TMath::Cos(randomphi);
3512
3513         if(uexE < 0.) uexE = -uexE;
3514
3515         fhMixXEUeCharged->Fill(ptTrig,uexE);
3516       }
3517
3518       if(bin < 0) continue ; // this pt bin was not considered
3519
3520       if(TMath::Abs(deltaEta) > 0.8)
3521         fhMixDeltaPhiChargedAssocPtBinDEta08  [bin]->Fill(ptTrig,   deltaPhi);
3522       if(TMath::Abs(deltaEta) < 0.01)
3523         fhMixDeltaPhiChargedAssocPtBinDEta0   [bin]->Fill(ptTrig,   deltaPhi);
3524
3525         fhMixDeltaPhiChargedAssocPtBin        [bin]->Fill(ptTrig,   deltaPhi);
3526         fhMixDeltaPhiDeltaEtaChargedAssocPtBin[bin]->Fill(deltaPhi, deltaEta);
3527
3528     } // track loop
3529   } // mixed event loop
3530 }
3531
3532
3533 //________________________________________________________________________________________________________________
3534 Bool_t  AliAnaParticleHadronCorrelation::MakeNeutralCorrelation(AliAODPWG4ParticleCorrelation * const aodParticle,
3535                                                                 const TObjArray* pi0list, Bool_t bFillHisto)
3536 {
3537   // Neutral Pion Correlation Analysis
3538   if(GetDebug() > 1) printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Make trigger particle - pi0 correlation, %d pi0's \n",
3539                             pi0list->GetEntriesFast());
3540
3541   Int_t evtIndex11 = 0 ;
3542   Int_t evtIndex12 = 0 ;
3543   if (GetMixedEvent())
3544   {
3545     evtIndex11 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(0)) ;
3546     evtIndex12 = GetMixedEvent()->EventIndexForCaloCluster(aodParticle->GetCaloLabel(1)) ;
3547   }
3548
3549   Float_t pt   = -100. ;
3550 //  Float_t zT   = -100. ;
3551   Float_t phi  = -100. ;
3552   Float_t eta  = -100. ;
3553   Float_t xE   = -100. ;
3554   Float_t hbpXE= -100. ;
3555   //Float_t hbpZT= -100. ;
3556
3557   Float_t ptTrig  = aodParticle->Pt();
3558   Float_t phiTrig = aodParticle->Phi();
3559   Float_t etaTrig = aodParticle->Eta();
3560   Float_t deltaPhi= -100. ;
3561
3562   TLorentzVector photonMom ;
3563
3564   // In case of pi0/eta trigger, we may want to check their decay correlation,
3565   // get their decay children
3566   TLorentzVector decayMom1;
3567   TLorentzVector decayMom2;
3568   Bool_t decayFound = kFALSE;
3569   if(fPi0Trigger && bFillHisto) decayFound = GetDecayPhotonMomentum(aodParticle,decayMom1, decayMom2);
3570
3571   TObjArray * refpi0 = 0x0;
3572   Int_t nrefs        = 0;
3573
3574   //Loop on stored AOD pi0
3575
3576   Int_t naod = pi0list->GetEntriesFast();
3577   if(GetDebug() > 0)
3578     printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelationFillHistograms() -  aod branch entries %d\n", naod);
3579
3580   for(Int_t iaod = 0; iaod < naod ; iaod++)
3581   {
3582     AliAODPWG4Particle* pi0 =  (AliAODPWG4Particle*) (pi0list->At(iaod));
3583
3584     Int_t evtIndex2 = 0 ;
3585     Int_t evtIndex3 = 0 ;
3586     if (GetMixedEvent())
3587     {
3588       evtIndex2 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(0)) ;
3589       evtIndex3 = GetMixedEvent()->EventIndexForCaloCluster(pi0->GetCaloLabel(1)) ;
3590
3591       if (evtIndex11 == evtIndex2 || evtIndex12 == evtIndex2 ||
3592           evtIndex11 == evtIndex3 || evtIndex12 == evtIndex3) // trigger and pi0 are not from different events
3593         continue ;
3594     }
3595
3596     pt  = pi0->Pt();
3597
3598     if(pt < fMinAssocPt || pt > fMaxAssocPt) continue ;
3599
3600     //jump out this event if near side associated particle pt larger than trigger
3601     if (fMakeNearSideLeading)
3602     {
3603       if(pt > ptTrig && TMath::Abs(phi-phiTrig) < TMath::PiOver2())  return kFALSE;
3604     }
3605     //jump out this event if there is any other particle with pt larger than trigger
3606     else if(fMakeAbsoluteLeading)
3607     {
3608       if(pt > ptTrig)  return kFALSE;
3609     }
3610
3611     if(bFillHisto)
3612     {
3613       phi = pi0->Phi() ;
3614       eta = pi0->Eta() ;
3615
3616       FillNeutralAngularCorrelationHistograms(pt, ptTrig, phi, phiTrig, deltaPhi, eta, etaTrig);
3617
3618       //zT  = pt/ptTrig ;
3619       xE  =-pt/ptTrig*TMath::Cos(deltaPhi); // -(px*pxTrig+py*pyTrig)/(ptTrig*ptTrig);
3620
3621       //if(xE <0.)xE =-xE;
3622
3623       hbpXE = -100;
3624       //hbpZT = -100;
3625
3626       if(xE > 0 ) hbpXE = TMath::Log(1./xE);
3627       //if(zT > 0 ) hbpZT = TMath::Log(1./zT);
3628
3629       if(fPi0Trigger && decayFound)
3630         FillDecayPhotonCorrelationHistograms(pt, phi, decayMom1,decayMom2,kFALSE) ;
3631
3632       //delta phi cut for correlation
3633       if( (deltaPhi > fDeltaPhiMinCut) && ( deltaPhi < fDeltaPhiMaxCut) )
3634       {
3635         fhDeltaPhiNeutralPt->Fill(pt,deltaPhi);
3636         fhXENeutral        ->Fill(ptTrig,xE);
3637         fhPtHbpXENeutral   ->Fill(ptTrig,hbpXE);
3638       }
3639       else if ( (deltaPhi > fUeDeltaPhiMinCut) && (deltaPhi < fUeDeltaPhiMaxCut) )
3640       {
3641         fhDeltaPhiUeNeutralPt->Fill(pt,deltaPhi);
3642         fhXEUeNeutral        ->Fill(ptTrig,xE);
3643         fhPtHbpXEUeNeutral   ->Fill(ptTrig,hbpXE);
3644       }
3645
3646       //several UE calculation
3647       if(fMakeSeveralUE) FillChargedUnderlyingEventSidesHistograms(ptTrig,pt,deltaPhi);
3648
3649           }
3650     else
3651     {
3652       nrefs++;
3653       if(nrefs==1)
3654       {
3655         refpi0 = new TObjArray(0);
3656         refpi0->SetName(GetAODObjArrayName()+"Pi0s");
3657         refpi0->SetOwner(kFALSE);
3658       }
3659       refpi0->Add(pi0);
3660     }//put references in trigger AOD
3661
3662     if(GetDebug() > 2 )
3663       printf("AliAnaParticleHadronCorrelation::MakeNeutralCorrelation() - Selected neutral for momentum imbalance: pt %2.2f, phi %2.2f, eta %2.2f \n",pt,phi,eta);
3664
3665   }//loop
3666
3667   return kTRUE;
3668 }
3669
3670 //_________________________________________________________________________________________________________
3671 void  AliAnaParticleHadronCorrelation::MakeMCChargedCorrelation(AliAODPWG4ParticleCorrelation *aodParticle)
3672 {
3673   // Charged Hadron Correlation Analysis with MC information
3674
3675   if ( GetDebug() > 1 )
3676     AliInfo("Make trigger particle - charged hadron correlation in AOD MC level");
3677
3678   AliStack         * stack        = 0x0 ;
3679   TParticle        * primary      = 0x0 ;
3680   TClonesArray     * mcparticles  = 0x0 ;
3681   AliAODMCParticle * aodprimary   = 0x0 ;
3682
3683   Double_t eprim   = 0 ;
3684   Double_t ptprim  = 0 ;
3685   Double_t phiprim = 0 ;
3686   Double_t etaprim = 0 ;
3687   Int_t    nTracks = 0 ;
3688   Int_t iParticle  = 0 ;
3689
3690   Bool_t lead = kFALSE;
3691
3692   Int_t label= aodParticle->GetLabel();
3693   if( label < 0 )
3694   {
3695     if( GetDebug() > 0 ) AliInfo(Form(" *** bad label ***:  label %d", label));
3696     return;
3697   }
3698
3699   if( GetReader()->ReadStack() )
3700   {
3701     stack =  GetMCStack() ;
3702     if(!stack)
3703     {
3704       AliFatal("Stack not available, is the MC handler called? STOP");
3705       return;
3706     }
3707
3708     //nTracks = stack->GetNtrack() ;
3709     nTracks = stack->GetNprimary();
3710     if( label >=  stack->GetNtrack() )
3711     {
3712       if(GetDebug() > 2)
3713         AliInfo(Form("*** large label ***:  label %d, n tracks %d", label, stack->GetNtrack()));
3714       return ;
3715     }
3716
3717     primary = stack->Particle(label);
3718     if ( !primary )
3719     {
3720       AliInfo(Form(" *** no primary ***:  label %d", label));
3721       return;
3722     }
3723
3724     eprim    = primary->Energy();
3725     ptprim   = primary->Pt();
3726     phiprim  = primary->Phi();
3727     etaprim  = primary->Eta();
3728
3729     if(ptprim < 0.01 || eprim < 0.01) return ;
3730
3731     for (iParticle = 0 ; iParticle <  nTracks ; iParticle++)
3732     {
3733       TParticle * particle = stack->Particle(iParticle);
3734       TLorentzVector momentum;
3735
3736       //keep only final state particles
3737       if( particle->GetStatusCode() != 1 ) continue ;
3738
3739       if ( particle->Pt() < GetReader()->GetCTSPtMin()) continue;
3740
3741       //---------- Charged particles ----------------------
3742       Int_t pdg    = particle->GetPdgCode();
3743       Int_t charge = (Int_t) TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
3744       if(charge == 0) continue;
3745
3746       particle->Momentum(momentum);
3747
3748       //Particles in CTS acceptance, make sure to use the same selection as in the reader
3749       Bool_t inCTS =  GetReader()->GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
3750       //printf("Accepted? %d; pt %2.2f, eta %2.2f, phi %2.2f\n",inCTS,momentum.Pt(),momentum.Eta(),momentum.Phi()*TMath::RadToDeg());
3751       if( !inCTS ) continue;
3752
3753       // Remove conversions
3754       if ( TMath::Abs(pdg) == 11 && stack->Particle(particle->GetFirstMother())->GetPdgCode() == 22 ) continue ;
3755
3756       if ( label == iParticle ) continue; // avoid trigger particle
3757
3758       lead = FillChargedMCCorrelationHistograms(particle->Pt(),particle->Phi(),particle->Eta(),ptprim,phiprim,etaprim);
3759       if ( !lead ) return;
3760
3761     } //track loop
3762
3763   } //ESD MC
3764
3765   else if( GetReader()->ReadAODMCParticles() )
3766   {
3767     //Get the list of MC particles
3768     mcparticles = GetReader()->GetAODMCParticles();
3769     if( !mcparticles ) return;
3770
3771     nTracks = mcparticles->GetEntriesFast() ;
3772
3773     if( label >= nTracks )
3774     {
3775       if(GetDebug() > 2)
3776         AliInfo(Form(" *** large label ***:  label %d, n tracks %d", label,nTracks));
3777       return;
3778     }
3779
3780     //Get the particle
3781     aodprimary = (AliAODMCParticle*) mcparticles->At(label);
3782     if( !aodprimary )
3783     {
3784       AliInfo(Form(" *** no AOD primary ***:  label %d", label));
3785       return;
3786     }
3787
3788     ptprim  = aodprimary->Pt();
3789     phiprim = aodprimary->Phi();
3790     etaprim = aodprimary->Eta();
3791     eprim   = aodprimary->E();
3792
3793     if(ptprim < 0.01 || eprim < 0.01) return ;
3794
3795     for (iParticle = 0; iParticle < nTracks; iParticle++)
3796     {
3797       AliAODMCParticle *part = (AliAODMCParticle*) mcparticles->At(iParticle);
3798
3799       if (!part->IsPhysicalPrimary() ) continue; // same as part->GetStatus() !=1
3800
3801       if ( part->Charge() == 0 ) continue;
3802
3803       if ( part->Pt() < GetReader()->GetCTSPtMin()) continue;
3804
3805       TLorentzVector momentum(part->Px(),part->Py(),part->Pz(),part->E());
3806
3807       //Particles in CTS acceptance, make sure to use the same selection as in the reader
3808       Bool_t inCTS =  GetReader()->GetFiducialCut()->IsInFiducialCut(momentum,"CTS");
3809       //printf("Accepted? %d; pt %2.2f, eta %2.2f, phi %2.2f\n",inCTS,momentum.Pt(),momentum.Eta(),momentum.Phi()*TMath::RadToDeg());
3810       if( !inCTS ) continue;
3811
3812       // Remove conversions
3813       Int_t indexmother = part->GetMother();
3814       if ( indexmother > -1 )
3815       {
3816         Int_t pdg = part->GetPdgCode();
3817         Int_t mPdg = ((AliAODMCParticle*) mcparticles->At(indexmother)) ->GetPdgCode();
3818         if (TMath::Abs(pdg) == 11 && mPdg == 22) continue;
3819       }
3820
3821       if ( label == iParticle ) continue; // avoid trigger particle
3822
3823       lead = FillChargedMCCorrelationHistograms(part->Pt(),part->Phi(),part->Eta(),ptprim,phiprim,etaprim);
3824       if ( !lead ) return;
3825
3826     }  //MC particle loop
3827   }// AOD MC
3828
3829   // Leading MC particle histograms
3830   if (lead)
3831   {
3832     fhMCPtLeading ->Fill(ptprim);
3833     fhMCPhiLeading->Fill(ptprim,phiprim);
3834     fhMCEtaLeading->Fill(ptprim,etaprim);
3835   }
3836 }
3837
3838 //_____________________________________________________________________
3839 void AliAnaParticleHadronCorrelation::Print(const Option_t * opt) const
3840 {
3841
3842   //Print some relevant parameters set for the analysis
3843   if(! opt)
3844     return;
3845
3846   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
3847   AliAnaCaloTrackCorrBaseClass::Print(" ");
3848   printf("Pt trigger           >  %3.2f\n", fMinTriggerPt) ;
3849   printf("Pt associated hadron <  %3.2f\n", fMaxAssocPt) ;
3850   printf("Pt associated hadron >  %3.2f\n", fMinAssocPt) ;
3851   printf("Phi trigger particle-Hadron   <  %3.2f\n", fDeltaPhiMaxCut) ;
3852   printf("Phi trigger particle-Hadron   >  %3.2f\n", fDeltaPhiMinCut) ;
3853   printf("Phi trigger particle-UeHadron <  %3.2f\n", fUeDeltaPhiMaxCut) ;
3854   printf("Phi trigger particle-UeHadron >  %3.2f\n", fUeDeltaPhiMinCut) ;
3855   printf("Isolated Trigger?  %d\n"     , fSelectIsolated) ;
3856   printf("Several UE?  %d\n"           , fMakeSeveralUE) ;
3857   printf("Name of AOD Pi0 Branch %s \n", fPi0AODBranchName.Data());
3858   printf("Do Decay-hadron correlation ?  %d\n", fPi0Trigger) ;
3859   printf("Select absolute leading for cluster triggers ?  %d or Near Side %d\n",
3860          fMakeAbsoluteLeading, fMakeNearSideLeading) ;
3861   printf("Trigger pt bins  %d\n", fNAssocPtBins) ;
3862   for (Int_t ibin = 0; ibin<fNAssocPtBins; ibin++) {
3863     printf("\t bin %d = [%2.1f,%2.1f]\n", ibin, fAssocPtBinLimit[ibin], fAssocPtBinLimit[ibin+1]) ;
3864   }
3865
3866 }
3867
3868 //____________________________________________________________
3869 void AliAnaParticleHadronCorrelation::SetNAssocPtBins(Int_t n)
3870 {
3871   // Set number of bins
3872
3873   fNAssocPtBins  = n ;
3874
3875
3876   if(n < 20 && n > 0)
3877   {
3878     fNAssocPtBins  = n ;
3879   }
3880   else
3881   {
3882     printf("n = larger than 19 or too small, set to 19 \n");
3883     fNAssocPtBins = 19;
3884   }
3885 }
3886
3887 //______________________________________________________________________________
3888 void AliAnaParticleHadronCorrelation::SetAssocPtBinLimit(Int_t ibin, Float_t pt)
3889 {
3890   // Set the list of limits for the trigger pt bins
3891
3892   if(ibin <= fNAssocPtBins || ibin >= 0)
3893   {
3894     fAssocPtBinLimit[ibin] = pt  ;
3895   }
3896   else
3897   {
3898     printf("AliAnaParticleHadronCorrelation::SetAssocPtBinLimit() - bin  number too large %d > %d or small, nothing done\n", ibin, fNAssocPtBins) ;
3899   }
3900 }
3901