]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG4/PartCorrDep/AliAnaCalorimeterQA.cxx
Setters added (Magnus)
[u/mrichter/AliRoot.git] / PWG4 / PartCorrDep / AliAnaCalorimeterQA.cxx
index 588470bd3088f5edef3c0da89ec7369665310a67..eff534682da1d529c60d2e1f62e5908aae808b45 100755 (executable)
@@ -42,7 +42,7 @@
 #include "AliStack.h"
 #include "AliVCaloCells.h"
 #include "AliFiducialCut.h"
-#include "AliVTrack.h"
+#include "AliAODTrack.h"
 #include "AliVCluster.h"
 #include "AliVEvent.h"
 #include "AliVEventHandler.h"
@@ -81,7 +81,8 @@ fh2E(0),fh2Pt(0),fh2Phi(0),fh2Eta(0),
 fhLambda(0), fhDispersion(0), 
 fhIM(0), fhIMCellCut(0),fhAsym(0), 
 fhNCellsPerCluster(0),fhNCellsPerClusterMIP(0), fhNCellsPerClusterMIPCharged(0), fhNClusters(0), 
-fhClusterTimeEnergy(0),fhCellTimeSpreadRespectToCellMax(0),fhCellIdCellLargeTimeSpread(0),
+fhClusterTimeEnergy(0),fhCellTimeSpreadRespectToCellMax(0),fhCellIdCellLargeTimeSpread(0), 
+fhBadClusterMaxCellTimeEnergy(0), fhBadClusterMaxCellCloseCellRatio(0),fhClusterMaxCellTimeEnergy(0), fhClusterMaxCellCloseCellRatio(0), 
 fhRNCells(0),fhXNCells(0),fhYNCells(0),fhZNCells(0),
 fhRE(0),     fhXE(0),     fhYE(0),     fhZE(0),    fhXYZ(0),
 fhRCellE(0), fhXCellE(0), fhYCellE(0), fhZCellE(0),fhXYZCell(0),
@@ -225,7 +226,30 @@ TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
   fhClusterTimeEnergy->SetXTitle("E (GeV) ");
   fhClusterTimeEnergy->SetYTitle("TOF (ns)");
   outputContainer->Add(fhClusterTimeEnergy);
