]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGGA/CaloTrackCorrelations/AliAnaPhoton.cxx
put different cluster parameters (time, n cells, n SM) in the AOD particle, recover...
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaPhoton.cxx
index 86d566a3a14c220574e04c6f71b9bef3c3ea337e..41def4aa60950a5d5db463349feea6dd8d7e95a9 100755 (executable)
@@ -27,7 +27,6 @@
 
 // --- ROOT system ---
 #include <TH2F.h>
-#include <TH3D.h>
 #include <TClonesArray.h>
 #include <TObjString.h>
 #include "TParticle.h"
@@ -54,7 +53,7 @@ ClassImp(AliAnaPhoton)
 
 //____________________________
 AliAnaPhoton::AliAnaPhoton() :
-AliAnaCaloTrackCorrBaseClass(), fCalorimeter(""),
+AliAnaCaloTrackCorrBaseClass(),
 fMinDist(0.),                 fMinDist2(0.),                fMinDist3(0.),
 fRejectTrackMatch(0),         fFillTMHisto(kFALSE),
 fTimeCutMin(-10000),          fTimeCutMax(10000),
@@ -62,51 +61,16 @@ fNCellsCut(0),
 fNLMCutMin(-1),               fNLMCutMax(10),
 fFillSSHistograms(kFALSE),    fFillOnlySimpleSSHisto(1),
 fNOriginHistograms(8),        fNPrimaryHistograms(4),
-fFillPileUpHistograms(0),     fFillEMCALBCHistograms(0),
+fMomentum(),                  fPrimaryMom(),
 // Histograms
-fhNCellsE(0),                 fhCellsE(0),   // Control histograms
-fhMaxCellDiffClusterE(0),     fhTimePt(0),   // Control histograms
-fhEtaPhi(0),                  fhEtaPhiEMCALBC0(0),
-fhEtaPhiEMCALBC1(0),          fhEtaPhiEMCALBCN(0),
-fhTimeTriggerEMCALBCCluster(0),
-fhTimeTriggerEMCALBCUMCluster(0),
-fhEtaPhiTriggerEMCALBCClusterOverTh(0),
-fhEtaPhiTriggerEMCALBCUMClusterOverTh(0),
-fhEtaPhiTriggerEMCALBCClusterBelowTh1(0),
-fhEtaPhiTriggerEMCALBCUMClusterBelowTh1(0),
-fhEtaPhiTriggerEMCALBCClusterBelowTh2(0),
-fhEtaPhiTriggerEMCALBCUMClusterBelowTh2(0),
-fhEtaPhiTriggerEMCALBCExotic(0),             fhTimeTriggerEMCALBCExotic(0),
-fhEtaPhiTriggerEMCALBCUMExotic(0),           fhTimeTriggerEMCALBCUMExotic(0),
-fhEtaPhiTriggerEMCALBCBad(0),                fhTimeTriggerEMCALBCBad(0),
-fhEtaPhiTriggerEMCALBCUMBad(0),              fhTimeTriggerEMCALBCUMBad(0),
-fhEtaPhiTriggerEMCALBCBadExotic(0),          fhTimeTriggerEMCALBCBadExotic(0),
-fhEtaPhiTriggerEMCALBCUMBadExotic(0),        fhTimeTriggerEMCALBCUMBadExotic(0),
-fhEtaPhiTriggerEMCALBCExoticCluster(0),      fhTimeTriggerEMCALBCExoticCluster(0),
-fhEtaPhiTriggerEMCALBCUMExoticCluster(0),    fhTimeTriggerEMCALBCUMExoticCluster(0),
-fhEtaPhiTriggerEMCALBCBadCluster(0),         fhTimeTriggerEMCALBCBadCluster(0),
-fhEtaPhiTriggerEMCALBCUMBadCluster(0),       fhTimeTriggerEMCALBCUMBadCluster(0),
-fhEtaPhiTriggerEMCALBCBadExoticCluster(0),   fhTimeTriggerEMCALBCBadExoticCluster(0),
-fhEtaPhiTriggerEMCALBCUMBadExoticCluster(0), fhTimeTriggerEMCALBCUMBadExoticCluster(0),
-fhTimeTriggerEMCALBCBadMaxCell(0),           fhTimeTriggerEMCALBCUMBadMaxCell(0),
-fhTimeTriggerEMCALBCBadMaxCellExotic(0),     fhTimeTriggerEMCALBCUMBadMaxCellExotic(0),
-fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster (0), fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster(0),
-fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster(0),fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster(0),
-fhEtaPhiTriggerEMCALBCUMReMatchBothCluster(0),      fhTimeTriggerEMCALBCUMReMatchBothCluster(0),
-fhTimeTriggerEMCALBC0UMReMatchOpenTime(0),
-fhTimeTriggerEMCALBC0UMReMatchCheckNeigh(0),
-fhTimeTriggerEMCALBC0UMReMatchBoth(0),
 
-fhEtaPhiNoTrigger(0),                        fhTimeNoTrigger(0),
+// Control histograms
+fhNCellsE(0),                 fhCellsE(0),
+fhMaxCellDiffClusterE(0),     fhTimePt(0),                  fhEtaPhi(0),
 
 fhEPhoton(0),                 fhPtPhoton(0),
 fhPhiPhoton(0),               fhEtaPhoton(0),
 fhEtaPhiPhoton(0),            fhEtaPhi05Photon(0),
-fhEtaPhiPhotonEMCALBC0(0),    fhEtaPhiPhotonEMCALBC1(0),   fhEtaPhiPhotonEMCALBCN(0),
-fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime(0),
-fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh(0),
-fhTimePhotonTriggerEMCALBC0UMReMatchBoth(0),
-
 fhPtCentralityPhoton(0),      fhPtEventPlanePhoton(0),
 
 // Shower shape histograms
@@ -137,22 +101,15 @@ fhEmbedPhotonELambda0FullSignal(0),   fhEmbedPhotonELambda0MostlySignal(0),
 fhEmbedPhotonELambda0MostlyBkg(0),    fhEmbedPhotonELambda0FullBkg(0),
 fhEmbedPi0ELambda0FullSignal(0),      fhEmbedPi0ELambda0MostlySignal(0),
 fhEmbedPi0ELambda0MostlyBkg(0),       fhEmbedPi0ELambda0FullBkg(0),
-// PileUp
-fhTimePtNoCut(0),                     fhTimePtSPD(0),
+
 fhTimePtPhotonNoCut(0),               fhTimePtPhotonSPD(0),
 fhTimeNPileUpVertSPD(0),              fhTimeNPileUpVertTrack(0),
-fhTimeNPileUpVertContributors(0),
-fhTimePileUpMainVertexZDistance(0),   fhTimePileUpMainVertexZDiamond(0),
-fhClusterMultSPDPileUp(),             fhClusterMultNoPileUp(),
-fhEtaPhiBC0(0),  fhEtaPhiBCPlus(0),   fhEtaPhiBCMinus(0),
-fhEtaPhiBC0PileUpSPD(0),
-fhEtaPhiBCPlusPileUpSPD(0),           fhEtaPhiBCMinusPileUpSPD(0),
-fhPtNPileUpSPDVtx(0),                 fhPtNPileUpTrkVtx(0),
-fhPtNPileUpSPDVtxTimeCut(0),          fhPtNPileUpTrkVtxTimeCut(0),
-fhPtNPileUpSPDVtxTimeCut2(0),         fhPtNPileUpTrkVtxTimeCut2(0),
 fhPtPhotonNPileUpSPDVtx(0),           fhPtPhotonNPileUpTrkVtx(0),
 fhPtPhotonNPileUpSPDVtxTimeCut(0),    fhPtPhotonNPileUpTrkVtxTimeCut(0),
-fhPtPhotonNPileUpSPDVtxTimeCut2(0),   fhPtPhotonNPileUpTrkVtxTimeCut2(0)
+fhPtPhotonNPileUpSPDVtxTimeCut2(0),   fhPtPhotonNPileUpTrkVtxTimeCut2(0),
+
+fhEClusterSM(0),                      fhEPhotonSM(0),
+fhPtClusterSM(0),                     fhPtPhotonSM(0)
 {
   //default ctor
   
@@ -185,21 +142,10 @@ fhPtPhotonNPileUpSPDVtxTimeCut2(0),   fhPtPhotonNPileUpTrkVtxTimeCut2(0)
     fhDispEtaDispPhi[i] = 0;
     fhLambda0DispPhi[i] = 0;
     fhLambda0DispEta[i] = 0;
-    
-    fhPtPileUp       [i] = 0;
-    fhPtChargedPileUp[i] = 0;
-    fhPtPhotonPileUp [i] = 0;
-    
-    fhLambda0PileUp       [i] = 0;
-    fhLambda0ChargedPileUp[i] = 0;
-    
-    fhClusterEFracLongTimePileUp  [i] = 0;
-    
-    fhClusterCellTimePileUp       [i] = 0;
-    fhClusterTimeDiffPileUp       [i] = 0;
-    fhClusterTimeDiffChargedPileUp[i] = 0;
+
+    fhPtPhotonPileUp[i] = 0;
     fhClusterTimeDiffPhotonPileUp [i] = 0;
-    
+
     for(Int_t j = 0; j < 6; j++)
     {
       fhMCDispEtaDispPhi[i][j] = 0;
@@ -234,7 +180,8 @@ fhPtPhotonNPileUpSPDVtxTimeCut2(0),   fhPtPhotonNPileUpTrkVtxTimeCut2(0)
   
   for(Int_t i = 0; i < 5; i++)
   {
-    fhClusterCuts[i] = 0;
+    fhClusterCutsE [i] = 0;
+    fhClusterCutsPt[i] = 0;
   }
   
   // Track matching residuals
@@ -252,104 +199,86 @@ fhPtPhotonNPileUpSPDVtxTimeCut2(0),   fhPtPhotonNPileUpTrkVtxTimeCut2(0)
     fhEOverPTRD[i] = 0;
   }
   
-  for(Int_t i = 0; i < 4; i++)
-  {
-    fhClusterMultSPDPileUp[i] = 0;
-    fhClusterMultNoPileUp [i] = 0;
-  }
-  
-  for(Int_t i = 0; i < 11; i++)
-  {
-    fhEtaPhiTriggerEMCALBC             [i] = 0 ;
-    fhTimeTriggerEMCALBC               [i] = 0 ;
-    fhEtaPhiTriggerEMCALBCUM           [i] = 0 ;
-    fhTimeTriggerEMCALBCUM             [i] = 0 ;
-    
-    fhEtaPhiPhotonTriggerEMCALBC       [i] = 0 ;
-    fhTimePhotonTriggerEMCALBC         [i] = 0 ;
-    fhEtaPhiPhotonTriggerEMCALBCUM     [i] = 0 ;
-    fhTimePhotonTriggerEMCALBCUM       [i] = 0 ;
-    
-    fhTimePhotonTriggerEMCALBCPileUpSPD[i] = 0 ;
-    fhTimeTriggerEMCALBCPileUpSPD      [i] = 0 ;
-    
-    fhEtaPhiTriggerEMCALBCCluster      [i] = 0 ;
-    fhEtaPhiTriggerEMCALBCUMCluster    [i] = 0 ;    
-  }
-  
   //Initialize parameters
   InitParameters();
   
 }
 
 //_________________________________________________________________________________________
-Bool_t  AliAnaPhoton::ClusterSelected(AliVCluster* calo, TLorentzVector mom, Int_t nMaxima)
+Bool_t  AliAnaPhoton::ClusterSelected(AliVCluster* calo, Int_t nMaxima)
 {
   //Select clusters if they pass different cuts
   
-  Float_t ptcluster  = mom.Pt();
-  Float_t ecluster   = mom.E();
-  Float_t l0cluster  = calo->GetM02();
-  Float_t etacluster = mom.Eta();
-  Float_t phicluster = mom.Phi();
-  if(phicluster<0) phicluster+=TMath::TwoPi();
-  Float_t tofcluster   = calo->GetTOF()*1.e9;
+  Float_t ptcluster  = fMomentum.Pt();
+  Float_t ecluster   = fMomentum.E();
+  Float_t etacluster = fMomentum.Eta();
+  Float_t phicluster = fMomentum.Phi();
+
+  if(phicluster < 0) phicluster+=TMath::TwoPi();
   
   Bool_t matched = IsTrackMatched(calo,GetReader()->GetInputEvent());
   
-  if(GetDebug() > 2)
-    printf("AliAnaPhoton::ClusterSelected() - Current Event %d; Before selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f\n",
+  AliDebug(2,Form("Current Event %d; Before selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f",
            GetReader()->GetEventNumber(),
-           ecluster,ptcluster, phicluster*TMath::RadToDeg(),etacluster);
+           ecluster,ptcluster, phicluster*TMath::RadToDeg(),etacluster));
   
-  fhClusterCuts[1]->Fill(ecluster);
+  fhClusterCutsE [1]->Fill( ecluster);
+  fhClusterCutsPt[1]->Fill(ptcluster);
   
   if(ecluster > 0.5) fhEtaPhi->Fill(etacluster, phicluster);
   
-  FillEMCALTriggerClusterBCHistograms(calo->GetID(),ecluster,tofcluster,etacluster,phicluster);
+  Int_t   nSM  = GetModuleNumber(calo);
+  if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
+  {
+    fhEClusterSM ->Fill(ecluster ,nSM);
+    fhPtClusterSM->Fill(ptcluster,nSM);
+  }
   
   //.......................................
   //If too small or big energy, skip it
   if(ecluster < GetMinEnergy() || ecluster > GetMaxEnergy() ) return kFALSE ;
   
-  if(GetDebug() > 2) printf("\t Cluster %d Pass E Cut \n",calo->GetID());
+  AliDebug(2,Form("\t Cluster %d Pass E Cut",calo->GetID()));
   
-  fhClusterCuts[2]->Fill(ecluster);
-  
-  FillClusterPileUpHistograms(calo,matched,ptcluster,etacluster,phicluster,l0cluster);
+  fhClusterCutsE [2]->Fill( ecluster);
+  fhClusterCutsPt[2]->Fill(ptcluster);
   
   //.......................................
   // TOF cut, BE CAREFUL WITH THIS CUT
   Double_t tof = calo->GetTOF()*1e9;
   if(tof < fTimeCutMin || tof > fTimeCutMax) return kFALSE;
   
-  if(GetDebug() > 2)  printf("\t Cluster %d Pass Time Cut \n",calo->GetID());
+  AliDebug(2,Form("\t Cluster %d Pass Time Cut",calo->GetID()));
   
-  fhClusterCuts[3]->Fill(ecluster);
+  fhClusterCutsE [3]->Fill( ecluster);
+  fhClusterCutsPt[3]->Fill(ptcluster);
   
   //.......................................
   if(calo->GetNCells() <= fNCellsCut && GetReader()->GetDataType() != AliCaloTrackReader::kMC) return kFALSE;
   
-  if(GetDebug() > 2) printf("\t Cluster %d Pass NCell Cut \n",calo->GetID());
+  AliDebug(2,Form("\t Cluster %d Pass NCell Cut",calo->GetID()));
   
-  fhClusterCuts[4]->Fill(ecluster);
+  fhClusterCutsE [4]->Fill( ecluster);
+  fhClusterCutsPt[4]->Fill(ptcluster);
   
   if(nMaxima < fNLMCutMin || nMaxima > fNLMCutMax) return kFALSE ;
-  if(GetDebug() > 2) printf(" \t Cluster %d pass NLM %d of out of range \n",calo->GetID(), nMaxima);
+  AliDebug(2,Form("\t Cluster %d pass NLM %d of out of range",calo->GetID(), nMaxima));
   
-  fhClusterCuts[5]->Fill(ecluster);
+  fhClusterCutsE [5]->Fill( ecluster);
+  fhClusterCutsPt[5]->Fill(ptcluster);
   
   //.......................................
   //Check acceptance selection
   if(IsFiducialCutOn())
   {
-    Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
+    Bool_t in = GetFiducialCut()->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),GetCalorimeter()) ;
     if(! in ) return kFALSE ;
   }
   
-  if(GetDebug() > 2) printf("Fiducial cut passed \n");
+  AliDebug(2,Form("\t Fiducial cut passed"));
   
-  fhClusterCuts[6]->Fill(ecluster);
+  fhClusterCutsE [6]->Fill( ecluster);
+  fhClusterCutsPt[6]->Fill(ptcluster);
   
   //.......................................
   //Skip matched clusters with tracks
@@ -361,25 +290,15 @@ Bool_t  AliAnaPhoton::ClusterSelected(AliVCluster* calo, TLorentzVector mom, Int
   {
     if(matched)
     {
-      if(GetDebug() > 2) printf("\t Reject track-matched clusters\n");
+      AliDebug(2,"\t Reject track-matched clusters");
       return kFALSE ;
     }
     else
-      if(GetDebug() > 2)  printf(" Track-matching cut passed \n");
+      AliDebug(2,"\t Track-matching cut passed");
   }// reject matched clusters
   
-  fhClusterCuts[7]->Fill(ecluster);
-  
-  if(fFillPileUpHistograms)
-  {
-    if(GetReader()->IsPileUpFromSPD())               {fhPtChargedPileUp[0]->Fill(ptcluster); fhLambda0ChargedPileUp[0]->Fill(ecluster,l0cluster); }
-    if(GetReader()->IsPileUpFromEMCal())             {fhPtChargedPileUp[1]->Fill(ptcluster); fhLambda0ChargedPileUp[1]->Fill(ecluster,l0cluster); }
-    if(GetReader()->IsPileUpFromSPDOrEMCal())        {fhPtChargedPileUp[2]->Fill(ptcluster); fhLambda0ChargedPileUp[2]->Fill(ecluster,l0cluster); }
-    if(GetReader()->IsPileUpFromSPDAndEMCal())       {fhPtChargedPileUp[3]->Fill(ptcluster); fhLambda0ChargedPileUp[3]->Fill(ecluster,l0cluster); }
-    if(GetReader()->IsPileUpFromSPDAndNotEMCal())    {fhPtChargedPileUp[4]->Fill(ptcluster); fhLambda0ChargedPileUp[4]->Fill(ecluster,l0cluster); }
-    if(GetReader()->IsPileUpFromEMCalAndNotSPD())    {fhPtChargedPileUp[5]->Fill(ptcluster); fhLambda0ChargedPileUp[5]->Fill(ecluster,l0cluster); }
-    if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtChargedPileUp[6]->Fill(ptcluster); fhLambda0ChargedPileUp[6]->Fill(ecluster,l0cluster); }
-  }
+  fhClusterCutsE [7]->Fill( ecluster);
+  fhClusterCutsPt[7]->Fill(ptcluster);
   
   //.......................................
   //Check Distance to Bad channel, set bit.