-  
+    
+  fhClusterMaxCellCloseCellRatio  = new TH2F ("hClusterMaxCellCloseCell","energy vs ratio of max cell / neighbour cell, reconstructed clusters",
+                                          nptbins,ptmin,ptmax, 100,0,1.); 
+  fhClusterMaxCellCloseCellRatio->SetXTitle("E_{cluster} (GeV) ");
+  fhClusterMaxCellCloseCellRatio->SetYTitle("ratio");
+  outputContainer->Add(fhClusterMaxCellCloseCellRatio);
+  
+  fhBadClusterMaxCellCloseCellRatio  = new TH2F ("hBadClusterMaxCellCloseCell","energy vs ratio of max cell / neighbour cell constributing cell, reconstructed bad clusters",
+                                             nptbins,ptmin,ptmax, 100,0,1.); 
+  fhBadClusterMaxCellCloseCellRatio->SetXTitle("E_{cluster} (GeV) ");
+  fhBadClusterMaxCellCloseCellRatio->SetYTitle("ratio");
+  outputContainer->Add(fhBadClusterMaxCellCloseCellRatio);
+  
+  fhClusterMaxCellTimeEnergy  = new TH2F ("hClusterMaxCellTimeEnergy","energy vs TOF of maximum constributing cell, reconstructed clusters",
+                                          nptbins,ptmin,ptmax, ntimebins,timemin,timemax); 
+  fhClusterMaxCellTimeEnergy->SetXTitle("E_{cluster} (GeV) ");
+  fhClusterMaxCellTimeEnergy->SetYTitle("TOF (ns)");
+  outputContainer->Add(fhClusterMaxCellTimeEnergy);
+  
+  fhBadClusterMaxCellTimeEnergy  = new TH2F ("hBadClusterMaxCellTimeEnergy","energy vs TOF of maximum constributing cell, reconstructed clusters",
+                                             nptbins,ptmin,ptmax, ntimebins,timemin,timemax); 
+  fhBadClusterMaxCellTimeEnergy->SetXTitle("E_{cluster} (GeV) ");
+  fhBadClusterMaxCellTimeEnergy->SetYTitle("TOF (ns)");
+  outputContainer->Add(fhBadClusterMaxCellTimeEnergy);    
   
   //Shower shape
   fhLambda  = new TH3F ("hLambda","#lambda_{0}^{2} vs #lambda_{1}^{2} vs energy, reconstructed clusters",
@@ -491,7 +515,7 @@ TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
     fhCellIdCellLargeTimeSpread= new TH1F ("hCellIdCellLargeTimeSpread","", colmax*rowmax*fNModules,0,colmax*rowmax*fNModules); 
     fhCellIdCellLargeTimeSpread->SetXTitle("Absolute Cell Id");
     outputContainer->Add(fhCellIdCellLargeTimeSpread);
-    
+
     fhTime  = new TH1F ("hTime","Cell Time",ntimebins,timemin,timemax); 
     fhTime->SetXTitle("Cell Time (ns)");
     outputContainer->Add(fhTime);
@@ -583,7 +607,7 @@ TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
     fhCaloV0MCorrECells  = new TH2F ("hCaloV0MECells",Form("summed energy of Cells in %s vs V0 signal",fCalorimeter.Data()), nv0mbins,nv0mmin,nv0mmax,nptbins,ptmin,ptmax); 
     fhCaloV0MCorrECells->SetXTitle("V0 signal");
     fhCaloV0MCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
-    outputContainer->Add(fhCaloV0SCorrECells);    
+    outputContainer->Add(fhCaloV0MCorrECells);    
     
     //Calorimeter VS Track multiplicity
     fhCaloTrackMCorrNClusters  = new TH2F ("hCaloTrackMNClusters",Form("# clusters in %s vs V0 signal",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nbins,nmin,nmax); 
@@ -604,7 +628,7 @@ TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
     fhCaloTrackMCorrECells  = new TH2F ("hCaloTrackMECells",Form("summed energy of Cells in %s vs V0 signal",fCalorimeter.Data()), ntrmbins,ntrmmin,ntrmmax,nptbins,ptmin,ptmax); 
     fhCaloTrackMCorrECells->SetXTitle("Track Multiplicity");
     fhCaloTrackMCorrECells->SetYTitle(Form("#Sigma E of Cells in %s (GeV)",fCalorimeter.Data()));
-    outputContainer->Add(fhCaloV0SCorrECells);    
+    outputContainer->Add(fhCaloTrackMCorrECells);    
     
     
   }//correlate calorimeters
@@ -672,43 +696,44 @@ TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
     
     if(fCalorimeter == "EMCAL"){
       for(Int_t ifrac = 0; ifrac < 3; ifrac++){
-       fhAmplitudeModFraction[imod*3+ifrac]  = new TH1F (Form("hAmplitude_Mod%d_Frac%d",imod,ifrac),Form("Cell reconstructed Energy in Module %d, Fraction %d ",imod,ifrac), nptbins,ptmin,ptmax); 
-       fhAmplitudeModFraction[imod*3+ifrac]->SetXTitle("E (GeV)");
-       outputContainer->Add(fhAmplitudeModFraction[imod*3+ifrac]);
+        fhAmplitudeModFraction[imod*3+ifrac]  = new TH1F (Form("hAmplitude_Mod%d_Frac%d",imod,ifrac),Form("Cell reconstructed Energy in Module %d, Fraction %d ",imod,ifrac), nptbins,ptmin,ptmax); 
+        fhAmplitudeModFraction[imod*3+ifrac]->SetXTitle("E (GeV)");
+        outputContainer->Add(fhAmplitudeModFraction[imod*3+ifrac]);
       }
       
     }
-    
-    for(Int_t ircu = 0; ircu < fNRCU; ircu++){
-      fhTimeAmpPerRCU[imod*fNRCU+ircu]  = new TH2F (Form("hTimeAmp_Mod%d_RCU%d",imod,ircu),
-                                                    Form("Cell Energy vs Cell Time in Module %d, RCU %d ",imod,ircu), 
-                                                    nptbins,ptmin,ptmax,ntimebins,timemin,timemax); 
-      fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
-      fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("time (ns)");
-      outputContainer->Add(fhTimeAmpPerRCU[imod*fNRCU+ircu]);
+    if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
       
-      //                               fhT0TimeAmpPerRCU[imod*fNRCU+ircu]  = new TH2F (Form("hT0TimeAmp_Mod%d_RCU%d",imod,ircu),
-      //                                                                                                                         Form("Cell Energy vs T0-Cell Time in Module %d, RCU %d ",imod,ircu), 
-      //                                                                                                                         nptbins,ptmin,ptmax,ntimebins,timemin,timemax); 
-      //                               fhT0TimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
-      //                               fhT0TimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("T_{0} - T_{EMCal} (ns)");
-      //                               outputContainer->Add(fhT0TimeAmpPerRCU[imod*fNRCU+ircu]);
-      //                       
-                       
-      //                               for(Int_t imod2 = 0; imod2 < fNModules; imod2++){
-      //                                               for(Int_t ircu2 = 0; ircu2 < fNModules; ircu2++){
-      //                                                       Int_t index =  (imod2*fNRCU+ircu2)+(fNModules*fNRCU)*(ircu+imod)+fNRCU*fNModules*imod; 
-      //                                                       fhTimeCorrRCU[index]  = new TH2F (Form("hTimeCorrRCU_Mod%d_RCU%d_CompareTo_Mod%d_RCU%d",imod, ircu,imod2, ircu2),
-      //                                                                                                                                                       Form("Cell Energy > 0.3, Correlate cell Time in Module %d, RCU %d to Module %d, RCU %d",imod,ircu,imod2, ircu2),
-      //                                                                                                                                                       ntimebins,timemin,timemax,ntimebins,timemin,timemax); 
-      //                                                       fhTimeCorrRCU[index]->SetXTitle("Trigger Cell Time (ns)");
-      //                                                       fhTimeCorrRCU[index]->SetYTitle("Cell Time (ns)");
-      //                                                       outputContainer->Add(fhTimeCorrRCU[index]);
-      //                                               }
-      //                               }
+      for(Int_t ircu = 0; ircu < fNRCU; ircu++){
+        fhTimeAmpPerRCU[imod*fNRCU+ircu]  = new TH2F (Form("hTimeAmp_Mod%d_RCU%d",imod,ircu),
+                                                      Form("Cell Energy vs Cell Time in Module %d, RCU %d ",imod,ircu), 
+                                                      nptbins,ptmin,ptmax,ntimebins,timemin,timemax); 
+        fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
+        fhTimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("time (ns)");
+        outputContainer->Add(fhTimeAmpPerRCU[imod*fNRCU+ircu]);
+        
+        //                             fhT0TimeAmpPerRCU[imod*fNRCU+ircu]  = new TH2F (Form("hT0TimeAmp_Mod%d_RCU%d",imod,ircu),
+        //                                                                                                                       Form("Cell Energy vs T0-Cell Time in Module %d, RCU %d ",imod,ircu), 
+        //                                                                                                                       nptbins,ptmin,ptmax,ntimebins,timemin,timemax); 
+        //                             fhT0TimeAmpPerRCU[imod*fNRCU+ircu]->SetXTitle("E (GeV)");
+        //                             fhT0TimeAmpPerRCU[imod*fNRCU+ircu]->SetYTitle("T_{0} - T_{EMCal} (ns)");
+        //                             outputContainer->Add(fhT0TimeAmpPerRCU[imod*fNRCU+ircu]);
+        //                     
+        
+        //                             for(Int_t imod2 = 0; imod2 < fNModules; imod2++){
+        //                                             for(Int_t ircu2 = 0; ircu2 < fNModules; ircu2++){
+        //                                                     Int_t index =  (imod2*fNRCU+ircu2)+(fNModules*fNRCU)*(ircu+imod)+fNRCU*fNModules*imod; 
+        //                                                     fhTimeCorrRCU[index]  = new TH2F (Form("hTimeCorrRCU_Mod%d_RCU%d_CompareTo_Mod%d_RCU%d",imod, ircu,imod2, ircu2),
+        //                                                                                                                                                     Form("Cell Energy > 0.3, Correlate cell Time in Module %d, RCU %d to Module %d, RCU %d",imod,ircu,imod2, ircu2),
+        //                                                                                                                                                     ntimebins,timemin,timemax,ntimebins,timemin,timemax); 
+        //                                                     fhTimeCorrRCU[index]->SetXTitle("Trigger Cell Time (ns)");
+        //                                                     fhTimeCorrRCU[index]->SetYTitle("Cell Time (ns)");
+        //                                                     outputContainer->Add(fhTimeCorrRCU[index]);
+        //                                             }
+        //                             }
+      }
     }
     
-    
     fhIMMod[imod]  = new TH2F (Form("hIM_Mod%d",imod),
                                Form("Cluster pairs Invariant mass vs reconstructed pair energy in Module %d",imod),
                                nptbins,ptmin,ptmax,nmassbins,massmin,massmax); 
@@ -1187,6 +1212,9 @@ TList *  AliAnaCalorimeterQA::GetCreateOutputObjects()
     outputContainer->Add(fhMCNeutral1pOverER02);
   }
   
+//  for(Int_t i = 0; i < outputContainer->GetEntries() ; i++)
+//    printf("i=%d, name= %s\n",i,outputContainer->At(i)->GetName());
+  
   return outputContainer;
 }
 
@@ -1205,15 +1233,12 @@ Int_t AliAnaCalorimeterQA::GetNewRebinForRePlotting(TH1D* histo, const Float_t n
 void AliAnaCalorimeterQA::Init()
 { 
   //Check if the data or settings are ok
-  if(fCalorimeter != "PHOS" && fCalorimeter !="EMCAL"){
-    printf("AliAnaCalorimeterQA::Init() - Wrong calorimeter name <%s>, END\n", fCalorimeter.Data());
-    abort();
-  }    
   
-  if(GetReader()->GetDataType()== AliCaloTrackReader::kMC){
-    printf("AliAnaCalorimeterQA::Init() - Analysis of reconstructed data, MC reader not aplicable\n");
-    abort();
-  }    
+  if(fCalorimeter != "PHOS" && fCalorimeter !="EMCAL")
+    AliFatal(Form("Wrong calorimeter name <%s>", fCalorimeter.Data()));
+
+  if(GetReader()->GetDataType()== AliCaloTrackReader::kMC)
+    AliFatal("Analysis of reconstructed data, MC reader not aplicable");
   
 }
 
@@ -1303,10 +1328,9 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
   if(IsDataMC()){
     if(GetReader()->ReadStack()){
       
-      if(!GetMCStack()) {
-        printf("AliAnaPhoton::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
-        abort();
-      }
+      if(!GetMCStack()) 
+        AliFatal("Stack not available, is the MC handler called?\n");
+        
       //Fill some pure MC histograms, only primaries.
       for(Int_t i=0 ; i<GetMCStack()->GetNprimary(); i++){//Only primary particles, for all MC transport put GetNtrack()
         TParticle *primary = GetMCStack()->Particle(i) ;
@@ -1318,10 +1342,9 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
     }
     else if(GetReader()->ReadAODMCParticles()){
       
-      if(!GetReader()->GetAODMCParticles(0))   {
-        printf("AliAnaPhoton::MakeAnalysisFillHistograms() -  AODMCParticles not available!\n");
-        abort();
-      }
+      if(!GetReader()->GetAODMCParticles(0))   
+        AliFatal("AODMCParticles not available!");
+      
       //Fill some pure MC histograms, only primaries.
       for(Int_t i=0 ; i < (GetReader()->GetAODMCParticles(0))->GetEntriesFast(); i++){
         AliAODMCParticle *aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(i) ;
@@ -1358,7 +1381,6 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
     //----------------------------------------------------------
     if(fCorrelate)     Correlate();
     
-    
     //----------------------------------------------------------
     // CALOCLUSTERS
     //----------------------------------------------------------
@@ -1381,38 +1403,38 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
       if(GetDebug() > 0) printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - cluster: %d/%d, data %d \n",
                                 iclus+1,nCaloClusters,GetReader()->GetDataType());
       
-      if(GetReader()->GetDataType()==AliCaloTrackReader::kESD){
-        AliVCluster* clus =  (AliVCluster*)caloClusters->At(iclus);
-        AliVCaloCells * cell = 0x0; 
-        if(fCalorimeter == "PHOS") cell =  GetPHOSCells();
-        else                                    cell =  GetEMCALCells();
-        
-        //Get cluster kinematics
-        clus->GetPosition(pos);
-        clus->GetMomentum(mom,v);
-        tof = clus->GetTOF()*1e9;
-        if(tof < fTimeCutMin || tof > fTimeCutMax) continue;
-        
-        //Check only certain regions
-        Bool_t in = kTRUE;
-        if(IsFiducialCutOn()) in =  GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
-        if(!in) continue;
-        
-        //Get module of cluster
-        nCaloClustersAccepted++;
-        nModule = GetModuleNumber(clus);
-        if(nModule >=0 && nModule < fNModules) nClustersInModule[nModule]++;
-        
-        //MC labels
-        nLabel = clus->GetNLabels();
-        labels = clus->GetLabels();
-        
-        //Cells per cluster
-        nCaloCellsPerCluster =  clus->GetNCells();
-        //if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster);
-        
-        //matched cluster with tracks
-        nTracksMatched = clus->GetNTracksMatched();
+      AliVCluster* clus =  (AliVCluster*)caloClusters->At(iclus);
+      AliVCaloCells * cell = 0x0; 
+      if(fCalorimeter == "PHOS") cell =  GetPHOSCells();
+      else                                      cell =  GetEMCALCells();
+      
+      //Get cluster kinematics
+      clus->GetPosition(pos);
+      clus->GetMomentum(mom,v);
+      tof = clus->GetTOF()*1e9;
+      if(tof < fTimeCutMin || tof > fTimeCutMax) continue;
+      
+      //Check only certain regions
+      Bool_t in = kTRUE;
+      if(IsFiducialCutOn()) in =  GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
+      if(!in) continue;
+      
+      //Get module of cluster
+      nCaloClustersAccepted++;
+      nModule = GetModuleNumber(clus);
+      if(nModule >=0 && nModule < fNModules) nClustersInModule[nModule]++;
+      
+      //MC labels
+      nLabel = clus->GetNLabels();
+      labels = clus->GetLabels();
+      
+      //Cells per cluster
+      nCaloCellsPerCluster =  clus->GetNCells();
+      //if(mom.E() > 10 && nCaloCellsPerCluster == 1 ) printf("%s:************** E = %f ********** ncells = %d\n",fCalorimeter.Data(), mom.E(),nCaloCellsPerCluster);
+      
+      //matched cluster with tracks
+      nTracksMatched = clus->GetNTracksMatched();
+      if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
         trackIndex     = clus->GetTrackMatchedIndex();
         if(trackIndex >= 0){
           track = (AliVTrack*)GetReader()->GetInputEvent()->GetTrack(trackIndex);
@@ -1421,172 +1443,207 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
           if(nTracksMatched == 1) nTracksMatched = 0;
           track = 0;
         }
+      }//kESD
+      else{//AODs
+        if(nTracksMatched > 0) track = (AliVTrack*)clus->GetTrackMatched(0);
+      }
+      
+      //Shower shape parameters
+      showerShape[0] = clus->GetM20();
+      showerShape[1] = clus->GetM02();
+      showerShape[2] = clus->GetDispersion();
+      
+      //======================
+      //Cells in cluster
+      //======================
+      
+      //Get list of contributors
+      UShort_t * indexList = clus->GetCellsAbsId() ;
+      // check time of cells respect to max energy cell
+      //Get maximum energy cell
+      Float_t emax  = -1;
+      Double_t tmax = -1;
+      Int_t imax    = -1;
+      Int_t absId   = -1 ;
+      //printf("nCaloCellsPerCluster %d\n",nCaloCellsPerCluster);
+      //Loop on cluster cells
+      for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
+        //     printf("Index %d\n",ipos);
+        absId  = indexList[ipos]; 
         
-        //Shower shape parameters
-        showerShape[0] = clus->GetM20();
-        showerShape[1] = clus->GetM02();
-        showerShape[2] = clus->GetDispersion();
+        //Get position of cell compare to cluster
+        if(fFillAllPosHisto){
+          if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
+            
+            Double_t cellpos[] = {0, 0, 0};
+            GetEMCALGeometry()->GetGlobal(absId, cellpos);
+            
+            fhDeltaCellClusterXNCells->Fill(pos[0]-cellpos[0],nCaloCellsPerCluster) ; 
+            fhDeltaCellClusterYNCells->Fill(pos[1]-cellpos[1],nCaloCellsPerCluster) ; 
+            fhDeltaCellClusterZNCells->Fill(pos[2]-cellpos[2],nCaloCellsPerCluster) ;
+            
+            fhDeltaCellClusterXE->Fill(pos[0]-cellpos[0],mom.E())  ; 
+            fhDeltaCellClusterYE->Fill(pos[1]-cellpos[1],mom.E())  ; 
+            fhDeltaCellClusterZE->Fill(pos[2]-cellpos[2],mom.E())  ; 
+            
+            Float_t r     = TMath::Sqrt(pos[0]*pos[0]        +pos[1]*pos[1]);//     +pos[2]*pos[2]);
+            Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
+            fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ; 
+            fhDeltaCellClusterRE     ->Fill(r-rcell, mom.E())  ; 
+            
+            //                                 Float_t celleta = 0, cellphi = 0;
+            //                                 GetEMCALGeometry()->EtaPhiFromIndex(absId, celleta, cellphi); 
+            //                                 Int_t imod = -1, iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1;
+            //                                 GetEMCALGeometry()->GetCellIndex(absId,imod,iTower,iIphi,iIeta); 
+            //                                 GetEMCALGeometry()->GetCellPhiEtaIndexInSModule(imod,iTower,
+            //                                                                                                                                                          iIphi, iIeta,iphi,ieta);
+            //                                 printf("AbsId %d, SM %d, Index eta %d, phi %d\n", absId, imod, ieta, iphi);
+            //                                 printf("Cluster E %f, eta %f, phi %f; Cell: Amp %f, eta %f, phi%f\n", mom.E(),mom.Eta(), mom.Phi()*TMath::RadToDeg(), cell->GetCellAmplitude(absId),celleta, cellphi*TMath::RadToDeg());
+            //                                 printf("x cluster %f, x cell %f, cluster-cell %f\n",pos[0], cellpos[0],pos[0]-cellpos[0]);
+            //                                 printf("y cluster %f, y cell %f, cluster-cell %f\n",pos[1], cellpos[1],pos[1]-cellpos[1]);
+            //                                 printf("z cluster %f, z cell %f, cluster-cell %f\n",pos[2], cellpos[2],pos[2]-cellpos[2]);
+            //                                 printf("r cluster %f, r cell %f, cluster-cell %f\n",r,      rcell,     r-rcell);
+            //                                 
+            
+          }//EMCAL and its matrices are available
+          else if(fCalorimeter=="PHOS" && GetCaloUtils()->IsPHOSGeoMatrixSet()){
+            TVector3 xyz;
+            Int_t relId[4], module;
+            Float_t xCell, zCell;
+            
+            GetPHOSGeometry()->AbsToRelNumbering(absId,relId);
+            module = relId[0];
+            GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
+            GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
+            
+            fhDeltaCellClusterXNCells->Fill(pos[0]-xyz.X(),nCaloCellsPerCluster) ; 
+            fhDeltaCellClusterYNCells->Fill(pos[1]-xyz.Y(),nCaloCellsPerCluster) ; 
+            fhDeltaCellClusterZNCells->Fill(pos[2]-xyz.Z(),nCaloCellsPerCluster) ;
+            
+            fhDeltaCellClusterXE->Fill(pos[0]-xyz.X(),mom.E())  ; 
+            fhDeltaCellClusterYE->Fill(pos[1]-xyz.Y(),mom.E())  ; 
+            fhDeltaCellClusterZE->Fill(pos[2]-xyz.Z(),mom.E())  ; 
+            
+            Float_t r     = TMath::Sqrt(pos[0]*pos[0]  +pos[1]*pos[1]);//     +pos[2]*pos[2]);
+            Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());//+xyz.Z()*xyz.Z());
+            fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ; 
+            fhDeltaCellClusterRE     ->Fill(r-rcell, mom.E())  ; 
+            
+            //                   printf("x cluster %f, x cell %f, cluster-cell %f\n",pos[0], cellpos[0],pos[0]-cellpos[0]);
+            //                 printf("y cluster %f, y cell %f, cluster-cell %f\n",pos[1], cellpos[1],pos[1]-cellpos[1]);
+            //                 printf("z cluster %f, z cell %f, cluster-cell %f\n",pos[2], cellpos[2],pos[2]-cellpos[2]);
+            //                                 printf("r cluster %f, r cell %f, cluster-cell %f\n",r,      rcell,     r-rcell);
+          }//PHOS and its matrices are available
+        }//Fill all position histograms
         
-        //======================
-        //Cells in cluster
-        //======================
+        //Find maximum energy cluster
+        if(cell->GetCellAmplitude(absId) > emax) {
+          imax = ipos;
+          emax = cell->GetCellAmplitude(absId);
+          tmax = cell->GetCellTime(absId);
+        } 
         
-        //Get list of contributors
-        UShort_t * indexList = clus->GetCellsAbsId() ;
-        // check time of cells respect to max energy cell
-        //Get maximum energy cell
-        Float_t emax  = -1;
-        Double_t tmax = -1;
-        Int_t imax    = -1;
-        Int_t absId   = -1 ;
-        //printf("nCaloCellsPerCluster %d\n",nCaloCellsPerCluster);
-        //Loop on cluster cells
+      }// cluster cell loop
+      //Bad clusters histograms
+      Float_t minNCells = 1+mom.E()/3;//-x*x*0.0033
+      if(nCaloCellsPerCluster < minNCells) {
+        if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) 
+          fhBadClusterMaxCellTimeEnergy->Fill(mom.E(),tmax);
+        else 
+          fhBadClusterMaxCellTimeEnergy->Fill(mom.E(),tof);
+        //printf("bad tof : %2.3f\n",tof);
+
         for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
           //   printf("Index %d\n",ipos);
-          absId  = indexList[ipos]; 
-          
-          //Get position of cell compare to cluster
-          if(fFillAllPosHisto){
-            if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
-              
-              Double_t cellpos[] = {0, 0, 0};
-              GetEMCALGeometry()->GetGlobal(absId, cellpos);
-              
-              fhDeltaCellClusterXNCells->Fill(pos[0]-cellpos[0],nCaloCellsPerCluster) ; 
-              fhDeltaCellClusterYNCells->Fill(pos[1]-cellpos[1],nCaloCellsPerCluster) ; 
-              fhDeltaCellClusterZNCells->Fill(pos[2]-cellpos[2],nCaloCellsPerCluster) ;
-              
-              fhDeltaCellClusterXE->Fill(pos[0]-cellpos[0],mom.E())  ; 
-              fhDeltaCellClusterYE->Fill(pos[1]-cellpos[1],mom.E())  ; 
-              fhDeltaCellClusterZE->Fill(pos[2]-cellpos[2],mom.E())  ; 
-              
-              Float_t r     = TMath::Sqrt(pos[0]*pos[0]        +pos[1]*pos[1]);//     +pos[2]*pos[2]);
-              Float_t rcell = TMath::Sqrt(cellpos[0]*cellpos[0]+cellpos[1]*cellpos[1]);//+cellpos[2]*cellpos[2]);
-              fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ; 
-              fhDeltaCellClusterRE     ->Fill(r-rcell, mom.E())  ; 
-              
-              //                                       Float_t celleta = 0, cellphi = 0;
-              //                                       GetEMCALGeometry()->EtaPhiFromIndex(absId, celleta, cellphi); 
-              //                                       Int_t imod = -1, iTower = -1, iIphi = -1, iIeta = -1, iphi = -1, ieta = -1;
-              //                                       GetEMCALGeometry()->GetCellIndex(absId,imod,iTower,iIphi,iIeta); 
-              //                                       GetEMCALGeometry()->GetCellPhiEtaIndexInSModule(imod,iTower,
-              //                                                                                                                                                                iIphi, iIeta,iphi,ieta);
-              //                                       printf("AbsId %d, SM %d, Index eta %d, phi %d\n", absId, imod, ieta, iphi);
-              //                                       printf("Cluster E %f, eta %f, phi %f; Cell: Amp %f, eta %f, phi%f\n", mom.E(),mom.Eta(), mom.Phi()*TMath::RadToDeg(), cell->GetCellAmplitude(absId),celleta, cellphi*TMath::RadToDeg());
-              //                                       printf("x cluster %f, x cell %f, cluster-cell %f\n",pos[0], cellpos[0],pos[0]-cellpos[0]);
-              //                                       printf("y cluster %f, y cell %f, cluster-cell %f\n",pos[1], cellpos[1],pos[1]-cellpos[1]);
-              //                                       printf("z cluster %f, z cell %f, cluster-cell %f\n",pos[2], cellpos[2],pos[2]-cellpos[2]);
-              //                                       printf("r cluster %f, r cell %f, cluster-cell %f\n",r,      rcell,     r-rcell);
-              //                                       
-              
-            }//EMCAL and its matrices are available
-            else if(fCalorimeter=="PHOS" && GetCaloUtils()->IsPHOSGeoMatrixSet()){
-              TVector3 xyz;
-              Int_t relId[4], module;
-              Float_t xCell, zCell;
-              
-              GetPHOSGeometry()->AbsToRelNumbering(absId,relId);
-              module = relId[0];
-              GetPHOSGeometry()->RelPosInModule(relId,xCell,zCell);
-              GetPHOSGeometry()->Local2Global(module,xCell,zCell,xyz);
-              
-              fhDeltaCellClusterXNCells->Fill(pos[0]-xyz.X(),nCaloCellsPerCluster) ; 
-              fhDeltaCellClusterYNCells->Fill(pos[1]-xyz.Y(),nCaloCellsPerCluster) ; 
-              fhDeltaCellClusterZNCells->Fill(pos[2]-xyz.Z(),nCaloCellsPerCluster) ;
-              
-              fhDeltaCellClusterXE->Fill(pos[0]-xyz.X(),mom.E())  ; 
-              fhDeltaCellClusterYE->Fill(pos[1]-xyz.Y(),mom.E())  ; 
-              fhDeltaCellClusterZE->Fill(pos[2]-xyz.Z(),mom.E())  ; 
-              
-              Float_t r     = TMath::Sqrt(pos[0]*pos[0]  +pos[1]*pos[1]);//     +pos[2]*pos[2]);
-              Float_t rcell = TMath::Sqrt(xyz.X()*xyz.X()+xyz.Y()*xyz.Y());//+xyz.Z()*xyz.Z());
-              fhDeltaCellClusterRNCells->Fill(r-rcell, nCaloCellsPerCluster) ; 
-              fhDeltaCellClusterRE     ->Fill(r-rcell, mom.E())  ; 
-              
-              //                         printf("x cluster %f, x cell %f, cluster-cell %f\n",pos[0], cellpos[0],pos[0]-cellpos[0]);
-              //                       printf("y cluster %f, y cell %f, cluster-cell %f\n",pos[1], cellpos[1],pos[1]-cellpos[1]);
-              //                       printf("z cluster %f, z cell %f, cluster-cell %f\n",pos[2], cellpos[2],pos[2]-cellpos[2]);
-              //                               printf("r cluster %f, r cell %f, cluster-cell %f\n",r,      rcell,     r-rcell);
-            }//PHOS and its matrices are available
-          }//Fill all position histograms
-          
-          //Find maximum energy cluster
-          if(cell->GetCellAmplitude(absId) > emax) {
-            imax = ipos;
-            emax = cell->GetCellAmplitude(absId);
-            tmax = cell->GetCellTime(absId);
-          } 
-          
-        }// cluster cell loop
-        
-        // check time of cells respect to max energy cell
-        if(nCaloCellsPerCluster > 1){
-          for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
-            if(imax == ipos) continue;
+          if(ipos!=imax){
+            absId  = indexList[ipos]; 
+            Float_t frac = cell->GetCellAmplitude(absId)/emax;
+            //printf("bad frac : %2.3f, e %2.2f, ncells %d, min %2.1f\n",frac,mom.E(),nCaloCellsPerCluster,minNCells);
+            fhBadClusterMaxCellCloseCellRatio->Fill(mom.E(),frac);
+          }
+        }
+      }//Bad cluster
+      else{
+        if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) 
+          fhClusterMaxCellTimeEnergy->Fill(mom.E(),tmax);
+        else 
+          fhClusterMaxCellTimeEnergy->Fill(mom.E(),tof);
+        for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
+          //   printf("Index %d\n",ipos);
+          if(ipos!=imax){
             absId  = indexList[ipos]; 
-            Float_t diff = (tmax-cell->GetCellTime(absId))*1e9;
-            fhCellTimeSpreadRespectToCellMax->Fill(diff);
-            if(TMath::Abs(TMath::Abs(diff) > 100)) fhCellIdCellLargeTimeSpread->Fill(absId);
-          }// fill cell-cluster histogram loop
+            Float_t frac = cell->GetCellAmplitude(absId)/emax;
+            //printf("good frac : %2.3f\n",frac);
+            fhClusterMaxCellCloseCellRatio->Fill(mom.E(),frac);
+          }
+        }
+      }//good cluster
           