@@ -389,14 +308,14 @@ Bool_t  AliAnaPhoton::ClusterSelected(AliVCluster* calo, TLorentzVector mom, Int
   {//In bad channel (PHOS cristal size 2.2x2.2 cm), EMCAL ( cell units )
     return kFALSE ;
   }
-  else if(GetDebug() > 2) printf("\t Bad channel cut passed %4.2f > %2.2f \n",distBad, fMinDist);
+  else AliDebug(2,Form("\t Bad channel cut passed %4.2f > %2.2f",distBad, fMinDist));
   
-  fhClusterCuts[8]->Fill(ecluster);
+  fhClusterCutsE [8]->Fill( ecluster);
+  fhClusterCutsPt[8]->Fill(ptcluster);
   
-  if(GetDebug() > 0)
-    printf("AliAnaPhoton::ClusterSelected() Current Event %d; After  selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f\n",
+  AliDebug(1,Form("Current Event %d; After  selection : E %2.2f, pT %2.2f, phi %2.2f, eta %2.2f",
            GetReader()->GetEventNumber(),
-           ecluster, ptcluster,mom.Phi()*TMath::RadToDeg(),mom.Eta());
+           ecluster, ptcluster,fMomentum.Phi()*TMath::RadToDeg(),fMomentum.Eta()));
   
   //All checks passed, cluster selected
   return kTRUE;
@@ -416,758 +335,215 @@ void AliAnaPhoton::FillAcceptanceHistograms()
   
   Int_t    pdg       =  0 ;
   Int_t    tag       =  0 ;
+  Int_t    status    =  0 ;
   Int_t    mcIndex   =  0 ;
-  Bool_t   inacceptance = kFALSE;
+  Int_t    nprim     =  0 ;
+  Bool_t   inacceptance = kFALSE ;
   
-  if(GetReader()->ReadStack())
-  {
-    AliStack * stack = GetMCStack();
-    if(stack)
-    {
-      for(Int_t i=0 ; i<stack->GetNtrack(); i++)
-      {
-        if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
-        
-        TParticle * prim = stack->Particle(i) ;
-        pdg = prim->GetPdgCode();
-        //printf("i %d, %s %d  %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
-        //                             prim->GetName(), prim->GetPdgCode());
-        
-        if(pdg == 22)
-        {
-          // Get tag of this particle photon from fragmentation, decay, prompt ...
-          tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader());
-          if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
-          {
-            //A conversion photon from a hadron, skip this kind of photon
-            // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
-            // GetMCAnalysisUtils()->PrintMCTag(tag);
-            
-            return;
-          }
-          
-          //Get photon kinematics
-          if(prim->Energy() == TMath::Abs(prim->Pz()))  continue ; //Protection against floating point exception
-          
-          photonY   = 0.5*TMath::Log((prim->Energy()-prim->Pz())/(prim->Energy()+prim->Pz())) ;
-          photonE   = prim->Energy() ;
-          photonPt  = prim->Pt() ;
-          photonPhi = TMath::RadToDeg()*prim->Phi() ;
-          if(photonPhi < 0) photonPhi+=TMath::TwoPi();
-          photonEta = prim->Eta() ;
-          
-          //Check if photons hit the Calorimeter
-          TLorentzVector lv;
-          prim->Momentum(lv);
-          inacceptance = kFALSE;
-          if     (fCalorimeter == "PHOS")
-          {
-            if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet())
-            {
-              Int_t mod ;
-              Double_t x,z ;
-              if(GetPHOSGeometry()->ImpactOnEmc(prim,mod,z,x))
-                inacceptance = kTRUE;
-              if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-            else
-            {
-              if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
-                inacceptance = kTRUE ;
-              if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-          }
-          else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet())
-          {
-            if(GetEMCALGeometry())
-            {
-              Int_t absID=0;
-              
-              GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
-              
-              if( absID >= 0)
-                inacceptance = kTRUE;
-              
-              //                  if(GetEMCALGeometry()->Impact(phot1) && GetEMCALGeometry()->Impact(phot2))
-              //                    inacceptance = kTRUE;
-              if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-            else
-            {
-              if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
-                inacceptance = kTRUE ;
-              if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-          }      //In EMCAL
-          
-          //Fill histograms
-          fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY) ;
-          if(TMath::Abs(photonY) < 1.0)
-          {
-            fhEPrimMC  [kmcPPhoton]->Fill(photonE ) ;
-            fhPtPrimMC [kmcPPhoton]->Fill(photonPt) ;
-            fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta) ;
-          }
-          if(inacceptance)
-          {
-            fhEPrimMCAcc  [kmcPPhoton]->Fill(photonE ) ;
-            fhPtPrimMCAcc [kmcPPhoton]->Fill(photonPt) ;
-            fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta) ;
-            fhYPrimMCAcc  [kmcPPhoton]->Fill(photonE , photonY) ;
-          }//Accepted
-          
-          //Origin of photon
-          if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhEPrimMC[kmcPPrompt])
-          {
-            mcIndex = kmcPPrompt;
-          }
-          else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) && fhEPrimMC[kmcPFragmentation])
-          {
-            mcIndex = kmcPFragmentation ;
-          }
-          else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) && fhEPrimMC[kmcPISR])
-          {
-            mcIndex = kmcPISR;
-          }
-          else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)&& fhEPrimMC[kmcPPi0Decay])
-          {
-            mcIndex = kmcPPi0Decay;
-          }
-          else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
-                    GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) && fhEPrimMC[kmcPOtherDecay])
-          {
-            mcIndex = kmcPOtherDecay;
-          }
-          else if(fhEPrimMC[kmcPOther])
-          {
-            mcIndex = kmcPOther;
-          }//Other origin
-          
-          fhYPrimMC[mcIndex]->Fill(photonPt, photonY) ;
-          if(TMath::Abs(photonY) < 1.0)
-          {
-            fhEPrimMC  [mcIndex]->Fill(photonE ) ;
-            fhPtPrimMC [mcIndex]->Fill(photonPt) ;
-            fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta) ;
-          }
-          
-          if(inacceptance)
-          {
-            fhEPrimMCAcc  [mcIndex]->Fill(photonE ) ;
-            fhPtPrimMCAcc [mcIndex]->Fill(photonPt) ;
-            fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta) ;
-            fhYPrimMCAcc  [mcIndex]->Fill(photonE , photonY) ;
-          }//Accepted
-          
-        }// Primary photon
-      }//loop on primaries
-    }//stack exists and data is MC
-  }//read stack
-  else if(GetReader()->ReadAODMCParticles())
-  {
-    TClonesArray * mcparticles = GetReader()->GetAODMCParticles();
-    if(mcparticles)
-    {
-      Int_t nprim = mcparticles->GetEntriesFast();
-      
-      for(Int_t i=0; i < nprim; i++)
-      {
-        if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
-        
-        AliAODMCParticle * prim = (AliAODMCParticle *) mcparticles->At(i);
-        
-        pdg = prim->GetPdgCode();
-        
-        if(pdg == 22)
-        {
-          // Get tag of this particle photon from fragmentation, decay, prompt ...
-          tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader());
-          if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
-          {
-            //A conversion photon from a hadron, skip this kind of photon
-            //            printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
-            //            GetMCAnalysisUtils()->PrintMCTag(tag);
-            
-            return;
-          }
-          
-          //Get photon kinematics
-          if(prim->E() == TMath::Abs(prim->Pz()))  continue ; //Protection against floating point exception
-          
-          photonY   = 0.5*TMath::Log((prim->E()-prim->Pz())/(prim->E()+prim->Pz())) ;
-          photonE   = prim->E() ;
-          photonPt  = prim->Pt() ;
-          photonPhi = prim->Phi() ;
-          if(photonPhi < 0) photonPhi+=TMath::TwoPi();
-          photonEta = prim->Eta() ;
-          
-          //Check if photons hit the Calorimeter
-          TLorentzVector lv;
-          lv.SetPxPyPzE(prim->Px(),prim->Py(),prim->Pz(),prim->E());
-          inacceptance = kFALSE;
-          if     (fCalorimeter == "PHOS")
-          {
-            if(GetPHOSGeometry() && GetCaloUtils()->IsPHOSGeoMatrixSet())
-            {
-              Int_t mod ;
-              Double_t x,z ;
-              Double_t vtx[]={prim->Xv(),prim->Yv(),prim->Zv()};
-              if(GetPHOSGeometry()->ImpactOnEmc(vtx, prim->Theta(),prim->Phi(),mod,z,x))
-                inacceptance = kTRUE;
-              if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-            else
-            {
-              if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
-                inacceptance = kTRUE ;
-              if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-          }
-          else if(fCalorimeter == "EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet())
-          {
-            if(GetEMCALGeometry())
-            {
-              Int_t absID=0;
-              
-              GetEMCALGeometry()->GetAbsCellIdFromEtaPhi(prim->Eta(),prim->Phi(),absID);
-              
-              if( absID >= 0)
-                inacceptance = kTRUE;
-              
-              if(GetDebug() > 2) printf("In %s Real acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-            else
-            {
-              if(GetFiducialCut()->IsInFiducialCut(lv,fCalorimeter))
-                inacceptance = kTRUE ;
-              if(GetDebug() > 2) printf("In %s fiducial cut acceptance? %d\n",fCalorimeter.Data(),inacceptance);
-            }
-          }      //In EMCAL
-          
-          //Fill histograms
-          
-          fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY) ;
-          if(TMath::Abs(photonY) < 1.0)
-          {
-            fhEPrimMC  [kmcPPhoton]->Fill(photonE ) ;
-            fhPtPrimMC [kmcPPhoton]->Fill(photonPt) ;
-            fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta) ;
-          }
-          
-          if(inacceptance)
-          {
-            fhEPrimMCAcc[kmcPPhoton]  ->Fill(photonE ) ;
-            fhPtPrimMCAcc[kmcPPhoton] ->Fill(photonPt) ;
-            fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta) ;
-            fhYPrimMCAcc[kmcPPhoton]  ->Fill(photonE , photonY) ;
-          }//Accepted
-          
-          //Origin of photon
-          if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt) && fhEPrimMC[kmcPPrompt])
-          {
-            mcIndex = kmcPPrompt;
-          }
-          else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation) && fhEPrimMC[kmcPFragmentation])
-          {
-            mcIndex = kmcPFragmentation ;
-          }
-          else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR) && fhEPrimMC[kmcPISR])
-          {
-            mcIndex = kmcPISR;
-          }
-          else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay)&& fhEPrimMC[kmcPPi0Decay])
-          {
-            mcIndex = kmcPPi0Decay;
-          }
-          else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
-                    GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)) && fhEPrimMC[kmcPOtherDecay])
-          {
-            mcIndex = kmcPOtherDecay;
-          }
-          else if(fhEPrimMC[kmcPOther])
-          {
-            mcIndex = kmcPOther;
-          }//Other origin
-          
-          fhYPrimMC[mcIndex]->Fill(photonPt, photonY) ;
-          if(TMath::Abs(photonY) < 1.0)
-          {
-            fhEPrimMC  [mcIndex]->Fill(photonE ) ;
-            fhPtPrimMC [mcIndex]->Fill(photonPt) ;
-            fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta) ;
-          }
-          if(inacceptance)
-          {
-            fhEPrimMCAcc  [mcIndex]->Fill(photonE ) ;
-            fhPtPrimMCAcc [mcIndex]->Fill(photonPt) ;
-            fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi) ;
-            fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta) ;
-            fhYPrimMCAcc  [mcIndex]->Fill(photonE , photonY) ;
-          }//Accepted
-          
-        }// Primary photon
-      }//loop on primaries
-      
-    }//kmc array exists and data is MC
-  }    // read AOD MC
-}
-
-//________________________________________________________________________________________________________________
-void  AliAnaPhoton::FillEMCALTriggerClusterBCHistograms(Int_t idcalo,       Float_t ecluster,  Float_t tofcluster,
-                                                        Float_t etacluster, Float_t phicluster)
-
-{
-  // Fill trigger related histograms
-  
-  if(!fFillEMCALBCHistograms || fCalorimeter!="EMCAL") return ;
+  TParticle        * primStack = 0;
+  AliAODMCParticle * primAOD   = 0;
   
-  Float_t tofclusterUS = TMath::Abs(tofcluster);
-  
-  if(ecluster > 2)
+  // Get the ESD MC particles container
+  AliStack * stack = 0;
+  if( GetReader()->ReadStack() )
   {
-    if      (tofclusterUS < 25) fhEtaPhiEMCALBC0->Fill(etacluster, phicluster);
-    else if (tofclusterUS < 75) fhEtaPhiEMCALBC1->Fill(etacluster, phicluster);
-    else                        fhEtaPhiEMCALBCN->Fill(etacluster, phicluster);
+    stack = GetMCStack();
+    if(!stack ) return;
+    nprim = stack->GetNtrack();
   }
   
-  Int_t  bc     = GetReader()->GetTriggerClusterBC();
-  Int_t  id     = GetReader()->GetTriggerClusterId();
-  Bool_t badMax = GetReader()->IsBadMaxCellTriggerEvent();
-  
-  Int_t histoBC = bc+5;
-  if(GetReader()->AreBadTriggerEventsRemoved()) histoBC=0; // histograms created only for one BC since the others where rejected
-
-  if(id==-2)
+  // Get the AOD MC particles container
+  TClonesArray * mcparticles = 0;
+  if( GetReader()->ReadAODMCParticles() )
   {
-    //printf("AliAnaPhoton::ClusterSelected() - No trigger found bc=%d\n",bc);
-    fhEtaPhiNoTrigger->Fill(etacluster, phicluster);
-    fhTimeNoTrigger  ->Fill(ecluster, tofcluster);
+    mcparticles = GetReader()->GetAODMCParticles();
+    if( !mcparticles ) return;
+    nprim = mcparticles->GetEntriesFast();
   }
-  else if(TMath::Abs(bc) < 6)
+  
+  for(Int_t i=0 ; i < nprim; i++)
   {
-    if(!GetReader()->IsBadCellTriggerEvent() && !GetReader()->IsExoticEvent() )
-    {
-      if(GetReader()->IsTriggerMatched())
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBC[histoBC]->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBC[histoBC]->Fill(ecluster, tofcluster);
-        if(GetReader()->IsPileUpFromSPD()) fhTimeTriggerEMCALBCPileUpSPD[histoBC]->Fill(ecluster, tofcluster);
-        
-        if(idcalo ==  GetReader()->GetTriggerClusterId())
-        {
-          fhEtaPhiTriggerEMCALBCCluster[histoBC]->Fill(etacluster, phicluster);
-          fhTimeTriggerEMCALBCCluster        ->Fill(ecluster, tofcluster);
-          
-          if(bc==0)
-          {
-            Float_t threshold = GetReader()->GetEventTriggerL1Threshold() ;
-            if(GetReader()->IsEventEMCALL0()) threshold = GetReader()->GetEventTriggerL0Threshold() ;
-            
-            if(ecluster > threshold)
-              fhEtaPhiTriggerEMCALBCClusterOverTh->Fill(etacluster, phicluster);
-            else if(ecluster > threshold-1)
-              fhEtaPhiTriggerEMCALBCClusterBelowTh1->Fill(etacluster, phicluster);
-            else
-              fhEtaPhiTriggerEMCALBCClusterBelowTh2->Fill(etacluster, phicluster);
-          }
-        }
-      }
-      else
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCUM[histoBC]->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCUM[histoBC]->Fill(ecluster, tofcluster);
-        
-        if(bc==0)
-        {
-          if(GetReader()->IsTriggerMatchedOpenCuts(0)) fhTimeTriggerEMCALBC0UMReMatchOpenTime   ->Fill(ecluster, tofcluster);
-          if(GetReader()->IsTriggerMatchedOpenCuts(1)) fhTimeTriggerEMCALBC0UMReMatchCheckNeigh ->Fill(ecluster, tofcluster);
-          if(GetReader()->IsTriggerMatchedOpenCuts(2)) fhTimeTriggerEMCALBC0UMReMatchBoth       ->Fill(ecluster, tofcluster);
-        }
-        
-        if(idcalo ==  GetReader()->GetTriggerClusterId())
-        {
-          fhEtaPhiTriggerEMCALBCUMCluster[histoBC]->Fill(etacluster, phicluster);
-          fhTimeTriggerEMCALBCUMCluster->Fill(ecluster, tofcluster);
-          if(bc==0)
-          {
-            Float_t threshold = GetReader()->GetEventTriggerL1Threshold() ;
-            if(GetReader()->IsEventEMCALL0()) threshold = GetReader()->GetEventTriggerL0Threshold() ;
-            
-            if(ecluster > threshold)
-              fhEtaPhiTriggerEMCALBCUMClusterOverTh->Fill(etacluster, phicluster);
-            else if(ecluster > threshold-1)
-              fhEtaPhiTriggerEMCALBCUMClusterBelowTh1->Fill(etacluster, phicluster);
-            else
-              fhEtaPhiTriggerEMCALBCUMClusterBelowTh2->Fill(etacluster, phicluster);
-            
-            if(GetReader()->IsTriggerMatchedOpenCuts(0))
-            {
-              fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster->Fill(etacluster, phicluster);
-              fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster  ->Fill(ecluster, tofcluster);
-            }
-            if(GetReader()->IsTriggerMatchedOpenCuts(1))
-            {
-              fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster->Fill(etacluster, phicluster);
-              fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster  ->Fill(ecluster, tofcluster);
-            }
-            if(GetReader()->IsTriggerMatchedOpenCuts(2))
-            {
-              fhEtaPhiTriggerEMCALBCUMReMatchBothCluster->Fill(etacluster, phicluster);
-              fhTimeTriggerEMCALBCUMReMatchBothCluster  ->Fill(ecluster, tofcluster);
-            }
-            
-          }
-        }
-      }
-    }// neither bad nor exotic
-    else if(GetReader()->IsBadCellTriggerEvent() && GetReader()->IsExoticEvent())
-    {
-      if(GetReader()->IsTriggerMatched())
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCBadExotic->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCBadExotic->Fill(ecluster, tofcluster);
-        if(badMax)  fhTimeTriggerEMCALBCBadMaxCellExotic->Fill(ecluster, tofcluster);
-      }
-      else
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCUMBadExotic->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCUMBadExotic->Fill(ecluster, tofcluster);
-        if(badMax)  fhTimeTriggerEMCALBCUMBadMaxCellExotic->Fill(ecluster, tofcluster);
-        
-      }
-    }// Bad and exotic cluster trigger
-    else if(GetReader()->IsBadCellTriggerEvent() )
-    {
-      if(GetReader()->IsTriggerMatched())
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCBad->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCBad->Fill(ecluster, tofcluster);
-        if(badMax)  fhTimeTriggerEMCALBCBadMaxCell->Fill(ecluster, tofcluster);
-      }
-      else
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCUMBad->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCUMBad->Fill(ecluster, tofcluster);
-        if(badMax)  fhTimeTriggerEMCALBCUMBadMaxCell->Fill(ecluster, tofcluster);
-      }
-    }// Bad cluster trigger
-    else if(GetReader()->IsExoticEvent() )
+    if(GetReader()->AcceptOnlyHIJINGLabels() && !GetReader()->IsHIJINGLabel(i)) continue ;
+    
+    if(GetReader()->ReadStack())
     {
-      if(GetReader()->IsTriggerMatched())
-      {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCExotic->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCExotic->Fill(ecluster, tofcluster);
-      }
-      else
+      primStack = stack->Particle(i) ;
+      if(!primStack)
       {
-        if(ecluster > 2) fhEtaPhiTriggerEMCALBCUMExotic->Fill(etacluster, phicluster);
-        fhTimeTriggerEMCALBCUMExotic->Fill(ecluster, tofcluster);
+        AliWarning("ESD primaries pointer not available!!");
+        continue;
       }
-    }
-  }
-  else if(TMath::Abs(bc) >= 6)
-    printf("AliAnaPhoton::ClusterSelected() - Trigger BC not expected = %d\n",bc);
-  
-}
-
-//_________________________________________________________________________________________________________
-void  AliAnaPhoton::FillClusterPileUpHistograms(AliVCluster * calo, Bool_t matched,     Float_t ptcluster,
-                                                Float_t etacluster, Float_t phicluster, Float_t l0cluster)
-{
-  // Fill some histograms related to pile up before any cluster cut is applied
-  
-  if(!fFillPileUpHistograms) return ;
-  
-  // Get the fraction of the cluster energy that carries the cell with highest energy and its absId
-  AliVCaloCells* cells = 0;
-  if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
-  else                        cells = GetPHOSCells();
-  
-  Float_t maxCellFraction = 0.;
-  Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, calo,maxCellFraction);
-  
-  Double_t tmax  = cells->GetCellTime(absIdMax);
-  GetCaloUtils()->RecalibrateCellTime(tmax, fCalorimeter, absIdMax,GetReader()->GetInputEvent()->GetBunchCrossNumber());
-  tmax*=1.e9;
-  
-  Bool_t okPhoton = kFALSE;
-  if( GetCaloPID()->GetIdentifiedParticleType(calo)== AliCaloPID::kPhoton) okPhoton = kTRUE;
-  
-  Float_t clusterLongTimePt = 0;
-  Float_t clusterOKTimePt   = 0;
-    
-  //Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
-  if(cells->GetCellAmplitude(absIdMax) > 0.1 && TMath::Abs(tmax) < 30)
-  {
-    for (Int_t ipos = 0; ipos < calo->GetNCells(); ipos++)
-    {
-      Int_t absId  = calo->GetCellsAbsId()[ipos];
-      
-      if( absId == absIdMax ) continue ;
-      
-      Double_t time  = cells->GetCellTime(absId);
-      Float_t  amp   = cells->GetCellAmplitude(absId);
-      Int_t    bc    = GetReader()->GetInputEvent()->GetBunchCrossNumber();
-      GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,time,cells);
-      time*=1e9;
       
-      Float_t diff = (tmax-time);
+      pdg    = primStack->GetPdgCode();
+      status = primStack->GetStatusCode();
       
-      if(GetReader()->IsInTimeWindow(time,amp)) clusterOKTimePt   += amp;
-      else                                      clusterLongTimePt += amp;
+      if(primStack->Energy() == TMath::Abs(primStack->Pz()))  continue ; //Protection against floating point exception
       
-      if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
+      //printf("i %d, %s %d  %s %d \n",i, stack->Particle(i)->GetName(), stack->Particle(i)->GetPdgCode(),
+      //       prim->GetName(), prim->GetPdgCode());
       
-      if(GetReader()->IsPileUpFromSPD())
-      {
-        fhClusterCellTimePileUp[0]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[0]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[0]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[0]->Fill(ptcluster, diff);
-        }
-      }
+      //Photon kinematics
+      primStack->Momentum(fMomentum);
       
-      if(GetReader()->IsPileUpFromEMCal())
-      {
-        fhClusterCellTimePileUp[1]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[1]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[1]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[1]->Fill(ptcluster, diff);
-        }
-      }
-      
-      if(GetReader()->IsPileUpFromSPDOrEMCal())
-      {
-        fhClusterCellTimePileUp[2]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[2]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[2]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[2]->Fill(ptcluster, diff);
-        }
-      }
-      
-      if(GetReader()->IsPileUpFromSPDAndEMCal())
+      photonY = 0.5*TMath::Log((primStack->Energy()+primStack->Pz())/(primStack->Energy()-primStack->Pz())) ;
+    }
+    else
+    {
+      primAOD = (AliAODMCParticle *) mcparticles->At(i);
+      if(!primAOD)
       {
-        fhClusterCellTimePileUp[3]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[3]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[3]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[3]->Fill(ptcluster, diff);
-        }
+        AliWarning("AOD primaries pointer not available!!");
+        continue;
       }
       
-      if(GetReader()->IsPileUpFromSPDAndNotEMCal())
-      {
-        fhClusterCellTimePileUp[4]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[4]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[4]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[4]->Fill(ptcluster, diff);
-        }
-      }
+      pdg    = primAOD->GetPdgCode();
+      status = primAOD->GetStatus();
       
-      if(GetReader()->IsPileUpFromEMCalAndNotSPD())
-      {
-        fhClusterCellTimePileUp[5]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[5]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[5]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[5]->Fill(ptcluster, diff);
-        }
-      }
+      if(primAOD->E() == TMath::Abs(primAOD->Pz()))  continue ; //Protection against floating point exception
       
-      if(GetReader()->IsPileUpFromNotSPDAndNotEMCal())
-      {
-        fhClusterCellTimePileUp[6]->Fill(ptcluster, time);
-        fhClusterTimeDiffPileUp[6]->Fill(ptcluster, diff);
-        if(!matched)
-        {
-          fhClusterTimeDiffChargedPileUp[6]->Fill(ptcluster, diff);
-          if(okPhoton)  fhClusterTimeDiffPhotonPileUp[6]->Fill(ptcluster, diff);
-        }
-      }
-    }//loop
+      //Photon kinematics
+      fMomentum.SetPxPyPzE(primAOD->Px(),primAOD->Py(),primAOD->Pz(),primAOD->E());
+
+      photonY = 0.5*TMath::Log((primAOD->E()+primAOD->Pz())/(primAOD->E()-primAOD->Pz())) ;
+    }
+
+    // Select only photons in the final state
+    if(pdg != 22 ) continue ;
     
+    // If too small or too large pt, skip, same cut as for data analysis
+    photonPt  = fMomentum.Pt () ;
     
-    Float_t frac = 0;
-    if(clusterLongTimePt+clusterOKTimePt > 0.001)
-      frac = clusterLongTimePt/(clusterLongTimePt+clusterOKTimePt);
-    //printf("E long %f, E OK %f, Fraction large time %f, E %f\n",clusterLongTimePt,clusterOKTimePt,frac,ptcluster);
+    if(photonPt < GetMinPt() || photonPt > GetMaxPt() ) continue ;
     
-    if(GetReader()->IsPileUpFromSPD())               {fhPtPileUp[0]->Fill(ptcluster); fhLambda0PileUp[0]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[0]->Fill(ptcluster,frac);}
-    if(GetReader()->IsPileUpFromEMCal())             {fhPtPileUp[1]->Fill(ptcluster); fhLambda0PileUp[1]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[1]->Fill(ptcluster,frac);}
-    if(GetReader()->IsPileUpFromSPDOrEMCal())        {fhPtPileUp[2]->Fill(ptcluster); fhLambda0PileUp[2]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[2]->Fill(ptcluster,frac);}
-    if(GetReader()->IsPileUpFromSPDAndEMCal())       {fhPtPileUp[3]->Fill(ptcluster); fhLambda0PileUp[3]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[3]->Fill(ptcluster,frac);}
-    if(GetReader()->IsPileUpFromSPDAndNotEMCal())    {fhPtPileUp[4]->Fill(ptcluster); fhLambda0PileUp[4]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[4]->Fill(ptcluster,frac);}
-    if(GetReader()->IsPileUpFromEMCalAndNotSPD())    {fhPtPileUp[5]->Fill(ptcluster); fhLambda0PileUp[5]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[5]->Fill(ptcluster,frac);}
-    if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) {fhPtPileUp[6]->Fill(ptcluster); fhLambda0PileUp[6]->Fill(ptcluster,l0cluster); fhClusterEFracLongTimePileUp[6]->Fill(ptcluster,frac);}
+    photonE   = fMomentum.E  () ;
+    photonEta = fMomentum.Eta() ;
+    photonPhi = fMomentum.Phi() ;
     
-    fhEtaPhiBC0->Fill(etacluster,phicluster);
-    if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBC0PileUpSPD    ->Fill(etacluster,phicluster);
-  }
-  
-  else if (tmax > 25)         {fhEtaPhiBCPlus ->Fill(etacluster,phicluster); if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBCPlusPileUpSPD ->Fill(etacluster,phicluster); }
-  else if (tmax <-25)         {fhEtaPhiBCMinus->Fill(etacluster,phicluster); if(GetReader()->IsPileUpFromSPD()) fhEtaPhiBCMinusPileUpSPD->Fill(etacluster,phicluster); }
-}
-
-//_______________________________________________
-void AliAnaPhoton::FillPileUpHistogramsPerEvent()
-{
-  // Fill some histograms per event to understand pile-up
-  // Open the time cut in the reader to be more meaningful
-  
-  if(!fFillPileUpHistograms) return;
-  
-  AliVEvent * event = GetReader()->GetInputEvent();
-       
-  AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
-  AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
-       
-  // N pile up vertices
-  Int_t nVtxSPD = -1;
-  Int_t nVtxTrk = -1;
-  TLorentzVector mom;
-       
-  if      (esdEv)
-  {
-               nVtxSPD = esdEv->GetNumberOfPileupVerticesSPD();
-               nVtxTrk = esdEv->GetNumberOfPileupVerticesTracks();
-  }//ESD
-  else if (aodEv)
-  {
-               nVtxSPD = aodEv->GetNumberOfPileupVerticesSPD();
-               nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
-  }//AOD
-       
-       
-       // Get the appropriate list of clusters
-       TClonesArray * clusterList = 0;
-       TString  clusterListName   = GetReader()->GetEMCALClusterListName();
-       if     (event->FindListObject(clusterListName))
-               clusterList = dynamic_cast<TClonesArray*> (event->FindListObject(clusterListName));
-       else if(GetReader()->GetOutputEvent())
-               clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetOutputEvent()->FindListObject(clusterListName));
-       
-       // Loop on clusters, get the maximum energy cluster as reference
-  Int_t nclusters = 0;
-       if(clusterList) nclusters = clusterList->GetEntriesFast();
-       else            nclusters = event->GetNumberOfCaloClusters();
-       
-  Int_t   idMax = 0;
-  Float_t  eMax = 0;
-  Float_t  tMax = 0;
-  for(Int_t iclus = 0; iclus < nclusters ; iclus++)
-  {
-               AliVCluster * clus = 0;
-               if(clusterList) clus = (AliVCluster*) (clusterList->At(iclus));
-    else            clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
-               
-               if(!clus)            continue;
-               
-               if(!clus->IsEMCAL()) continue;
-               
-               Float_t tof = clus->GetTOF()*1e9;
-               if(clus->E() > eMax && TMath::Abs(tof) < 30)
+    if(photonPhi < 0) photonPhi+=TMath::TwoPi();
+    
+    // Check if photons hit desired acceptance
+    inacceptance = kTRUE;
+    
+    // Check same fidutial borders as in data analysis on top of real acceptance if real was requested.
+    if( IsFiducialCutOn() && !GetFiducialCut()->IsInFiducialCut(fMomentum.Eta(),fMomentum.Phi(),GetCalorimeter())) inacceptance = kFALSE ;
+    
+    // Check if photons hit the Calorimeter acceptance
+    if(IsRealCaloAcceptanceOn()) // defined on base class
     {
-      eMax  = clus->E();
-                       tMax  = tof;
-      idMax = iclus;
+      if(GetReader()->ReadStack()          &&
+         !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primStack)) inacceptance = kFALSE ;
+      if(GetReader()->ReadAODMCParticles() &&
+         !GetCaloUtils()->IsMCParticleInCalorimeterAcceptance(GetCalorimeter(), primAOD  )) inacceptance = kFALSE ;
     }