-        }//check time of cells respect to max energy cell
-        
-        //-----------------------------------------------------------
-        //Fill histograms related to single cluster or track matching
-        //-----------------------------------------------------------
-        
-        ClusterHistograms(mom, tof, pos, showerShape, nCaloCellsPerCluster, nModule, nTracksMatched, track, labels, nLabel);   
-        
-        
-        //-----------------------------------------------------------
-        //Invariant mass
-        //-----------------------------------------------------------
-        if(GetDebug()>1) printf("Invariant mass \n");
+      // check time of cells respect to max energy cell
+      if(nCaloCellsPerCluster > 1 &&  GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
+        for (Int_t ipos = 0; ipos < nCaloCellsPerCluster; ipos++) {
+          if(imax == ipos) continue;
+          absId  = indexList[ipos]; 
+          Float_t diff = (tmax-cell->GetCellTime(absId))*1e9;
+          fhCellTimeSpreadRespectToCellMax->Fill(diff);
+          if(TMath::Abs(TMath::Abs(diff) > 100)) fhCellIdCellLargeTimeSpread->Fill(absId);
+        }// fill cell-cluster histogram loop
+      }//check time of cells respect to max energy cell
+      
+      //-----------------------------------------------------------
+      //Fill histograms related to single cluster or track matching
+      //-----------------------------------------------------------
+      ClusterHistograms(mom, tof, pos, showerShape, nCaloCellsPerCluster, nModule, nTracksMatched, track, labels, nLabel);     
+      
+      
+      //-----------------------------------------------------------
+      //Invariant mass
+      //-----------------------------------------------------------
+      if(GetDebug()>1) printf("Invariant mass \n");
+      
+      //do not do for bad vertex
+      // Float_t fZvtxCut = 40. ;      
+      if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
+      
+      Int_t nModule2 = -1;
+      Int_t nCaloCellsPerCluster2=0;
+      if (nCaloClusters > 1 ) {
+        for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++) {
+          AliVCluster* clus2 =  (AliVCluster*)caloClusters->At(jclus);
+          
+          //Get cluster kinematics
+          clus2->GetMomentum(mom2,v);
+          //Check only certain regions
+          Bool_t in2 = kTRUE;
+          if(IsFiducialCutOn()) in2 =  GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;
+          if(!in2) continue;   
+          //Get module of cluster
+          nModule2 = GetModuleNumber(clus2);
+          //Cells per cluster
+          nCaloCellsPerCluster2 = clus2->GetNCells();
+        }
+        //Fill invariant mass histograms
+        //All modules
         
-        //do not do for bad vertex
-        // Float_t fZvtxCut = 40. ;    
-        if(v[2]<-GetZvertexCut() || v[2]> GetZvertexCut()) continue ; //Event can not be used (vertex, centrality,... cuts not fulfilled)
+        //printf("QA : Fill inv mass histo: pt1 %f, pt2 %f, pt12 %f, mass %f, calo %s \n",mom.Pt(),mom2.Pt(),(mom+mom2).Pt(),(mom+mom2).M(), fCalorimeter.Data());
+        fhIM  ->Fill((mom+mom2).Pt(),(mom+mom2).M());
+        //Single module
+        if(nModule == nModule2 && nModule >=0 && nModule < fNModules)
+          fhIMMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
         
-        Int_t nModule2 = -1;
-        Int_t nCaloCellsPerCluster2=0;
-        if (nCaloClusters > 1 ) {
-          for(Int_t jclus = iclus + 1 ; jclus < nCaloClusters ; jclus++) {
-            AliVCluster* clus2 =  (AliVCluster*)caloClusters->At(jclus);
-            
-            //Get cluster kinematics
-            clus2->GetMomentum(mom2,v);
-            //Check only certain regions
-            Bool_t in2 = kTRUE;
-            if(IsFiducialCutOn()) in2 =  GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;
-            if(!in2) continue; 
-            //Get module of cluster
-            nModule2 = GetModuleNumber(clus2);
-            //Cells per cluster
-            nCaloCellsPerCluster2 = clus2->GetNCells();
-          }
-          //Fill invariant mass histograms
+        //Select only clusters with at least 2 cells
+        if(nCaloCellsPerCluster > 1 && nCaloCellsPerCluster2 > 1) {
           //All modules
-          
-          //printf("QA : Fill inv mass histo: pt1 %f, pt2 %f, pt12 %f, mass %f, calo %s \n",mom.Pt(),mom2.Pt(),(mom+mom2).Pt(),(mom+mom2).M(), fCalorimeter.Data());
-          fhIM  ->Fill((mom+mom2).Pt(),(mom+mom2).M());
-          //Single module
+          fhIMCellCut  ->Fill((mom+mom2).Pt(),(mom+mom2).M());
+          //Single modules
           if(nModule == nModule2 && nModule >=0 && nModule < fNModules)
-            fhIMMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
-          
-          //Select only clusters with at least 2 cells
-          if(nCaloCellsPerCluster > 1 && nCaloCellsPerCluster2 > 1) {
-            //All modules
-            fhIMCellCut  ->Fill((mom+mom2).Pt(),(mom+mom2).M());
-            //Single modules
-            if(nModule == nModule2 && nModule >=0 && nModule < fNModules)
-              fhIMCellCutMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
-          }
-          
-          //Asymetry histograms
-          fhAsym->Fill((mom+mom2).Pt(),TMath::Abs((mom.E()-mom2.E())/(mom.E()+mom2.E())));
-          
-        }// 2nd cluster loop
-      }////more than 1 cluster in calorimeter          
+            fhIMCellCutMod[nModule]->Fill((mom+mom2).Pt(),(mom+mom2).M());
+        }
+        
+        //Asymetry histograms
+        fhAsym->Fill((mom+mom2).Pt(),TMath::Abs((mom.E()-mom2.E())/(mom.E()+mom2.E())));
+        
+      }// 2nd cluster loop
     }//cluster loop
     
     //Number of clusters histograms