-         
-               clus->GetMomentum(mom,GetVertex(0));
-               Float_t pt = mom.Pt();
-         
-               fhPtNPileUpSPDVtx->Fill(pt,nVtxSPD);
-               fhPtNPileUpTrkVtx->Fill(pt,nVtxTrk);
     
-               if(TMath::Abs(tof) < 30)
-               {
-                       fhPtNPileUpSPDVtxTimeCut->Fill(pt,nVtxSPD);
-                       fhPtNPileUpTrkVtxTimeCut->Fill(pt,nVtxTrk);
-               }
+    // Get tag of this particle photon from fragmentation, decay, prompt ...
+    // Set the origin of the photon.
+    tag = GetMCAnalysisUtils()->CheckOrigin(i,GetReader(),GetCalorimeter());
     
-    if(tof < 75 && tof > -30)
+    if(!GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton))
     {
-      fhPtNPileUpSPDVtxTimeCut2->Fill(pt,nVtxSPD);
-      fhPtNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
+      // A conversion photon from a hadron, skip this kind of photon
+      // printf("AliAnaPhoton::FillAcceptanceHistograms() - not a photon, weird!\n ");
+      // GetMCAnalysisUtils()->PrintMCTag(tag);
+      
+      continue;
     }
     
-    fhTimePtNoCut->Fill(pt,tof);
-    if(GetReader()->IsPileUpFromSPD()) fhTimePtSPD->Fill(pt,tof);
-
-  }
-       
-  if(eMax < 5) return;
-  
-  // Loop again on clusters to compare this max cluster t and the rest of the clusters, if E > 0.3
-  Int_t n20  = 0;
-  Int_t n40  = 0;
-  Int_t n    = 0;
-  Int_t nOK  = 0;
-  
-  for(Int_t iclus = 0; iclus < nclusters ; iclus++)
-  {
-               AliVCluster * clus = 0;
-               if(clusterList) clus = (AliVCluster*) (clusterList->At(iclus));
-    else            clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
-               
-               if(!clus)            continue;
-               
-               if(!clus->IsEMCAL()) continue;
+    // Consider only final state particles, but this depends on generator,
+    // status 1 is the usual one, in case of not being ok, leave the possibility
+    // to not consider this.
+    if(status > 1) continue ; // Avoid "partonic" photons
     
-    if(clus->E() < 0.3 || iclus==idMax) continue;
+    Bool_t takeIt  = kFALSE ;
+    if(status == 1 && GetMCAnalysisUtils()->GetMCGenerator() != AliMCAnalysisUtils::kBoxLike ) takeIt = kTRUE ;
+
+    if      (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)) continue;
     
-    Float_t tdiff = TMath::Abs(tMax-clus->GetTOF()*1e9);
-    n++;
-    if(tdiff < 25) nOK++;
+    //Origin of photon
+    if     (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPrompt))
+    {
+      mcIndex = kmcPPrompt;
+    }
+    else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCFragmentation))
+    {
+      mcIndex = kmcPFragmentation ;
+    }
+    else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
+    {
+      mcIndex = kmcPISR;
+    }
+    else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0Decay))
+    {
+      mcIndex = kmcPPi0Decay;
+    }
+    else if( (GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) ||
+              GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay)))
+    {
+      mcIndex = kmcPOtherDecay;
+    }
     else
     {
-      n20++;
-      if(tdiff > 40 ) n40++;
+      // Other decay but from non final state particle
+      mcIndex = kmcPOtherDecay;
+    }//Other origin
+    
+    if(!takeIt &&  (mcIndex == kmcPPi0Decay || mcIndex == kmcPOtherDecay)) takeIt = kTRUE ;
+
+    if(!takeIt) continue ;
+      
+    //Fill histograms for all photons
+    fhYPrimMC[kmcPPhoton]->Fill(photonPt, photonY) ;
+    if(TMath::Abs(photonY) < 1.0)
+    {
+      fhEPrimMC  [kmcPPhoton]->Fill(photonE ) ;
+      fhPtPrimMC [kmcPPhoton]->Fill(photonPt) ;
+      fhPhiPrimMC[kmcPPhoton]->Fill(photonE , photonPhi) ;
+      fhEtaPrimMC[kmcPPhoton]->Fill(photonE , photonEta) ;
     }
-  }
-  
-  // Check pile-up and fill histograms depending on the different cluster multiplicities
-  if(GetReader()->IsPileUpFromSPD())
-  {
-    fhClusterMultSPDPileUp[0]->Fill(eMax,n  );
-    fhClusterMultSPDPileUp[1]->Fill(eMax,nOK);
-    fhClusterMultSPDPileUp[2]->Fill(eMax,n20);
-    fhClusterMultSPDPileUp[3]->Fill(eMax,n40);
-  }
-  else
-  {
-    fhClusterMultNoPileUp[0]->Fill(eMax,n  );
-    fhClusterMultNoPileUp[1]->Fill(eMax,nOK);
-    fhClusterMultNoPileUp[2]->Fill(eMax,n20);
-    fhClusterMultNoPileUp[3]->Fill(eMax,n40);
-  }
-  
-}
+    
+    if(inacceptance)
+    {
+      fhEPrimMCAcc  [kmcPPhoton]->Fill(photonE ) ;
+      fhPtPrimMCAcc [kmcPPhoton]->Fill(photonPt) ;
+      fhPhiPrimMCAcc[kmcPPhoton]->Fill(photonE , photonPhi) ;
+      fhEtaPrimMCAcc[kmcPPhoton]->Fill(photonE , photonEta) ;
+      fhYPrimMCAcc  [kmcPPhoton]->Fill(photonE , photonY) ;
+    }//Accepted
+    
+    //Fill histograms for photons origin
+    if(mcIndex < fNPrimaryHistograms)
+    {
+      fhYPrimMC[mcIndex]->Fill(photonPt, photonY) ;
+      if(TMath::Abs(photonY) < 1.0)
+      {
+        fhEPrimMC  [mcIndex]->Fill(photonE ) ;
+        fhPtPrimMC [mcIndex]->Fill(photonPt) ;
+        fhPhiPrimMC[mcIndex]->Fill(photonE , photonPhi) ;
+        fhEtaPrimMC[mcIndex]->Fill(photonE , photonEta) ;
+      }
+      
+      if(inacceptance)
+      {
+        fhEPrimMCAcc  [mcIndex]->Fill(photonE ) ;
+        fhPtPrimMCAcc [mcIndex]->Fill(photonPt) ;
+        fhPhiPrimMCAcc[mcIndex]->Fill(photonE , photonPhi) ;
+        fhEtaPrimMCAcc[mcIndex]->Fill(photonE , photonEta) ;
+        fhYPrimMCAcc  [mcIndex]->Fill(photonE , photonY) ;
+      }//Accepted
+    }
+  }//loop on primaries
 
+}
 
-//_________________________________________________________________________________________________
-void AliAnaPhoton::FillPileUpHistograms(Float_t energy, Float_t pt, Float_t time)
+//________________________________________________________________________________
+void AliAnaPhoton::FillPileUpHistograms(AliVCluster* cluster, AliVCaloCells *cells,
+                                        Int_t absIdMax)
 {
   // Fill some histograms to understand pile-up
-  if(!fFillPileUpHistograms) return;
   
-  //printf("E %f, time %f\n",energy,time);
+  Float_t pt   = fMomentum.Pt();
+  Float_t time = cluster->GetTOF()*1.e9;
+  
   AliVEvent * event = GetReader()->GetInputEvent();
   
   if(GetReader()->IsPileUpFromSPD())               fhPtPhotonPileUp[0]->Fill(pt);
@@ -1181,7 +557,38 @@ void AliAnaPhoton::FillPileUpHistograms(Float_t energy, Float_t pt, Float_t time
   fhTimePtPhotonNoCut->Fill(pt,time);
   if(GetReader()->IsPileUpFromSPD()) fhTimePtPhotonSPD->Fill(pt,time);
   
-  if(energy < 8) return; // Fill time figures for high energy clusters not too close to trigger threshold
+  // cells inside the cluster
+
+  //Loop on cells inside cluster, max cell must be over 100 MeV and time in BC=0
+  if(cells->GetCellAmplitude(absIdMax) > 0.1 && TMath::Abs(time) < 30)
+  {
+    for (Int_t ipos = 0; ipos < cluster->GetNCells(); ipos++)
+    {
+      Int_t absId  = cluster->GetCellsAbsId()[ipos];
+      
+      if( absId == absIdMax ) continue ;
+      
+      Double_t tcell = cells->GetCellTime(absId);
+      Float_t  amp   = cells->GetCellAmplitude(absId);
+      Int_t    bc    = GetReader()->GetInputEvent()->GetBunchCrossNumber();
+      
+      GetCaloUtils()->GetEMCALRecoUtils()->AcceptCalibrateCell(absId,bc,amp,tcell,cells);
+      tcell*=1e9;
+      
+      Float_t diff = (time-tcell);
+      
+      if( cells->GetCellAmplitude(absIdMax) < 0.1 ) continue ;
+      
+      if(GetReader()->IsPileUpFromSPD())               fhClusterTimeDiffPhotonPileUp[0]->Fill(pt, diff);
+      if(GetReader()->IsPileUpFromEMCal())             fhClusterTimeDiffPhotonPileUp[1]->Fill(pt, diff);
+      if(GetReader()->IsPileUpFromSPDOrEMCal())        fhClusterTimeDiffPhotonPileUp[2]->Fill(pt, diff);
+      if(GetReader()->IsPileUpFromSPDAndEMCal())       fhClusterTimeDiffPhotonPileUp[3]->Fill(pt, diff);
+      if(GetReader()->IsPileUpFromSPDAndNotEMCal())    fhClusterTimeDiffPhotonPileUp[4]->Fill(pt, diff);
+      if(GetReader()->IsPileUpFromEMCalAndNotSPD())    fhClusterTimeDiffPhotonPileUp[5]->Fill(pt, diff);
+      if(GetReader()->IsPileUpFromNotSPDAndNotEMCal()) fhClusterTimeDiffPhotonPileUp[6]->Fill(pt, diff);
+
+    }//loop
+  }
   
   AliESDEvent* esdEv = dynamic_cast<AliESDEvent*> (event);
   AliAODEvent* aodEv = dynamic_cast<AliAODEvent*> (event);
@@ -1202,8 +609,11 @@ void AliAnaPhoton::FillPileUpHistograms(Float_t energy, Float_t pt, Float_t time
     nVtxTrk = aodEv->GetNumberOfPileupVerticesTracks();
   }//AOD
   
-  fhTimeNPileUpVertSPD  ->Fill(time,nVtxSPD);
-  fhTimeNPileUpVertTrack->Fill(time,nVtxTrk);
+  if(pt < 8)
+  {
+    fhTimeNPileUpVertSPD  ->Fill(time,nVtxSPD);
+    fhTimeNPileUpVertTrack->Fill(time,nVtxTrk);
+  }
   
   fhPtPhotonNPileUpSPDVtx->Fill(pt,nVtxSPD);
   fhPtPhotonNPileUpTrkVtx->Fill(pt,nVtxTrk);
@@ -1217,47 +627,14 @@ void AliAnaPhoton::FillPileUpHistograms(Float_t energy, Float_t pt, Float_t time
   if(time < 75 && time > -25)
   {
     fhPtPhotonNPileUpSPDVtxTimeCut2->Fill(pt,nVtxSPD);
-    fhPtPhotonNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
-  }
-  
-  //printf("Is SPD %d, Is SPD Multi %d, n spd %d, n track %d\n",
-  //       GetReader()->IsPileUpFromSPD(),event->IsPileupFromSPDInMultBins(),nVtxSPD,nVtxTrk);
-  
-  Int_t ncont = -1;
-  Float_t z1 = -1, z2 = -1;
-  Float_t diamZ = -1;
-  for(Int_t iVert=0; iVert<nVtxSPD;iVert++)
-  {
-    if      (esdEv)
-    {
-      const AliESDVertex* pv=esdEv->GetPileupVertexSPD(iVert);
-      ncont=pv->GetNContributors();
-      z1 = esdEv->GetPrimaryVertexSPD()->GetZ();
-      z2 = pv->GetZ();
-      diamZ = esdEv->GetDiamondZ();
-    }//ESD
-    else if (aodEv)
-    {
-      AliAODVertex *pv=aodEv->GetVertex(iVert);
-      if(pv->GetType()!=AliAODVertex::kPileupSPD) continue;
-      ncont=pv->GetNContributors();
-      z1=aodEv->GetPrimaryVertexSPD()->GetZ();
-      z2=pv->GetZ();
-      diamZ = aodEv->GetDiamondZ();
-    }// AOD
-    
-    Double_t distZ  = TMath::Abs(z2-z1);
-    diamZ  = TMath::Abs(z2-diamZ);
-    
-    fhTimeNPileUpVertContributors  ->Fill(time,ncont);
-    fhTimePileUpMainVertexZDistance->Fill(time,distZ);
-    fhTimePileUpMainVertexZDiamond ->Fill(time,diamZ);
-    
-  }// loop
+    fhPtPhotonNPileUpTrkVtxTimeCut2->Fill(pt,nVtxTrk);
+  }
+  
 }
 
-//____________________________________________________________________________________
-void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
+//_______________________________________________________________________________
+void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster,
+                                              Int_t mcTag, Int_t maxCellFraction)
 {
   //Fill cluster Shower Shape histograms
   
@@ -1269,26 +646,16 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
   Float_t lambda1 = cluster->GetM20();
   Float_t disp    = cluster->GetDispersion()*cluster->GetDispersion();
   
-  TLorentzVector mom;
-  if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
-  {
-    cluster->GetMomentum(mom,GetVertex(0)) ;
-  }//Assume that come from vertex in straight line
-  else
-  {
-    Double_t vertex[]={0,0,0};
-    cluster->GetMomentum(mom,vertex) ;
-  }
-  
-  Float_t eta = mom.Eta();
-  Float_t phi = mom.Phi();
+  Float_t eta = fMomentum.Eta();
+  Float_t phi = fMomentum.Phi();
   if(phi < 0) phi+=TMath::TwoPi();
   
   fhLam0E ->Fill(energy,lambda0);
   fhLam1E ->Fill(energy,lambda1);
   fhDispE ->Fill(energy,disp);
   
-  if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5)
+  if(GetCalorimeter() == kEMCAL &&  GetFirstSMCoveredByTRD() >= 0 &&
+     GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD()  )
   {
     fhLam0ETRD->Fill(energy,lambda0);
     fhLam1ETRD->Fill(energy,lambda1);
@@ -1298,7 +665,7 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
   Float_t l0   = 0., l1   = 0.;
   Float_t dispp= 0., dEta = 0., dPhi    = 0.;
   Float_t sEta = 0., sPhi = 0., sEtaPhi = 0.;
-  if(fCalorimeter == "EMCAL" && !fFillOnlySimpleSSHisto)
+  if(GetCalorimeter() == kEMCAL && !fFillOnlySimpleSSHisto)
   {
     GetCaloUtils()->GetEMCALRecoUtils()->RecalculateClusterShowerShapeParameters(GetEMCALGeometry(), GetReader()->GetInputEvent()->GetEMCALCells(), cluster,
                                                                                  l0, l1, dispp, dEta, dPhi, sEta, sPhi, sEtaPhi);
@@ -1349,7 +716,8 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
       fhLam1ETM ->Fill(energy,lambda1);
       fhDispETM ->Fill(energy,disp);
       
-      if(fCalorimeter == "EMCAL" && GetModuleNumber(cluster) > 5)
+      if(GetCalorimeter() == kEMCAL &&  GetFirstSMCoveredByTRD()   >= 0 &&
+         GetModuleNumber(cluster) >= GetFirstSMCoveredByTRD()  )
       {
         fhLam0ETMTRD->Fill(energy,lambda0);
         fhLam1ETMTRD->Fill(energy,lambda1);
@@ -1359,7 +727,8 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
   }// if track-matching was of, check effect of matching residual cut
   
   
-  if(!fFillOnlySimpleSSHisto){
+  if(!fFillOnlySimpleSSHisto)
+  {
     if(energy < 2)
     {
       fhNCellsLam0LowE ->Fill(ncells,lambda0);
@@ -1385,10 +754,11 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
       fhPhiLam0HighE   ->Fill(phi, lambda0);
     }
   }
+  
   if(IsDataMC())
   {
     AliVCaloCells* cells = 0;
-    if(fCalorimeter == "EMCAL") cells = GetEMCALCells();
+    if(GetCalorimeter() == kEMCAL) cells = GetEMCALCells();
     else                        cells = GetPHOSCells();
     
     //Fill histograms to check shape of embedded clusters
@@ -1410,19 +780,12 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
       //Fraction of total energy due to the embedded signal
       fraction/=clusterE;
       
-      if(GetDebug() > 1 )
-        printf("AliAnaPhoton::FillShowerShapeHistogram() - Energy fraction of embedded signal %2.3f, Energy %2.3f\n",fraction, clusterE);
+      AliDebug(1,Form("Energy fraction of embedded signal %2.3f, Energy %2.3f",fraction, clusterE));
       
       fhEmbeddedSignalFractionEnergy->Fill(clusterE,fraction);
       
     }  // embedded fraction
-    
-    // Get the fraction of the cluster energy that carries the cell with highest energy
-    Int_t   absID           =-1 ;
-    Float_t maxCellFraction = 0.;
-    
-    absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
-    
+        
     // Check the origin and fill histograms
     
     Int_t mcIndex = -1;
@@ -1460,7 +823,7 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
         }
         else
         {
-          printf("AliAnaPhoton::FillShowerShapeHistogram() - n overlaps = %d!!\n", noverlaps);
+          AliWarning(Form("n overlaps = %d!!", noverlaps));
         }
       }//No embedding
       
@@ -1557,7 +920,7 @@ void  AliAnaPhoton::FillShowerShapeHistograms(AliVCluster* cluster, Int_t mcTag)
         fhMCNCellsvsClusterMaxCellDiffE6 [mcIndex]->Fill(ncells,  maxCellFraction);
       }
       
-      if(fCalorimeter == "EMCAL")
+      if(GetCalorimeter() == kEMCAL)
       {
         fhMCEDispEta        [mcIndex]-> Fill(energy,dEta);
         fhMCEDispPhi        [mcIndex]-> Fill(energy,dPhi);
@@ -1628,7 +991,8 @@ void AliAnaPhoton::FillTrackMatchingResidualHistograms(AliVCluster* cluster,
     
     Int_t nSMod = GetModuleNumber(cluster);
     
-    if(fCalorimeter=="EMCAL" &&  nSMod > 5)
+    if(GetCalorimeter()==kEMCAL &&   GetFirstSMCoveredByTRD() >= 0 &&
+       nSMod >= GetFirstSMCoveredByTRD()   )
     {
       fhTrackMatchedDEtaTRD[cut]->Fill(cluster->E(),dZ);
       fhTrackMatchedDPhiTRD[cut]->Fill(cluster->E(),dR);
@@ -1646,20 +1010,19 @@ void AliAnaPhoton::FillTrackMatchingResidualHistograms(AliVCluster* cluster,
         fhdEdx[cut]  ->Fill(cluster->E(), dEdx);
         fhEOverP[cut]->Fill(cluster->E(), eOverp);
         
-        if(fCalorimeter=="EMCAL" && nSMod > 5)
+        if(GetCalorimeter()==kEMCAL &&  GetFirstSMCoveredByTRD() >= 0 &&
+           nSMod >= GetFirstSMCoveredByTRD()  )
           fhEOverPTRD[cut]->Fill(cluster->E(), eOverp);
         
         
       }
       else
-        printf("AliAnaPhoton::FillTrackMatchingResidualHistograms() - Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT? \n", dR,dZ);
-      
-      
+        AliWarning(Form("Residual OK but (dR, dZ)= (%2.4f,%2.4f) no track associated WHAT?", dR,dZ));
       
       if(IsDataMC())
       {
         
-        Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(),GetReader());
+        Int_t tag = GetMCAnalysisUtils()->CheckOrigin(cluster->GetLabels(),cluster->GetNLabels(),GetReader(),GetCalorimeter());
         
         if  ( !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion)  )
         {
@@ -1722,17 +1085,17 @@ TObjString *  AliAnaPhoton::GetAnalysisCuts()
   const Int_t buffersize = 255;
   char onePar[buffersize] ;
   
-  snprintf(onePar,buffersize,"--- AliAnaPhoton ---\n") ;
+  snprintf(onePar,buffersize,"--- AliAnaPhoton ---:") ;
   parList+=onePar ;
-  snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
+  snprintf(onePar,buffersize,"Calorimeter: %s;",GetCalorimeterString().Data()) ;
   parList+=onePar ;
-  snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
+  snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster);",fMinDist) ;
   parList+=onePar ;
-  snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
+  snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation);",fMinDist2) ;
   parList+=onePar ;
-  snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
+  snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study);",fMinDist3) ;
   parList+=onePar ;
-  snprintf(onePar,buffersize,"fRejectTrackMatch: %d\n",fRejectTrackMatch) ;
+  snprintf(onePar,buffersize,"fRejectTrackMatch: %d",fRejectTrackMatch) ;
   parList+=onePar ;
   
   //Get parameters set in base class.
@@ -1778,534 +1141,112 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
   
   Int_t bin[] = {0,2,4,6,10,15,20,100}; // energy bins for SS studies
   