@@ -1622,18 +1679,18 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
   else                       
     cell = GetEMCALCells();
   
-  if(!cell) {
-    printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - STOP: No %s ESD CELLS available for analysis\n",fCalorimeter.Data());
-    abort();
-  }
+  if(!cell) 
+    AliFatal(Form("No %s CELLS available for analysis",fCalorimeter.Data()));
   
   if(GetDebug() > 0) 
-    printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - In ESD %s cell entries %d\n", fCalorimeter.Data(), cell->GetNumberOfCells());    
+    printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - %s cell entries %d\n", fCalorimeter.Data(), cell->GetNumberOfCells());    
   
   for (Int_t iCell = 0; iCell < cell->GetNumberOfCells(); iCell++) {      
-    if(GetDebug() > 2)  printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Cell : amp %f, absId %d \n", cell->GetAmplitude(iCell), cell->GetCellNumber(iCell));
+    if(GetDebug() > 2)  
+      printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - Cell : amp %f, absId %d \n", cell->GetAmplitude(iCell), cell->GetCellNumber(iCell));
     nModule = GetModuleNumberCellIndexes(cell->GetCellNumber(iCell),fCalorimeter, icol, irow, iRCU);
-    if(GetDebug() > 2) printf("\t module %d, column %d, row %d \n", nModule,icol,irow);
+    if(GetDebug() > 2) 
+      printf("\t module %d, column %d, row %d \n", nModule,icol,irow);
     
     if(nModule < fNModules) {  
       //Check if the cell is a bad channel
@@ -1648,7 +1705,11 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
           }
         }
       }
-      
+      else {
+        delete [] nCellsInModule;
+        return;
+      }
+
       //Get Recalibration factor if set
       if (GetCaloUtils()->IsRecalibrationOn()) {
         if(fCalorimeter == "PHOS") recalF = GetCaloUtils()->GetPHOSChannelRecalibrationFactor(nModule,icol,irow);
@@ -1658,26 +1719,17 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
       
       amp     = cell->GetAmplitude(iCell)*recalF;
       time    = cell->GetTime(iCell)*1e9;//transform time to ns
-      if(time < fTimeCutMin || time > fTimeCutMax) continue;
       
+      //Remove noisy channels, only possible in ESDs
+      if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
+        if(time < fTimeCutMin || time > fTimeCutMax) continue;
+      }
       //if(amp > 3 && fCalorimeter=="EMCAL") printf("Amp = %f, time = %f, (mod, col, row)= (%d,%d,%d)\n",
       //                                                                                  amp,time,nModule,icol,irow);
       
-      //printf("%s: time %g\n",fCalorimeter.Data(), time);
       id      = cell->GetCellNumber(iCell);
       fhAmplitude->Fill(amp);
       fhAmpId    ->Fill(amp,id);
-      fhTime     ->Fill(time);
-      fhTimeId   ->Fill(time,id);
-      fhTimeAmp  ->Fill(amp,time);
-      //Double_t t0 = GetReader()->GetInputEvent()->GetT0();
-      //printf("---->>> Time EMCal %e, T0 %e, T0 vertex %e, T0 clock %e, T0 trig %d \n",time,t0, 
-      //          GetReader()->GetInputEvent()->GetT0zVertex(),
-      //          GetReader()->GetInputEvent()->GetT0clock(),
-      //          GetReader()->GetInputEvent()->GetT0Trig());
-      //fhT0Time     ->Fill(time-t0);
-      //fhT0TimeId   ->Fill(time-t0,id);
-      //fhT0TimeAmp  ->Fill(amp,time-t0);
       
       fhAmplitudeMod[nModule]->Fill(amp);
       if(fCalorimeter=="EMCAL"){
@@ -1687,48 +1739,69 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
         fhAmplitudeModFraction[nModule*3+ifrac]->Fill(amp);
       }
       
-      fhTimeAmpPerRCU  [nModule*fNRCU+iRCU]->Fill(amp, time);
-      //printf("id %d, nModule %d, iRCU %d: Histo Name %s\n",id, nModule,iRCU, fhTimeAmpPerRCU[nModule*fNRCU+iRCU]->GetName());
-      //fhT0TimeAmpPerRCU[nModule*fNRCU+iRCU]->Fill(amp, time-t0);
       nCellsInModule[nModule]++;
       fhGridCellsMod[nModule]    ->Fill(icol,irow);
       fhGridCellsEMod[nModule]   ->Fill(icol,irow,amp);
-      if(amp > 0.3){
-        fhGridCellsTimeMod[nModule]->Fill(icol,irow,time);
+      
+      if(GetReader()->GetDataType() == AliCaloTrackReader::kESD){
+        //printf("%s: time %g\n",fCalorimeter.Data(), time);
+        fhTime     ->Fill(time);
+        fhTimeId   ->Fill(time,id);
+        fhTimeAmp  ->Fill(amp,time);
         
-        //                                     AliESDCaloCells * cell2 = 0x0; 
-        //                                     if(fCalorimeter == "PHOS") cell2 =  GetReader()->GetInputEvent()->GetPHOSCells();
-        //                                     else                       cell2 = GetReader()->GetInputEvent()->GetEMCALCells();
-        //                                     Int_t icol2    = -1;
-        //                                     Int_t irow2    = -1;
-        //                                     Int_t iRCU2    = -1;
-        //                                     Float_t amp2   =  0.;
-        //                                     Float_t time2  =  0.;
-        //                                     Int_t id2      = -1;
-        //                                     Int_t nModule2 = -1;
-        //                                     for (Int_t iCell2 = 0; iCell2 < ncells; iCell2++) {  
-        //                                             amp2    = cell2->GetAmplitude(iCell2);
-        //                                             if(amp2 < 0.3) continue;
-        //                                             if(iCell2 == iCell) continue;
-        //                                             time2    = cell2->GetTime(iCell2)*1e9;//transform time to ns
-        //                                             //printf("%s: time %g\n",fCalorimeter.Data(), time);
-        //                                             id2      = cell2->GetCellNumber(iCell2);
-        //                                             nModule2 = GetModuleNumberCellIndexes(cell2->GetCellNumber(iCell2), fCalorimeter, icol2, irow2, iRCU2);
-        //                                             Int_t index = (nModule2*fNRCU+iRCU2)+(fNModules*fNRCU)*(iRCU+fNRCU*nModule); 
-        //                                             //printf("id %d, nModule %d, iRCU %d, id2 %d, nModule2 %d, iRCU2 %d, index %d: Histo Name %s\n",id, nModule,iRCU,cell2->GetCellNumber(iCell2),nModule2,iRCU2,index, fhTimeCorrRCU[index]->GetName());
-        //                                             fhTimeCorrRCU[index]->Fill(time,time2); 
-        //                                             
-        //                                     }// second cell loop
-      }// amplitude cut
-    }//nmodules
+        //Double_t t0 = GetReader()->GetInputEvent()->GetT0();
+        //printf("---->>> Time EMCal %e, T0 %e, T0 vertex %e, T0 clock %e, T0 trig %d \n",time,t0, 
+        //        GetReader()->GetInputEvent()->GetT0zVertex(),
+        //        GetReader()->GetInputEvent()->GetT0clock(),
+        //        GetReader()->GetInputEvent()->GetT0Trig());
+        //fhT0Time     ->Fill(time-t0);
+        //fhT0TimeId   ->Fill(time-t0,id);
+        //fhT0TimeAmp  ->Fill(amp,time-t0);
+        
+        //printf("id %d, nModule %d, iRCU %d: Histo Name %s\n",id, nModule,iRCU, fhTimeAmpPerRCU[nModule*fNRCU+iRCU]->GetName());
+        //fhT0TimeAmpPerRCU[nModule*fNRCU+iRCU]->Fill(amp, time-t0);
+        
+        fhTimeAmpPerRCU  [nModule*fNRCU+iRCU]->Fill(amp, time);
+        
+        if(amp > 0.3){
+          fhGridCellsTimeMod[nModule]->Fill(icol,irow,time);
+          
+          //                                   AliESDCaloCells * cell2 = 0x0; 
+          //                                   if(fCalorimeter == "PHOS") cell2 =  GetReader()->GetInputEvent()->GetPHOSCells();
+          //                                   else                       cell2 = GetReader()->GetInputEvent()->GetEMCALCells();
+          //                                   Int_t icol2    = -1;
+          //                                   Int_t irow2    = -1;
+          //                                   Int_t iRCU2    = -1;
+          //                                   Float_t amp2   =  0.;
+          //                                   Float_t time2  =  0.;
+          //                                   Int_t id2      = -1;
+          //                                   Int_t nModule2 = -1;
+          //                                   for (Int_t iCell2 = 0; iCell2 < ncells; iCell2++) {  
+          //                                           amp2    = cell2->GetAmplitude(iCell2);
+          //                                           if(amp2 < 0.3) continue;
+          //                                           if(iCell2 == iCell) continue;
+          //                                           time2    = cell2->GetTime(iCell2)*1e9;//transform time to ns
+          //                                           //printf("%s: time %g\n",fCalorimeter.Data(), time);
+          //                                           id2      = cell2->GetCellNumber(iCell2);
+          //                                           nModule2 = GetModuleNumberCellIndexes(cell2->GetCellNumber(iCell2), fCalorimeter, icol2, irow2, iRCU2);
+          //                                           Int_t index = (nModule2*fNRCU+iRCU2)+(fNModules*fNRCU)*(iRCU+fNRCU*nModule); 
+          //                                           //printf("id %d, nModule %d, iRCU %d, id2 %d, nModule2 %d, iRCU2 %d, index %d: Histo Name %s\n",id, nModule,iRCU,cell2->GetCellNumber(iCell2),nModule2,iRCU2,index, fhTimeCorrRCU[index]->GetName());
+          //                                           fhTimeCorrRCU[index]->Fill(time,time2); 
+          //                                           
+          //                                   }// second cell loop
+          
+        }// amplitude cut
+      }
+    
     
     //Get Eta-Phi position of Cell
-    //if(fFillAllPosHisto)
+    if(fFillAllPosHisto)
     {
       if(fCalorimeter=="EMCAL" && GetCaloUtils()->IsEMCALGeoMatrixSet()){
         Float_t celleta = 0.;
         Float_t cellphi = 0.;
         GetEMCALGeometry()->EtaPhiFromIndex(id, celleta, cellphi); 
+        
         fhEtaPhiAmp->Fill(celleta,cellphi,amp);
         Double_t cellpos[] = {0, 0, 0};
         GetEMCALGeometry()->GetGlobal(id, cellpos);
@@ -1756,10 +1829,12 @@ void  AliAnaCalorimeterQA::MakeAnalysisFillHistograms()
         fhXYZCell->Fill(xyz.X(),xyz.Y(),xyz.Z())  ;
       }//PHOS cells
     }//fill cell position histograms
+    
     if     (fCalorimeter=="EMCAL" && amp > fEMCALCellAmpMin) ncells ++ ;
     else if(fCalorimeter=="PHOS"  && amp > fPHOSCellAmpMin)  ncells ++ ;
     //else  
-    //  printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - no %s CELLS passed the analysis cut\n",fCalorimeter.Data());       
+    //  printf("AliAnaCalorimeterQA::MakeAnalysisFillHistograms() - no %s CELLS passed the analysis cut\n",fCalorimeter.Data());    
+    }//nmodules
   }//cell loop
   if(ncells > 0 )fhNCells->Fill(ncells) ; //fill the cells after the cut 
   
@@ -1816,20 +1891,17 @@ void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom, const Doub
   fhNCellsPerClusterMIP->Fill(e, nCaloCellsPerCluster,eta);
   
   //Position