-  Int_t nTrigBC  = 1;
-  Int_t iBCShift = 0;
-  if(!GetReader()->AreBadTriggerEventsRemoved())
-  {
-    nTrigBC = 11;
-    iBCShift = 5;
-  }
-  
   TString cut[] = {"Open","Reader","E","Time","NCells","NLM","Fidutial","Matching","Bad","PID"};
   for (Int_t i = 0; i < 10 ;  i++)
   {
-    fhClusterCuts[i] = new TH1F(Form("hCut_%d_%s", i, cut[i].Data()),
+    fhClusterCutsE[i] = new TH1F(Form("hE_Cut_%d_%s", i, cut[i].Data()),
+                                Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
+                                nptbins,ptmin,ptmax);
+    fhClusterCutsE[i]->SetYTitle("d#it{N}/d#it{E} ");
+    fhClusterCutsE[i]->SetXTitle("#it{E} (GeV)");
+    outputContainer->Add(fhClusterCutsE[i]) ;
+    
+    fhClusterCutsPt[i] = new TH1F(Form("hPt_Cut_%d_%s", i, cut[i].Data()),
                                 Form("Number of clusters that pass cuts <= %d, %s", i, cut[i].Data()),
                                 nptbins,ptmin,ptmax);
-    fhClusterCuts[i]->SetYTitle("dN/dE ");
-    fhClusterCuts[i]->SetXTitle("E (GeV)");
-    outputContainer->Add(fhClusterCuts[i]) ;
+    fhClusterCutsPt[i]->SetYTitle("d#it{N}/d#it{E} ");
+    fhClusterCutsPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
+    outputContainer->Add(fhClusterCutsPt[i]) ;
   }
   
+  fhEClusterSM = new TH2F("hEClusterSM","Raw clusters E and super-module number",
+                          nptbins,ptmin,ptmax,
+                          GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
+  fhEClusterSM->SetYTitle("SuperModule ");
+  fhEClusterSM->SetXTitle("#it{E} (GeV)");
+  outputContainer->Add(fhEClusterSM) ;
+
+  fhPtClusterSM = new TH2F("hPtClusterSM","Raw clusters #it{p}_{T} and super-module number",
+                          nptbins,ptmin,ptmax,
+                          GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
+  fhPtClusterSM->SetYTitle("SuperModule ");
+  fhPtClusterSM->SetXTitle("#it{E} (GeV)");
+  outputContainer->Add(fhPtClusterSM) ;
+  
+  fhEPhotonSM = new TH2F("hEPhotonSM","Selected clusters E and super-module number",
+                          nptbins,ptmin,ptmax,
+                          GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
+  fhEPhotonSM->SetYTitle("SuperModule ");
+  fhEPhotonSM->SetXTitle("#it{E} (GeV)");
+  outputContainer->Add(fhEPhotonSM) ;
+  
+  fhPtPhotonSM = new TH2F("hPtPhotonSM","Selected clusters #it{p}_{T} and super-module number",
+                           nptbins,ptmin,ptmax,
+                           GetCaloUtils()->GetNumberOfSuperModulesUsed(),0,GetCaloUtils()->GetNumberOfSuperModulesUsed());
+  fhPtPhotonSM->SetYTitle("SuperModule ");
+  fhPtPhotonSM->SetXTitle("#it{E} (GeV)");
+  outputContainer->Add(fhPtPhotonSM) ;
+  
   fhNCellsE  = new TH2F ("hNCellsE","# of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nbins,nmin,nmax);
-  fhNCellsE->SetXTitle("E (GeV)");
+  fhNCellsE->SetXTitle("#it{E} (GeV)");
   fhNCellsE->SetYTitle("# of cells in cluster");
   outputContainer->Add(fhNCellsE);
   
   fhCellsE  = new TH2F ("hCellsE","energy of cells in cluster vs E of clusters", nptbins,ptmin,ptmax, nptbins*2,ptmin,ptmax);
-  fhCellsE->SetXTitle("E_{cluster} (GeV)");
-  fhCellsE->SetYTitle("E_{cell} (GeV)");
+  fhCellsE->SetXTitle("#it{E}_{cluster} (GeV)");
+  fhCellsE->SetYTitle("#it{E}_{cell} (GeV)");
   outputContainer->Add(fhCellsE);
   
   fhTimePt  = new TH2F ("hTimePt","time of cluster vs pT of clusters", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-  fhTimePt->SetXTitle("p_{T} (GeV/c)");
-  fhTimePt->SetYTitle("time (ns)");
+  fhTimePt->SetXTitle("#it{p}_{T} (GeV/#it{c})");
+  fhTimePt->SetYTitle("#it{time} (ns)");
   outputContainer->Add(fhTimePt);
   
   fhMaxCellDiffClusterE  = new TH2F ("hMaxCellDiffClusterE","energy vs difference of cluster energy - max cell energy / cluster energy, good clusters",
                                      nptbins,ptmin,ptmax, 500,0,1.);
-  fhMaxCellDiffClusterE->SetXTitle("E_{cluster} (GeV) ");
-  fhMaxCellDiffClusterE->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+  fhMaxCellDiffClusterE->SetXTitle("#it{E}_{cluster} (GeV) ");
+  fhMaxCellDiffClusterE->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
   outputContainer->Add(fhMaxCellDiffClusterE);
   
   fhEPhoton  = new TH1F("hEPhoton","Number of #gamma over calorimeter vs energy",nptbins,ptmin,ptmax);
-  fhEPhoton->SetYTitle("N");
-  fhEPhoton->SetXTitle("E_{#gamma}(GeV)");
+  fhEPhoton->SetYTitle("#it{counts}");
+  fhEPhoton->SetXTitle("#it{E}_{#gamma}(GeV)");
   outputContainer->Add(fhEPhoton) ;
   
-  fhPtPhoton  = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs p_{T}",nptbins,ptmin,ptmax);
-  fhPtPhoton->SetYTitle("N");
-  fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/c)");
+  fhPtPhoton  = new TH1F("hPtPhoton","Number of #gamma over calorimeter vs #it{p}_{T}",nptbins,ptmin,ptmax);
+  fhPtPhoton->SetYTitle("#it{counts}");
+  fhPtPhoton->SetXTitle("p_{T #gamma}(GeV/#it{c})");
   outputContainer->Add(fhPtPhoton) ;
   
-  fhPtCentralityPhoton  = new TH2F("hPtCentralityPhoton","centrality vs p_{T}",nptbins,ptmin,ptmax, 100,0,100);
-  fhPtCentralityPhoton->SetYTitle("Centrality");
-  fhPtCentralityPhoton->SetXTitle("p_{T}(GeV/c)");
-  outputContainer->Add(fhPtCentralityPhoton) ;
-  
-  fhPtEventPlanePhoton  = new TH2F("hPtEventPlanePhoton","centrality vs p_{T}",nptbins,ptmin,ptmax, 100,0,TMath::Pi());
-  fhPtEventPlanePhoton->SetYTitle("Event plane angle (rad)");
-  fhPtEventPlanePhoton->SetXTitle("p_{T} (GeV/c)");
-  outputContainer->Add(fhPtEventPlanePhoton) ;
+  if(IsHighMultiplicityAnalysisOn())
+  {
+    fhPtCentralityPhoton  = new TH2F("hPtCentralityPhoton","centrality vs #it{p}_{T}",nptbins,ptmin,ptmax, 100,0,100);
+    fhPtCentralityPhoton->SetYTitle("Centrality");
+    fhPtCentralityPhoton->SetXTitle("#it{p}_{T}(GeV/#it{c})");
+    outputContainer->Add(fhPtCentralityPhoton) ;
+    
+    fhPtEventPlanePhoton  = new TH2F("hPtEventPlanePhoton","centrality vs #it{p}_{T}",nptbins,ptmin,ptmax, 100,0,TMath::Pi());
+    fhPtEventPlanePhoton->SetYTitle("Event plane angle (rad)");
+    fhPtEventPlanePhoton->SetXTitle("#it{p}_{T} (GeV/#it{c})");
+    outputContainer->Add(fhPtEventPlanePhoton) ;
+  }
   
   fhEtaPhi  = new TH2F
-  ("hEtaPhi","cluster,E > 0.5 GeV, #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
+  ("hEtaPhi","cluster,#it{E} > 0.5 GeV, #eta vs #phi",netabins,etamin,etamax,nphibins,phimin,phimax);
   fhEtaPhi->SetYTitle("#phi (rad)");
   fhEtaPhi->SetXTitle("#eta");
   outputContainer->Add(fhEtaPhi) ;
   
-  if(fCalorimeter=="EMCAL" && fFillEMCALBCHistograms)
-  {
-    fhEtaPhiEMCALBC0  = new TH2F
-    ("hEtaPhiEMCALBC0","cluster,E > 2 GeV, #eta vs #phi, for clusters with |time| < 25 ns, EMCAL-BC=0",netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiEMCALBC0->SetYTitle("#phi (rad)");
-    fhEtaPhiEMCALBC0->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiEMCALBC0) ;
-    
-    fhEtaPhiEMCALBC1  = new TH2F
-    ("hEtaPhiEMCALBC1","cluster,E > 2 GeV, #eta vs #phi, for clusters with 25 < |time| < 75 ns, EMCAL-BC=1",netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiEMCALBC1->SetYTitle("#phi (rad)");
-    fhEtaPhiEMCALBC1->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiEMCALBC1) ;
-    
-    fhEtaPhiEMCALBCN  = new TH2F
-    ("hEtaPhiEMCALBCN","cluster,E > 2 GeV, #eta vs #phi, for clusters with |time| > 75 ns, EMCAL-BC>1",netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiEMCALBCN->SetYTitle("#phi (rad)");
-    fhEtaPhiEMCALBCN->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiEMCALBCN) ;
-    
-    for(Int_t i = 0; i < nTrigBC; i++)
-    {
-      fhEtaPhiTriggerEMCALBC[i] = new TH2F
-      (Form("hEtaPhiTriggerEMCALBC%d",i-iBCShift),
-       Form("cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=%d",i-iBCShift),
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBC[i]->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBC[i]->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBC[i]) ;
-      
-      fhTimeTriggerEMCALBC[i] = new TH2F
-      (Form("hTimeTriggerEMCALBC%d",i-iBCShift),
-       Form("cluster time vs E of clusters, Trigger EMCAL-BC=%d",i-iBCShift),
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBC[i]->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBC[i]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBC[i]);
-      
-      fhTimeTriggerEMCALBCPileUpSPD[i] = new TH2F
-      (Form("hTimeTriggerEMCALBC%dPileUpSPD",i-iBCShift),
-       Form("cluster time vs E of clusters, Trigger EMCAL-BC=%d",i-iBCShift),
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCPileUpSPD[i]->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCPileUpSPD[i]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCPileUpSPD[i]);
-      
-      fhEtaPhiTriggerEMCALBCUM[i] = new TH2F
-      (Form("hEtaPhiTriggerEMCALBC%d_UnMatch",i-iBCShift),
-       Form("cluster E > 2 GeV, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-iBCShift),
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUM[i]->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUM[i]->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUM[i]) ;
-      
-      fhTimeTriggerEMCALBCUM[i] = new TH2F
-      (Form("hTimeTriggerEMCALBC%d_UnMatch",i-iBCShift),
-       Form("cluster time vs E of clusters, unmatched trigger EMCAL-BC=%d",i-iBCShift),
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUM[i]->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUM[i]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUM[i]);
-      
-      fhEtaPhiTriggerEMCALBCCluster[i] = new TH2F
-      (Form("hEtaPhiTriggerEMCALBC%d_OnlyTrigger",i-iBCShift),
-       Form("trigger cluster, #eta vs #phi, Trigger EMCAL-BC=%d",i-iBCShift),
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCCluster[i]->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCCluster[i]->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCCluster[i]) ;
-            
-      fhEtaPhiTriggerEMCALBCUMCluster[i] = new TH2F
-      (Form("hEtaPhiTriggerEMCALBC%d_OnlyTrigger_UnMatch",i-iBCShift),
-       Form("trigger cluster, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-iBCShift),
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMCluster[i]->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMCluster[i]->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMCluster[i]) ;
-    }
-    
-    fhTimeTriggerEMCALBCCluster = new TH2F("hTimeTriggerEMCALBC_OnlyTrigger",
-                                           "trigger cluster time vs E of clusters",
-                                           nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBCCluster->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBCCluster->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBCCluster);
-    
-    fhTimeTriggerEMCALBCUMCluster = new TH2F("hTimeTriggerEMCALBC_OnlyTrigger_UnMatch",
-                                             "trigger cluster time vs E of clusters, unmatched trigger",
-                                             nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBCUMCluster->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBCUMCluster->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBCUMCluster);
-    
-    fhEtaPhiTriggerEMCALBCClusterOverTh = new TH2F
-    ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_OverThreshold",
-     "trigger cluster E > trigger threshold, #eta vs #phi, Trigger EMCAL-BC=0",
-     netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCClusterOverTh->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCClusterOverTh->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCClusterOverTh) ;
-    
-    fhEtaPhiTriggerEMCALBCUMClusterOverTh = new TH2F
-    ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_OverThreshold_UnMatch",
-     "trigger cluster E > trigger threshold, #eta vs #phi, unmatched trigger EMCAL-BC=0",
-     netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCUMClusterOverTh->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCUMClusterOverTh->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCUMClusterOverTh) ;
-    
-    fhEtaPhiTriggerEMCALBCClusterBelowTh1 = new TH2F
-    ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold1",
-     "trigger cluster thresh-1 < E < thres, #eta vs #phi, Trigger EMCAL-BC=0",
-     netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCClusterBelowTh1->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCClusterBelowTh1->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCClusterBelowTh1) ;
-    
-    fhEtaPhiTriggerEMCALBCUMClusterBelowTh1 = new TH2F
-    ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold1_UnMatch",
-     "trigger cluster thresh-1 < E < thres, #eta vs #phi, unmatched trigger EMCAL-BC=0",
-     netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCUMClusterBelowTh1->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCUMClusterBelowTh1->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCUMClusterBelowTh1) ;
-    
-    fhEtaPhiTriggerEMCALBCClusterBelowTh2 = new TH2F
-    ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold2",
-     "trigger cluster thresh-2 < E < thres, #eta vs #phi, Trigger EMCAL-BC=0",
-     netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCClusterBelowTh2->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCClusterBelowTh2->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCClusterBelowTh2) ;
-    
-    fhEtaPhiTriggerEMCALBCUMClusterBelowTh2 = new TH2F
-    ("hEtaPhiTriggerEMCALBC0_OnlyTrigger_BelowThreshold2_UnMatch",
-     "trigger cluster thresh-2 < E < thres, #eta vs #phi, unmatched trigger EMCAL-BC=0",
-     netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCUMClusterBelowTh2->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCUMClusterBelowTh2->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCUMClusterBelowTh2) ;
-    
-    if(!GetReader()->AreBadTriggerEventsRemoved())
-    {
-      fhEtaPhiTriggerEMCALBCExotic = new TH2F
-      ("hEtaPhiTriggerExotic",
-       "cluster E > 2 GeV, #eta vs #phi, Trigger Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCExotic->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCExotic->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCExotic) ;
-      
-      fhTimeTriggerEMCALBCExotic = new TH2F
-      ("hTimeTriggerExotic",
-       "cluster time vs E of clusters, Trigger Exotic ",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCExotic->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCExotic->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCExotic);
-      
-      fhEtaPhiTriggerEMCALBCUMExotic = new TH2F
-      ("hEtaPhiTriggerExotic_UnMatch",
-       "cluster E > 2 GeV, #eta vs #phi, unmatched trigger Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMExotic->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMExotic->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMExotic) ;
-      
-      fhTimeTriggerEMCALBCUMExotic = new TH2F
-      ("hTimeTriggerExotic_UnMatch",
-       "cluster time vs E of clusters, unmatched trigger Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMExotic->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMExotic->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMExotic);
-      
-      fhEtaPhiTriggerEMCALBCExoticCluster = new TH2F
-      ("hEtaPhiTriggerExotic_OnlyTrigger",
-       "trigger cluster E > 2 GeV, #eta vs #phi, Trigger Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCExoticCluster->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCExoticCluster->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCExoticCluster) ;
-      
-      fhTimeTriggerEMCALBCExoticCluster = new TH2F
-      ("hTimeTriggerExotic_OnlyTrigger",
-       "trigger cluster time vs E of clusters, Trigger Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCExoticCluster->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCExoticCluster->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCExoticCluster);
-      
-      fhEtaPhiTriggerEMCALBCUMExoticCluster = new TH2F
-      ("hEtaPhiTriggerExotic_OnlyTrigger_UnMatch",
-       "trigger cluster E > 2 GeV, #eta vs #phi, unmatched trigger Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMExoticCluster->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMExoticCluster->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMExoticCluster) ;
-      
-      fhTimeTriggerEMCALBCUMExoticCluster = new TH2F
-      ("hTimeTriggerExotic_OnlyTrigger_UnMatch",
-       "trigger cluster time vs E of clusters, unmatched trigger Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMExoticCluster->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMExoticCluster->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMExoticCluster);
-      
-      fhEtaPhiTriggerEMCALBCBad = new TH2F
-      ("hEtaPhiTriggerBad",
-       "cluster E > 2 GeV, #eta vs #phi, Trigger Bad",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCBad->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCBad->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCBad) ;
-      
-      fhTimeTriggerEMCALBCBad = new TH2F
-      ("hTimeTriggerBad",
-       "cluster time vs E of clusters, Trigger Bad ",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCBad->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCBad->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCBad);
-      
-      fhEtaPhiTriggerEMCALBCUMBad = new TH2F
-      ("hEtaPhiTriggerBad_UnMatch",
-       "cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMBad->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMBad->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBad) ;
-      
-      fhTimeTriggerEMCALBCUMBad = new TH2F
-      ("hTimeTriggerBad_UnMatch",
-       "cluster time vs E of clusters, unmatched trigger Bad",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMBad->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMBad->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMBad);
-      
-      fhEtaPhiTriggerEMCALBCBadCluster = new TH2F
-      ("hEtaPhiTriggerBad_OnlyTrigger",
-       "trigger cluster E > 2 GeV, #eta vs #phi, Trigger Bad",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCBadCluster->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCBadCluster->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCBadCluster) ;
-      
-      fhTimeTriggerEMCALBCBadCluster = new TH2F
-      ("hTimeTriggerBad_OnlyTrigger",
-       "trigger cluster time vs E of clusters, Trigger Bad",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCBadCluster->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCBadCluster->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCBadCluster);
-      
-      fhEtaPhiTriggerEMCALBCUMBadCluster = new TH2F
-      ("hEtaPhiTriggerBad_OnlyTrigger_UnMatch",
-       "trigger cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMBadCluster->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMBadCluster->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBadCluster) ;
-      
-      fhTimeTriggerEMCALBCUMBadCluster = new TH2F
-      ("hTimeTriggerBad_OnlyTrigger_UnMatch",
-       "trigger cluster time vs E of clusters, unmatched trigger Bad",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMBadCluster->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMBadCluster->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMBadCluster);
-      
-      fhEtaPhiTriggerEMCALBCBadExotic = new TH2F
-      ("hEtaPhiTriggerBadExotic",
-       "cluster E > 2 GeV, #eta vs #phi, Trigger Bad&Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCBadExotic->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCBadExotic->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCBadExotic) ;
-      
-      fhTimeTriggerEMCALBCBadExotic = new TH2F
-      ("hTimeTriggerBadExotic",
-       "cluster time vs E of clusters, Trigger Bad&Exotic ",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCBadExotic->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCBadExotic->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCBadExotic);
-      
-      fhEtaPhiTriggerEMCALBCUMBadExotic = new TH2F
-      ("hEtaPhiTriggerBadExotic_UnMatch",
-       "cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad&Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMBadExotic->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMBadExotic->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBadExotic) ;
-      
-      fhTimeTriggerEMCALBCUMBadExotic = new TH2F
-      ("hTimeTriggerBadExotic_UnMatch",
-       "cluster time vs E of clusters, unmatched trigger Bad&Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMBadExotic->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMBadExotic->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMBadExotic);
-      
-      fhEtaPhiTriggerEMCALBCBadExoticCluster = new TH2F
-      ("hEtaPhiTriggerBadExotic_OnlyTrigger",
-       "trigger cluster E > 2 GeV, #eta vs #phi, Trigger Bad&Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCBadExoticCluster->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCBadExoticCluster->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCBadExoticCluster) ;
-      
-      fhTimeTriggerEMCALBCBadExoticCluster = new TH2F
-      ("hTimeTriggerBadExotic_OnlyTrigger",
-       "trigger cluster time vs E of clusters, Trigger Bad&Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCBadExoticCluster->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCBadExoticCluster->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCBadExoticCluster);
-      
-      fhEtaPhiTriggerEMCALBCUMBadExoticCluster = new TH2F
-      ("hEtaPhiTriggerBadExotic_OnlyTrigger_UnMatch",
-       "trigger cluster E > 2 GeV, #eta vs #phi, unmatched trigger Bad&Exotic",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiTriggerEMCALBCUMBadExoticCluster->SetYTitle("#phi (rad)");
-      fhEtaPhiTriggerEMCALBCUMBadExoticCluster->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiTriggerEMCALBCUMBadExoticCluster) ;
-      
-      fhTimeTriggerEMCALBCUMBadExoticCluster = new TH2F
-      ("hTimeTriggerBadExotic_OnlyTrigger_UnMatch",
-       "trigger cluster time vs E of clusters, unmatched trigger Bad&Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMBadExoticCluster->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMBadExoticCluster->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMBadExoticCluster);
-      
-      fhTimeTriggerEMCALBCBadMaxCell = new TH2F
-      ("hTimeTriggerBadMaxCell",
-       "cluster time vs E of clusters, Trigger BadMaxCell",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCBadMaxCell->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCBadMaxCell->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCBadMaxCell);
-      
-      fhTimeTriggerEMCALBCUMBadMaxCell = new TH2F
-      ("hTimeTriggerBadMaxCell_UnMatch",
-       "cluster time vs E of clusters, unmatched trigger BadMaxCell",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMBadMaxCell->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMBadMaxCell->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMBadMaxCell);
-      
-      
-      fhTimeTriggerEMCALBCBadMaxCellExotic = new TH2F
-      ("hTimeTriggerBadMaxCellExotic",
-       "cluster time vs E of clusters, Trigger BadMaxCell&Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCBadMaxCellExotic->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCBadMaxCellExotic->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCBadMaxCellExotic);
-      
-      fhTimeTriggerEMCALBCUMBadMaxCellExotic = new TH2F
-      ("hTimeTriggerBadMaxCellExotic_UnMatch",
-       "cluster time vs E of clusters, unmatched trigger BadMaxCell&Exotic",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeTriggerEMCALBCUMBadMaxCellExotic->SetXTitle("E (GeV)");
-      fhTimeTriggerEMCALBCUMBadMaxCellExotic->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeTriggerEMCALBCUMBadMaxCellExotic);
-      
-      fhTimeNoTrigger = new TH2F
-      ("hTimeNoTrigger",
-       "events with no foundable trigger, time vs e of clusters",
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimeNoTrigger->SetXTitle("E (GeV)");
-      fhTimeNoTrigger->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeNoTrigger);
-      
-      fhEtaPhiNoTrigger = new TH2F
-      ("hEtaPhiNoTrigger",
-       "events with no foundable trigger, eta vs phi of clusters",
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiNoTrigger->SetYTitle("#phi (rad)");
-      fhEtaPhiNoTrigger->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiNoTrigger) ;
-    }
-    
-    fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster = new TH2F("hEtaPhiTriggerEMCALBC0_OnlyTrigger_UnMatch_ReMatch_OpenTime",
-                                                       "cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=0, un match, rematch open time",
-                                                       netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCUMReMatchOpenTimeCluster) ;
-    
-    fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster = new TH2F("hTimeTrigger_OnlyTrigger_UnMatch_ReMatch_OpenTime",
-                                                     "cluster time vs E of clusters, no match, rematch open time",
-                                                     nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBCUMReMatchOpenTimeCluster);
-    
-    
-    fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster = new TH2F("hEtaPhiTriggerEMCALBC0_OnlyTrigger_UnMatch_ReMatch_CheckNeighbours",
-                                                         "cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=0, un match, rematch with neighbour patches",
-                                                         netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCUMReMatchCheckNeighCluster) ;
-    
-    fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster = new TH2F("hTimeTrigger_OnlyTrigger_UnMatch_ReMatch_CheckNeighbours",
-                                                       "cluster time vs E of clusters, no match, rematch with neigbour parches",
-                                                       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBCUMReMatchCheckNeighCluster);
-    
-    fhEtaPhiTriggerEMCALBCUMReMatchBothCluster = new TH2F("hEtaPhiTriggerEMCALBC0_OnlyTrigger_UnMatch_ReMatch_Both",
-                                                   "cluster E > 2 GeV, #eta vs #phi, Trigger EMCAL-BC=0, un match, rematch open time and neighbour",
-                                                   netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiTriggerEMCALBCUMReMatchBothCluster->SetYTitle("#phi (rad)");
-    fhEtaPhiTriggerEMCALBCUMReMatchBothCluster->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiTriggerEMCALBCUMReMatchBothCluster) ;
-    
-    fhTimeTriggerEMCALBCUMReMatchBothCluster = new TH2F("hTimeTrigger_OnlyTrigger_UnMatch_ReMatch_Both",
-                                                 "cluster time vs E of clusters, no match, rematch open time and neigbour",
-                                                 nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBCUMReMatchBothCluster->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBCUMReMatchBothCluster->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBCUMReMatchBothCluster);
-    
-    fhTimeTriggerEMCALBC0UMReMatchOpenTime = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_OpenTime",
-                                                            "cluster time vs E of clusters, no match, rematch open time",
-                                                            nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBC0UMReMatchOpenTime->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBC0UMReMatchOpenTime->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchOpenTime);
-    
-        
-    fhTimeTriggerEMCALBC0UMReMatchCheckNeigh = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_CheckNeighbours",
-                                                              "cluster time vs E of clusters, no match, rematch with neigbour parches",
-                                                              nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBC0UMReMatchCheckNeigh->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBC0UMReMatchCheckNeigh->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchCheckNeigh);
-        
-    fhTimeTriggerEMCALBC0UMReMatchBoth = new TH2F("hTimeTriggerBC0_UnMatch_ReMatch_Both",
-                                                        "cluster time vs E of clusters, no match, rematch open time and neigbour",
-                                                        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimeTriggerEMCALBC0UMReMatchBoth->SetXTitle("E (GeV)");
-    fhTimeTriggerEMCALBC0UMReMatchBoth->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimeTriggerEMCALBC0UMReMatchBoth);
-    
-  }
-  
   fhPhiPhoton  = new TH2F
-  ("hPhiPhoton","#phi_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
+  ("hPhiPhoton","#phi_{#gamma} vs #it{p}_{T}",nptbins,ptmin,ptmax,nphibins,phimin,phimax);
   fhPhiPhoton->SetYTitle("#phi (rad)");
-  fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
+  fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/#it{c})");
   outputContainer->Add(fhPhiPhoton) ;
   
   fhEtaPhoton  = new TH2F
-  ("hEtaPhoton","#eta_{#gamma} vs p_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
+  ("hEtaPhoton","#eta_{#gamma} vs #it{p}_{T}",nptbins,ptmin,ptmax,netabins,etamin,etamax);
   fhEtaPhoton->SetYTitle("#eta");
-  fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
+  fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/#it{c})");
   outputContainer->Add(fhEtaPhoton) ;
   
   fhEtaPhiPhoton  = new TH2F
@@ -2322,98 +1263,10 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
     outputContainer->Add(fhEtaPhi05Photon) ;
   }
   
-  if(fCalorimeter=="EMCAL" && fFillEMCALBCHistograms)
-  {
-    fhEtaPhiPhotonEMCALBC0  = new TH2F
-    ("hEtaPhiPhotonEMCALBC0","identified photon, E > 2 GeV, #eta vs #phi, for clusters with |time| < 25 ns, EMCAL-BC=0",netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiPhotonEMCALBC0->SetYTitle("#phi (rad)");
-    fhEtaPhiPhotonEMCALBC0->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiPhotonEMCALBC0) ;
-    
-    fhEtaPhiPhotonEMCALBC1  = new TH2F
-    ("hEtaPhiPhotonEMCALBC1","identified photon, E > 2 GeV, #eta vs #phi, for clusters with 25 < |time| < 75 ns, EMCAL-BC=1",netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiPhotonEMCALBC1->SetYTitle("#phi (rad)");
-    fhEtaPhiPhotonEMCALBC1->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiPhotonEMCALBC1) ;
-    
-    fhEtaPhiPhotonEMCALBCN  = new TH2F
-    ("hEtaPhiPhotonEMCALBCN","identified photon, E > 2 GeV, #eta vs #phi, for clusters with |time| > 75 ns, EMCAL-BC>1",netabins,etamin,etamax,nphibins,phimin,phimax);
-    fhEtaPhiPhotonEMCALBCN->SetYTitle("#phi (rad)");
-    fhEtaPhiPhotonEMCALBCN->SetXTitle("#eta");
-    outputContainer->Add(fhEtaPhiPhotonEMCALBCN) ;
-    
-    for(Int_t i = 0; i < nTrigBC; i++)
-    {
-      fhEtaPhiPhotonTriggerEMCALBC[i] = new TH2F
-      (Form("hEtaPhiPhotonTriggerEMCALBC%d",i-iBCShift),
-       Form("photon E > 2 GeV, #eta vs #phi, PhotonTrigger EMCAL-BC=%d",i-iBCShift),
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiPhotonTriggerEMCALBC[i]->SetYTitle("#phi (rad)");
-      fhEtaPhiPhotonTriggerEMCALBC[i]->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiPhotonTriggerEMCALBC[i]) ;
-      
-      fhTimePhotonTriggerEMCALBC[i] = new TH2F
-      (Form("hTimePhotonTriggerEMCALBC%d",i-iBCShift),
-       Form("photon time vs E of clusters, PhotonTrigger EMCAL-BC=%d",i-iBCShift),
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimePhotonTriggerEMCALBC[i]->SetXTitle("E (GeV)");
-      fhTimePhotonTriggerEMCALBC[i]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimePhotonTriggerEMCALBC[i]);
-      
-      fhTimePhotonTriggerEMCALBCPileUpSPD[i] = new TH2F
-      (Form("hTimePhotonTriggerEMCALBC%dPileUpSPD",i-iBCShift),
-       Form("photon time vs E, PhotonTrigger EMCAL-BC=%d",i-iBCShift),
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimePhotonTriggerEMCALBCPileUpSPD[i]->SetXTitle("E (GeV)");
-      fhTimePhotonTriggerEMCALBCPileUpSPD[i]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimePhotonTriggerEMCALBCPileUpSPD[i]);
-      
-      fhEtaPhiPhotonTriggerEMCALBCUM[i] = new TH2F
-      (Form("hEtaPhiPhotonTriggerEMCALBC%d_UnMatch",i-iBCShift),
-       Form("photon E > 2 GeV, #eta vs #phi, unmatched trigger EMCAL-BC=%d",i-iBCShift),
-       netabins,etamin,etamax,nphibins,phimin,phimax);
-      fhEtaPhiPhotonTriggerEMCALBCUM[i]->SetYTitle("#phi (rad)");
-      fhEtaPhiPhotonTriggerEMCALBCUM[i]->SetXTitle("#eta");
-      outputContainer->Add(fhEtaPhiPhotonTriggerEMCALBCUM[i]) ;
-      
-      fhTimePhotonTriggerEMCALBCUM[i] = new TH2F
-      (Form("hTimePhotonTriggerEMCALBC%d_UnMatch",i-iBCShift),
-       Form("photon time vs E, unmatched trigger EMCAL-BC=%d",i-iBCShift),
-       nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-      fhTimePhotonTriggerEMCALBCUM[i]->SetXTitle("E (GeV)");
-      fhTimePhotonTriggerEMCALBCUM[i]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimePhotonTriggerEMCALBCUM[i]);
-      
-    }
-    
-    fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime = new TH2F("hTimePhotonTriggerBC0_UnMatch_ReMatch_OpenTime",
-                                                      "cluster time vs E of photons, no match, rematch open time",
-                                                      nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime->SetXTitle("E (GeV)");
-    fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime);
-    
-    
-    fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh = new TH2F("hTimePhotonTriggerBC0_UnMatch_ReMatch_CheckNeighbours",
-                                                        "cluster time vs E of photons, no match, rematch with neigbour parches",
-                                                        nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh->SetXTitle("E (GeV)");
-    fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh);
-    
-    fhTimePhotonTriggerEMCALBC0UMReMatchBoth = new TH2F("hTimePhotonTriggerBC0_UnMatch_ReMatch_Both",
-                                                  "cluster time vs E of photons, no match, rematch open time and neigbour",
-                                                  nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePhotonTriggerEMCALBC0UMReMatchBoth->SetXTitle("E (GeV)");
-    fhTimePhotonTriggerEMCALBC0UMReMatchBoth->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimePhotonTriggerEMCALBC0UMReMatchBoth);
-
-  }
-  
   fhNLocMax = new TH2F("hNLocMax","Number of local maxima in cluster",
                        nptbins,ptmin,ptmax,10,0,10);
   fhNLocMax ->SetYTitle("N maxima");
-  fhNLocMax ->SetXTitle("E (GeV)");
+  fhNLocMax ->SetXTitle("#it{E} (GeV)");
   outputContainer->Add(fhNLocMax) ;
   
   //Shower shape
@@ -2421,69 +1274,69 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
   {
     fhLam0E  = new TH2F ("hLam0E","#lambda_{0}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
     fhLam0E->SetYTitle("#lambda_{0}^{2}");
-    fhLam0E->SetXTitle("E (GeV)");
+    fhLam0E->SetXTitle("#it{E} (GeV)");
     outputContainer->Add(fhLam0E);
     
     fhLam1E  = new TH2F ("hLam1E","#lambda_{1}^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
     fhLam1E->SetYTitle("#lambda_{1}^{2}");
-    fhLam1E->SetXTitle("E (GeV)");
+    fhLam1E->SetXTitle("#it{E} (GeV)");
     outputContainer->Add(fhLam1E);
     
     fhDispE  = new TH2F ("hDispE"," dispersion^{2} vs E", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
     fhDispE->SetYTitle("D^{2}");
-    fhDispE->SetXTitle("E (GeV) ");
+    fhDispE->SetXTitle("#it{E} (GeV) ");
     outputContainer->Add(fhDispE);
     
     if(!fRejectTrackMatch)
     {
       fhLam0ETM  = new TH2F ("hLam0ETM","#lambda_{0}^{2} vs E, cut on track-matching residual |#Delta #eta| < 0.05,  |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
       fhLam0ETM->SetYTitle("#lambda_{0}^{2}");
-      fhLam0ETM->SetXTitle("E (GeV)");
+      fhLam0ETM->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhLam0ETM);
       
       fhLam1ETM  = new TH2F ("hLam1ETM","#lambda_{1}^{2} vs E, cut on track-matching residual |#Delta #eta| < 0.05,  |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
       fhLam1ETM->SetYTitle("#lambda_{1}^{2}");
-      fhLam1ETM->SetXTitle("E (GeV)");
+      fhLam1ETM->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhLam1ETM);
       
       fhDispETM  = new TH2F ("hDispETM"," dispersion^{2} vs E, cut on track-matching residual |#Delta #eta| < 0.05,  |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
       fhDispETM->SetYTitle("D^{2}");
-      fhDispETM->SetXTitle("E (GeV) ");
+      fhDispETM->SetXTitle("#it{E} (GeV) ");
       outputContainer->Add(fhDispETM);
     }
     
-    if(fCalorimeter == "EMCAL")
+    if(GetCalorimeter() == kEMCAL &&  GetFirstSMCoveredByTRD() >= 0)
     {
       fhLam0ETRD  = new TH2F ("hLam0ETRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
       fhLam0ETRD->SetYTitle("#lambda_{0}^{2}");
-      fhLam0ETRD->SetXTitle("E (GeV)");
+      fhLam0ETRD->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhLam0ETRD);
       
       fhLam1ETRD  = new TH2F ("hLam1ETRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
       fhLam1ETRD->SetYTitle("#lambda_{1}^{2}");
-      fhLam1ETRD->SetXTitle("E (GeV)");
+      fhLam1ETRD->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhLam1ETRD);
       
       fhDispETRD  = new TH2F ("hDispETRD"," dispersion^{2} vs E, EMCAL SM covered by TRD", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
       fhDispETRD->SetYTitle("Dispersion^{2}");
-      fhDispETRD->SetXTitle("E (GeV) ");
+      fhDispETRD->SetXTitle("#it{E} (GeV) ");
       outputContainer->Add(fhDispETRD);
       
-      if(!fRejectTrackMatch)
+      if(!fRejectTrackMatch &&  GetFirstSMCoveredByTRD() >=0 )
       {
         fhLam0ETMTRD  = new TH2F ("hLam0ETMTRD","#lambda_{0}^{2} vs E, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05,  |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhLam0ETMTRD->SetYTitle("#lambda_{0}^{2}");
-        fhLam0ETMTRD->SetXTitle("E (GeV)");
+        fhLam0ETMTRD->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhLam0ETMTRD);
         
         fhLam1ETMTRD  = new TH2F ("hLam1ETMTRD","#lambda_{1}^{2} vs E, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05,  |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhLam1ETMTRD->SetYTitle("#lambda_{1}^{2}");
-        fhLam1ETMTRD->SetXTitle("E (GeV)");
+        fhLam1ETMTRD->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhLam1ETMTRD);
         
         fhDispETMTRD  = new TH2F ("hDispETMTRD"," dispersion^{2} vs E, EMCAL SM covered by TRD, cut on track-matching residual |#Delta #eta| < 0.05,  |#Delta #phi| < 0.05", nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhDispETMTRD->SetYTitle("Dispersion^{2}");
-        fhDispETMTRD->SetXTitle("E (GeV) ");
+        fhDispETMTRD->SetXTitle("#it{E} (GeV) ");
         outputContainer->Add(fhDispETMTRD);
       }
     }
@@ -2495,7 +1348,7 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhNCellsLam0LowE->SetYTitle("#lambda_{0}^{2}");
       outputContainer->Add(fhNCellsLam0LowE);
       
-      fhNCellsLam0HighE  = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, E > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
+      fhNCellsLam0HighE  = new TH2F ("hNCellsLam0HighE","N_{cells} in cluster vs #lambda_{0}^{2}, #it{E} > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
       fhNCellsLam0HighE->SetXTitle("N_{cells}");
       fhNCellsLam0HighE->SetYTitle("#lambda_{0}^{2}");
       outputContainer->Add(fhNCellsLam0HighE);
@@ -2505,7 +1358,7 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhNCellsLam1LowE->SetYTitle("#lambda_{0}^{2}");
       outputContainer->Add(fhNCellsLam1LowE);
       
-      fhNCellsLam1HighE  = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, E > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
+      fhNCellsLam1HighE  = new TH2F ("hNCellsLam1HighE","N_{cells} in cluster vs #lambda_{1}^{2}, #it{E} > 2 GeV", nbins,nmin, nmax, ssbins,ssmin,ssmax);
       fhNCellsLam1HighE->SetXTitle("N_{cells}");
       fhNCellsLam1HighE->SetYTitle("#lambda_{0}^{2}");
       outputContainer->Add(fhNCellsLam1HighE);
@@ -2530,12 +1383,12 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhPhiLam0LowE->SetXTitle("#phi");
       outputContainer->Add(fhPhiLam0LowE);
       
-      fhEtaLam0HighE  = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, E > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
+      fhEtaLam0HighE  = new TH2F ("hEtaLam0HighE","#eta vs #lambda_{0}^{2}, #it{E} > 2 GeV", netabins,etamin,etamax, ssbins,ssmin,ssmax);
       fhEtaLam0HighE->SetYTitle("#lambda_{0}^{2}");
       fhEtaLam0HighE->SetXTitle("#eta");
       outputContainer->Add(fhEtaLam0HighE);
       
-      fhPhiLam0HighE  = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, E > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
+      fhPhiLam0HighE  = new TH2F ("hPhiLam0HighE","#phi vs #lambda_{0}^{2}, #it{E} > 2 GeV", nphibins,phimin,phimax, ssbins,ssmin,ssmax);
       fhPhiLam0HighE->SetYTitle("#lambda_{0}^{2}");
       fhPhiLam0HighE->SetXTitle("#phi");
       outputContainer->Add(fhPhiLam0HighE);
@@ -2545,7 +1398,7 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhLam1Lam0LowE->SetXTitle("#lambda_{1}^{2}");
       outputContainer->Add(fhLam1Lam0LowE);
       
-      fhLam1Lam0HighE  = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of E > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
+      fhLam1Lam0HighE  = new TH2F ("hLam1Lam0HighE","#lambda_{0}^{2} vs #lambda_{1}^{2} in cluster of #it{E} > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
       fhLam1Lam0HighE->SetYTitle("#lambda_{0}^{2}");
       fhLam1Lam0HighE->SetXTitle("#lambda_{1}^{2}");
       outputContainer->Add(fhLam1Lam0HighE);
@@ -2555,7 +1408,7 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhLam0DispLowE->SetYTitle("D^{2}");
       outputContainer->Add(fhLam0DispLowE);
       
-      fhLam0DispHighE  = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of E > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
+      fhLam0DispHighE  = new TH2F ("hLam0DispHighE","#lambda_{0}^{2} vs dispersion^{2} in cluster of #it{E} > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
       fhLam0DispHighE->SetXTitle("#lambda_{0}^{2}");
       fhLam0DispHighE->SetYTitle("D^{2}");
       outputContainer->Add(fhLam0DispHighE);
@@ -2565,59 +1418,59 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhDispLam1LowE->SetYTitle("#lambda_{1}^{2}");
       outputContainer->Add(fhDispLam1LowE);
       
-      fhDispLam1HighE  = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of E > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
+      fhDispLam1HighE  = new TH2F ("hDispLam1HighE","Dispersion^{2} vs #lambda_{1^{2}} in cluster of #it{E} > 2 GeV",  ssbins,ssmin,ssmax, ssbins,ssmin,ssmax);
       fhDispLam1HighE->SetXTitle("D^{2}");
       fhDispLam1HighE->SetYTitle("#lambda_{1}^{2}");
       outputContainer->Add(fhDispLam1HighE);
       
-      if(fCalorimeter == "EMCAL")
+      if(GetCalorimeter() == kEMCAL)
       {
         fhDispEtaE  = new TH2F ("hDispEtaE","#sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",  nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
-        fhDispEtaE->SetXTitle("E (GeV)");
+        fhDispEtaE->SetXTitle("#it{E} (GeV)");
         fhDispEtaE->SetYTitle("#sigma^{2}_{#eta #eta}");
         outputContainer->Add(fhDispEtaE);
         
         fhDispPhiE  = new TH2F ("hDispPhiE","#sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",  nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
-        fhDispPhiE->SetXTitle("E (GeV)");
+        fhDispPhiE->SetXTitle("#it{E} (GeV)");
         fhDispPhiE->SetYTitle("#sigma^{2}_{#phi #phi}");
         outputContainer->Add(fhDispPhiE);
         
         fhSumEtaE  = new TH2F ("hSumEtaE","#delta^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i})^{2} / #Sigma w_{i} - <#eta>^{2} vs E",  nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
-        fhSumEtaE->SetXTitle("E (GeV)");
+        fhSumEtaE->SetXTitle("#it{E} (GeV)");
         fhSumEtaE->SetYTitle("#delta^{2}_{#eta #eta}");
         outputContainer->Add(fhSumEtaE);
         
         fhSumPhiE  = new TH2F ("hSumPhiE","#delta^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i})^{2}/ #Sigma w_{i} - <#phi>^{2} vs E",
                                nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
-        fhSumPhiE->SetXTitle("E (GeV)");
+        fhSumPhiE->SetXTitle("#it{E} (GeV)");
         fhSumPhiE->SetYTitle("#delta^{2}_{#phi #phi}");
         outputContainer->Add(fhSumPhiE);
         
         fhSumEtaPhiE  = new TH2F ("hSumEtaPhiE","#delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",
                                   nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
-        fhSumEtaPhiE->SetXTitle("E (GeV)");
+        fhSumEtaPhiE->SetXTitle("#it{E} (GeV)");
         fhSumEtaPhiE->SetYTitle("#delta^{2}_{#eta #phi}");
         outputContainer->Add(fhSumEtaPhiE);
         
         fhDispEtaPhiDiffE  = new TH2F ("hDispEtaPhiDiffE","#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",
                                        nptbins,ptmin,ptmax,200, -10,10);
-        fhDispEtaPhiDiffE->SetXTitle("E (GeV)");
+        fhDispEtaPhiDiffE->SetXTitle("#it{E} (GeV)");
         fhDispEtaPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
         outputContainer->Add(fhDispEtaPhiDiffE);
         
         fhSphericityE  = new TH2F ("hSphericityE","(#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",
                                    nptbins,ptmin,ptmax, 200, -1,1);
-        fhSphericityE->SetXTitle("E (GeV)");
+        fhSphericityE->SetXTitle("#it{E} (GeV)");
         fhSphericityE->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
         outputContainer->Add(fhSphericityE);
         
         fhDispSumEtaDiffE  = new TH2F ("hDispSumEtaDiffE","#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average vs E",  nptbins,ptmin,ptmax, 200,-0.01,0.01);
-        fhDispSumEtaDiffE->SetXTitle("E (GeV)");
+        fhDispSumEtaDiffE->SetXTitle("#it{E} (GeV)");
         fhDispSumEtaDiffE->SetYTitle("#sigma^{2}_{#eta #eta} - #delta^{2}_{#eta #eta} / average");
         outputContainer->Add(fhDispSumEtaDiffE);
         
         fhDispSumPhiDiffE  = new TH2F ("hDispSumPhiDiffE","#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average vs E",  nptbins,ptmin,ptmax, 200,-0.01,0.01);
-        fhDispSumPhiDiffE->SetXTitle("E (GeV)");
+        fhDispSumPhiDiffE->SetXTitle("#it{E} (GeV)");
         fhDispSumPhiDiffE->SetYTitle("#sigma^{2}_{#phi #phi} - #delta^{2}_{#phi #phi} / average");
         outputContainer->Add(fhDispSumPhiDiffE);
         
@@ -2658,14 +1511,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
        Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
        nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
       fhTrackMatchedDEta[i]->SetYTitle("d#eta");
-      fhTrackMatchedDEta[i]->SetXTitle("E_{cluster} (GeV)");
+      fhTrackMatchedDEta[i]->SetXTitle("#it{E}_{cluster} (GeV)");
       
       fhTrackMatchedDPhi[i]  = new TH2F
       (Form("hTrackMatchedDPhi%s",cutTM[i].Data()),
        Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
        nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
       fhTrackMatchedDPhi[i]->SetYTitle("d#phi (rad)");
-      fhTrackMatchedDPhi[i]->SetXTitle("E_{cluster} (GeV)");
+      fhTrackMatchedDPhi[i]->SetXTitle("#it{E}_{cluster} (GeV)");
       
       fhTrackMatchedDEtaDPhi[i]  = new TH2F
       (Form("hTrackMatchedDEtaDPhi%s",cutTM[i].Data()),
@@ -2679,14 +1532,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
        Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
        nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
       fhTrackMatchedDEtaPos[i]->SetYTitle("d#eta");
-      fhTrackMatchedDEtaPos[i]->SetXTitle("E_{cluster} (GeV)");
+      fhTrackMatchedDEtaPos[i]->SetXTitle("#it{E}_{cluster} (GeV)");
       
       fhTrackMatchedDPhiPos[i]  = new TH2F
       (Form("hTrackMatchedDPhiPos%s",cutTM[i].Data()),
        Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
        nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
       fhTrackMatchedDPhiPos[i]->SetYTitle("d#phi (rad)");
-      fhTrackMatchedDPhiPos[i]->SetXTitle("E_{cluster} (GeV)");
+      fhTrackMatchedDPhiPos[i]->SetXTitle("#it{E}_{cluster} (GeV)");
       
       fhTrackMatchedDEtaDPhiPos[i]  = new TH2F
       (Form("hTrackMatchedDEtaDPhiPos%s",cutTM[i].Data()),
@@ -2700,14 +1553,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
        Form("d#eta of cluster-track vs cluster energy, %s",cutTM[i].Data()),
        nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
       fhTrackMatchedDEtaNeg[i]->SetYTitle("d#eta");
-      fhTrackMatchedDEtaNeg[i]->SetXTitle("E_{cluster} (GeV)");
+      fhTrackMatchedDEtaNeg[i]->SetXTitle("#it{E}_{cluster} (GeV)");
       
       fhTrackMatchedDPhiNeg[i]  = new TH2F
       (Form("hTrackMatchedDPhiNeg%s",cutTM[i].Data()),
        Form("d#phi of cluster-track vs cluster energy, %s",cutTM[i].Data()),
        nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
       fhTrackMatchedDPhiNeg[i]->SetYTitle("d#phi (rad)");
-      fhTrackMatchedDPhiNeg[i]->SetXTitle("E_{cluster} (GeV)");
+      fhTrackMatchedDPhiNeg[i]->SetXTitle("#it{E}_{cluster} (GeV)");
       
       fhTrackMatchedDEtaDPhiNeg[i]  = new TH2F
       (Form("hTrackMatchedDEtaDPhiNeg%s",cutTM[i].Data()),
@@ -2718,12 +1571,12 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       
       fhdEdx[i]  = new TH2F (Form("hdEdx%s",cutTM[i].Data()),Form("matched track <dE/dx> vs cluster E, %s",cutTM[i].Data()),
                              nptbins,ptmin,ptmax,ndedxbins, dedxmin, dedxmax);
-      fhdEdx[i]->SetXTitle("E (GeV)");
+      fhdEdx[i]->SetXTitle("#it{E} (GeV)");
       fhdEdx[i]->SetYTitle("<dE/dx>");
       
       fhEOverP[i]  = new TH2F (Form("hEOverP%s",cutTM[i].Data()),Form("matched track E/p vs cluster E, %s",cutTM[i].Data()),
                                nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
-      fhEOverP[i]->SetXTitle("E (GeV)");
+      fhEOverP[i]->SetXTitle("#it{E} (GeV)");
       fhEOverP[i]->SetYTitle("E/p");
       
       outputContainer->Add(fhTrackMatchedDEta[i]) ;
@@ -2738,27 +1591,27 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       outputContainer->Add(fhdEdx[i]);
       outputContainer->Add(fhEOverP[i]);
       
-      if(fCalorimeter=="EMCAL")
+      if(GetCalorimeter()==kEMCAL &&  GetFirstSMCoveredByTRD() >=0 )
       {
         fhTrackMatchedDEtaTRD[i]  = new TH2F
         (Form("hTrackMatchedDEtaTRD%s",cutTM[i].Data()),
          Form("d#eta of cluster-track vs cluster energy, SM behind TRD, %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
         fhTrackMatchedDEtaTRD[i]->SetYTitle("d#eta");
-        fhTrackMatchedDEtaTRD[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDEtaTRD[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         fhTrackMatchedDPhiTRD[i]  = new TH2F
         (Form("hTrackMatchedDPhiTRD%s",cutTM[i].Data()),
          Form("d#phi of cluster-track vs cluster energy, SM behing TRD, %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
         fhTrackMatchedDPhiTRD[i]->SetYTitle("d#phi (rad)");
-        fhTrackMatchedDPhiTRD[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDPhiTRD[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         fhEOverPTRD[i]  = new TH2F
         (Form("hEOverPTRD%s",cutTM[i].Data()),
          Form("matched track E/p vs cluster E, behind TRD, %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nPoverEbins,pOverEmin,pOverEmax);
-        fhEOverPTRD[i]->SetXTitle("E (GeV)");
+        fhEOverPTRD[i]->SetXTitle("#it{E} (GeV)");
         fhEOverPTRD[i]->SetYTitle("E/p");
         
         outputContainer->Add(fhTrackMatchedDEtaTRD[i]) ;
@@ -2773,14 +1626,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
          Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
         fhTrackMatchedDEtaMCNoOverlap[i]->SetYTitle("d#eta");
-        fhTrackMatchedDEtaMCNoOverlap[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDEtaMCNoOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         fhTrackMatchedDPhiMCNoOverlap[i]  = new TH2F
         (Form("hTrackMatchedDPhiMCNoOverlap%s",cutTM[i].Data()),
          Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
         fhTrackMatchedDPhiMCNoOverlap[i]->SetYTitle("d#phi (rad)");
-        fhTrackMatchedDPhiMCNoOverlap[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDPhiMCNoOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         outputContainer->Add(fhTrackMatchedDEtaMCNoOverlap[i]) ;
         outputContainer->Add(fhTrackMatchedDPhiMCNoOverlap[i]) ;
@@ -2789,14 +1642,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
          Form("d#eta of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
         fhTrackMatchedDEtaMCOverlap[i]->SetYTitle("d#eta");
-        fhTrackMatchedDEtaMCOverlap[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDEtaMCOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         fhTrackMatchedDPhiMCOverlap[i]  = new TH2F
         (Form("hTrackMatchedDPhiMCOverlap%s",cutTM[i].Data()),
          Form("d#phi of cluster-track vs cluster energy, several MC particles overlap %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
         fhTrackMatchedDPhiMCOverlap[i]->SetYTitle("d#phi (rad)");
-        fhTrackMatchedDPhiMCOverlap[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDPhiMCOverlap[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         outputContainer->Add(fhTrackMatchedDEtaMCOverlap[i]) ;
         outputContainer->Add(fhTrackMatchedDPhiMCOverlap[i]) ;
@@ -2806,14 +1659,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
          Form("d#eta of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresetabins,resetamin,resetamax);
         fhTrackMatchedDEtaMCConversion[i]->SetYTitle("d#eta");
-        fhTrackMatchedDEtaMCConversion[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDEtaMCConversion[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         fhTrackMatchedDPhiMCConversion[i]  = new TH2F
         (Form("hTrackMatchedDPhiMCConversion%s",cutTM[i].Data()),
          Form("d#phi of cluster-track vs cluster energy, no other MC particles overlap appart from conversions %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,nresphibins,resphimin,resphimax);
         fhTrackMatchedDPhiMCConversion[i]->SetYTitle("d#phi (rad)");
-        fhTrackMatchedDPhiMCConversion[i]->SetXTitle("E_{cluster} (GeV)");
+        fhTrackMatchedDPhiMCConversion[i]->SetXTitle("#it{E}_{cluster} (GeV)");
         
         outputContainer->Add(fhTrackMatchedDEtaMCConversion[i]) ;
         outputContainer->Add(fhTrackMatchedDPhiMCConversion[i]) ;
@@ -2822,7 +1675,7 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
         (Form("hTrackMatchedMCParticle%s",cutTM[i].Data()),
          Form("Origin of particle vs energy %s",cutTM[i].Data()),
          nptbins,ptmin,ptmax,8,0,8);
-        fhTrackMatchedMCParticle[i]->SetXTitle("E (GeV)");
+        fhTrackMatchedMCParticle[i]->SetXTitle("#it{E} (GeV)");
         //fhTrackMatchedMCParticle[i]->SetYTitle("Particle type");
         
         fhTrackMatchedMCParticle[i]->GetYaxis()->SetBinLabel(1 ,"Photon");
@@ -2839,246 +1692,84 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
     }
   }
   
-  if(fFillPileUpHistograms)
+  if(IsPileUpAnalysisOn())
   {
-    
     TString pileUpName[] = {"SPD","EMCAL","SPDOrEMCAL","SPDAndEMCAL","SPDAndNotEMCAL","EMCALAndNotSPD","NotSPDAndNotEMCAL"} ;
     
     for(Int_t i = 0 ; i < 7 ; i++)
     {
-      fhPtPileUp[i]  = new TH1F(Form("hPtPileUp%s",pileUpName[i].Data()),
-                                Form("Cluster  p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
-      fhPtPileUp[i]->SetXTitle("p_{T} (GeV/c)");
-      outputContainer->Add(fhPtPileUp[i]);
-      
-      fhPtChargedPileUp[i]  = new TH1F(Form("hPtChargedPileUp%s",pileUpName[i].Data()),
-                                       Form("Charged clusters p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
-      fhPtChargedPileUp[i]->SetXTitle("p_{T} (GeV/c)");
-      outputContainer->Add(fhPtChargedPileUp[i]);
-      
       fhPtPhotonPileUp[i]  = new TH1F(Form("hPtPhotonPileUp%s",pileUpName[i].Data()),
-                                      Form("Selected photon p_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
-      fhPtPhotonPileUp[i]->SetXTitle("p_{T} (GeV/c)");
+                                      Form("Selected photon #it{p}_{T} distribution, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax);
+      fhPtPhotonPileUp[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
       outputContainer->Add(fhPtPhotonPileUp[i]);
       
-      
-      fhClusterEFracLongTimePileUp[i]  = new TH2F(Form("hClusterEFracLongTimePileUp%s",pileUpName[i].Data()),
-                                                  Form("Cluster E vs fraction of cluster energy from large T cells, %s Pile-Up event",pileUpName[i].Data()),
-                                                  nptbins,ptmin,ptmax,200,0,1);
-      fhClusterEFracLongTimePileUp[i]->SetXTitle("E (GeV)");
-      fhClusterEFracLongTimePileUp[i]->SetYTitle("E(large time) / E");
-      outputContainer->Add(fhClusterEFracLongTimePileUp[i]);
-      
-      fhClusterCellTimePileUp[i]  = new TH2F(Form("hClusterCellTimePileUp%s",pileUpName[i].Data()),
-                                             Form("Cluster E vs cell time in cluster, %s Pile-Up event",pileUpName[i].Data()),
-                                             nptbins,ptmin,ptmax,ntimebins,timemin,timemax);
-      fhClusterCellTimePileUp[i]->SetXTitle("E (GeV)");
-      fhClusterCellTimePileUp[i]->SetYTitle("t_{cell} (ns)");
-      outputContainer->Add(fhClusterCellTimePileUp[i]);
-      
-      fhClusterTimeDiffPileUp[i]  = new TH2F(Form("hClusterTimeDiffPileUp%s",pileUpName[i].Data()),
-                                             Form("Cluster E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
-                                             nptbins,ptmin,ptmax,400,-200,200);
-      fhClusterTimeDiffPileUp[i]->SetXTitle("E (GeV)");
-      fhClusterTimeDiffPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
-      outputContainer->Add(fhClusterTimeDiffPileUp[i]);
-      
-      fhClusterTimeDiffChargedPileUp[i]  = new TH2F(Form("hClusterTimeDiffChargedPileUp%s",pileUpName[i].Data()),
-                                                    Form("Charged clusters E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
-                                                    nptbins,ptmin,ptmax,400,-200,200);
-      fhClusterTimeDiffChargedPileUp[i]->SetXTitle("E (GeV)");
-      fhClusterTimeDiffChargedPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
-      outputContainer->Add(fhClusterTimeDiffChargedPileUp[i]);
-      
       fhClusterTimeDiffPhotonPileUp[i]  = new TH2F(Form("hClusterTimeDiffPhotonPileUp%s",pileUpName[i].Data()),
-                                                   Form("Selected photon E vs t_{max}-t_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
-                                                   nptbins,ptmin,ptmax,400,-200,200);
-      fhClusterTimeDiffPhotonPileUp[i]->SetXTitle("E (GeV)");
-      fhClusterTimeDiffPhotonPileUp[i]->SetYTitle("t_{max}-t_{cell} (ns)");
+                                             Form("Photon cluster E vs #it{t}_{max}-#it{t}_{cell} in cluster, %s Pile-Up event",pileUpName[i].Data()),
+                                             nptbins,ptmin,ptmax,400,-200,200);
+      fhClusterTimeDiffPhotonPileUp[i]->SetXTitle("#it{E} (GeV)");
+      fhClusterTimeDiffPhotonPileUp[i]->SetYTitle("#it{t}_{max}-#it{t}_{cell} (ns)");
       outputContainer->Add(fhClusterTimeDiffPhotonPileUp[i]);
-      
-      fhLambda0PileUp[i]  = new TH2F(Form("hLambda0PileUp%s",pileUpName[i].Data()),
-                                     Form("Cluster E vs #lambda^{2}_{0} in cluster, %s Pile-Up event",pileUpName[i].Data()),
-                                     nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
-      fhLambda0PileUp[i]->SetXTitle("E (GeV)");
-      fhLambda0PileUp[i]->SetYTitle("#lambda^{2}_{0}");
-      outputContainer->Add(fhLambda0PileUp[i]);
-      
-      fhLambda0ChargedPileUp[i]  = new TH2F(Form("hLambda0ChargedPileUp%s",pileUpName[i].Data()),
-                                            Form("Charged clusters E vs #lambda^{2}_{0}in cluster, %s Pile-Up event",pileUpName[i].Data()), nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
-      fhLambda0ChargedPileUp[i]->SetXTitle("E (GeV)");
-      fhLambda0ChargedPileUp[i]->SetYTitle("#lambda^{2}_{0}");
-      outputContainer->Add(fhLambda0ChargedPileUp[i]);
-      
     }
     
-    fhEtaPhiBC0  = new TH2F ("hEtaPhiBC0","eta-phi for clusters tof corresponding to BC=0",netabins,etamin,etamax, nphibins,phimin,phimax);
-    fhEtaPhiBC0->SetXTitle("#eta ");
-    fhEtaPhiBC0->SetYTitle("#phi (rad)");
-    outputContainer->Add(fhEtaPhiBC0);
-    
-    fhEtaPhiBCPlus  = new TH2F ("hEtaPhiBCPlus","eta-phi for clusters tof corresponding to BC>0",netabins,etamin,etamax, nphibins,phimin,phimax);
-    fhEtaPhiBCPlus->SetXTitle("#eta ");
-    fhEtaPhiBCPlus->SetYTitle("#phi (rad)");
-    outputContainer->Add(fhEtaPhiBCPlus);
-    
-    fhEtaPhiBCMinus  = new TH2F ("hEtaPhiBCMinus","eta-phi for clusters tof corresponding to BC<0",netabins,etamin,etamax, nphibins,phimin,phimax);
-    fhEtaPhiBCMinus->SetXTitle("#eta ");
-    fhEtaPhiBCMinus->SetYTitle("#phi (rad)");
-    outputContainer->Add(fhEtaPhiBCMinus);
-    
-    fhEtaPhiBC0PileUpSPD  = new TH2F ("hEtaPhiBC0PileUpSPD","eta-phi for clusters tof corresponding to BC=0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
-    fhEtaPhiBC0PileUpSPD->SetXTitle("#eta ");
-    fhEtaPhiBC0PileUpSPD->SetYTitle("#phi (rad)");
-    outputContainer->Add(fhEtaPhiBC0PileUpSPD);
-    
-    fhEtaPhiBCPlusPileUpSPD  = new TH2F ("hEtaPhiBCPlusPileUpSPD","eta-phi for clusters tof corresponding to BC>0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
-    fhEtaPhiBCPlusPileUpSPD->SetXTitle("#eta ");
-    fhEtaPhiBCPlusPileUpSPD->SetYTitle("#phi (rad)");
-    outputContainer->Add(fhEtaPhiBCPlusPileUpSPD);
-    
-    fhEtaPhiBCMinusPileUpSPD  = new TH2F ("hEtaPhiBCMinusPileUpSPD","eta-phi for clusters tof corresponding to BC<0, SPD pile-up",netabins,etamin,etamax, nphibins,phimin,phimax);
-    fhEtaPhiBCMinusPileUpSPD->SetXTitle("#eta ");
-    fhEtaPhiBCMinusPileUpSPD->SetYTitle("#phi (rad)");
-    outputContainer->Add(fhEtaPhiBCMinusPileUpSPD);
-    
-    fhTimePtNoCut  = new TH2F ("hTimePt_NoCut","time of cluster vs pT of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePtNoCut->SetXTitle("p_{T} (GeV/c)");
-    fhTimePtNoCut->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimePtNoCut);
-    
-    fhTimePtSPD  = new TH2F ("hTimePt_SPD","time of cluster vs pT of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePtSPD->SetXTitle("p_{T} (GeV/c)");
-    fhTimePtSPD->SetYTitle("time (ns)");
-    outputContainer->Add(fhTimePtSPD);
-    
     fhTimePtPhotonNoCut  = new TH2F ("hTimePtPhoton_NoCut","time of photon cluster vs pT of clusters, no cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePtPhotonNoCut->SetXTitle("p_{T} (GeV/c)");
-    fhTimePtPhotonNoCut->SetYTitle("time (ns)");
+    fhTimePtPhotonNoCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
+    fhTimePtPhotonNoCut->SetYTitle("#it{time} (ns)");
     outputContainer->Add(fhTimePtPhotonNoCut);
     
     fhTimePtPhotonSPD  = new TH2F ("hTimePtPhoton_SPD","time of  photon cluster vs pT of clusters, SPD cut", nptbins,ptmin,ptmax, ntimebins,timemin,timemax);
-    fhTimePtPhotonSPD->SetXTitle("p_{T} (GeV/c)");
-    fhTimePtPhotonSPD->SetYTitle("time (ns)");
+    fhTimePtPhotonSPD->SetXTitle("#it{p}_{T} (GeV/#it{c})");
+    fhTimePtPhotonSPD->SetYTitle("#it{time} (ns)");
     outputContainer->Add(fhTimePtPhotonSPD);
-        
+    
     fhTimeNPileUpVertSPD  = new TH2F ("hTime_NPileUpVertSPD","time of cluster vs N pile-up SPD vertex", ntimebins,timemin,timemax,20,0,20);
     fhTimeNPileUpVertSPD->SetYTitle("# vertex ");
-    fhTimeNPileUpVertSPD->SetXTitle("time (ns)");
+    fhTimeNPileUpVertSPD->SetXTitle("#it{time} (ns)");
     outputContainer->Add(fhTimeNPileUpVertSPD);
     
     fhTimeNPileUpVertTrack  = new TH2F ("hTime_NPileUpVertTracks","time of cluster vs N pile-up Tracks vertex", ntimebins,timemin,timemax, 20,0,20 );
     fhTimeNPileUpVertTrack->SetYTitle("# vertex ");
-    fhTimeNPileUpVertTrack->SetXTitle("time (ns)");
+    fhTimeNPileUpVertTrack->SetXTitle("#it{time} (ns)");
     outputContainer->Add(fhTimeNPileUpVertTrack);
     
-    fhTimeNPileUpVertContributors  = new TH2F ("hTime_NPileUpVertContributors","time of cluster vs N constributors to pile-up SPD vertex", ntimebins,timemin,timemax,50,0,50);
-    fhTimeNPileUpVertContributors->SetYTitle("# vertex ");
-    fhTimeNPileUpVertContributors->SetXTitle("time (ns)");
-    outputContainer->Add(fhTimeNPileUpVertContributors);
-    
-    fhTimePileUpMainVertexZDistance  = new TH2F ("hTime_PileUpMainVertexZDistance","time of cluster vs distance in Z pile-up SPD vertex - main SPD vertex",ntimebins,timemin,timemax,100,0,50);
-    fhTimePileUpMainVertexZDistance->SetYTitle("distance Z (cm) ");
-    fhTimePileUpMainVertexZDistance->SetXTitle("time (ns)");
-    outputContainer->Add(fhTimePileUpMainVertexZDistance);
-    
-    fhTimePileUpMainVertexZDiamond  = new TH2F ("hTime_PileUpMainVertexZDiamond","time of cluster vs distance in Z pile-up SPD vertex - z diamond",ntimebins,timemin,timemax,100,0,50);
-    fhTimePileUpMainVertexZDiamond->SetYTitle("diamond distance Z (cm) ");
-    fhTimePileUpMainVertexZDiamond->SetXTitle("time (ns)");
-    outputContainer->Add(fhTimePileUpMainVertexZDiamond);
-    
-    TString title[] = {"no |t diff| cut","|t diff|<20 ns","|t diff|>20 ns","|t diff|>40 ns"};
-    TString name [] = {"TDiffNoCut","TDiffSmaller20ns","TDiffLarger20ns","TDiffLarger40ns"};
-    for(Int_t i = 0; i < 4; i++)
-    {
-      fhClusterMultSPDPileUp[i] = new TH2F(Form("fhClusterMultSPDPileUp_%s", name[i].Data()),
-                                           Form("Number of clusters per pile up event with E > 0.5 and %s respect cluster max vs cluster max E ",title[i].Data()),
-                                           nptbins,ptmin,ptmax,100,0,100);
-      fhClusterMultSPDPileUp[i]->SetYTitle("n clusters ");
-      fhClusterMultSPDPileUp[i]->SetXTitle("E_{cluster max} (GeV)");
-      outputContainer->Add(fhClusterMultSPDPileUp[i]) ;
-      
-      fhClusterMultNoPileUp[i] = new TH2F(Form("fhClusterMultNoPileUp_%s", name[i].Data()),
-                                          Form("Number of clusters per non pile up event with E > 0.5 and %s respect cluster max vs cluster max E ",title[i].Data()),
-                                          nptbins,ptmin,ptmax,100,0,100);
-      fhClusterMultNoPileUp[i]->SetYTitle("n clusters ");
-      fhClusterMultNoPileUp[i]->SetXTitle("E_{cluster max} (GeV)");
-      outputContainer->Add(fhClusterMultNoPileUp[i]) ;
-    }
-    
-    fhPtNPileUpSPDVtx  = new TH2F ("hPt_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
-                                   nptbins,ptmin,ptmax,20,0,20);
-    fhPtNPileUpSPDVtx->SetYTitle("# vertex ");
-    fhPtNPileUpSPDVtx->SetXTitle("p_{T} (GeV/c)");
-    outputContainer->Add(fhPtNPileUpSPDVtx);
-         
-    fhPtNPileUpTrkVtx  = new TH2F ("hPt_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
-                                   nptbins,ptmin,ptmax, 20,0,20 );
-    fhPtNPileUpTrkVtx->SetYTitle("# vertex ");
-    fhPtNPileUpTrkVtx->SetXTitle("p_{T} (GeV/c)");
-    outputContainer->Add(fhPtNPileUpTrkVtx);
-    
-    fhPtNPileUpSPDVtxTimeCut  = new TH2F ("hPt_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
-                                          nptbins,ptmin,ptmax,20,0,20);
-    fhPtNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
-    fhPtNPileUpSPDVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
-    outputContainer->Add(fhPtNPileUpSPDVtxTimeCut);
-         
-    fhPtNPileUpTrkVtxTimeCut  = new TH2F ("hPt_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
-                                          nptbins,ptmin,ptmax, 20,0,20 );
-    fhPtNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
-    fhPtNPileUpTrkVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
-    outputContainer->Add(fhPtNPileUpTrkVtxTimeCut);
-    
-               fhPtNPileUpSPDVtxTimeCut2  = new TH2F ("hPt_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
-                                           nptbins,ptmin,ptmax,20,0,20);
-    fhPtNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
-    fhPtNPileUpSPDVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
-    outputContainer->Add(fhPtNPileUpSPDVtxTimeCut2);
-         
-    fhPtNPileUpTrkVtxTimeCut2  = new TH2F ("hPt_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex, -25 < tof < 75 ns",
-                                           nptbins,ptmin,ptmax, 20,0,20 );
-    fhPtNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
-    fhPtNPileUpTrkVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
-    outputContainer->Add(fhPtNPileUpTrkVtxTimeCut2);
-    
     fhPtPhotonNPileUpSPDVtx  = new TH2F ("hPtPhoton_NPileUpVertSPD","pT of cluster vs N pile-up SPD vertex",
                                          nptbins,ptmin,ptmax,20,0,20);
     fhPtPhotonNPileUpSPDVtx->SetYTitle("# vertex ");
-    fhPtPhotonNPileUpSPDVtx->SetXTitle("p_{T} (GeV/c)");
+    fhPtPhotonNPileUpSPDVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
     outputContainer->Add(fhPtPhotonNPileUpSPDVtx);
          
     fhPtPhotonNPileUpTrkVtx  = new TH2F ("hPtPhoton_NPileUpVertTracks","pT of cluster vs N pile-up Tracks vertex",
                                          nptbins,ptmin,ptmax, 20,0,20 );
     fhPtPhotonNPileUpTrkVtx->SetYTitle("# vertex ");
-    fhPtPhotonNPileUpTrkVtx->SetXTitle("p_{T} (GeV/c)");
+    fhPtPhotonNPileUpTrkVtx->SetXTitle("#it{p}_{T} (GeV/#it{c})");
     outputContainer->Add(fhPtPhotonNPileUpTrkVtx);
          
     fhPtPhotonNPileUpSPDVtxTimeCut  = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut","pT of cluster vs N pile-up SPD vertex, |tof| < 25 ns",
                                                 nptbins,ptmin,ptmax,20,0,20);
     fhPtPhotonNPileUpSPDVtxTimeCut->SetYTitle("# vertex ");
-    fhPtPhotonNPileUpSPDVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
+    fhPtPhotonNPileUpSPDVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
     outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut);
          
     fhPtPhotonNPileUpTrkVtxTimeCut  = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut","pT of cluster vs N pile-up Tracks vertex, |tof| < 25 ns",
                                                 nptbins,ptmin,ptmax, 20,0,20 );
     fhPtPhotonNPileUpTrkVtxTimeCut->SetYTitle("# vertex ");
-    fhPtPhotonNPileUpTrkVtxTimeCut->SetXTitle("p_{T} (GeV/c)");
+    fhPtPhotonNPileUpTrkVtxTimeCut->SetXTitle("#it{p}_{T} (GeV/#it{c})");
     outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut);
          
     fhPtPhotonNPileUpSPDVtxTimeCut2  = new TH2F ("hPtPhoton_NPileUpVertSPD_TimeCut2","pT of cluster vs N pile-up SPD vertex, -25 < tof < 75 ns",
                                                  nptbins,ptmin,ptmax,20,0,20);
     fhPtPhotonNPileUpSPDVtxTimeCut2->SetYTitle("# vertex ");
-    fhPtPhotonNPileUpSPDVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
+    fhPtPhotonNPileUpSPDVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
     outputContainer->Add(fhPtPhotonNPileUpSPDVtxTimeCut2);
          
     fhPtPhotonNPileUpTrkVtxTimeCut2  = new TH2F ("hPtPhoton_NPileUpVertTracks_TimeCut2","pT of cluster vs N pile-up Tracks vertex,  -25 < tof < 75 ns",
                                                  nptbins,ptmin,ptmax, 20,0,20 );
     fhPtPhotonNPileUpTrkVtxTimeCut2->SetYTitle("# vertex ");
-    fhPtPhotonNPileUpTrkVtxTimeCut2->SetXTitle("p_{T} (GeV/c)");
+    fhPtPhotonNPileUpTrkVtxTimeCut2->SetXTitle("#it{p}_{T} (GeV/#it{c})");
     outputContainer->Add(fhPtPhotonNPileUpTrkVtxTimeCut2);
     
   }
+
+  
   
   if(IsDataMC())
   {
@@ -3095,65 +1786,65 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhMCE[i]  = new TH1F(Form("hE_MC%s",pname[i].Data()),
                            Form("cluster from %s : E ",ptype[i].Data()),
                            nptbins,ptmin,ptmax);
-      fhMCE[i]->SetXTitle("E (GeV)");
+      fhMCE[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhMCE[i]) ;
       
       fhMCPt[i]  = new TH1F(Form("hPt_MC%s",pname[i].Data()),
-                            Form("cluster from %s : p_{T} ",ptype[i].Data()),
+                            Form("cluster from %s : #it{p}_{T} ",ptype[i].Data()),
                             nptbins,ptmin,ptmax);
-      fhMCPt[i]->SetXTitle("p_{T} (GeV/c)");
+      fhMCPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
       outputContainer->Add(fhMCPt[i]) ;
       
       fhMCEta[i]  = new TH2F(Form("hEta_MC%s",pname[i].Data()),
                              Form("cluster from %s : #eta ",ptype[i].Data()),
                              nptbins,ptmin,ptmax,netabins,etamin,etamax);
       fhMCEta[i]->SetYTitle("#eta");
-      fhMCEta[i]->SetXTitle("E (GeV)");
+      fhMCEta[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhMCEta[i]) ;
       
       fhMCPhi[i]  = new TH2F(Form("hPhi_MC%s",pname[i].Data()),
                              Form("cluster from %s : #phi ",ptype[i].Data()),
                              nptbins,ptmin,ptmax,nphibins,phimin,phimax);
       fhMCPhi[i]->SetYTitle("#phi (rad)");
-      fhMCPhi[i]->SetXTitle("E (GeV)");
+      fhMCPhi[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhMCPhi[i]) ;
       
       
       fhMCDeltaE[i]  = new TH2F (Form("hDeltaE_MC%s",pname[i].Data()),
                                  Form("MC - Reco E from %s",pname[i].Data()),
                                  nptbins,ptmin,ptmax, 200,-50,50);
-      fhMCDeltaE[i]->SetYTitle("#Delta E (GeV)");
-      fhMCDeltaE[i]->SetXTitle("E (GeV)");
+      fhMCDeltaE[i]->SetYTitle("#Delta #it{E} (GeV)");
+      fhMCDeltaE[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhMCDeltaE[i]);
       
       fhMCDeltaPt[i]  = new TH2F (Form("hDeltaPt_MC%s",pname[i].Data()),
-                                  Form("MC - Reco p_{T} from %s",pname[i].Data()),
+                                  Form("MC - Reco #it{p}_{T} from %s",pname[i].Data()),
                                   nptbins,ptmin,ptmax, 200,-50,50);
-      fhMCDeltaPt[i]->SetXTitle("p_{T,rec} (GeV/c)");
-      fhMCDeltaPt[i]->SetYTitle("#Delta p_{T} (GeV/c)");
+      fhMCDeltaPt[i]->SetXTitle("p_{T,rec} (GeV/#it{c})");
+      fhMCDeltaPt[i]->SetYTitle("#Delta #it{p}_{T} (GeV/#it{c})");
       outputContainer->Add(fhMCDeltaPt[i]);
       
       fhMC2E[i]  = new TH2F (Form("h2E_MC%s",pname[i].Data()),
                              Form("E distribution, reconstructed vs generated from %s",pname[i].Data()),
                              nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
-      fhMC2E[i]->SetXTitle("E_{rec} (GeV)");
-      fhMC2E[i]->SetYTitle("E_{gen} (GeV)");
+      fhMC2E[i]->SetXTitle("#it{E}_{rec} (GeV)");
+      fhMC2E[i]->SetYTitle("#it{E}_{gen} (GeV)");
       outputContainer->Add(fhMC2E[i]);
       
       fhMC2Pt[i]  = new TH2F (Form("h2Pt_MC%s",pname[i].Data()),
                               Form("p_T distribution, reconstructed vs generated from %s",pname[i].Data()),
                               nptbins,ptmin,ptmax,nptbins,ptmin,ptmax);
-      fhMC2Pt[i]->SetXTitle("p_{T,rec} (GeV/c)");
-      fhMC2Pt[i]->SetYTitle("p_{T,gen} (GeV/c)");
+      fhMC2Pt[i]->SetXTitle("p_{T,rec} (GeV/#it{c})");
+      fhMC2Pt[i]->SetYTitle("p_{T,gen} (GeV/#it{c})");
       outputContainer->Add(fhMC2Pt[i]);
       
       
     }
     
-    TString pptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}","hadron?",
+    TString pptype[] = { "#gamma", "#gamma_{#pi decay}","#gamma_{other decay}",
       "#gamma_{prompt}","#gamma_{fragmentation}","#gamma_{ISR}"} ;
     
-    TString ppname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay","Hadron",
+    TString ppname[] = { "Photon","PhotonPi0Decay","PhotonOtherDecay",
       "PhotonPrompt","PhotonFragmentation","PhotonISR"} ;
     
     for(Int_t i = 0; i < fNPrimaryHistograms; i++)
@@ -3161,68 +1852,68 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
       fhEPrimMC[i]  = new TH1F(Form("hEPrim_MC%s",ppname[i].Data()),
                                Form("primary photon %s : E ",pptype[i].Data()),
                                nptbins,ptmin,ptmax);
-      fhEPrimMC[i]->SetXTitle("E (GeV)");
+      fhEPrimMC[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhEPrimMC[i]) ;
       
       fhPtPrimMC[i]  = new TH1F(Form("hPtPrim_MC%s",ppname[i].Data()),
-                                Form("primary photon %s : p_{T} ",pptype[i].Data()),
+                                Form("primary photon %s : #it{p}_{T} ",pptype[i].Data()),
                                 nptbins,ptmin,ptmax);
-      fhPtPrimMC[i]->SetXTitle("p_{T} (GeV/c)");
+      fhPtPrimMC[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
       outputContainer->Add(fhPtPrimMC[i]) ;
       
       fhYPrimMC[i]  = new TH2F(Form("hYPrim_MC%s",ppname[i].Data()),
                                Form("primary photon %s : Rapidity ",pptype[i].Data()),
                                nptbins,ptmin,ptmax,200,-2,2);
       fhYPrimMC[i]->SetYTitle("Rapidity");
-      fhYPrimMC[i]->SetXTitle("E (GeV)");
+      fhYPrimMC[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhYPrimMC[i]) ;
       
       fhEtaPrimMC[i]  = new TH2F(Form("hEtaPrim_MC%s",ppname[i].Data()),
                                Form("primary photon %s : #eta",pptype[i].Data()),
                                nptbins,ptmin,ptmax,200,-2,2);
       fhEtaPrimMC[i]->SetYTitle("#eta");
-      fhEtaPrimMC[i]->SetXTitle("E (GeV)");
+      fhEtaPrimMC[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhEtaPrimMC[i]) ;
       
       fhPhiPrimMC[i]  = new TH2F(Form("hPhiPrim_MC%s",ppname[i].Data()),
                                  Form("primary photon %s : #phi ",pptype[i].Data()),
                                  nptbins,ptmin,ptmax,nphibins,0,TMath::TwoPi());
       fhPhiPrimMC[i]->SetYTitle("#phi (rad)");
-      fhPhiPrimMC[i]->SetXTitle("E (GeV)");
+      fhPhiPrimMC[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhPhiPrimMC[i]) ;
       
       
       fhEPrimMCAcc[i]  = new TH1F(Form("hEPrimAcc_MC%s",ppname[i].Data()),
                                   Form("primary photon %s in acceptance: E ",pptype[i].Data()),
                                   nptbins,ptmin,ptmax);
-      fhEPrimMCAcc[i]->SetXTitle("E (GeV)");
+      fhEPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhEPrimMCAcc[i]) ;
       
       fhPtPrimMCAcc[i]  = new TH1F(Form("hPtPrimAcc_MC%s",ppname[i].Data()),
-                                   Form("primary photon %s in acceptance: p_{T} ",pptype[i].Data()),
+                                   Form("primary photon %s in acceptance: #it{p}_{T} ",pptype[i].Data()),
                                    nptbins,ptmin,ptmax);
-      fhPtPrimMCAcc[i]->SetXTitle("p_{T} (GeV/c)");
+      fhPtPrimMCAcc[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
       outputContainer->Add(fhPtPrimMCAcc[i]) ;
       
       fhYPrimMCAcc[i]  = new TH2F(Form("hYPrimAcc_MC%s",ppname[i].Data()),
                                   Form("primary photon %s in acceptance: Rapidity ",pptype[i].Data()),
                                   nptbins,ptmin,ptmax,100,-1,1);
       fhYPrimMCAcc[i]->SetYTitle("Rapidity");
-      fhYPrimMCAcc[i]->SetXTitle("E (GeV)");
+      fhYPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhYPrimMCAcc[i]) ;
 
       fhEtaPrimMCAcc[i]  = new TH2F(Form("hEtaPrimAcc_MC%s",ppname[i].Data()),
                                   Form("primary photon %s in acceptance: #eta ",pptype[i].Data()),
                                   nptbins,ptmin,ptmax,netabins,etamin,etamax);
       fhEtaPrimMCAcc[i]->SetYTitle("#eta");
-      fhEtaPrimMCAcc[i]->SetXTitle("E (GeV)");
+      fhEtaPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhEtaPrimMCAcc[i]) ;
       
       fhPhiPrimMCAcc[i]  = new TH2F(Form("hPhiPrimAcc_MC%s",ppname[i].Data()),
                                     Form("primary photon %s in acceptance: #phi ",pptype[i].Data()),
                                     nptbins,ptmin,ptmax,nphibins,phimin,phimax);
       fhPhiPrimMCAcc[i]->SetYTitle("#phi (rad)");
-      fhPhiPrimMCAcc[i]->SetXTitle("E (GeV)");
+      fhPhiPrimMCAcc[i]->SetXTitle("#it{E} (GeV)");
       outputContainer->Add(fhPhiPrimMCAcc[i]) ;
       
     }
@@ -3239,35 +1930,35 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
                                     Form("cluster from %s : E vs #lambda_{0}^{2}",ptypess[i].Data()),
                                     nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhMCELambda0[i]->SetYTitle("#lambda_{0}^{2}");
-        fhMCELambda0[i]->SetXTitle("E (GeV)");
+        fhMCELambda0[i]->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhMCELambda0[i]) ;
         
         fhMCELambda1[i]  = new TH2F(Form("hELambda1_MC%s",pnamess[i].Data()),
                                     Form("cluster from %s : E vs #lambda_{1}^{2}",ptypess[i].Data()),
                                     nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhMCELambda1[i]->SetYTitle("#lambda_{1}^{2}");
-        fhMCELambda1[i]->SetXTitle("E (GeV)");
+        fhMCELambda1[i]->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhMCELambda1[i]) ;
         
         fhMCEDispersion[i]  = new TH2F(Form("hEDispersion_MC%s",pnamess[i].Data()),
                                        Form("cluster from %s : E vs dispersion^{2}",ptypess[i].Data()),
                                        nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhMCEDispersion[i]->SetYTitle("D^{2}");
-        fhMCEDispersion[i]->SetXTitle("E (GeV)");
+        fhMCEDispersion[i]->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhMCEDispersion[i]) ;
         
         fhMCNCellsE[i]  = new TH2F (Form("hNCellsE_MC%s",pnamess[i].Data()),
                                     Form("# of cells in cluster from %s vs E of clusters",ptypess[i].Data()),
                                     nptbins,ptmin,ptmax, nbins,nmin,nmax);
-        fhMCNCellsE[i]->SetXTitle("E (GeV)");
+        fhMCNCellsE[i]->SetXTitle("#it{E} (GeV)");
         fhMCNCellsE[i]->SetYTitle("# of cells in cluster");
         outputContainer->Add(fhMCNCellsE[i]);
         
         fhMCMaxCellDiffClusterE[i]  = new TH2F (Form("hMaxCellDiffClusterE_MC%s",pnamess[i].Data()),
                                                 Form("energy vs difference of cluster energy from %s - max cell energy / cluster energy, good clusters",ptypess[i].Data()),
                                                 nptbins,ptmin,ptmax, 500,0,1.);
-        fhMCMaxCellDiffClusterE[i]->SetXTitle("E_{cluster} (GeV) ");
-        fhMCMaxCellDiffClusterE[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+        fhMCMaxCellDiffClusterE[i]->SetXTitle("#it{E}_{cluster} (GeV) ");
+        fhMCMaxCellDiffClusterE[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
         outputContainer->Add(fhMCMaxCellDiffClusterE[i]);
         
         if(!fFillOnlySimpleSSHisto)
@@ -3276,78 +1967,78 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
                                                            Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
                                                            ssbins,ssmin,ssmax,500,0,1.);
           fhMCLambda0vsClusterMaxCellDiffE0[i]->SetXTitle("#lambda_{0}^{2}");
-          fhMCLambda0vsClusterMaxCellDiffE0[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+          fhMCLambda0vsClusterMaxCellDiffE0[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
           outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE0[i]) ;
           
           fhMCLambda0vsClusterMaxCellDiffE2[i]  = new TH2F(Form("hLambda0vsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
                                                            Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
                                                            ssbins,ssmin,ssmax,500,0,1.);
           fhMCLambda0vsClusterMaxCellDiffE2[i]->SetXTitle("#lambda_{0}^{2}");
-          fhMCLambda0vsClusterMaxCellDiffE2[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+          fhMCLambda0vsClusterMaxCellDiffE2[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
           outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE2[i]) ;
           
           fhMCLambda0vsClusterMaxCellDiffE6[i]  = new TH2F(Form("hLambda0vsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
-                                                           Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, E > 6 GeV",ptypess[i].Data()),
+                                                           Form("cluster from %s : #lambda^{2}_{0} vs fraction of energy carried by max cell, #it{E} > 6 GeV",ptypess[i].Data()),
                                                            ssbins,ssmin,ssmax,500,0,1.);
           fhMCLambda0vsClusterMaxCellDiffE6[i]->SetXTitle("#lambda_{0}^{2}");
-          fhMCLambda0vsClusterMaxCellDiffE6[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+          fhMCLambda0vsClusterMaxCellDiffE6[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
           outputContainer->Add(fhMCLambda0vsClusterMaxCellDiffE6[i]) ;
           
           fhMCNCellsvsClusterMaxCellDiffE0[i]  = new TH2F(Form("hNCellsvsClusterMaxCellDiffE0_MC%s",pnamess[i].Data()),
                                                           Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E < 2 GeV",ptypess[i].Data()),
                                                           nbins/5,nmin,nmax/5,500,0,1.);
           fhMCNCellsvsClusterMaxCellDiffE0[i]->SetXTitle("N cells in cluster");
-          fhMCNCellsvsClusterMaxCellDiffE0[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+          fhMCNCellsvsClusterMaxCellDiffE0[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
           outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE0[i]) ;
           
           fhMCNCellsvsClusterMaxCellDiffE2[i]  = new TH2F(Form("hNCellsvsClusterMaxCellDiffE2_MC%s",pnamess[i].Data()),
                                                           Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, 2< E < 6 GeV",ptypess[i].Data()),
                                                           nbins/5,nmin,nmax/5,500,0,1.);
           fhMCNCellsvsClusterMaxCellDiffE2[i]->SetXTitle("N cells in cluster");
-          fhMCNCellsvsClusterMaxCellDiffE2[i]->SetYTitle("(E_{cluster} - E_{cell max})/ E_{cluster}");
+          fhMCNCellsvsClusterMaxCellDiffE2[i]->SetYTitle("(#it{E}_{cluster} - #it{E}_{cell max})/ #it{E}_{cluster}");
           outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE2[i]) ;
           
           fhMCNCellsvsClusterMaxCellDiffE6[i]  = new TH2F(Form("hNCellsvsClusterMaxCellDiffE6_MC%s",pnamess[i].Data()),
-                                                          Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, E > 6 GeV",ptypess[i].Data()),
+                                                          Form("cluster from %s : N cells in cluster vs fraction of energy carried by max cell, #it{E} > 6 GeV",ptypess[i].Data()),
                                                           nbins/5,nmin,nmax/5,500,0,1.);
           fhMCNCellsvsClusterMaxCellDiffE6[i]->SetXTitle("N cells in cluster");
-          fhMCNCellsvsClusterMaxCellDiffE6[i]->SetYTitle("E (GeV)");
+          fhMCNCellsvsClusterMaxCellDiffE6[i]->SetYTitle("#it{E} (GeV)");
           outputContainer->Add(fhMCNCellsvsClusterMaxCellDiffE6[i]) ;
           
-          if(fCalorimeter=="EMCAL")
+          if(GetCalorimeter()==kEMCAL)
           {
             fhMCEDispEta[i]  = new TH2F (Form("hEDispEtaE_MC%s",pnamess[i].Data()),
                                          Form("cluster from %s : #sigma^{2}_{#eta #eta} = #Sigma w_{i}(#eta_{i} - <#eta>)^{2}/ #Sigma w_{i} vs E",ptypess[i].Data()),
                                          nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
-            fhMCEDispEta[i]->SetXTitle("E (GeV)");
+            fhMCEDispEta[i]->SetXTitle("#it{E} (GeV)");
             fhMCEDispEta[i]->SetYTitle("#sigma^{2}_{#eta #eta}");
             outputContainer->Add(fhMCEDispEta[i]);
             
             fhMCEDispPhi[i]  = new TH2F (Form("hEDispPhiE_MC%s",pnamess[i].Data()),
                                          Form("cluster from %s : #sigma^{2}_{#phi #phi} = #Sigma w_{i}(#phi_{i} - <#phi>)^{2} / #Sigma w_{i} vs E",ptypess[i].Data()),
                                          nptbins,ptmin,ptmax, ssbins,ssmin,ssmax);
-            fhMCEDispPhi[i]->SetXTitle("E (GeV)");
+            fhMCEDispPhi[i]->SetXTitle("#it{E} (GeV)");
             fhMCEDispPhi[i]->SetYTitle("#sigma^{2}_{#phi #phi}");
             outputContainer->Add(fhMCEDispPhi[i]);
             
             fhMCESumEtaPhi[i]  = new TH2F (Form("hESumEtaPhiE_MC%s",pnamess[i].Data()),
                                            Form("cluster from %s : #delta^{2}_{#eta #phi} = #Sigma w_{i}(#phi_{i} #eta_{i} ) / #Sigma w_{i} - <#phi><#eta> vs E",ptypess[i].Data()),
                                            nptbins,ptmin,ptmax, 2*ssbins,-ssmax,ssmax);
-            fhMCESumEtaPhi[i]->SetXTitle("E (GeV)");
+            fhMCESumEtaPhi[i]->SetXTitle("#it{E} (GeV)");
             fhMCESumEtaPhi[i]->SetYTitle("#delta^{2}_{#eta #phi}");
             outputContainer->Add(fhMCESumEtaPhi[i]);
             
             fhMCEDispEtaPhiDiff[i]  = new TH2F (Form("hEDispEtaPhiDiffE_MC%s",pnamess[i].Data()),
                                                 Form("cluster from %s : #sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta} vs E",ptypess[i].Data()),
                                                 nptbins,ptmin,ptmax,200,-10,10);
-            fhMCEDispEtaPhiDiff[i]->SetXTitle("E (GeV)");
+            fhMCEDispEtaPhiDiff[i]->SetXTitle("#it{E} (GeV)");
             fhMCEDispEtaPhiDiff[i]->SetYTitle("#sigma^{2}_{#phi #phi}-#sigma^{2}_{#eta #eta}");
             outputContainer->Add(fhMCEDispEtaPhiDiff[i]);
             
             fhMCESphericity[i]  = new TH2F (Form("hESphericity_MC%s",pnamess[i].Data()),
                                             Form("cluster from %s : (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi}) vs E",ptypess[i].Data()),
                                             nptbins,ptmin,ptmax, 200,-1,1);
-            fhMCESphericity[i]->SetXTitle("E (GeV)");
+            fhMCESphericity[i]->SetXTitle("#it{E} (GeV)");
             fhMCESphericity[i]->SetYTitle("s = (#sigma^{2}_{#phi #phi} - #sigma^{2}_{#eta #eta}) / (#sigma^{2}_{#eta #eta} + #sigma^{2}_{#phi #phi})");
             outputContainer->Add(fhMCESphericity[i]);
             
@@ -3384,21 +2075,21 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
                                                 "cluster from Photon : E vs #lambda_{0}^{2}",
                                                 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhMCPhotonELambda0NoOverlap->SetYTitle("#lambda_{0}^{2}");
-        fhMCPhotonELambda0NoOverlap->SetXTitle("E (GeV)");
+        fhMCPhotonELambda0NoOverlap->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhMCPhotonELambda0NoOverlap) ;
         
         fhMCPhotonELambda0TwoOverlap  = new TH2F("hELambda0_MCPhoton_TwoOverlap",
                                                  "cluster from Photon : E vs #lambda_{0}^{2}",
                                                  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhMCPhotonELambda0TwoOverlap->SetYTitle("#lambda_{0}^{2}");
-        fhMCPhotonELambda0TwoOverlap->SetXTitle("E (GeV)");
+        fhMCPhotonELambda0TwoOverlap->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhMCPhotonELambda0TwoOverlap) ;
         
         fhMCPhotonELambda0NOverlap  = new TH2F("hELambda0_MCPhoton_NOverlap",
                                                "cluster from Photon : E vs #lambda_{0}^{2}",
                                                nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhMCPhotonELambda0NOverlap->SetYTitle("#lambda_{0}^{2}");
-        fhMCPhotonELambda0NOverlap->SetXTitle("E (GeV)");
+        fhMCPhotonELambda0NOverlap->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhMCPhotonELambda0NOverlap) ;
         
       } //No embedding
@@ -3410,63 +2101,63 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
                                                    "Energy Fraction of embedded signal versus cluster energy",
                                                    nptbins,ptmin,ptmax,100,0.,1.);
         fhEmbeddedSignalFractionEnergy->SetYTitle("Fraction");
-        fhEmbeddedSignalFractionEnergy->SetXTitle("E (GeV)");
+        fhEmbeddedSignalFractionEnergy->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbeddedSignalFractionEnergy) ;
         
         fhEmbedPhotonELambda0FullSignal  = new TH2F("hELambda0_EmbedPhoton_FullSignal",
                                                     "cluster from Photon embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
                                                     nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPhotonELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPhotonELambda0FullSignal->SetXTitle("E (GeV)");
+        fhEmbedPhotonELambda0FullSignal->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPhotonELambda0FullSignal) ;
         
         fhEmbedPhotonELambda0MostlySignal  = new TH2F("hELambda0_EmbedPhoton_MostlySignal",
                                                       "cluster from Photon embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
                                                       nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPhotonELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPhotonELambda0MostlySignal->SetXTitle("E (GeV)");
+        fhEmbedPhotonELambda0MostlySignal->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPhotonELambda0MostlySignal) ;
         
         fhEmbedPhotonELambda0MostlyBkg  = new TH2F("hELambda0_EmbedPhoton_MostlyBkg",
                                                    "cluster from Photon embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
                                                    nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPhotonELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPhotonELambda0MostlyBkg->SetXTitle("E (GeV)");
+        fhEmbedPhotonELambda0MostlyBkg->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPhotonELambda0MostlyBkg) ;
         
         fhEmbedPhotonELambda0FullBkg  = new TH2F("hELambda0_EmbedPhoton_FullBkg",
                                                  "cluster from Photonm embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
                                                  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPhotonELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPhotonELambda0FullBkg->SetXTitle("E (GeV)");
+        fhEmbedPhotonELambda0FullBkg->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPhotonELambda0FullBkg) ;
         
         fhEmbedPi0ELambda0FullSignal  = new TH2F("hELambda0_EmbedPi0_FullSignal",
                                                  "cluster from Pi0 embedded with more than 90% energy in cluster : E vs #lambda_{0}^{2}",
                                                  nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPi0ELambda0FullSignal->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPi0ELambda0FullSignal->SetXTitle("E (GeV)");
+        fhEmbedPi0ELambda0FullSignal->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPi0ELambda0FullSignal) ;
         
         fhEmbedPi0ELambda0MostlySignal  = new TH2F("hELambda0_EmbedPi0_MostlySignal",
                                                    "cluster from Pi0 embedded with 50% to 90% energy in cluster : E vs #lambda_{0}^{2}",
                                                    nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPi0ELambda0MostlySignal->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPi0ELambda0MostlySignal->SetXTitle("E (GeV)");
+        fhEmbedPi0ELambda0MostlySignal->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPi0ELambda0MostlySignal) ;
         
         fhEmbedPi0ELambda0MostlyBkg  = new TH2F("hELambda0_EmbedPi0_MostlyBkg",
                                                 "cluster from Pi0 embedded with 10% to 50% energy in cluster : E vs #lambda_{0}^{2}",
                                                 nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPi0ELambda0MostlyBkg->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPi0ELambda0MostlyBkg->SetXTitle("E (GeV)");
+        fhEmbedPi0ELambda0MostlyBkg->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPi0ELambda0MostlyBkg) ;
         
         fhEmbedPi0ELambda0FullBkg  = new TH2F("hELambda0_EmbedPi0_FullBkg",
                                               "cluster from Pi0 embedded with 0% to 10% energy in cluster : E vs #lambda_{0}^{2}",
                                               nptbins,ptmin,ptmax,ssbins,ssmin,ssmax);
         fhEmbedPi0ELambda0FullBkg->SetYTitle("#lambda_{0}^{2}");
-        fhEmbedPi0ELambda0FullBkg->SetXTitle("E (GeV)");
+        fhEmbedPi0ELambda0FullBkg->SetXTitle("#it{E} (GeV)");
         outputContainer->Add(fhEmbedPi0ELambda0FullBkg) ;
         
       }// embedded histograms
@@ -3483,19 +2174,14 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
 //_______________________
 void AliAnaPhoton::Init()
 {
-  
   //Init
   //Do some checks
-  if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD())
-  {
-    printf("AliAnaPhoton::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");
-    abort();
-  }
-  else  if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD())
-  {
-    printf("AliAnaPhoton::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
-    abort();
-  }
+  
+  if      ( GetCalorimeter() == kPHOS  && !GetReader()->IsPHOSSwitchedOn()  && NewOutputAOD() )
+    AliFatal("!!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!");
+  else  if( GetCalorimeter() == kEMCAL && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD() )
+    AliFatal("!!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!");
   
   if(GetReader()->GetDataType() == AliCaloTrackReader::kMC) GetCaloPID()->SwitchOnBayesian();
   
@@ -3508,7 +2194,6 @@ void AliAnaPhoton::InitParameters()
   //Initialize the parameters of the analysis.
   AddToHistogramsName("AnaPhoton_");
   
-  fCalorimeter = "EMCAL" ;
   fMinDist     = 2.;
   fMinDist2    = 4.;
   fMinDist3    = 5.;
@@ -3533,12 +2218,12 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
   //Select the Calorimeter of the photon
   TObjArray * pl = 0x0;
   AliVCaloCells* cells    = 0;
-  if      (fCalorimeter == "PHOS" )
+  if      (GetCalorimeter() == kPHOS )
   {
     pl    = GetPHOSClusters();
     cells = GetPHOSCells();
   }
-  else if (fCalorimeter == "EMCAL")
+  else if (GetCalorimeter() == kEMCAL)
   {
     pl    = GetEMCALClusters();
     cells = GetEMCALCells();
@@ -3546,20 +2231,30 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
   
   if(!pl)
   {
-    Info("MakeAnalysisFillAOD","TObjArray with %s clusters is NULL!\n",fCalorimeter.Data());
+    AliWarning(Form("TObjArray with %s clusters is NULL!",GetCalorimeterString().Data()));
     return;
   }
   
-  FillPileUpHistogramsPerEvent();
-  
   // Loop on raw clusters before filtering in the reader and fill control histogram
-  if((GetReader()->GetEMCALClusterListName()=="" && fCalorimeter=="EMCAL") || fCalorimeter=="PHOS")
+  if((GetReader()->GetEMCALClusterListName()=="" && GetCalorimeter()==kEMCAL) || GetCalorimeter()==kPHOS)
   {
     for(Int_t iclus = 0; iclus < GetReader()->GetInputEvent()->GetNumberOfCaloClusters(); iclus++ )
     {
       AliVCluster * clus = GetReader()->GetInputEvent()->GetCaloCluster(iclus);
-      if     (fCalorimeter == "PHOS"  && clus->IsPHOS()  && clus->E() > GetReader()->GetPHOSPtMin() ) fhClusterCuts[0]->Fill(clus->E());
-      else if(fCalorimeter == "EMCAL" && clus->IsEMCAL() && clus->E() > GetReader()->GetEMCALPtMin()) fhClusterCuts[0]->Fill(clus->E());
+      if     (GetCalorimeter() == kPHOS  && clus->IsPHOS()  && clus->E() > GetReader()->GetPHOSPtMin() )
+      {
+        fhClusterCutsE [0]->Fill(clus->E());
+        
+        clus->GetMomentum(fMomentum,GetVertex(0)) ;
+        fhClusterCutsPt[0]->Fill(fMomentum.Pt());
+      }
+      else if(GetCalorimeter() == kEMCAL && clus->IsEMCAL() && clus->E() > GetReader()->GetEMCALPtMin())
+      {
+        fhClusterCutsE [0]->Fill(clus->E());
+        
+        clus->GetMomentum(fMomentum,GetVertex(0)) ;
+        fhClusterCutsPt[0]->Fill(fMomentum.Pt());
+      }
     }
   }
   else
@@ -3577,97 +2272,21 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
       for (Int_t iclus =  0; iclus <  nclusters; iclus++)
       {
         AliVCluster * clus = dynamic_cast<AliVCluster*> (clusterList->At(iclus));
-        if(clus)fhClusterCuts[0]->Fill(clus->E());
-      }
-    }
-  }
-  
-  // Fill some trigger related histograms
-  Int_t  idTrig = GetReader()->GetTriggerClusterIndex();
-  Bool_t exotic = GetReader()->IsExoticEvent();
-  Bool_t bad    = GetReader()->IsBadCellTriggerEvent();
-  
-  if( fFillEMCALBCHistograms && fCalorimeter=="EMCAL" &&
-     ( bad || exotic )  && idTrig >= 0 && !GetReader()->AreBadTriggerEventsRemoved())
-  {
-    //    printf("Index %d, Id %d,  bad %d, exo %d\n",
-    //           GetReader()->GetTriggerClusterIndex(),
-    //           GetReader()->GetTriggerClusterId(),
-    //           GetReader()->IsBadCellTriggerEvent(),
-    //           GetReader()->IsExoticEvent() );
-    
-    TClonesArray * clusterList = 0;
-    TString  clusterListName   = GetReader()->GetEMCALClusterListName();
-    if     (GetReader()->GetInputEvent()->FindListObject(clusterListName))
-      clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetInputEvent() ->FindListObject(clusterListName));
-    else if(GetReader()->GetOutputEvent())
-      clusterList = dynamic_cast<TClonesArray*> (GetReader()->GetOutputEvent()->FindListObject(clusterListName));
-    
-    AliVCluster  *  badClusTrig = 0;
-    if(clusterList) badClusTrig = (AliVCluster*) clusterList->At(idTrig);
-    else            badClusTrig = GetReader()->GetInputEvent()->GetCaloCluster(idTrig);
-    
-    if(!badClusTrig)
-      printf("AliAnaPhoton::MakeAnalysisFillAOD() - No cluster (bad-exotic trigger) found with requested index %d \n",idTrig);
-    else
-    {
-      TLorentzVector momBadClus;
-      
-      badClusTrig->GetMomentum(momBadClus,GetVertex(0));
-      
-      Float_t etaclusterBad = momBadClus.Eta();
-      Float_t phiclusterBad = momBadClus.Phi();
-      if( phiclusterBad < 0 ) phiclusterBad+=TMath::TwoPi();
-      Float_t tofclusterBad = badClusTrig->GetTOF()*1.e9;
-      Float_t eclusterBad   = badClusTrig->E();
-      
-      if( bad && exotic )
-      {
-        if(GetReader()->IsTriggerMatched())
-        {
-          fhEtaPhiTriggerEMCALBCBadExoticCluster->Fill(etaclusterBad, phiclusterBad);
-          fhTimeTriggerEMCALBCBadExoticCluster  ->Fill(eclusterBad,   tofclusterBad);
-        }
-        else
-        {
-          fhEtaPhiTriggerEMCALBCUMBadExoticCluster->Fill(etaclusterBad, phiclusterBad);
-          fhTimeTriggerEMCALBCUMBadExoticCluster  ->Fill(eclusterBad,   tofclusterBad);
-        }
-      }
-      else if( bad && !exotic )
-      {
-        if(GetReader()->IsTriggerMatched())
-        {
-          fhEtaPhiTriggerEMCALBCBadCluster->Fill(etaclusterBad, phiclusterBad);
-          fhTimeTriggerEMCALBCBadCluster  ->Fill(eclusterBad,   tofclusterBad);
-        }
-        else
-        {
-          fhEtaPhiTriggerEMCALBCUMBadCluster->Fill(etaclusterBad, phiclusterBad);
-          fhTimeTriggerEMCALBCUMBadCluster  ->Fill(eclusterBad,   tofclusterBad);
-        }
-      }// Bad cluster trigger
-      else if( !bad && exotic )
-      {
-        if(GetReader()->IsTriggerMatched())
-        {
-          fhEtaPhiTriggerEMCALBCExoticCluster->Fill(etaclusterBad, phiclusterBad);
-          fhTimeTriggerEMCALBCExoticCluster  ->Fill(eclusterBad, tofclusterBad);
-        }
-        else
+        if(clus)
         {
-          fhEtaPhiTriggerEMCALBCUMExoticCluster->Fill(etaclusterBad, phiclusterBad);
-          fhTimeTriggerEMCALBCUMExoticCluster  ->Fill(eclusterBad, tofclusterBad);
+          fhClusterCutsE [0]->Fill(clus->E());
+          
+          clus->GetMomentum(fMomentum,GetVertex(0)) ;
+          fhClusterCutsPt[0]->Fill(fMomentum.Pt());
         }
       }
-    }// cluster exists
-  } // study bad/exotic trigger BC
+    }
+  }
   
-  //Init arrays, variables, get number of clusters
-  TLorentzVector mom, mom2 ;
+  // Init arrays, variables, get number of clusters
   Int_t nCaloClusters = pl->GetEntriesFast();
   
-  if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillAOD() - input %s cluster entries %d\n", fCalorimeter.Data(), nCaloClusters);
+  AliDebug(1,Form("Input %s cluster entries %d", GetCalorimeterString().Data(), nCaloClusters));
   
   //----------------------------------------------------
   // Fill AOD with PHOS/EMCAL AliAODPWG4Particle objects
@@ -3690,31 +2309,31 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
     //Cluster selection, not charged, with photon id and in fiducial cut
     if(GetReader()->GetDataType() != AliCaloTrackReader::kMC)
     {
-      calo->GetMomentum(mom,GetVertex(evtIndex)) ;
+      calo->GetMomentum(fMomentum,GetVertex(evtIndex)) ;
     }//Assume that come from vertex in straight line
     else
     {
       Double_t vertex[]={0,0,0};
-      calo->GetMomentum(mom,vertex) ;
+      calo->GetMomentum(fMomentum,vertex) ;
     }
     
-    //--------------------------------------
+    //-----------------------------
     // Cluster selection
-    //--------------------------------------
+    //-----------------------------
     Int_t nMaxima = GetCaloUtils()->GetNumberOfLocalMaxima(calo, cells); // NLM
-    if(!ClusterSelected(calo,mom,nMaxima)) continue;
+    if(!ClusterSelected(calo,nMaxima)) continue;
     
     //----------------------------
-    //Create AOD for analysis
+    // Create AOD for analysis
     //----------------------------
-    AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
+    AliAODPWG4Particle aodph = AliAODPWG4Particle(fMomentum);
     
     //...............................................
     //Set the indeces of the original caloclusters (MC, ID), and calorimeter
     Int_t label = calo->GetLabel();
     aodph.SetLabel(label);
     aodph.SetCaloLabel(calo->GetID(),-1);
-    aodph.SetDetector(fCalorimeter);
+    aodph.SetDetectorTag(GetCalorimeter());
     //printf("Index %d, Id %d, iaod %d\n",icalo, calo->GetID(),GetOutputAODBranch()->GetEntriesFast());
     
     //...............................................
@@ -3725,31 +2344,40 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
     else                         aodph.SetDistToBad(0) ;
     //printf("DistBad %f Bit %d\n",distBad, aodph.DistToBad());
     
-    //--------------------------------------------------------------------------------------
+    //-------------------------------------
     // Play with the MC stack if available
-    //--------------------------------------------------------------------------------------
+    //-------------------------------------
     
     //Check origin of the candidates
     Int_t tag = -1;
     
     if(IsDataMC())
     {
-      tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader());
+      tag = GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(),GetCalorimeter());
       aodph.SetTag(tag);
       
-      if(GetDebug() > 0)
-        printf("AliAnaPhoton::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());
+      AliDebug(1,Form("Origin of candidate, bit map %d",aodph.GetTag()));
     }//Work with stack also
     
+    //--------------------------------------------------------
+    // Fill some shower shape histograms before PID is applied
+    //--------------------------------------------------------
     
-    //--------------------------------------------------------------------------------------
-    //Fill some shower shape histograms before PID is applied
-    //--------------------------------------------------------------------------------------
+    Float_t maxCellFraction = 0;
+    Int_t absIdMax = GetCaloUtils()->GetMaxEnergyCell(cells, calo, maxCellFraction);
+    if( absIdMax < 0 ) AliFatal("Wrong absID");
     
-    FillShowerShapeHistograms(calo,tag);
+    FillShowerShapeHistograms(calo,tag,maxCellFraction);
     
+    aodph.SetM02(calo->GetM02());
+    aodph.SetNLM(nMaxima);
+    aodph.SetTime(calo->GetTOF()*1e9);
+    aodph.SetNCells(calo->GetNCells());
+    Int_t nSM = GetModuleNumber(calo);
+    aodph.SetSModNumber(nSM);
+
     //-------------------------------------
-    //PID selection or bit setting
+    // PID selection or bit setting
     //-------------------------------------
     
     //...............................................
@@ -3761,7 +2389,7 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
       
       aodph.SetIdentifiedParticleType(GetCaloPID()->GetIdentifiedParticleType(calo));
       
-      if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PDG of identified particle %d\n",aodph.GetIdentifiedParticleType());
+      AliDebug(1,Form("PDG of identified particle %d",aodph.GetIdentifiedParticleType()));
       
       //If cluster does not pass pid, not photon, skip it.
       if(aodph.GetIdentifiedParticleType() != AliCaloPID::kPhoton) continue ;
@@ -3777,78 +2405,50 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
       
       GetCaloPID()->SetPIDBits(calo,&aodph, GetCaloUtils(),GetReader()->GetInputEvent());
       
-      if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - PID Bits set \n");
+      AliDebug(1,"PID Bits set");
     }
     
-    if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD() - Photon selection cuts passed: pT %3.2f, pdg %d\n",
-                              aodph.Pt(), aodph.GetIdentifiedParticleType());
+    AliDebug(1,Form("Photon selection cuts passed: pT %3.2f, pdg %d",aodph.Pt(),aodph.GetIdentifiedParticleType()));
+    
+    fhClusterCutsE [9]->Fill(calo->E());
+    fhClusterCutsPt[9]->Fill(fMomentum.Pt());
     
-    fhClusterCuts[9]->Fill(calo->E());
+    if(nSM < GetCaloUtils()->GetNumberOfSuperModulesUsed() && nSM >=0)
+    {
+      fhEPhotonSM ->Fill(fMomentum.E (),nSM);
+      fhPtPhotonSM->Fill(fMomentum.Pt(),nSM);
+    }
     
     fhNLocMax->Fill(calo->E(),nMaxima);
     
+    // Few more control histograms for selected clusters
+    fhMaxCellDiffClusterE->Fill(calo->E() ,maxCellFraction);
+    fhNCellsE            ->Fill(calo->E() ,calo->GetNCells());
+    fhTimePt             ->Fill(fMomentum.Pt()  ,calo->GetTOF()*1.e9);
+    
+    if(cells)
+    {
+      for(Int_t icell = 0; icell <  calo->GetNCells(); icell++)
+        fhCellsE->Fill(calo->E(),cells->GetCellAmplitude(calo->GetCellsAbsId()[icell]));
+    }
+    
     // Matching after cuts
-    if(fFillTMHisto) FillTrackMatchingResidualHistograms(calo,1);
+    if( fFillTMHisto )          FillTrackMatchingResidualHistograms(calo,1);
     
     // Fill histograms to undertand pile-up before other cuts applied
     // Remember to relax time cuts in the reader
-    FillPileUpHistograms(calo->E(),mom.Pt(),calo->GetTOF()*1e9);
-    
-    // Add number of local maxima to AOD, method name in AOD to be FIXED
-    aodph.SetFiducialArea(nMaxima);
-    
-    if(fFillEMCALBCHistograms && fCalorimeter=="EMCAL" && idTrig >= 0)
-    {
-      Double_t calotof   = calo->GetTOF()*1e9;
-      Float_t  calotofUS = TMath::Abs(calotof);
-      Float_t phicluster = aodph.Phi();
-      if(phicluster < 0) phicluster+=TMath::TwoPi();
-      
-      if(calo->E() > 2)
-      {
-        if      (calotofUS < 25) fhEtaPhiPhotonEMCALBC0->Fill(aodph.Eta(), phicluster);
-        else if (calotofUS < 75) fhEtaPhiPhotonEMCALBC1->Fill(aodph.Eta(), phicluster);
-        else                     fhEtaPhiPhotonEMCALBCN->Fill(aodph.Eta(), phicluster);
-      }
-      
-      Int_t bc = GetReader()->GetTriggerClusterBC();
-      Int_t histoBC = bc-5;
-      if(GetReader()->AreBadTriggerEventsRemoved()) histoBC=0; // histograms created only for one BC since the others where rejected
-      if(TMath::Abs(bc) < 6 && !GetReader()->IsBadCellTriggerEvent() && !GetReader()->IsExoticEvent())
-      {
-        if(GetReader()->IsTriggerMatched())
-        {
-          if(calo->E() > 2) fhEtaPhiPhotonTriggerEMCALBC[histoBC]->Fill(aodph.Eta(), phicluster);
-          fhTimePhotonTriggerEMCALBC[histoBC]->Fill(calo->E(), calotof);
-          if(GetReader()->IsPileUpFromSPD()) fhTimePhotonTriggerEMCALBCPileUpSPD[histoBC]->Fill(calo->E(), calotof);
-        }
-        else
-        {
-          if(calo->E() > 2) fhEtaPhiPhotonTriggerEMCALBCUM[histoBC]->Fill(aodph.Eta(), phicluster);
-          fhTimePhotonTriggerEMCALBCUM[histoBC]->Fill(calo->E(), calotof);
-          
-          if(bc==0)
-          {
-            if(GetReader()->IsTriggerMatchedOpenCuts(0)) fhTimePhotonTriggerEMCALBC0UMReMatchOpenTime   ->Fill(calo->E(), calotof);
-            if(GetReader()->IsTriggerMatchedOpenCuts(1)) fhTimePhotonTriggerEMCALBC0UMReMatchCheckNeigh ->Fill(calo->E(), calotof);
-            if(GetReader()->IsTriggerMatchedOpenCuts(2)) fhTimePhotonTriggerEMCALBC0UMReMatchBoth       ->Fill(calo->E(), calotof);
-          }
-        }
-      }
-      else if(TMath::Abs(bc) >= 6)
-        printf("AliAnaPhoton::MakeAnalysisFillAOD() - Trigger BC not expected = %d\n",bc);
-    }
+    if( IsPileUpAnalysisOn() ) FillPileUpHistograms(calo,cells, absIdMax);
     
-    //Add AOD with photon object to aod branch
+    // Add AOD with photon object to aod branch
     AddAODParticle(aodph);
     
   }//loop
   
-  if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillAOD()  End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());
+  AliDebug(1,Form("End fill AODs, with %d entries",GetOutputAODBranch()->GetEntriesFast()));
   
 }
 
-//__________________________________________________________________
+//______________________________________________
 void  AliAnaPhoton::MakeAnalysisFillHistograms()
 {
   //Fill histograms
@@ -3865,7 +2465,7 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
   //----------------------------------
   //Loop on stored AOD photons
   Int_t naod = GetOutputAODBranch()->GetEntriesFast();
-  if(GetDebug() > 0) printf("AliAnaPhoton::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
+  AliDebug(1,Form("AOD branch entries %d", naod));
   
   Float_t cen = GetEventCentrality();
   // printf("++++++++++ GetEventCentrality() %f\n",cen);
@@ -3877,16 +2477,14 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
     AliAODPWG4Particle* ph =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
     Int_t pdg = ph->GetIdentifiedParticleType();
     
-    if(GetDebug() > 3)
-      printf("AliAnaPhoton::MakeAnalysisFillHistograms() - PDG %d, MC TAG %d, Calorimeter %s\n",
-             ph->GetIdentifiedParticleType(),ph->GetTag(), (ph->GetDetector()).Data()) ;
+    AliDebug(2,Form("PDG %d, MC TAG %d, Calorimeter <%d>",ph->GetIdentifiedParticleType(),ph->GetTag(), ph->GetDetectorTag())) ;
     
-    //If PID used, fill histos with photons in Calorimeter fCalorimeter
+    //If PID used, fill histos with photons in Calorimeter GetCalorimeter()
     if(IsCaloPIDOn() && pdg != AliCaloPID::kPhoton) continue;
-    if(ph->GetDetector() != fCalorimeter) continue;
     
-    if(GetDebug() > 2)
-      printf("AliAnaPhoton::MakeAnalysisFillHistograms() - ID Photon: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;
+    if(((Int_t) ph->GetDetectorTag()) != GetCalorimeter()) continue;
+    
+    AliDebug(2,Form("ID Photon: pt %f, phi %f, eta %f", ph->Pt(),ph->Phi(),ph->Eta())) ;
     
     //................................
     //Fill photon histograms
@@ -3902,76 +2500,54 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
     if     (ecluster   > 0.5) fhEtaPhiPhoton  ->Fill(etacluster, phicluster);
     else if(GetMinPt() < 0.5) fhEtaPhi05Photon->Fill(etacluster, phicluster);
     
-    fhPtCentralityPhoton ->Fill(ptcluster,cen) ;
-    fhPtEventPlanePhoton ->Fill(ptcluster,ep ) ;
-    
-    //Get original cluster, to recover some information
-    Int_t absID             = 0;
-    Float_t maxCellFraction = 0;
-    AliVCaloCells* cells    = 0;
-    TObjArray * clusters    = 0;
-    if(fCalorimeter == "EMCAL")
-    {
-      cells    = GetEMCALCells();
-      clusters = GetEMCALClusters();
-    }
-    else
-    {
-      cells    = GetPHOSCells();
-      clusters = GetPHOSClusters();
-    }
-    
-    Int_t iclus = -1;
-    AliVCluster *cluster = FindCluster(clusters,ph->GetCaloLabel(0),iclus);
-    if(cluster)
+    if(IsHighMultiplicityAnalysisOn())
     {
-      absID = GetCaloUtils()->GetMaxEnergyCell(cells, cluster,maxCellFraction);
-      
-      // Control histograms
-      fhMaxCellDiffClusterE->Fill(ph->E() ,maxCellFraction);
-      fhNCellsE            ->Fill(ph->E() ,cluster->GetNCells());
-      fhTimePt             ->Fill(ph->Pt(),cluster->GetTOF()*1.e9);
-      if(cells)
-      {
-        for(Int_t icell = 0; icell <  cluster->GetNCells(); icell++)
-          fhCellsE->Fill(ph->E(),cells->GetCellAmplitude(cluster->GetCellsAbsId()[icell]));
-      }
+      fhPtCentralityPhoton ->Fill(ptcluster,cen) ;
+      fhPtEventPlanePhoton ->Fill(ptcluster,ep ) ;
     }
+  
+//    Comment this part, not needed but in case to know how to do it in the future
+//    //Get original cluster, to recover some information
+//    AliVCaloCells* cells    = 0;
+//    TObjArray * clusters    = 0;
+//    if(GetCalorimeter() == kEMCAL)
+//    {
+//      cells    = GetEMCALCells();
+//      clusters = GetEMCALClusters();
+//    }
+//    else
+//    {
+//      cells    = GetPHOSCells();
+//      clusters = GetPHOSClusters();
+//    }
+    
+//    Int_t iclus = -1;
+//    AliVCluster *cluster = FindCluster(clusters,ph->GetCaloLabel(0),iclus);
+//    if(cluster)
     
     //.......................................
     //Play with the MC data if available
     if(IsDataMC())
     {
-      if(GetDebug()>0)
-      {
-        if(GetReader()->ReadStack() && !GetMCStack())
-        {
-          printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called?\n");
-        }
-        else if(GetReader()->ReadAODMCParticles() && !GetReader()->GetAODMCParticles())
-        {
-          printf("AliAnaPhoton::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n");
-        }
-      }
-            
       //....................................................................
       // Access MC information in stack if requested, check that it exists.
-      Int_t label =ph->GetLabel();
+      Int_t label = ph->GetLabel();
       
       if(label < 0)
       {
-        if(GetDebug() > 1) printf("AliAnaPhoton::MakeAnalysisFillHistograms() *** bad label ***:  label %d \n", label);
+        AliDebug(1,Form("*** bad label ***:  label %d", label));
         continue;
       }
       
       Float_t eprim   = 0;
       Float_t ptprim  = 0;
       Bool_t ok = kFALSE;
-      TLorentzVector primary = GetMCAnalysisUtils()->GetMother(label,GetReader(),ok);
+      fPrimaryMom = GetMCAnalysisUtils()->GetMother(label,GetReader(),ok);
+                 
       if(ok)
       {
-        eprim   = primary.Energy();
-        ptprim  = primary.Pt();
+        eprim   = fPrimaryMom.Energy();
+        ptprim  = fPrimaryMom.Pt();
       }
       
       Int_t tag =ph->GetTag();
@@ -3989,7 +2565,6 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
         fhMCDeltaPt[kmcPhoton] ->Fill(ptcluster,ptprim-ptcluster);
         
         if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCConversion) &&
-           GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)     &&
            fhMCE[kmcConversion])
         {
           fhMCE  [kmcConversion] ->Fill(ecluster);
@@ -4011,7 +2586,7 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
         {
           mcParticleTag = kmcFragmentation;
         }
-        else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR)&& fhMCE[kmcISR])
+        else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCISR))
         {
           mcParticleTag = kmcISR;
         }
@@ -4020,20 +2595,19 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
         {
           mcParticleTag = kmcPi0Decay;
         }
-        else if((( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay) &&
-                  !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta)        ) ||
-                 GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCOtherDecay) ))
-        {
-          mcParticleTag = kmcOtherDecay;
-        }
-        else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))
+        else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))
         {
           mcParticleTag = kmcPi0;
         }
-        else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))
+        else if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta) &&
+                !GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEtaDecay))
         {
           mcParticleTag = kmcEta;
         }
+        else
+        {
+          mcParticleTag = kmcOtherDecay;
+        }
       }
       else if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiNeutron))
       {
@@ -4060,7 +2634,7 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
         
       }
       
-      if(mcParticleTag >= 0 && fhMCE  [mcParticleTag])
+      if(mcParticleTag >= 0 && fhMCE[mcParticleTag])
       {
         fhMCE  [mcParticleTag]->Fill(ecluster);
         fhMCPt [mcParticleTag]->Fill(ptcluster);
@@ -4090,7 +2664,7 @@ void AliAnaPhoton::Print(const Option_t * opt) const
   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
   AliAnaCaloTrackCorrBaseClass::Print(" ");
   
-  printf("Calorimeter            =     %s\n", fCalorimeter.Data()) ;
+  printf("Calorimeter            =     %s\n", GetCalorimeterString().Data()) ;
   printf("Min Distance to Bad Channel   = %2.1f\n",fMinDist);
   printf("Min Distance to Bad Channel 2 = %2.1f\n",fMinDist2);
   printf("Min Distance to Bad Channel 3 = %2.1f\n",fMinDist3);