-  //if(fFillAllPosHisto)
-  {
-    fhXE     ->Fill(pos[0],e);
-    fhYE     ->Fill(pos[1],e);
-    fhZE     ->Fill(pos[2],e);
-    fhXYZ    ->Fill(pos[0], pos[1],pos[2]);
-    
-    fhXNCells->Fill(pos[0],nCaloCellsPerCluster);
-    fhYNCells->Fill(pos[1],nCaloCellsPerCluster);
-    fhZNCells->Fill(pos[2],nCaloCellsPerCluster);
-    Float_t rxyz = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);//+pos[2]*pos[2]);
-    fhRE     ->Fill(rxyz,e);
-    fhRNCells->Fill(rxyz  ,nCaloCellsPerCluster);
-  }
+  fhXE     ->Fill(pos[0],e);
+  fhYE     ->Fill(pos[1],e);
+  fhZE     ->Fill(pos[2],e);
+  fhXYZ    ->Fill(pos[0], pos[1],pos[2]);
+    
+  fhXNCells->Fill(pos[0],nCaloCellsPerCluster);
+  fhYNCells->Fill(pos[1],nCaloCellsPerCluster);
+  fhZNCells->Fill(pos[2],nCaloCellsPerCluster);
+  Float_t rxyz = TMath::Sqrt(pos[0]*pos[0]+pos[1]*pos[1]);//+pos[2]*pos[2]);
+  fhRE     ->Fill(rxyz,e);
+  fhRNCells->Fill(rxyz  ,nCaloCellsPerCluster);
   
   fhClusterTimeEnergy->Fill(e,tof);
        
@@ -1932,10 +2004,8 @@ void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom, const Doub
     }
     else if(GetReader()->ReadAODMCParticles() && !GetMCAnalysisUtils()->CheckTagBit(tag, AliMCAnalysisUtils::kMCUnknown)){//it MC AOD and known tag
       //Get the list of MC particles
-      if(!GetReader()->GetAODMCParticles(0))   {
-        printf("AliAnaCalorimeterQA::ClusterHistograms() -  MCParticles not available!\n");
-        abort();
-      }                
+      if(!GetReader()->GetAODMCParticles(0)) 
+        AliFatal("MCParticles not available!");
       
       aodprimary = (AliAODMCParticle*) (GetReader()->GetAODMCParticles(0))->At(label);
       iMother = label;
@@ -2119,9 +2189,7 @@ void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom, const Doub
        
   //Match tracks and clusters
   //To be Modified in case of AODs
-       
-  //if(ntracksmatched==1 && trackIndex==-1) ntracksmatched=0;
-       
+               
   if( nTracksMatched > 0){
     if(fFillAllTH12){
       fhECharged      ->Fill(e);       
@@ -2152,54 +2220,53 @@ void AliAnaCalorimeterQA::ClusterHistograms(const TLorentzVector mom, const Doub
     Int_t nTPC   = 0;
     
     //In case of ESDs get the parameters in this way
-    //         if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
-    if (track->GetOuterParam() ) {
-      okout = kTRUE;
-      
-      bfield = GetReader()->GetInputEvent()->GetMagneticField();
-      okpos = track->GetOuterParam()->GetXYZAt(radius,bfield,emcpos);
-      okmom = track->GetOuterParam()->GetPxPyPzAt(radius,bfield,emcmom);
-      if(!(okpos && okmom)) return;
-      
-      TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
-      TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
-      tphi = position.Phi();
-      teta = position.Eta();
-      tmom = momentum.Mag();
-      
-      //Double_t tphi  = track->GetOuterParam()->Phi();
-      //Double_t teta  = track->GetOuterParam()->Eta();
-      //Double_t tmom  = track->GetOuterParam()->P();
-      tpt       = track->Pt();
-      tmom2     = track->P();
-      tpcSignal = track->GetTPCsignal();
-      
-      nITS = track->GetNcls(0);
-      nTPC = track->GetNcls(1);
-    }//Outer param available 
-    //}// ESDs
-    //                 else if(GetReader()->GetDataType()==AliCaloTrackReader::kAOD) {
-    //                         AliAODPid* pid = (AliAODPid*) ((AliAODTrack *) track)->GetDetPid();
-    //                         if (pid) {
-    //                                 okout = kTRUE;
-    //                                 pid->GetEMCALPosition(emcpos);
-    //                                 pid->GetEMCALMomentum(emcmom);  
-    //                                 
-    //                                 TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
-    //                                 TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
-    //                                 tphi = position.Phi();
-    //                                 teta = position.Eta();
-    //                                 tmom = momentum.Mag();
-    //                                 
-    //                                 tpt       = ((AliAODTrack*)track)->Pt();
-    //                                 tmom2     = ((AliAODTrack*)track)->P();
-    //                                 tpcSignal = pid->GetTPCsignal();
-    //                         
-    //                                 //nITS = ((AliAODTrack*)track)->GetNcls(0);
-    //                                 //nTPC = ((AliAODTrack*)track)->GetNcls(1);
-    //                         }//Outer param available 
-    //                 }//AODs
-    //                 else return; //Do nothing case not implemented.
+    if(GetReader()->GetDataType()==AliCaloTrackReader::kESD) {
+      if (track->GetOuterParam() ) {
+        okout = kTRUE;
+        
+        bfield = GetReader()->GetInputEvent()->GetMagneticField();
+        okpos = track->GetOuterParam()->GetXYZAt(radius,bfield,emcpos);
+        okmom = track->GetOuterParam()->GetPxPyPzAt(radius,bfield,emcmom);
+        if(!(okpos && okmom)) return;
+        
+        TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
+        TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
+        tphi = position.Phi();
+        teta = position.Eta();
+        tmom = momentum.Mag();
+        
+        //Double_t tphi  = track->GetOuterParam()->Phi();
+        //Double_t teta  = track->GetOuterParam()->Eta();
+        //Double_t tmom  = track->GetOuterParam()->P();
+        tpt       = track->Pt();
+        tmom2     = track->P();
+        tpcSignal = track->GetTPCsignal();
+        
+        nITS = track->GetNcls(0);
+        nTPC = track->GetNcls(1);
+      }//Outer param available 
+    }// ESDs
+    else if(GetReader()->GetDataType()==AliCaloTrackReader::kAOD) {
+      AliAODPid* pid = (AliAODPid*) ((AliAODTrack *) track)->GetDetPid();
+      if (pid) {
+        okout = kTRUE;
+        pid->GetEMCALPosition(emcpos);
+        pid->GetEMCALMomentum(emcmom); 
+        
+        TVector3 position(emcpos[0],emcpos[1],emcpos[2]);
+        TVector3 momentum(emcmom[0],emcmom[1],emcmom[2]);
+        tphi = position.Phi();
+        teta = position.Eta();
+        tmom = momentum.Mag();
+        
+        tpt       = track->Pt();
+        tmom2     = track->P();
+        tpcSignal = pid->GetTPCsignal();
+        
+        //nITS = ((AliAODTrack*)track)->GetNcls(0);
+        //nTPC = ((AliAODTrack*)track)->GetNcls(1);
+      }//pid 
+    }//AODs
                
     if(okout){
       Double_t deta = teta - eta;
@@ -3177,7 +3244,7 @@ void  AliAnaCalorimeterQA::Terminate(TList* outputList)
   Int_t rbX = 1;
   Int_t rbY = 1;
   Int_t rbZ = 1;
-  //if(fFillAllPosHisto)
+  if(fFillAllPosHisto)
   {
     snprintf(cname,buffersize,"%s_QA_ClusterXY",fCalorimeter.Data());
     TCanvas  * cxyz = new TCanvas(cname, "Cluster XY distributions", 1200, 400) ;