adapted to part corr fw
authorrvernet <rvernet@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 25 Mar 2011 18:43:05 +0000 (18:43 +0000)
committerrvernet <rvernet@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 25 Mar 2011 18:43:05 +0000 (18:43 +0000)
PWG4/PartCorrDep/AliAnaShowerParameter.cxx
PWG4/PartCorrDep/AliAnaShowerParameter.h

index 486ed4f..98dba8f 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes hereby granted      *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-/* $Id: AliAnaPhoton.cxx 28688 2008-09-11 15:04:07Z gconesab $ */
-
-//_________________________________________________________________________
-//
-// Class cloned from AliAnaPhoton, main aim is shower shape studies
-// 
-// 
-//
-//-- Author: Jocelyn Mlynarz (WSU) and Gustavo Conesa (LPSC-CNRS)
-//////////////////////////////////////////////////////////////////////////////
-
-
-// --- ROOT system --- 
-#include <TH3F.h>
-#include <TH2F.h>
-#include <TClonesArray.h>
-//#include <TObjString.h>
-#include <Riostream.h>
-#include "TParticle.h"
-//#include <fstream>
-
-// --- Analysis system --- 
-#include "AliAnaShowerParameter.h" 
-#include "AliCaloTrackReader.h"
-#include "AliStack.h"
-#include "AliCaloPID.h"
-#include "AliMCAnalysisUtils.h"
-#include "AliFiducialCut.h"
-#include "AliAODCaloCluster.h"
-#include "AliAODMCParticle.h"
-#include "AliAnaShowerParameter.h" 
-#include "AliEMCALGeoUtils.h"
-#include "AliAODCaloCells.h"
-#include "AliAODEvent.h"
-
-
-ClassImp(AliAnaShowerParameter)
-
-//____________________________________________________________________________
-AliAnaShowerParameter::AliAnaShowerParameter() : 
-AliAnaPartCorrBaseClass(), fCalorimeter(""), 
-fMinDist(0.),fMinDist2(0.),fMinDist3(0.),fRejectTrackMatch(0),
-fCheckConversion(kFALSE),fAddConvertedPairsToAOD(kFALSE), fMassCut(0), fNCellsCut(0),
-fTimeCutMin(-1), fTimeCutMax(9999999),fEMCALGeo(0),fNumClusters(0),
-fhNClusters(0),fhNCellCluster(0),fhPtCluster(0),fhPhiCluster(0),fhEtaCluster(0),fhDeltaPhiClusters(0),fhDeltaEtaClusters(0),
-fhLambdaCluster(0),fhDispersionCluster(0),fhELambdaCluster(0),fhELambdaCellCluster(0),
-fhNCellPhoton(0),fhLambdaPhoton(0),fhDispersionPhoton(0),fhELambdaPhoton(0),fhELambdaCellPhoton(0),
-fhNCellPi0(0),fhLambdaPi0(0),fhDispersionPi0(0),fhELambdaPi0(0),fhELambdaCellPi0(0),
-fhNCellChargedHadron(0),fhLambdaChargedHadron(0),fhDispersionChargedHadron(0),fhELambdaChargedHadron(0),fhELambdaCellChargedHadron(0),  
-//MC
-fhDeltaE(0), fhDeltaPt(0),fhRatioE(0), fhRatioPt(0),fh2E(0),fh2Pt(0),fhMCPdg(0),
-fhEMCPhoton(0),fhPhiMCPhoton(0),fhEtaMCPhoton(0),fhLambdaMCPhoton(0),fhPhotTrueE(0),fhRatioEPhoton(0),
-fhEMCPi0(0),fhPhiMCPi0(0),fhEtaMCPi0(0),fhLambdaMCPi0(0),fhPi0TrueE(0),
-fhProductionDistance(0),fhProductionRadius(0),fhDecayAngle(0),fhRatioEPi0(0),
-fhEMCPion(0),fhPhiMCPion(0),fhEtaMCPion(0),fhLambdaMCPion(0),fhPionTrueE(0),fhRatioEPion(0), 
-fhEMCProton(0),fhPhiMCProton(0),fhEtaMCProton(0),fhLambdaMCProton(0),fhProtonTrueE(0),fhRatioEProton(0),
-fhEMCAntiProton(0),fhPhiMCAntiProton(0),fhEtaMCAntiProton(0),fhLambdaMCAntiProton(0),fhAntiProtonTrueE(0),fhRatioEAntiProton(0),
-fhEMCNeutron(0),fhPhiMCNeutron(0),fhEtaMCNeutron(0),fhLambdaMCNeutron(0),fhNeutronTrueE(0),fhRatioENeutron(0),
-fhEMCEta(0),fhPhiMCEta(0),fhEtaMCEta(0),fhLambdaMCEta(0),
-fhPrimPt(0x0)
-{
-  //default ctor
-  
-  //Initialize parameters
-  InitParameters();
-  
-}
-
-//____________________________________________________________________________
-AliAnaShowerParameter::~AliAnaShowerParameter() 
-{
-  //dtor
-  
-}
-
-//________________________________________________________________________
-TObjString *  AliAnaShowerParameter::GetAnalysisCuts()
-{      
-  //Save parameters used for analysis
-  TString parList ; //this will be list of parameters used for this analysis.
-  const Int_t buffersize = 255;
-  char onePar[buffersize] ;
-  
-  snprintf(onePar,buffersize,"--- AliAnaShowerParameter ---\n") ;
-  parList+=onePar ;    
-  snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;
-  parList+=onePar ;
-  snprintf(onePar,buffersize,"fMinDist =%2.2f (Minimal distance to bad channel to accept cluster) \n",fMinDist) ;
-  parList+=onePar ;
-  snprintf(onePar,buffersize,"fMinDist2=%2.2f (Cuts on Minimal distance to study acceptance evaluation) \n",fMinDist2) ;
-  parList+=onePar ;
-  snprintf(onePar,buffersize,"fMinDist3=%2.2f (One more cut on distance used for acceptance-efficiency study) \n",fMinDist3) ;
-  parList+=onePar ;
-  snprintf(onePar,buffersize,"fRejectTrackMatch: %d\n",fRejectTrackMatch) ;
-  parList+=onePar ;  
-  
-  //Get parameters set in base class.
-  parList += GetBaseParametersList() ;
-  
-  //Get parameters set in PID class.
-  parList += GetCaloPID()->GetPIDParametersList() ;
-  
-  //Get parameters set in FiducialCut class (not available yet)
-  //parlist += GetFidCut()->GetFidCutParametersList() 
-  
-  return new TObjString(parList) ;
-}
-
-
-//________________________________________________________________________
-TList *  AliAnaShowerParameter::GetCreateOutputObjects()
-{  
-  // Create histograms to be saved in output file and 
-  // store them in outputContainer
-  TList * outputContainer = new TList() ; 
-  outputContainer->SetName("PhotonHistos") ; 
-       
-  Int_t nptbins  = GetHistoPtBins();
-  Int_t nphibins = GetHistoPhiBins();
-  Int_t netabins = GetHistoEtaBins();
-  Float_t ptmax  = GetHistoPtMax();
-  Float_t phimax = GetHistoPhiMax();
-  Float_t etamax = GetHistoEtaMax();
-  Float_t ptmin  = GetHistoPtMin();
-  Float_t phimin = GetHistoPhiMin();
-  Float_t etamin = GetHistoEtaMin();   
-  
-  //General non-MC Cluster histograms
-  fhNClusters = new TH1F ("hNClusters","NClusters",21,-0.5,20.5); 
-  fhNClusters->SetXTitle("N_{Clusters}");
-  outputContainer->Add(fhNClusters) ;
-  
-  fhNCellCluster = new TH2F ("hNCellCluster","Number of cell per cluster",nptbins,ptmin,ptmax,21,-0.5,20.5); 
-  fhNCellCluster->SetXTitle("p_{T}");
-  fhNCellCluster->SetYTitle("N_{Cells}");  
-  outputContainer->Add(fhNCellCluster) ;
-  
-  fhPtCluster  = new TH1F("hPtCluster","Number of #gamma over calorimeter",nptbins,ptmin,ptmax); 
-  fhPtCluster->SetXTitle("p_{Cluster}(GeV/c)");
-  outputContainer->Add(fhPtCluster) ; 
-  
-  fhPhiCluster  = new TH2F
-  ("hPhiCluster","#phi_{Cluster}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-  fhPhiCluster->SetYTitle("#phi");
-  fhPhiCluster->SetXTitle("E_{Cluster} (GeV/c)");
-  outputContainer->Add(fhPhiCluster) ; 
-  
-  fhEtaCluster  = new TH2F
-  ("hEtaCluster","#phi_{Cluster}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-  fhEtaCluster->SetYTitle("#eta");  
-  fhEtaCluster->SetXTitle("E_{Cluster} (GeV/c)");
-  outputContainer->Add(fhEtaCluster) ;
-
-  fhDeltaEtaClusters  = new TH1F("hDeltaEtaClusters","#Delta #eta of clusters over calorimeter",netabins,etamin,etamax); 
-  fhDeltaEtaClusters->SetXTitle("#Delta #eta_{Clusters}");
-  outputContainer->Add(fhDeltaEtaClusters) ; 
-
-  fhDeltaPhiClusters  = new TH1F("hDeltaPhiClusters","#Delta #phi of clusters over calorimeter",nphibins,phimin,phimax); 
-  fhDeltaPhiClusters->SetXTitle("#Delta #phi_{Clusters}");
-  outputContainer->Add(fhDeltaPhiClusters) ; 
-  
-  fhLambdaCluster  = new TH3F
-  ("hLambdaCluster","#lambda_{Cluster}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhLambdaCluster->SetZTitle("#lambda_{1}^{2}"); 
-  fhLambdaCluster->SetYTitle("#lambda_{0}^{2}");
-  fhLambdaCluster->SetXTitle("p_{T Cluster} (GeV/c)");
-  outputContainer->Add(fhLambdaCluster) ;
-  
-  fhELambdaCluster  = new TH3F
-  ("hELambdaCluster","#lambda_{Cluster}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhELambdaCluster->SetZTitle("#lambda_{1}^{2}"); 
-  fhELambdaCluster->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaCluster->SetXTitle("p_{T Cluster} (GeV)");
-  outputContainer->Add(fhELambdaCluster) ;
-  
-  fhDispersionCluster  = new TH2F
-  ("hDispersionCluster","Dispersion_{Cluster}",nptbins,ptmin,ptmax,200,0,2); 
-  fhDispersionCluster->SetYTitle("Dispersion");
-  fhDispersionCluster->SetXTitle("p_{T Cluster} (GeV/c)");
-  outputContainer->Add(fhDispersionCluster) ;
-  
-  fhELambdaCellCluster  = new TH3F
-  ("hELambdaCellCluster","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-  fhELambdaCellCluster->SetZTitle("N Cells"); 
-  fhELambdaCellCluster->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaCellCluster->SetXTitle("E_{Cluster} (GeV)");
-  outputContainer->Add(fhELambdaCellCluster) ;
-  
-  //Identified photon histograms
-  fhNCellPhoton = new TH2F ("hNCellPhoton","Number of cell per cluster",nptbins,ptmin,ptmax,21,-0.5,20.5); 
-  fhNCellPhoton->SetXTitle("E");
-  fhNCellPhoton->SetYTitle("N_{Cells}");  
-  outputContainer->Add(fhNCellPhoton) ;
-  
-  fhLambdaPhoton  = new TH3F
-  ("hLambdaPhoton","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhLambdaPhoton->SetZTitle("#lambda_{1}^{2}"); 
-  fhLambdaPhoton->SetYTitle("#lambda_{0}^{2}");
-  fhLambdaPhoton->SetXTitle("p_{T Photon} (GeV/c)");
-  outputContainer->Add(fhLambdaPhoton) ;
-  
-  fhELambdaPhoton  = new TH3F
-  ("hELambdaPhoton","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhELambdaPhoton->SetZTitle("#lambda_{1}^{2}"); 
-  fhELambdaPhoton->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaPhoton->SetXTitle("E_{Photon} (GeV)");
-  outputContainer->Add(fhELambdaPhoton) ;
-  
-  fhDispersionPhoton  = new TH2F
-  ("hLambdaDispersion","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2); 
-  fhDispersionPhoton->SetYTitle("Dispersion");
-  fhDispersionPhoton->SetXTitle("p_{T Photon} (GeV/c)");
-  outputContainer->Add(fhDispersionPhoton) ;
-  
-  fhELambdaCellPhoton  = new TH3F
-  ("hELambdaCellPhoton","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-  fhELambdaCellPhoton->SetZTitle("N Cells"); 
-  fhELambdaCellPhoton->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaCellPhoton->SetXTitle("E_{Photon} (GeV)");
-  outputContainer->Add(fhELambdaCellPhoton) ;
-  
-  //Identified Pi0 histograms
-  fhNCellPi0 = new TH2F ("hNCellPi0","Number of cell per cluster",nptbins,ptmin,ptmax,21,-0.5,20.5); 
-  fhNCellPi0->SetXTitle("E");
-  fhNCellPi0->SetYTitle("N_{Cells}");  
-  outputContainer->Add(fhNCellPi0) ;
-  
-  fhLambdaPi0  = new TH3F
-  ("hLambdaPi0","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhLambdaPi0->SetZTitle("#lambda_{1}^{2}"); 
-  fhLambdaPi0->SetYTitle("#lambda_{0}^{2}");
-  fhLambdaPi0->SetXTitle("p_{T Pi0} (GeV/c)");
-  outputContainer->Add(fhLambdaPi0) ;
-  
-  fhELambdaPi0  = new TH3F
-  ("hELambdaPi0","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhELambdaPi0->SetZTitle("#lambda_{1}^{2}"); 
-  fhELambdaPi0->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaPi0->SetXTitle("E_{Pi0} (GeV)");
-  outputContainer->Add(fhELambdaPi0) ;
-  
-  fhDispersionPi0  = new TH2F
-  ("hLambdaDispersion","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2); 
-  fhDispersionPi0->SetYTitle("Dispersion");
-  fhDispersionPi0->SetXTitle("p_{T Pi0} (GeV/c)");
-  outputContainer->Add(fhDispersionPi0) ;
-  
-  fhELambdaCellPi0  = new TH3F
-  ("hELambdaCellPi0","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-  fhELambdaCellPi0->SetZTitle("N Cells"); 
-  fhELambdaCellPi0->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaCellPi0->SetXTitle("E_{Pi0} (GeV)");
-  outputContainer->Add(fhELambdaCellPi0) ;
-  
-  //Identified charged hadron histograms
-  fhNCellChargedHadron = new TH2F ("hNCellChargedHadron","Number of cell per cluster",nptbins,ptmin,ptmax,21,-0.5,20.5); 
-  fhNCellChargedHadron->SetXTitle("p_{T}");
-  fhNCellChargedHadron->SetYTitle("N_{Cells}");  
-  outputContainer->Add(fhNCellChargedHadron) ;
-  
-  fhLambdaChargedHadron  = new TH3F
-  ("hLambdaChargedHadron","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhLambdaChargedHadron->SetZTitle("#lambda_{1}^{2}"); 
-  fhLambdaChargedHadron->SetYTitle("#lambda_{0}^{2}");
-  fhLambdaChargedHadron->SetXTitle("p_{T ChargedHadron} (GeV/c)");
-  outputContainer->Add(fhLambdaChargedHadron) ;
-  
-  fhELambdaChargedHadron  = new TH3F
-  ("hELambdaChargedHadron","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,200,0,2); 
-  fhELambdaChargedHadron->SetZTitle("#lambda_{1}^{2}"); 
-  fhELambdaChargedHadron->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaChargedHadron->SetXTitle("E_{ChargedHadron} (GeV)");
-  outputContainer->Add(fhELambdaChargedHadron) ;
-  
-  fhDispersionChargedHadron  = new TH2F
-  ("hLambdaDispersion","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2); 
-  fhDispersionChargedHadron->SetYTitle("Dispersion");
-  fhDispersionChargedHadron->SetXTitle("p_{T ChargedHadron} (GeV/c)");
-  outputContainer->Add(fhDispersionChargedHadron) ;
-  
-  fhELambdaCellChargedHadron  = new TH3F
-  ("hELambdaCellChargedHadron","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-  fhELambdaCellChargedHadron->SetZTitle("N Cells"); 
-  fhELambdaCellChargedHadron->SetYTitle("#lambda_{0}^{2}");
-  fhELambdaCellChargedHadron->SetXTitle("E_{ChargedHadron} (GeV)");
-  outputContainer->Add(fhELambdaCellChargedHadron) ;
-  
-  if(IsDataMC()){
-    fhDeltaE  = new TH1F ("hDeltaE","MC - Reco E ", 200,-50,50); 
-    fhDeltaE->SetXTitle("#Delta E (GeV)");
-    outputContainer->Add(fhDeltaE);
-    
-    fhDeltaPt  = new TH1F ("hDeltaPt","MC - Reco p_{T} ", 200,-50,50); 
-    fhDeltaPt->SetXTitle("#Delta p_{T} (GeV/c)");
-    outputContainer->Add(fhDeltaPt);
-    
-    fhRatioE  = new TH1F ("hRatioE","Reco/MC E ", 1000,0,2); 
-    fhRatioE->SetXTitle("E_{reco}/E_{gen}");
-    outputContainer->Add(fhRatioE);
-    
-    fhRatioPt  = new TH1F ("hRatioPt","Reco/MC p_{T} ", 1000,0,2); 
-    fhRatioPt->SetXTitle("p_{T, reco}/p_{T, gen}");
-    outputContainer->Add(fhRatioPt);    
-    
-    fh2E  = new TH2F ("h2E","E distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax); 
-    fh2E->SetXTitle("E_{rec} (GeV)");
-    fh2E->SetYTitle("E_{gen} (GeV)");
-    outputContainer->Add(fh2E);          
-    
-    fh2Pt  = new TH2F ("h2Pt","p_T distribution, reconstructed vs generated", nptbins,ptmin,ptmax,nptbins,ptmin,ptmax); 
-    fh2Pt->SetXTitle("p_{T,rec} (GeV/c)");
-    fh2Pt->SetYTitle("p_{T,gen} (GeV/c)");
-    outputContainer->Add(fh2Pt);
-
-    fhMCPdg  = new TH2F ("hMCPdg","PDG Code distribution",nptbins,ptmin,ptmax,6001,-3000.5,3000.5); 
-    fhMCPdg->SetXTitle("E_{Reco Cluster} (GeV)");
-    fhMCPdg->SetYTitle("PDG Code");
-    outputContainer->Add(fhMCPdg);
-    
-    fhEMCPhoton  = new TH1F("hEMCPhoton","Number of #gamma over calorimeter",nptbins,ptmin,ptmax); 
-    //fhEMCPhoton->SetYTitle("N");
-    fhEMCPhoton->SetXTitle("E_{#gamma}(GeV)");
-    outputContainer->Add(fhEMCPhoton) ; 
-    
-    fhPhiMCPhoton  = new TH2F
-    ("hPhiMCPhoton","#phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCPhoton->SetYTitle("#phi");
-    fhPhiMCPhoton->SetXTitle("E_{ #gamma} (GeV)");
-    outputContainer->Add(fhPhiMCPhoton) ; 
-    
-    fhEtaMCPhoton  = new TH2F
-    ("hEtaMCPhoton","#phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCPhoton->SetYTitle("#eta");
-    fhEtaMCPhoton->SetXTitle("E_{ #gamma} (GeV)");
-    outputContainer->Add(fhEtaMCPhoton) ;
-    
-    fhLambdaMCPhoton  = new TH3F
-    ("hLambdaMCPhoton","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCPhoton->SetZTitle("N_{Cells}");
-    fhLambdaMCPhoton->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCPhoton->SetXTitle("E_{#gamma} (GeV)");
-    outputContainer->Add(fhLambdaMCPhoton) ;
-    
-    fhPhotTrueE  = new TH3F
-    ("hPhotTrueE","#lambda_{#gamma}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax,200,0,2); 
-    fhPhotTrueE->SetZTitle("#lambda_{0}^{2}");
-    fhPhotTrueE->SetYTitle("Recons. E_{#gamma}");
-    fhPhotTrueE->SetXTitle("MC Truth E_{#gamma} (GeV)");
-    outputContainer->Add(fhPhotTrueE) ;
-    
-    fhRatioEPhoton  = new TH1F ("hRatioEPhoton","E_{Reco #gamma}/E_{MC truth #gamma}", 1000,0,2); 
-    fhRatioEPhoton->SetXTitle("E_{reco #gamma}/E_{gen #gamma}");
-    outputContainer->Add(fhRatioEPhoton);
-    
-    fhEMCPi0  = new TH1F("hEMCPi0","Number of #pi^0 over calorimeter",nptbins,ptmin,ptmax); 
-    fhEMCPi0->SetYTitle("N");
-    fhEMCPi0->SetXTitle("E_{ #pi^0}(GeV)");
-    outputContainer->Add(fhEMCPi0) ; 
-    
-    fhPhiMCPi0  = new TH2F
-    ("hPhiMCPi0","#phi_{#pi^0}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCPi0->SetYTitle("#phi");
-    fhPhiMCPi0->SetXTitle("E_{ #pi^0} (GeV)");
-    outputContainer->Add(fhPhiMCPi0) ; 
-    
-    fhEtaMCPi0  = new TH2F
-    ("hEtaMCPi0","#phi_{#pi^0}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCPi0->SetYTitle("#eta");
-    fhEtaMCPi0->SetXTitle("E_{ #pi^0} (GeV)");
-    outputContainer->Add(fhEtaMCPi0) ;
-    
-    fhLambdaMCPi0  = new TH3F
-    ("hLambdaMCPi0","#lambda_{#pi^0}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCPi0->SetZTitle("N_{Cells}");
-    fhLambdaMCPi0->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCPi0->SetXTitle("E_{ #pi^0} (GeV)");
-    outputContainer->Add(fhLambdaMCPi0) ;
-    
-    fhProductionDistance  = new TH2F
-    ("hProductionDistance","#lambda_{#pi^0}",nptbins,ptmin,ptmax,160,0,800); 
-    fhProductionDistance->SetYTitle("Distance of production of Pi0 from vertex (cm)");
-    fhProductionDistance->SetXTitle("E_{ #pi^0} (GeV)");
-    outputContainer->Add(fhProductionDistance) ;
-    
-    fhProductionRadius  = new TH2F
-    ("hProductionRadius","#lambda_{#pi^0}",nptbins,ptmin,ptmax,160,0,800); 
-    fhProductionRadius->SetYTitle("Radius of production of Pi0 from vertex (cm)");
-    fhProductionRadius->SetXTitle("E_{ #pi^0} (GeV)");
-    outputContainer->Add(fhProductionRadius) ;
-    
-    fhDecayAngle = new TH2F
-    ("hDecayAngle","#lambda_{#pi^0}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhDecayAngle->SetYTitle("Opening angle of the #pi^{0} decay (radians)");
-    fhDecayAngle->SetXTitle("E_{ #pi^0} (GeV)");
-    outputContainer->Add(fhDecayAngle) ;
-    
-    fhPi0TrueE  = new TH3F
-    ("hPi0TrueE","#lambda_{#pi^0}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax,200,0,2); 
-    fhPi0TrueE->SetZTitle("#lambda_{0}^{2}");
-    fhPi0TrueE->SetYTitle("Recons. E_{#pi^0}");
-    fhPi0TrueE->SetXTitle("MC Truth E_{#pi^0} (GeV)");
-    outputContainer->Add(fhPi0TrueE) ;
-    
-    fhRatioEPi0  = new TH1F ("hRatioEPi0","E_{Reco #pi^{0}}/E_{MC truth #pi^{0}}", 1000,0,2); 
-    fhRatioEPi0->SetXTitle("E_{reco #pi^{0}}/E_{gen #pi^{0}}");
-    outputContainer->Add(fhRatioEPi0);
-    
-    fhEMCPion  = new TH1F("hEMCPion","Number of #pi over calorimeter",nptbins,ptmin,ptmax); 
-    fhEMCPion->SetYTitle("N");
-    fhEMCPion->SetXTitle("E_{ #pi}(GeV)");
-    outputContainer->Add(fhEMCPion) ; 
-    
-    fhPhiMCPion  = new TH2F
-    ("hPhiMCPion","#phi_{#pi}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCPion->SetYTitle("#phi");
-    fhPhiMCPion->SetXTitle("E_{ #pi} (GeV)");
-    outputContainer->Add(fhPhiMCPion) ; 
-    
-    fhEtaMCPion  = new TH2F
-    ("hEtaMCPion","#phi_{#pi}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCPion->SetYTitle("#eta");
-    fhEtaMCPion->SetXTitle("E_{ #pi} (GeV)");
-    outputContainer->Add(fhEtaMCPion) ;
-    
-    fhLambdaMCPion  = new TH3F
-    ("hLambdaMCPion","#lambda_{#pi}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCPion->SetZTitle("N_{Cells}");
-    fhLambdaMCPion->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCPion->SetXTitle("E_{ #pi} (GeV)");
-    outputContainer->Add(fhLambdaMCPion) ;
-    
-    fhPionTrueE  = new TH3F
-    ("hPionTrueE","#lambda_{#pi}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax,200,0,2); 
-    fhPionTrueE->SetZTitle("#lambda_{0}^{2}");
-    fhPionTrueE->SetYTitle("Recons. E_{#pi}");
-    fhPionTrueE->SetXTitle("MC Truth E_{#pi} (GeV)");
-    outputContainer->Add(fhPionTrueE) ;
-    
-    fhRatioEPion  = new TH1F ("hRatioEPion","E_{Reco #pi^{#pm}}/E_{MC truth #pi^{#pm}}", 1000,0,2); 
-    fhRatioEPion->SetXTitle("E_{reco #pi^{#pm}}/E_{gen #pi^{#pm}}");
-    outputContainer->Add(fhRatioEPion);
-    
-    fhEMCProton  = new TH1F("hEMCProton","Number of p over calorimeter",nptbins,ptmin,ptmax); 
-    fhEMCProton->SetYTitle("N");
-    fhEMCProton->SetXTitle("E_{ p}(GeV)");
-    outputContainer->Add(fhEMCProton) ; 
-    
-    fhPhiMCProton  = new TH2F
-    ("hPhiMCProton","#phi_{p}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCProton->SetYTitle("#phi");
-    fhPhiMCProton->SetXTitle("E_{ p} (GeV)");
-    outputContainer->Add(fhPhiMCProton) ; 
-    
-    fhEtaMCProton  = new TH2F
-    ("hEtaMCProton","#phi_{p}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCProton->SetYTitle("#eta");
-    fhEtaMCProton->SetXTitle("E_{ p} (GeV)");
-    outputContainer->Add(fhEtaMCProton) ;
-    
-    fhLambdaMCProton  = new TH3F
-    ("hLambdaMCProton","#lambda_{p}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCProton->SetZTitle("N_{Cells}");
-    fhLambdaMCProton->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCProton->SetXTitle("E_{ p} (GeV)");
-    outputContainer->Add(fhLambdaMCProton) ;
-    
-    fhProtonTrueE  = new TH3F
-    ("hProtonTrueE","#lambda_{p}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax,200,0,2); 
-    fhProtonTrueE->SetZTitle("#lambda_{0}^{2}");
-    fhProtonTrueE->SetYTitle("Recons. E_{p}");
-    fhProtonTrueE->SetXTitle("MC Truth E_{p} (GeV)");
-    outputContainer->Add(fhProtonTrueE) ;
-    
-    fhRatioEProton  = new TH1F ("hRatioEProton","E_{Reco p}/E_{MC truth p}", 1000,0,2); 
-    fhRatioEProton->SetXTitle("E_{reco p}/E_{gen p}");
-    outputContainer->Add(fhRatioEProton);
-    
-    fhEMCAntiProton  = new TH1F("hEMCAntiProton","Number of #bar{p} over calorimeter",nptbins,ptmin,ptmax); 
-    fhEMCAntiProton->SetYTitle("N");
-    fhEMCAntiProton->SetXTitle("E_{#bar{p}}(GeV)");
-    outputContainer->Add(fhEMCAntiProton) ; 
-    
-    fhPhiMCAntiProton  = new TH2F
-    ("hPhiMCAntiProton","#phi_{#bar{p}}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCAntiProton->SetYTitle("#phi");
-    fhPhiMCAntiProton->SetXTitle("E_{#bar{p}} (GeV)");
-    outputContainer->Add(fhPhiMCAntiProton) ; 
-    
-    fhEtaMCAntiProton  = new TH2F
-    ("hEtaMCAntiProton","#phi_{#bar{p}}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCAntiProton->SetYTitle("#eta");
-    fhEtaMCAntiProton->SetXTitle("E_{#bar{p}} (GeV)");
-    outputContainer->Add(fhEtaMCAntiProton) ;
-    
-    fhLambdaMCAntiProton  = new TH3F
-    ("hLambdaMCAntiProton","#lambda_{#bar{p}}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCAntiProton->SetZTitle("N_{Cells}");
-    fhLambdaMCAntiProton->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCAntiProton->SetXTitle("E_{#bar{p}} (GeV)");
-    outputContainer->Add(fhLambdaMCAntiProton) ;
-    
-    fhAntiProtonTrueE  = new TH3F
-    ("hAntiProtonTrueE","#lambda_{#bar{p}}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax,200,0,2); 
-    fhAntiProtonTrueE->SetZTitle("#lambda_{0}^{2}");
-    fhAntiProtonTrueE->SetYTitle("Recons. E_{#bar{p}}");
-    fhAntiProtonTrueE->SetXTitle("MC Truth E_{#bar{p}} (GeV)");
-    outputContainer->Add(fhAntiProtonTrueE) ;
-    
-    fhRatioEAntiProton  = new TH1F ("hRatioEAntiProton","E_{Reco #bar{p}}/E_{MC truth #bar{p}}", 1000,0,2); 
-    fhRatioEAntiProton->SetXTitle("E_{reco #bar{p}}/E_{gen #bar{p}}");
-    outputContainer->Add(fhRatioEAntiProton);
-    
-    fhEMCNeutron  = new TH1F("hEMCNeutron","Number of p over calorimeter",nptbins,ptmin,ptmax); 
-    fhEMCNeutron->SetYTitle("N");
-    fhEMCNeutron->SetXTitle("E_{ p}(GeV)");
-    outputContainer->Add(fhEMCNeutron) ; 
-    
-    fhPhiMCNeutron  = new TH2F
-    ("hPhiMCNeutron","#phi_{p}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCNeutron->SetYTitle("#phi");
-    fhPhiMCNeutron->SetXTitle("E_{ p} (GeV)");
-    outputContainer->Add(fhPhiMCNeutron) ; 
-    
-    fhEtaMCNeutron  = new TH2F
-    ("hEtaMCNeutron","#phi_{p}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCNeutron->SetYTitle("#eta");
-    fhEtaMCNeutron->SetXTitle("E_{ p} (GeV)");
-    outputContainer->Add(fhEtaMCNeutron) ;
-    
-    fhLambdaMCNeutron  = new TH3F
-    ("hLambdaMCNeutron","#lambda_{p}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCNeutron->SetZTitle("N_{Cells}");
-    fhLambdaMCNeutron->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCNeutron->SetXTitle("E_{ p} (GeV)");
-    outputContainer->Add(fhLambdaMCNeutron) ;
-    
-    fhNeutronTrueE  = new TH3F
-    ("hNeutronTrueE","#lambda_{n}",nptbins,ptmin,ptmax,nptbins,ptmin,ptmax,200,0,2); 
-    fhNeutronTrueE->SetZTitle("#lambda_{0}^{2}");
-    fhNeutronTrueE->SetYTitle("Recons. E_{n}");
-    fhNeutronTrueE->SetXTitle("MC Truth E_{n} (GeV)");
-    outputContainer->Add(fhNeutronTrueE) ;
-    
-    fhRatioENeutron  = new TH1F ("hRatioENeutron","E_{Reco n}/E_{MC truth n}", 1000,0,2); 
-    fhRatioENeutron->SetXTitle("E_{reco n}/E_{gen n}");
-    outputContainer->Add(fhRatioENeutron);
-    
-    fhEMCEta  = new TH1F("hEMCEta","Number of #eta over calorimeter",nptbins,ptmin,ptmax); 
-    fhEMCEta->SetYTitle("N");
-    fhEMCEta->SetXTitle("E_{ #eta}(GeV)");
-    outputContainer->Add(fhEMCEta) ; 
-    
-    fhPhiMCEta  = new TH2F
-    ("hPhiMCEta","#phi_{#eta}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
-    fhPhiMCEta->SetYTitle("#phi");
-    fhPhiMCEta->SetXTitle("E_{ #eta} (GeV)");
-    outputContainer->Add(fhPhiMCEta) ; 
-    
-    fhEtaMCEta  = new TH2F
-    ("hEtaMCEta","#phi_{#eta}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
-    fhEtaMCEta->SetYTitle("#eta");
-    fhEtaMCEta->SetXTitle("E_{ #eta} (GeV)");
-    outputContainer->Add(fhEtaMCEta) ;
-    
-    fhLambdaMCEta  = new TH3F
-    ("hLambdaMCEta","#lambda_{#eta}",nptbins,ptmin,ptmax,200,0,2,21,-0.5,20.5); 
-    fhLambdaMCEta->SetZTitle("N_{Cells}");
-    fhLambdaMCEta->SetYTitle("#lambda_{0}^{2}");
-    fhLambdaMCEta->SetXTitle("E_{ #eta} (GeV)");
-    outputContainer->Add(fhLambdaMCEta) ;
-    
-    fhPrimPt     = new TH1D("hPrimPt","Primary pi0 pt",nptbins,ptmin,ptmax) ;
-    outputContainer->Add(fhPrimPt) ;
-    
-  }//Histos with MC
-  
-  return outputContainer ;
-  
-}
-
-//____________________________________________________________________________
-void AliAnaShowerParameter::Init()
-{
-  
-  //Init
-  //Do some checks
-  if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){
-    printf("AliAnaShowerParameter::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("AliAnaShowerParameter::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");
-    abort();
-  }
-  
-}
-
-//____________________________________________________________________________
-void AliAnaShowerParameter::InitParameters()
-{
-  
-  //Initialize the parameters of the analysis.
-  AddToHistogramsName("AnaPhoton_");
-  
-  fCalorimeter = "PHOS" ;
-  fMinDist  = 2.;
-  fMinDist2 = 4.;
-  fMinDist3 = 5.;
-  fMassCut  = 0.03; //30 MeV
-       
-  fTimeCutMin  = -1;
-  fTimeCutMax  = 9999999;
-  fNCellsCut = 0;
-       
-  fRejectTrackMatch       = kTRUE ;
-  fCheckConversion        = kFALSE;
-  fAddConvertedPairsToAOD = kFALSE;
-  
-  fNumClusters = -1; // By Default, select all events.
-       
-}
-
-//__________________________________________________________________
-void  AliAnaShowerParameter::MakeAnalysisFillAOD() 
-{
-  //Do analysis and fill aods
-  //Search for photons in fCalorimeter 
-  
-  //Get vertex for photon momentum calculation
-  
-  for (Int_t iev = 0; iev < GetNMixedEvent(); iev++) {
-    if (!GetMixedEvent()) 
-      GetReader()->GetVertex(GetVertex(iev));
-    else 
-      GetMixedEvent()->GetVertexOfEvent(iev)->GetXYZ(GetVertex(iev)); 
-  } 
-  
-  //Select the Calorimeter of the photon
-  TObjArray * pl = 0x0; 
-  if(fCalorimeter == "PHOS")
-    pl = GetPHOSClusters();
-  else if (fCalorimeter == "EMCAL")
-    pl = GetEMCALClusters();
-  
-  if(!pl){
-    printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Careful cluster array NULL!!\n");
-    return;
-  }
-  
-  //Fill AODCaloClusters and AODParticle with PHOS/EMCAL aods
-  TLorentzVector mom, mom2 ; 
-  Int_t nCaloClusters = pl->GetEntriesFast();   
-  //Cut on the number of clusters in the event.
-  if ((fNumClusters !=-1) && (nCaloClusters != fNumClusters)) return;
-  Bool_t * indexConverted = new Bool_t[nCaloClusters];
-  for (Int_t i = 0; i < nCaloClusters; i++) 
-    indexConverted[i] = kFALSE;
-  
-  for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){    
-    
-    AliAODCaloCluster * calo =  (AliAODCaloCluster*) (pl->At(icalo));  
-    Int_t evtIndex = 0 ; 
-    if (GetMixedEvent()) {
-      evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ; 
-    }
-    //Cluster selection, not charged, with photon id and in fiducial cut
-         
-    //Input from second AOD?
-    Int_t input = 0;
-    
-    //Get Momentum vector, 
-    if (input == 0) 
-      calo->GetMomentum(mom,GetVertex(evtIndex)) ;//Assume that come from vertex in straight line
-    
-    //Skip the cluster if it doesn't fit inside the cuts.
-    if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ; 
-    Double_t tof = calo->GetTOF()*1e9;    
-    if(tof < fTimeCutMin || tof > fTimeCutMax) continue;         
-    if(calo->GetNCells() <= fNCellsCut) continue;
-    
-    //Check acceptance selection
-    if(IsFiducialCutOn()){
-      Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;
-      if(! in ) continue ;
-    }
-    
-    //Create AOD for analysis
-    AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);
-    Int_t label = calo->GetLabel();
-    aodph.SetLabel(label);
-    aodph.SetInputFileIndex(input);
-    
-    //Set the indices of the original caloclusters  
-    aodph.SetCaloLabel(calo->GetID(),-1);
-    aodph.SetDetector(fCalorimeter);
-    if(GetDebug() > 1) 
-      printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Min pt cut and fiducial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",aodph.Pt(),aodph.Phi(),aodph.Eta());        
-    
-    //Check Distance to Bad channel, set bit.
-    Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel
-    if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;
-    if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)
-      continue ;
-    
-    if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Bad channel cut passed %4.2f\n",distBad);
-    
-    if     (distBad > fMinDist3) aodph.SetDistToBad(2) ;
-    else if(distBad > fMinDist2) aodph.SetDistToBad(1) ; 
-    else                         aodph.SetDistToBad(0) ;
-    
-    //Skip matched clusters with tracks
-    if(fRejectTrackMatch && IsTrackMatched(calo)) continue ;
-    if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - TrackMatching cut passed \n");
-    
-    //Set PID bits for later selection (AliAnaPi0 for example)
-    //GetPDG already called in SetPIDBits.
-    GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph, GetCaloUtils());
-    if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - PID Bits set \n");                   
-    
-    //Play with the MC stack if available
-    //Check origin of the candidates
-    if(IsDataMC()){
-      aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(), aodph.GetInputFileIndex()));
-      if(GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());
-    }
-    
-    // Check if cluster comes from a conversion in the material in front of the calorimeter
-    // Do invariant mass of all pairs, if mass is close to 0, then it is conversion.
-    
-    if(fCheckConversion && nCaloClusters > 1){
-      Bool_t bConverted = kFALSE;
-      Int_t id2 = -1;
-                 
-      //Check if set previously as converted couple, if so skip its use.
-      if (indexConverted[icalo]) continue;
-                 
-      for(Int_t jcalo = icalo + 1 ; jcalo < nCaloClusters ; jcalo++) {
-        //Check if set previously as converted couple, if so skip its use.
-        if (indexConverted[jcalo]) continue;
-        //printf("Check Conversion indeces %d and %d\n",icalo,jcalo);
-        AliAODCaloCluster * calo2 =  (AliAODCaloCluster*) (pl->At(jcalo));              //Get cluster kinematics
-        Int_t evtIndex2 = 0 ; 
-        if (GetMixedEvent()) {
-          evtIndex2=GetMixedEvent()->EventIndexForCaloCluster(calo2->GetID()) ; 
-        }        
-        calo2->GetMomentum(mom2,GetVertex(evtIndex2));
-        //Check only certain regions
-        Bool_t in2 = kTRUE;
-        if(IsFiducialCutOn()) in2 =  GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;
-        if(!in2) continue;      
-        
-        //Get mass of pair, if small, take this pair.
-        //printf("\t both in calo, mass %f, cut %f\n",(mom+mom2).M(),fMassCut);
-        if((mom+mom2).M() < fMassCut){  
-          bConverted = kTRUE;
-          id2 = calo2->GetID();
-          indexConverted[jcalo]=kTRUE;
-          break;
-        }
-                         
-      }//Mass loop
-                 
-      if(bConverted){ 
-        if(fAddConvertedPairsToAOD){
-          //Create AOD of pair analysis
-          TLorentzVector mpair = mom+mom2;
-          AliAODPWG4Particle aodpair = AliAODPWG4Particle(mpair);
-          aodpair.SetLabel(aodph.GetLabel());
-          aodpair.SetInputFileIndex(input);
-          
-          //printf("Index %d, Id %d\n",icalo, calo->GetID());
-          //Set the indeces of the original caloclusters  
-          aodpair.SetCaloLabel(calo->GetID(),id2);
-          aodpair.SetDetector(fCalorimeter);
-          aodpair.SetPdg(aodph.GetPdg());
-          aodpair.SetTag(aodph.GetTag());
-          
-          //Add AOD with pair object to aod branch
-          AddAODParticle(aodpair);
-          //printf("\t \t both added pair\n");
-        }
-        
-        //Do not add the current calocluster
-        continue;
-      }//converted pair
-    }//check conversion
-         
-    //Add AOD with photon object to aod branch
-    AddAODParticle(aodph);
-    
-  }//loop;
-  delete [] indexConverted;
-       
-  if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD()  End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());  
-  
-}
-
-//__________________________________________________________________
-void  AliAnaShowerParameter::MakeAnalysisFillHistograms() 
-{
-  
-  //Do analysis and fill histograms
-  
-  // Access MC information in stack if requested, check that it exists.        
-  AliStack * stack = 0x0;
-  TParticle * primary = 0x0;   
-  TClonesArray * mcparticles0 = 0x0;
-  //TClonesArray * mcparticles1 = 0x0;
-  AliAODMCParticle * aodprimary = 0x0; 
-  TObjArray * pl = 0x0;
-  Int_t iNumCell=0;
-  
-  //Check if the stack is available when analysing MC data.
-  if(IsDataMC()){
-    
-    if(GetReader()->ReadStack()){
-      stack =  GetMCStack() ;
-      if(!stack) {
-        printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");
-                               abort();
-      }
-      
-    }
-    else if(GetReader()->ReadAODMCParticles()){
-      
-      //Get the list of MC particles
-      mcparticles0 = GetReader()->GetAODMCParticles(0);
-      if(!mcparticles0 && GetDebug() > 0)      {
-        printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n"); 
-      }
-    }// is data and MC
-  }    
-  //Loop on stored AOD photons
-  Int_t naod = GetOutputAODBranch()->GetEntriesFast();
-  fhNClusters->Fill(naod);
-  if(GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);
-  
-  for(Int_t iaod = 0; iaod < naod ; iaod++){
-    AliAODPWG4Particle* ph =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
-    
-    if(ph->GetDetector() != fCalorimeter) continue;
-    
-    if(GetDebug() > 2) 
-      printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - ID Photon: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;
-    
-    //Fill Cluster histograms 
-    Float_t ptcluster  = ph->Pt();
-    Float_t phicluster = ph->Phi();
-    Float_t etacluster = ph->Eta();
-    Float_t ecluster   = ph->E();
-    Float_t lambdaMainCluster = -1; //Can't set their values here, but need a default -1 value in case something is wrong
-    Float_t lambdaSecondCluster = -1;
-    Float_t dispcluster = -1;
-    
-    //Get the list of clusters from the AOD and loop over them to find the onces corresponding to the reconstructed 'photons'.
-    if(fCalorimeter == "PHOS")
-      pl = GetPHOSClusters();
-    else if (fCalorimeter == "EMCAL")
-      pl = GetEMCALClusters();
-    if(pl){
-      //Some values are stored in AliAODCaloCluster objects only; we need to fetch them.
-      for(Int_t icalo = 0; icalo < pl->GetEntriesFast(); icalo++){
-        AliAODCaloCluster * calo =  (AliAODCaloCluster*) (pl->At(icalo));
-        if (calo->GetLabel()==ph->GetLabel()) {  //The Cluster is the right one for this particle
-          lambdaMainCluster = calo->GetM02();   //lambda_0
-          lambdaSecondCluster = calo->GetM20();     //lambda_1
-          dispcluster = calo->GetDispersion();    //Dispersion
-          iNumCell = calo->GetNCells();        
-          if(GetDebug() > 2) 
-            printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Cluster Lambda0: %3.2f, Lambda1: %3.2f, Dispersion: %3.2f, NCells: %d \n",lambdaMainCluster,lambdaSecondCluster,dispcluster,iNumCell) ;
-        }
-      }
-    }
-    
-    //Fill the basic non-MC information about the cluster.
-    fhPtCluster  ->Fill(ecluster);
-    fhPhiCluster ->Fill(ecluster,phicluster);
-    fhEtaCluster ->Fill(ecluster,etacluster);
-    fhDispersionCluster->Fill(ptcluster,dispcluster);
-    fhLambdaCluster ->Fill(ptcluster,lambdaMainCluster,lambdaSecondCluster);
-    fhELambdaCluster ->Fill(ecluster,lambdaMainCluster,lambdaSecondCluster);
-    fhELambdaCellCluster ->Fill(ecluster,lambdaMainCluster,iNumCell);
-    fhNCellCluster->Fill(ecluster,iNumCell);
-    
-    //In the case of an event with several clusters, calculate DeltaEta and DeltaPhi for the cluster pairs.
-    for (Int_t jaod=iaod+1;jaod<naod;jaod++){
-      AliAODPWG4Particle* phSecond = (AliAODPWG4Particle*) (GetOutputAODBranch()->At(jaod));   
-      if(phSecond->GetDetector() != fCalorimeter) continue;
-      fhDeltaPhiClusters->Fill(phicluster-(phSecond->Phi()));
-      fhDeltaEtaClusters->Fill(etacluster-(phSecond->Eta()));  
-    }
-    
-    //Fill the lambda distribution for identified particles
-    if(ph->GetPdg() ==  AliCaloPID::kPhoton){
-      fhDispersionPhoton->Fill(ptcluster,dispcluster);
-      fhLambdaPhoton ->Fill(ptcluster,lambdaMainCluster,lambdaSecondCluster);
-      fhELambdaPhoton ->Fill(ecluster,lambdaMainCluster,lambdaSecondCluster);
-      fhELambdaCellPhoton ->Fill(ecluster,lambdaMainCluster,iNumCell);
-      fhNCellPhoton->Fill(ecluster,iNumCell);
-    }
-    if(ph->GetPdg() ==  AliCaloPID::kPi0){
-      fhDispersionPi0->Fill(ptcluster,dispcluster);
-      fhLambdaPi0 ->Fill(ptcluster,lambdaMainCluster,lambdaSecondCluster);
-      fhELambdaPi0 ->Fill(ecluster,lambdaMainCluster,lambdaSecondCluster);
-      fhELambdaCellPi0 ->Fill(ecluster,lambdaMainCluster,iNumCell);
-      fhNCellPi0->Fill(ecluster,iNumCell);
-    }
-    if(ph->GetPdg() ==  AliCaloPID::kChargedHadron){
-      fhDispersionChargedHadron->Fill(ptcluster,dispcluster);
-      fhLambdaChargedHadron ->Fill(ptcluster,lambdaMainCluster,lambdaSecondCluster);
-      fhELambdaChargedHadron ->Fill(ecluster,lambdaMainCluster,lambdaSecondCluster);
-      fhELambdaCellChargedHadron ->Fill(ecluster,lambdaMainCluster,iNumCell);
-      fhNCellChargedHadron->Fill(ecluster,iNumCell);
-    }
-    
-    //Play with the MC data if available
-    if(IsDataMC()){
-      if(GetReader()->ReadStack() && !stack) return;
-      if(GetReader()->ReadAODMCParticles() && !mcparticles0) return;
-
-      //Get the tag from AliMCAnalysisUtils for PID
-      Int_t tag =ph->GetTag();
-      if ( ph->GetLabel() < 0){
-        if(GetDebug() > 0) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Label is -1; problem in the MC ESD? ");
-        continue;
-      }
-      
-      
-      //Check if the tag matches on of the different particle types and fill the corresponding histograms
-      if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)&&!(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))) //kMCPhoton; making sure that this is not a Pi0 cluster (it should not happen).
-      {
-        fhEMCPhoton  ->Fill(ecluster);
-        fhPhiMCPhoton ->Fill(ecluster,phicluster);
-        fhEtaMCPhoton ->Fill(ecluster,etacluster);
-        fhLambdaMCPhoton ->Fill(ecluster,lambdaMainCluster,iNumCell);
-        Double_t ePhot;
-        //Get the true energy of the photon
-        Int_t iCurrent = ph->GetLabel();
-        TParticle* pCurrent = stack->Particle(iCurrent);
-        ePhot = pCurrent->Energy();      
-        fhPhotTrueE->Fill(ePhot,ecluster,lambdaMainCluster);
-        fhRatioEPhoton->Fill(ecluster/ePhot);
-        fhMCPdg->Fill(ecluster,22);
-        if(GetDebug() > 3) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - MC Cluster identified as kMCPhoton: True E: %3.2f, Reco E: %3.2f, Lambda0:  %3.2f\n",ePhot,ecluster,lambdaMainCluster);
-      }//kMCPhoton
-      
-      if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) //kMCPi0 : single cluster created by 2 photons from the same Pi0 decay
-      {     
-        
-        //Fill the basic information about the Pi0 cluster
-        fhEMCPi0  ->Fill(ecluster);
-        fhPhiMCPi0 ->Fill(ecluster,phicluster);
-        fhEtaMCPi0 ->Fill(ecluster,etacluster);
-        fhLambdaMCPi0 ->Fill(ecluster,lambdaMainCluster,iNumCell);
-        
-        //Check the position of the Pi0: does it come from the vertex or was it created by hadron-material collision?
-        Int_t iCurrent = ph->GetLabel();
-        Double_t ePi0 = 0;
-        TParticle* pCurrent = stack->Particle(iCurrent);
-        while (pCurrent->GetPdgCode()!=111)
-        {
-          iCurrent = pCurrent->GetFirstMother();
-          pCurrent = stack->Particle(iCurrent);
-        }
-        Float_t fDistance = pCurrent->Rho();
-        Float_t fRadius = pCurrent->R();
-        ePi0 = pCurrent->Energy();
-        fhProductionDistance->Fill(ecluster,fDistance);
-        fhProductionRadius->Fill(ecluster,fRadius);
-        //Compare the cluster energy and true energy
-        fhPi0TrueE->Fill(ePi0,ecluster,lambdaMainCluster);
-        fhRatioEPi0->Fill(ecluster/ePi0);
-        fhMCPdg->Fill(ecluster,111);
-        if(GetDebug() > 3) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - MC Cluster identified as kMCPi0: True E: %3.2f, Reco E: %3.2f, Lambda0:  %3.2f\n",ePi0,ecluster,lambdaMainCluster);
-      }
-      if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion))//kMCPion
-      {
-        fhEMCPion  ->Fill(ecluster);
-        fhPhiMCPion ->Fill(ecluster,phicluster);
-        fhEtaMCPion ->Fill(ecluster,etacluster);
-        fhLambdaMCPion ->Fill(ecluster,lambdaMainCluster,iNumCell);
-        
-        Int_t iCurrent = ph->GetLabel();
-        Double_t ePion = 0;
-        TParticle* pCurrent = stack->Particle(iCurrent);
-        while ((TMath::Abs(pCurrent->GetPdgCode())!=211)&&(iCurrent>=0))
-        {
-          if (iCurrent!=pCurrent->GetFirstMother()) iCurrent = pCurrent->GetFirstMother(); //Avoiding infinite loops
-          if (iCurrent>=0) pCurrent = stack->Particle(iCurrent);
-        } 
-        if ((TMath::Abs(pCurrent->GetPdgCode())==211)&&(iCurrent>=0))
-        {
-          ePion = pCurrent->Energy();
-          fhPionTrueE->Fill(ePion,ecluster,lambdaMainCluster);
-          fhRatioEPion->Fill(ecluster/ePion);
-        }
-        fhMCPdg->Fill(ecluster,211);
-        if(GetDebug() > 3) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - MC Cluster identified as kMCPion: True E: %3.2f, Reco E: %3.2f, Lambda0:  %3.2f\n",ePion,ecluster,lambdaMainCluster);
-      }        
-      if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCProton))//kMCProton
-      {
-        fhEMCProton  ->Fill(ecluster);
-        fhPhiMCProton ->Fill(ecluster,phicluster);
-        fhEtaMCProton ->Fill(ecluster,etacluster);
-        fhLambdaMCProton ->Fill(ecluster,lambdaMainCluster,iNumCell);
-        
-        Int_t iCurrent = ph->GetLabel();
-        Double_t eProton = 0;
-        TParticle* pCurrent = stack->Particle(iCurrent);
-        while ((TMath::Abs(pCurrent->GetPdgCode())!=2212)&&(iCurrent>=0))
-        {
-          if (iCurrent!=pCurrent->GetFirstMother()) iCurrent = pCurrent->GetFirstMother(); //Avoiding infinite loops
-          if (iCurrent>=0) pCurrent = stack->Particle(iCurrent);
-        } 
-        if ((TMath::Abs(pCurrent->GetPdgCode())==2212)&&(iCurrent>=0))
-        {
-          eProton = pCurrent->Energy();
-          fhProtonTrueE->Fill(eProton,ecluster,lambdaMainCluster);
-          fhRatioEProton->Fill(ecluster/eProton);
-        }
-        fhMCPdg->Fill(ecluster,2212);
-        if(GetDebug() > 3) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - MC Cluster identified as kMCProton: True E: %3.2f, Reco E: %3.2f, Lambda0:  %3.2f\n",eProton,ecluster,lambdaMainCluster);
-      } 
-      if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCAntiProton))//kMCAntiProton
-      {
-        fhEMCAntiProton  ->Fill(ecluster);
-        fhPhiMCAntiProton ->Fill(ecluster,phicluster);
-        fhEtaMCAntiProton ->Fill(ecluster,etacluster);
-        fhLambdaMCAntiProton ->Fill(ecluster,lambdaMainCluster,iNumCell);
-        
-        Int_t iCurrent = ph->GetLabel();
-        Double_t eAntiProton = 0;
-        TParticle* pCurrent = stack->Particle(iCurrent);
-        while ((TMath::Abs(pCurrent->GetPdgCode())!=-2212)&&(iCurrent>=0))
-        {
-          if (iCurrent!=pCurrent->GetFirstMother()) iCurrent = pCurrent->GetFirstMother(); //Avoiding infinite loops
-          if (iCurrent>=0) pCurrent = stack->Particle(iCurrent);
-        } 
-        if ((TMath::Abs(pCurrent->GetPdgCode())==-2212)&&(iCurrent>=0))
-        {
-          eAntiProton = pCurrent->Energy();
-          fhAntiProtonTrueE->Fill(eAntiProton,ecluster,lambdaMainCluster);
-          fhRatioEAntiProton->Fill(ecluster/eAntiProton);
-        }      
-        fhMCPdg->Fill(ecluster,-2212);
-        if(GetDebug() > 3) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - MC Cluster identified as kMCAntiProton: True E: %3.2f, Reco E: %3.2f, Lambda0:  %3.2f\n",eAntiProton,ecluster,lambdaMainCluster);
-      }
-      if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCNeutron))//kMCNeutron
-      {
-        fhEMCNeutron  ->Fill(ecluster);
-        fhPhiMCNeutron ->Fill(ecluster,phicluster);
-        fhEtaMCNeutron ->Fill(ecluster,etacluster);
-        fhLambdaMCNeutron ->Fill(ecluster,lambdaMainCluster,iNumCell);
-        
-        Int_t iCurrent = ph->GetLabel();
-        Double_t eNeutron = 0;
-        TParticle* pCurrent = stack->Particle(iCurrent);
-        while ((TMath::Abs(pCurrent->GetPdgCode())!=2112)&&(iCurrent>=0))
-        {
-          if (iCurrent!=pCurrent->GetFirstMother()) iCurrent = pCurrent->GetFirstMother(); //Avoiding infinite loops
-          if (iCurrent>=0) pCurrent = stack->Particle(iCurrent);
-        } 
-        if ((TMath::Abs(pCurrent->GetPdgCode())==2112)&&(iCurrent>=0))
-        {
-          eNeutron = pCurrent->Energy();
-          fhNeutronTrueE->Fill(eNeutron,ecluster,lambdaMainCluster);
-          fhRatioENeutron->Fill(ecluster/eNeutron);
-        }
-        fhMCPdg->Fill(ecluster,2112);
-                 if(GetDebug() > 3) 
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - MC Cluster identified as kMCNeutron: True E: %3.2f, Reco E: %3.2f, Lambda0:  %3.2f\n",eNeutron,ecluster,lambdaMainCluster);
-      }
-      if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCEta))//kMCEta
-      {
-        fhEMCEta  ->Fill(ecluster);
-        fhPhiMCEta ->Fill(ecluster,phicluster);
-        fhEtaMCEta ->Fill(ecluster,etacluster);
-        fhLambdaMCEta ->Fill(ecluster,lambdaMainCluster,iNumCell);
-      }      
-      
-      // Access MC information in stack if requested, check that it exists.
-      Int_t label =ph->GetLabel();
-      if(label < 0) {
-        printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** bad label ***:  label %d \n", label);
-        continue;
-      }
-      
-      //Calculate Pi0 decay angles
-      if(stack && (IsDataMC() || (GetReader()->GetDataType() == AliCaloTrackReader::kMC)) ){
-        for(Int_t i=0 ; i<stack->GetNprimary(); i++){
-          TParticle * prim = stack->Particle(i) ;
-          if(prim->GetPdgCode() == 111){
-            TLorentzVector mom1 ;
-            (stack->Particle((prim->GetFirstDaughter())))->Momentum(mom1);
-            TLorentzVector mom2 ;
-            (stack->Particle((prim->GetLastDaughter())))->Momentum(mom2);
-            Float_t fDecayAngle = mom1.Angle(mom2.Vect());
-            fhDecayAngle->Fill(ecluster,fDecayAngle);
-          }
-        }
-      }
-      
-      
-      Float_t eprim   = 0;
-      Float_t ptprim  = 0;
-      if(GetReader()->ReadStack()){
-        
-        if(label >=  stack->GetNtrack()) {
-          if(GetDebug() > 2)  printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** large label ***:  label %d, n tracks %d \n", label, stack->GetNtrack());
-          continue ;
-        }
-        
-        primary = stack->Particle(label);
-        if(!primary){
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** no primary ***:  label %d \n", label);
-          continue;
-        }
-        eprim   = primary->Energy();
-        ptprim  = primary->Pt();
-      } 
-      else if(GetReader()->ReadAODMCParticles()){
-        //Check which is the input
-        if(ph->GetInputFileIndex() == 0){
-          if(!mcparticles0) continue;
-          if(label >=  mcparticles0->GetEntriesFast()) {
-            if(GetDebug() > 2)  printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** large label ***:  label %d, n tracks %d \n",label, mcparticles0->GetEntriesFast());
-            continue ;
-          }
-          //Get the particle
-          aodprimary = (AliAODMCParticle*) mcparticles0->At(label);
-          
-        }
-        //        else {//Second input
-        //          if(!mcparticles1) continue;
-        //          if(label >=  mcparticles1->GetEntriesFast()) {
-        //            if(GetDebug() > 2)  printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** large label ***:  label %d, n tracks %d \n",label, mcparticles1->GetEntriesFast());
-        //            continue ;
-        //          }
-        //          //Get the particle
-        //          aodprimary = (AliAODMCParticle*) mcparticles1->At(label); 
-        //        }//second input
-        
-        if(!aodprimary){
-          printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** no primary ***:  label %d \n", label);
-          continue;
-        }
-        
-        eprim   = aodprimary->E();
-        ptprim  = aodprimary->Pt();
-      }
-      //Write down MC truth information concerning all clusters
-      fh2E     ->Fill(ecluster, eprim);
-      fh2Pt    ->Fill(ptcluster, ptprim);     
-      fhDeltaE ->Fill(eprim-ecluster);
-      fhDeltaPt->Fill(ptprim-ptcluster);     
-      if(eprim > 0)  fhRatioE  ->Fill(ecluster/eprim);
-      if(ptprim > 0) fhRatioPt ->Fill(ptcluster/ptprim);               
-    }//Histograms with MC
-  }// aod loop 
-}
-
-
-//__________________________________________________________________
-void AliAnaShowerParameter::Print(const Option_t * opt) const
-{
-  //Print some relevant parameters set for the analysis
-  
-  if(! opt)
-    return;
-  
-  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;
-  AliAnaPartCorrBaseClass::Print(" ");
-  
-  printf("Calorimeter            =     %s\n", fCalorimeter.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);
-  printf("Reject clusters with a track matched = %d\n",fRejectTrackMatch);
-  printf("Check Pair Conversion                = %d\n",fCheckConversion);
-  printf("Add conversion pair to AOD           = %d\n",fAddConvertedPairsToAOD);
-  printf("Conversion pair mass cut             = %f\n",fMassCut);
-  printf("Time Cut: %3.1f < TOF  < %3.1f\n", fTimeCutMin, fTimeCutMax);
-  printf("Number of cells in cluster is        > %f \n", fNCellsCut);
-  printf("Number of clusters in envent is        > %d \n", fNumClusters);
-  printf("    \n") ;
-  
-} 
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes hereby granted      *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+/* $Id: AliAnaPhoton.cxx 28688 2008-09-11 15:04:07Z gconesab $ */\r
+\r
+//_________________________________________________________________________\r
+//\r
+// Class cloned from AliAnaPhoton, main aim is shower shape studies\r
+// \r
+// \r
+//\r
+//-- Author: Jocelyn Mlynarz (WSU) and Gustavo Conesa (LPSC-CNRS)\r
+//////////////////////////////////////////////////////////////////////////////\r
+\r
+\r
+// --- ROOT system --- \r
+#include <TH3F.h>\r
+#include <TH2F.h>\r
+#include <TClonesArray.h>\r
+//#include <TObjString.h>\r
+#include <Riostream.h>\r
+#include "TParticle.h"\r
+//#include <fstream>\r
+\r
+// --- Analysis system --- \r
+#include "AliAnaShowerParameter.h" \r
+#include "AliCaloTrackReader.h"\r
+#include "AliStack.h"\r
+#include "AliCaloPID.h"\r
+#include "AliMCAnalysisUtils.h"\r
+#include "AliFiducialCut.h"\r
+#include "AliAODCaloCluster.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliAnalysisTaskParticleCorrelation.h"\r
+#include "AliEMCALGeoUtils.h"\r
+#include "AliAODEvent.h"\r
+\r
+\r
+ClassImp(AliAnaShowerParameter)\r
+\r
+//____________________________________________________________________________\r
+AliAnaShowerParameter::AliAnaShowerParameter() : \r
+AliAnaPartCorrBaseClass(), fCalorimeter(""), fNCellsCutMin(0),\r
+fNCellsCutMax(0), fLambdaCut(0), fTimeCutMin(-1), fTimeCutMax(9999999),\r
+fhNClusters(0), fhNCellCluster(0), fhEtaPhiPtCluster(0), \r
+fhLambdaPtCluster(0), \r
+\r
+//MC\r
+\r
+fhLambdaPtPhoton(0), fhLambdaPtPi0(0), fhLambdaPtPion(0), fhPtTruthPi0(0)\r
+\r
+{\r
+  //default ctor\r
+  \r
+  //Initialize parameters\r
+  InitParameters();\r
+  \r
+}\r
+\r
+//____________________________________________________________________________\r
+AliAnaShowerParameter::~AliAnaShowerParameter() \r
+{\r
+  //dtor\r
+  \r
+}\r
+\r
+//________________________________________________________________________\r
+TObjString *  AliAnaShowerParameter::GetAnalysisCuts()\r
+{      \r
+  //Save parameters used for analysis\r
+  TString parList ; //this will be list of parameters used for this analysis.\r
+  const Int_t buffersize = 255;\r
+  char onePar[buffersize] ;\r
+  \r
+  snprintf(onePar,buffersize,"--- AliAnaShowerParameter ---\n") ;\r
+  parList+=onePar ;    \r
+  snprintf(onePar,buffersize,"Calorimeter: %s\n",fCalorimeter.Data()) ;\r
+  parList+=onePar ;\r
+  snprintf(onePar,buffersize,"fNCellsCutMin: (Cut in the minimum number of cells) %e\n",fNCellsCutMin) ;\r
+  parList+=onePar ;  \r
+  snprintf(onePar,buffersize,"fNCellsCutMax: (Cut in the maximum number of cells) %e\n",fNCellsCutMax) ;\r
+  parList+=onePar ;\r
+  snprintf(onePar,buffersize,"fLambdaCut: (Cut in the minimum lambda) %e\n",fLambdaCut) ;\r
+  parList+=onePar ;\r
+  \r
+  //Get parameters set in base class.\r
+  parList += GetBaseParametersList() ;\r
+  \r
+  //Get parameters set in PID class.\r
+  parList += GetCaloPID()->GetPIDParametersList() ;\r
+  \r
+  //Get parameters set in FiducialCut class (not available yet)\r
+  //parlist += GetFidCut()->GetFidCutParametersList() \r
+  \r
+  return new TObjString(parList) ;\r
+}\r
+\r
+\r
+//________________________________________________________________________\r
+TList *  AliAnaShowerParameter::GetCreateOutputObjects()\r
+{  \r
+  // Create histograms to be saved in output file and \r
+  // store them in outputContainer\r
+  TList * outputContainer = new TList() ; \r
+  outputContainer->SetName("PhotonHistos") ; \r
+       \r
+  Int_t nptbins  = GetHistoPtBins();\r
+  Int_t nphibins = GetHistoPhiBins();\r
+  Int_t netabins = GetHistoEtaBins();\r
+  Float_t ptmax  = GetHistoPtMax();\r
+  Float_t phimax = GetHistoPhiMax();\r
+  Float_t etamax = GetHistoEtaMax();\r
+  Float_t ptmin  = GetHistoPtMin();\r
+  Float_t phimin = GetHistoPhiMin();\r
+  Float_t etamin = GetHistoEtaMin();   \r
+  \r
+  //General non-MC Cluster histograms\r
+  fhNClusters = new TH1F ("hNClusters","NClusters",21,-0.5,20.5); \r
+  fhNClusters->SetXTitle("N_{Clusters}");\r
+  outputContainer->Add(fhNClusters) ;\r
+  \r
+  fhNCellCluster = new TH2F ("hNCellCluster","Number of cell per cluster",nptbins,ptmin,ptmax,21,-0.5,20.5); \r
+  fhNCellCluster->SetYTitle("N_{Cells}");  \r
+  fhNCellCluster->SetXTitle("p_{T}");\r
+  outputContainer->Add(fhNCellCluster) ;\r
+  \r
+  fhEtaPhiPtCluster = new TH3F\r
+  ("hEtaPhiPtCluster","#phi_{Cluster} and #eta_{Cluster}",nptbins,ptmin,ptmax,nphibins,phimin,phimax,netabins,etamin,etamax); \r
+  fhEtaPhiPtCluster->SetZTitle("#eta");  \r
+  fhEtaPhiPtCluster->SetYTitle("#phi");\r
+  fhEtaPhiPtCluster->SetXTitle("pT_{Cluster} (GeV/c)");\r
+  outputContainer->Add(fhEtaPhiPtCluster) ; \r
+  \r
+  fhLambdaPtCluster  = new TH2F\r
+  ("hLambdaCluster","#lambda_{Cluster}",nptbins,ptmin,ptmax,300,0,3); \r
+  fhLambdaPtCluster->SetYTitle("#lambda_{0}^{2}");\r
+  fhLambdaPtCluster->SetXTitle("p_{T Cluster} (GeV/c)");\r
+  outputContainer->Add(fhLambdaPtCluster) ;\r
+  \r
+  if(IsDataMC()){\r
+    \r
+    fhLambdaPtPhoton  = new TH2F\r
+    ("hLambdaPtPhoton","#lambda_{#gamma}",nptbins,ptmin,ptmax,200,0,2); \r
+    fhLambdaPtPhoton->SetYTitle("#lambda_{0}^{2}");\r
+    fhLambdaPtPhoton->SetXTitle("pT_{#gamma, Reco} (GeV)");\r
+    outputContainer->Add(fhLambdaPtPhoton) ;\r
+    \r
+    fhLambdaPtPi0  = new TH2F\r
+    ("hLambdaPtPi0","#lambda_{#pi^{0}}",nptbins,ptmin,ptmax,200,0,2); \r
+    fhLambdaPtPi0->SetYTitle("#lambda_{0}^{2}");\r
+    fhLambdaPtPi0->SetXTitle("pT_{#pi^{0}, Reco} (GeV)");\r
+    outputContainer->Add(fhLambdaPtPi0) ;\r
+    \r
+    fhLambdaPtPion  = new TH2F\r
+    ("hLambdaPtPion","#lambda_{#pi^{+}}",nptbins,ptmin,ptmax,200,0,2); \r
+    fhLambdaPtPion->SetYTitle("#lambda_{0}^{2}");\r
+    fhLambdaPtPion->SetXTitle("pT_{#pi^{+}, Reco} (GeV)");\r
+    outputContainer->Add(fhLambdaPtPion) ;\r
+   \r
+    fhPtTruthPi0  = new TH1D("hPtTruthPi0","#pi^{0} MC truth pT",nptbins,ptmin,ptmax) ;\r
+    fhPtTruthPi0->SetXTitle("pT_{#pi^{0}, Reco} (GeV)");\r
+    outputContainer->Add(fhPtTruthPi0) ;\r
+    \r
+  }//Histos with MC\r
+  \r
+  return outputContainer ;\r
+  \r
+}\r
+\r
+//____________________________________________________________________________\r
+void AliAnaShowerParameter::Init()\r
+{\r
+  \r
+  //Init\r
+  //Do some checks\r
+  if(fCalorimeter == "PHOS" && !GetReader()->IsPHOSSwitchedOn() && NewOutputAOD()){\r
+    printf("AliAnaShowerParameter::Init() - !!STOP: You want to use PHOS in analysis but it is not read!! \n!!Check the configuration file!!\n");\r
+    abort();\r
+  }\r
+  else  if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn() && NewOutputAOD()){\r
+    printf("AliAnaShowerParameter::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!! \n!!Check the configuration file!!\n");\r
+    abort();\r
+  }\r
+  \r
+}\r
+\r
+//____________________________________________________________________________\r
+void AliAnaShowerParameter::InitParameters()\r
+{\r
+  \r
+  //Initialize the parameters of the analysis.\r
+  AddToHistogramsName("AnaLambda_");\r
+  \r
+  fCalorimeter = "EMCAL" ;\r
+       \r
+  fTimeCutMin  = -1;\r
+  fTimeCutMax  = 9999999;\r
+  fNCellsCutMin = 0 ;\r
+  fNCellsCutMax = 0 ;\r
+  fLambdaCut = 0.01 ;\r
+       \r
+}\r
+\r
+//__________________________________________________________________\r
+void  AliAnaShowerParameter::MakeAnalysisFillAOD() \r
+{\r
+  /*//Do analysis and fill aods\r
+  //Search for photons in fCalorimeter \r
+  \r
+  //Get vertex for photon momentum calculation\r
+  \r
+  for (Int_t iev = 0; iev < GetNMixedEvent(); iev++) {\r
+    if (!GetMixedEvent()) \r
+      GetReader()->GetVertex(GetVertex(iev));\r
+    else \r
+      GetMixedEvent()->GetVertexOfEvent(iev)->GetXYZ(GetVertex(iev)); \r
+  } \r
+  \r
+  //Select the Calorimeter of the photon\r
+  TObjArray * pl = 0x0; \r
+  if(fCalorimeter == "PHOS")\r
+    pl = GetPHOSClusters();\r
+  else if (fCalorimeter == "EMCAL")\r
+    pl = GetEMCALClusters();\r
+  \r
+  if(!pl){\r
+    printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Careful cluster array NULL!!\n");\r
+    return;\r
+  }\r
+  \r
+  //Fill AODCaloClusters and AODParticle with PHOS/EMCAL aods\r
+  TLorentzVector mom, mom2 ; \r
+  Int_t nCaloClusters = pl->GetEntriesFast();   \r
+  //Cut on the number of clusters in the event.\r
+  if ((fNumClusters !=-1) && (nCaloClusters != fNumClusters)) return;\r
+  Bool_t * indexConverted = new Bool_t[nCaloClusters];\r
+  for (Int_t i = 0; i < nCaloClusters; i++) \r
+    indexConverted[i] = kFALSE;\r
+  \r
+  for(Int_t icalo = 0; icalo < nCaloClusters; icalo++){    \r
+    \r
+    AliAODCaloCluster * calo =  (AliAODCaloCluster*) (pl->At(icalo));  \r
+    Int_t evtIndex = 0 ; \r
+    if (GetMixedEvent()) {\r
+      evtIndex=GetMixedEvent()->EventIndexForCaloCluster(calo->GetID()) ; \r
+    }\r
+    //Cluster selection, not charged, with photon id and in fiducial cut\r
+         \r
+    //Input from second AOD?\r
+    Int_t input = 0;\r
+    \r
+    //Get Momentum vector, \r
+    if (input == 0) \r
+      calo->GetMomentum(mom,GetVertex(evtIndex)) ;//Assume that come from vertex in straight line\r
+    \r
+    //Skip the cluster if it doesn't fit inside the cuts.\r
+    if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ; \r
+    Double_t tof = calo->GetTOF()*1e9;    \r
+    if(tof < fTimeCutMin || tof > fTimeCutMax) continue;         \r
+    if(calo->GetNCells() <= fNCellsCut) continue;\r
+    \r
+    //Check acceptance selection\r
+    if(IsFiducialCutOn()){\r
+      Bool_t in = GetFiducialCut()->IsInFiducialCut(mom,fCalorimeter) ;\r
+      if(! in ) continue ;\r
+    }\r
+    \r
+    //Create AOD for analysis\r
+    AliAODPWG4Particle aodph = AliAODPWG4Particle(mom);\r
+    Int_t label = calo->GetLabel();\r
+    aodph.SetLabel(label);\r
+    aodph.SetInputFileIndex(input);\r
+    \r
+    //Set the indices of the original caloclusters  \r
+    aodph.SetCaloLabel(calo->GetID(),-1);\r
+    aodph.SetDetector(fCalorimeter);\r
+    if(GetDebug() > 1) \r
+      printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Min pt cut and fiducial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",aodph.Pt(),aodph.Phi(),aodph.Eta());        \r
+    \r
+    //Check Distance to Bad channel, set bit.\r
+    Double_t distBad=calo->GetDistanceToBadChannel() ; //Distance to bad channel\r
+    if(distBad < 0.) distBad=9999. ; //workout strange convension dist = -1. ;\r
+    if(distBad < fMinDist) //In bad channel (PHOS cristal size 2.2x2.2 cm)\r
+      continue ;\r
+    \r
+    if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Bad channel cut passed %4.2f\n",distBad);\r
+    \r
+    if     (distBad > fMinDist3) aodph.SetDistToBad(2) ;\r
+    else if(distBad > fMinDist2) aodph.SetDistToBad(1) ; \r
+    else                         aodph.SetDistToBad(0) ;\r
+    \r
+    //Skip matched clusters with tracks\r
+    if(fRejectTrackMatch && IsTrackMatched(calo)) continue ;\r
+    if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - TrackMatching cut passed \n");\r
+    \r
+    //Set PID bits for later selection (AliAnaPi0 for example)\r
+    //GetPDG already called in SetPIDBits.\r
+    GetCaloPID()->SetPIDBits(fCalorimeter,calo,&aodph, GetCaloUtils());\r
+    if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - PID Bits set \n");                   \r
+    \r
+    //Play with the MC stack if available\r
+    //Check origin of the candidates\r
+    if(IsDataMC()){\r
+      aodph.SetTag(GetMCAnalysisUtils()->CheckOrigin(calo->GetLabels(),calo->GetNLabels(),GetReader(), aodph.GetInputFileIndex()));\r
+      if(GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillAOD() - Origin of candidate, bit map %d\n",aodph.GetTag());\r
+    }\r
+    \r
+    // Check if cluster comes from a conversion in the material in front of the calorimeter\r
+    // Do invariant mass of all pairs, if mass is close to 0, then it is conversion.\r
+    \r
+    if(fCheckConversion && nCaloClusters > 1){\r
+      Bool_t bConverted = kFALSE;\r
+      Int_t id2 = -1;\r
+                 \r
+      //Check if set previously as converted couple, if so skip its use.\r
+      if (indexConverted[icalo]) continue;\r
+                 \r
+      for(Int_t jcalo = icalo + 1 ; jcalo < nCaloClusters ; jcalo++) {\r
+        //Check if set previously as converted couple, if so skip its use.\r
+        if (indexConverted[jcalo]) continue;\r
+        //printf("Check Conversion indeces %d and %d\n",icalo,jcalo);\r
+        AliAODCaloCluster * calo2 =  (AliAODCaloCluster*) (pl->At(jcalo));              //Get cluster kinematics\r
+        Int_t evtIndex2 = 0 ; \r
+        if (GetMixedEvent()) {\r
+          evtIndex2=GetMixedEvent()->EventIndexForCaloCluster(calo2->GetID()) ; \r
+        }        \r
+        calo2->GetMomentum(mom2,GetVertex(evtIndex2));\r
+        //Check only certain regions\r
+        Bool_t in2 = kTRUE;\r
+        if(IsFiducialCutOn()) in2 =  GetFiducialCut()->IsInFiducialCut(mom2,fCalorimeter) ;\r
+        if(!in2) continue;      \r
+        \r
+        //Get mass of pair, if small, take this pair.\r
+        //printf("\t both in calo, mass %f, cut %f\n",(mom+mom2).M(),fMassCut);\r
+        if((mom+mom2).M() < fMassCut){  \r
+          bConverted = kTRUE;\r
+          id2 = calo2->GetID();\r
+          indexConverted[jcalo]=kTRUE;\r
+          break;\r
+        }\r
+                         \r
+      }//Mass loop\r
+                 \r
+      if(bConverted){ \r
+        if(fAddConvertedPairsToAOD){\r
+          //Create AOD of pair analysis\r
+          TLorentzVector mpair = mom+mom2;\r
+          AliAODPWG4Particle aodpair = AliAODPWG4Particle(mpair);\r
+          aodpair.SetLabel(aodph.GetLabel());\r
+          aodpair.SetInputFileIndex(input);\r
+          \r
+          //printf("Index %d, Id %d\n",icalo, calo->GetID());\r
+          //Set the indeces of the original caloclusters  \r
+          aodpair.SetCaloLabel(calo->GetID(),id2);\r
+          aodpair.SetDetector(fCalorimeter);\r
+          aodpair.SetPdg(aodph.GetPdg());\r
+          aodpair.SetTag(aodph.GetTag());\r
+          \r
+          //Add AOD with pair object to aod branch\r
+          AddAODParticle(aodpair);\r
+          //printf("\t \t both added pair\n");\r
+        }\r
+        \r
+        //Do not add the current calocluster\r
+        continue;\r
+      }//converted pair\r
+    }//check conversion\r
+         \r
+    //Add AOD with photon object to aod branch\r
+    AddAODParticle(aodph);\r
+    \r
+  }//loop;\r
+  delete [] indexConverted;\r
+       \r
+  if(GetDebug() > 1) printf("AliAnaShowerParameter::MakeAnalysisFillAOD()  End fill AODs, with %d entries \n",GetOutputAODBranch()->GetEntriesFast());  \r
+  \r
+*/\r
+}\r
+//__________________________________________________________________\r
+void  AliAnaShowerParameter::MakeAnalysisFillHistograms() \r
+{\r
+  \r
+  //Do analysis and fill histograms\r
+  if ( GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Starting analysis.\n");\r
+    \r
+  // Access MC information in stack if requested, check that it exists.        \r
+  AliStack * stack = 0x0;\r
+  //TParticle * primary = 0x0;   \r
+  TClonesArray * mcparticles0 = 0x0;\r
+  //AliAODMCParticle * aodprimary = 0x0; \r
+  //TObjArray * pl = 0x0;\r
+  Int_t NClusters = 0 ;\r
+    TLorentzVector momCluster ;\r
+\r
+  \r
+  //Check if the stack is available when analysing MC data.\r
+  if(IsDataMC()){\r
+    \r
+    if(GetReader()->ReadStack()){\r
+      stack =  GetMCStack() ;\r
+      if(!stack) {\r
+        printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Stack not available, is the MC handler called? STOP\n");\r
+                               abort();\r
+      }\r
+      \r
+    }\r
+    else if(GetReader()->ReadAODMCParticles()){\r
+      \r
+      //Get the list of MC particles\r
+      mcparticles0 = GetReader()->GetAODMCParticles(0);\r
+      if(!mcparticles0 && GetDebug() > 0)      {\r
+        printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() -  Standard MCParticles not available !\n"); \r
+      }\r
+    }// is data and MC\r
+  }    \r
+  //Loop on stored AOD photons\r
+\r
+  TClonesArray *  clustArray = 0x0 ; \r
+  clustArray = GetAODCaloClusters() ;\r
+  NClusters = clustArray->GetEntriesFast() ;    \r
+  fhNClusters->Fill(NClusters) ;\r
+\r
+  if(GetDebug() > 0) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - aod branch entries %d\n", NClusters);\r
+  \r
+  for(Int_t iCluster = 0 ; iCluster < NClusters ; iCluster++){    \r
+    Int_t input = 0 ;\r
+\r
+    AliAODCaloCluster * caloCluster =  (AliAODCaloCluster*) (clustArray->At(iCluster)) ;\r
+    caloCluster->GetMomentum(momCluster,GetVertex(0)) ;\r
+    AliAODPWG4Particle * aodCluster = new AliAODPWG4Particle(momCluster) ;\r
+    Int_t labelCluster = caloCluster->GetLabel() ;\r
+    aodCluster->SetLabel(labelCluster) ;;\r
+    aodCluster->SetInputFileIndex(input) ;\r
+    aodCluster->SetCaloLabel(caloCluster->GetID(),-1) ;\r
+    aodCluster->SetDetector(fCalorimeter) ;\r
+    aodCluster->SetTag(GetMCAnalysisUtils()->CheckOrigin(caloCluster->GetLabels(),caloCluster->GetNLabels(),GetReader(), aodCluster->GetInputFileIndex())) ;\r
+\r
+    if (GetDebug() > 3) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Cluster particle label = %d\n",labelCluster) ;\r
+    \r
+    //Fill Cluster histograms \r
+    Float_t ptcluster  = aodCluster->Pt() ;\r
+    Float_t phicluster = aodCluster->Phi() ;\r
+    Float_t etacluster = aodCluster->Eta() ;\r
+    Float_t lambdaMainCluster = caloCluster->GetM02() ;\r
+    Int_t iNumCell = caloCluster->GetNCells() ;        \r
+    \r
+    if (iNumCell>=fNCellsCutMin&&iNumCell<=fNCellsCutMax&&lambdaMainCluster>=fLambdaCut){\r
+      \r
+      //Fill the basic non-MC information about the cluster.\r
+      fhNCellCluster->Fill(ptcluster,iNumCell) ;\r
+      fhEtaPhiPtCluster->Fill(ptcluster,phicluster,etacluster) ;\r
+      fhLambdaPtCluster->Fill(ptcluster,lambdaMainCluster) ;\r
+      \r
+    //Play with the MC data if available\r
+      if(IsDataMC()){\r
+       \r
+       if(GetReader()->ReadStack() && !stack) return;\r
+       if(GetReader()->ReadAODMCParticles() && !mcparticles0) return;\r
+       \r
+       //Get the tag from AliMCAnalysisUtils for PID\r
+       Int_t tag = aodCluster->GetTag();\r
+       if (GetDebug() > 3) printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Cluster particle tag= %d\n",tag) ;\r
+       if ( aodCluster->GetLabel() < 0){\r
+         if(GetDebug() > 0) \r
+           printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() - Label is -1; problem in the MC ESD? ");\r
+         continue;\r
+       }\r
+       \r
+       //Check if the tag matches one of the different particle types and fill the corresponding histograms\r
+       if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPhoton)&&!(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0))) {\r
+         fhLambdaPtPhoton->Fill(ptcluster,lambdaMainCluster) ;\r
+       }//kMCPhoton\r
+       \r
+       if(GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPi0)) {     \r
+         fhLambdaPtPi0->Fill(ptcluster,lambdaMainCluster) ;        \r
+       }\r
+       if( GetMCAnalysisUtils()->CheckTagBit(tag,AliMCAnalysisUtils::kMCPion)) {\r
+         fhLambdaPtPion->Fill(ptcluster,lambdaMainCluster) ;\r
+       }\r
+       \r
+       // Access MC information in stack if requested, check that it exists.\r
+       Int_t label = aodCluster->GetLabel();\r
+       if(label < 0) {\r
+         printf("AliAnaShowerParameter::MakeAnalysisFillHistograms() *** bad label ***:  label %d \n", label);\r
+         continue;\r
+       }\r
+      }\r
+    }\r
+  }\r
+  \r
+  Float_t fPtGener = 0 ;\r
+  if(IsDataMC()){\r
+    TParticle* pGener = stack->Particle(0) ;\r
+    fPtGener = pGener->Pt() ;\r
+  }\r
+\r
+  \r
+}\r
+\r
+\r
+//__________________________________________________________________\r
+void AliAnaShowerParameter::Print(const Option_t * opt) const\r
+{\r
+  //Print some relevant parameters set for the analysis\r
+  \r
+  if(! opt)\r
+    return;\r
+  \r
+  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
+  AliAnaPartCorrBaseClass::Print(" ");\r
+  printf("Calorimeter            =     %s\n", fCalorimeter.Data()) ;\r
+  printf("Min number of cells in cluster is        > %f \n", fNCellsCutMin);\r
+  printf("Max number of cells in cluster is        > %f \n", fNCellsCutMax);\r
+  printf("Min lambda of cluster is        > %f \n", fLambdaCut);\r
+  printf("    \n") ;\r
+  \r
+} \r
index 7a2604d..ecb15da 100644 (file)
-#ifndef ALIANASHOWERPARAMETER_H
-#define ALIANASHOWERPARAMETER_H
-/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice     */
-/* $Id: AliAnaShowerParameter.h 27413 2008-07-18 13:28:12Z gconesab $ */
-
-//_________________________________________________________________________
-//
-// Class cloned from AliAnaPhoton, main aim is shower shape studies
-// 
-// 
-//
-//-- Author: Jocelyn Mlynarz (WSU) and Gustavo Conesa (LPSC)
-
-// --- ROOT system ---
-class TH3F;
-class TH2F ;
-class TH1F;
-class TString ;
-
-// --- ANALYSIS system ---
-#include "AliAnaPartCorrBaseClass.h"
-//#include "AliStack.h"
-//#include "TParticle.h"
-class AliStack;
-class TParticle;
-
-class TList ;
-class AliEMCALGeoUtils;
-class AliAnaShowerParameter : public AliAnaPartCorrBaseClass {
-
-public: 
-
-  AliAnaShowerParameter() ; // default ctor
-  AliAnaShowerParameter(const AliAnaShowerParameter & g) ; // cpy ctor
-  AliAnaShowerParameter & operator = (const AliAnaShowerParameter & g) ;//cpy assignment
-  virtual ~AliAnaShowerParameter() ; //virtual dtor
-  
-  TList *  GetCreateOutputObjects();
-
-  void Init();
-
-  TObjString* GetAnalysisCuts();
-
-  void MakeAnalysisFillAOD()  ;
-    
-  void MakeAnalysisFillHistograms() ; 
-  
-  void Print(const Option_t * opt)const;
-  
-  TString GetCalorimeter()   const {return fCalorimeter ; }
-  void SetCalorimeter(TString det)    {fCalorimeter = det ; }
-
-  Bool_t IsTrackMatchRejectionOn()   const {return fRejectTrackMatch ; }
-  void SwitchOnTrackMatchRejection()  {fRejectTrackMatch = kTRUE ; }
-  void SwitchOffTrackMatchRejection() {fRejectTrackMatch = kFALSE ; }  
-
-  Bool_t IsCheckConversionOn()   const {return fCheckConversion ; }
-  void SwitchOnConversionChecker()  {fCheckConversion = kTRUE ; }
-  void SwitchOffConversionChecker() {fCheckConversion = kFALSE ; }  
-       
-  Bool_t AreConvertedPairsInAOD()   const {return fAddConvertedPairsToAOD ; }
-  void SwitchOnAdditionConvertedPairsToAOD()  {fAddConvertedPairsToAOD = kTRUE ; }
-  void SwitchOffAdditionConvertedPairsToAOD() {fAddConvertedPairsToAOD = kFALSE ; }  
-       
-  void InitParameters();
-  void SetMinDistanceToBadChannel(Float_t m1, Float_t m2, Float_t m3) {
-    fMinDist = m1;
-    fMinDist2 = m2;
-    fMinDist3 = m3;
-  }
-       
-  Float_t GetMassCut()    const {return fMassCut ; }
-  void SetMassCut(Float_t m)    {fMassCut = m ; }
-  void SetNClusterCut(Int_t nc)   {fNumClusters = nc;}
-  Int_t GetNClusterCut()   {return fNumClusters;}
-  void SetTimeCut(Double_t min, Double_t max) {fTimeCutMin = min; fTimeCutMax = max;}
-  Double_t GetTimeCutMin() const {return fTimeCutMin;}
-  Double_t GetTimeCutMax() const {return fTimeCutMax;} 
-       
-  private:
-  TString fCalorimeter ; // Calorimeter where the gamma is searched;
-  Float_t fMinDist ;     // Minimal distance to bad channel to accept cluster
-  Float_t fMinDist2;     // Cuts on Minimal distance to study acceptance evaluation
-  Float_t fMinDist3;     // One more cut on distance used for acceptance-efficiency study
-  Bool_t  fRejectTrackMatch ;      //If PID on, reject clusters which have an associated TPC track
-  Bool_t  fCheckConversion;        // Combine pairs of clusters with mass close to 0
-  Bool_t  fAddConvertedPairsToAOD; // Put Converted pairs in AOD
-  Float_t fMassCut;                // Mass cut for the conversion pairs selection
-  Float_t fNCellsCut;
-  Double_t fTimeCutMin  ;    // Remove clusters/cells with time smaller than this value, in ns
-  Double_t fTimeCutMax  ;    // Remove clusters/cells with time larger than this value, in ns
-  AliEMCALGeoUtils * fEMCALGeo;
-  Int_t fNumClusters;        //Cut that selects events with a specific number of clusters
-
-  //Histograms   
-  TH1F * fhNClusters  ; //! Number of clusters per event.
-  TH2F * fhNCellCluster;//! Number of cells per cluster
-  TH1F * fhPtCluster   ; //! Number of clusters vs transerse momentum 
-  TH2F * fhPhiCluster  ; //! Azimuthal angle of clusters vs transerse momentum 
-  TH2F * fhEtaCluster  ; //! Pseudorapidity of clusters vs transerse momentum
-  TH1F * fhDeltaPhiClusters  ; //! Delta phi of cluster pairs
-  TH1F * fhDeltaEtaClusters  ; //! Delta eta of cluster pairs 
-  TH3F * fhLambdaCluster  ; //! Shower parameters of clusters vs transerse momentum
-  TH2F * fhDispersionCluster  ; //! Dispersion of the clusters
-  TH3F * fhELambdaCluster  ; //! Shower parameters of clusters vs Energy
-  TH3F * fhELambdaCellCluster  ; //! Shower parameters of clusters vs Energy
-
-  TH2F * fhNCellPhoton;//! Number of cells per photon
-  TH3F * fhLambdaPhoton  ; //! Shower parameters of photons vs transerse momentum
-  TH2F * fhDispersionPhoton  ; //! Dispersion of the photons
-  TH3F * fhELambdaPhoton  ; //! Shower parameters of photons vs Energy
-  TH3F * fhELambdaCellPhoton  ; //! Shower parameters of photons vs Energy
-
-  TH2F * fhNCellPi0;//! Number of cells per neutral pion
-  TH3F * fhLambdaPi0  ; //! Shower parameters of neutral pions vs transerse momentum
-  TH2F * fhDispersionPi0  ; //! Dispersion of the neutral pions
-  TH3F * fhELambdaPi0  ; //! Shower parameters of neutral pions vs Energy
-  TH3F * fhELambdaCellPi0  ; //! Shower parameters of neutral pions vs Energy
-
-  TH2F * fhNCellChargedHadron;//! Number of cells per charged hadron
-  TH3F * fhLambdaChargedHadron  ; //! Shower parameters of charged hadrons vs transerse momentum
-  TH2F * fhDispersionChargedHadron  ; //! Dispersion of the charged hadrons
-  TH3F * fhELambdaChargedHadron  ; //! Shower parameters of charged hadrons vs Energy
-  TH3F * fhELambdaCellChargedHadron  ; //! Shower parameters of charged hadrons vs Energy
-
-  //MC
-  TH1F * fhDeltaE  ; //! MC-Reco E distribution      
-  TH1F * fhDeltaPt ; //! MC-Reco pT distribution
-  TH1F * fhRatioE  ; //! Reco/MC E distribution      
-  TH1F * fhRatioPt ; //! Reco/MC pT distribution
-  TH2F * fh2E  ; //! E distribution, Reco vs MC
-  TH2F * fh2Pt ; //! pT distribution, Reco vs MC
-  TH2F * fhMCPdg; //! Complete list of PDG Codes.
-  
-  TH1F * fhEMCPhoton;   //! Number of identified gamma 
-  TH2F * fhPhiMCPhoton;  //! Phi of identified gamma
-  TH2F * fhEtaMCPhoton;  //! eta of identified gamma   
-  TH3F * fhLambdaMCPhoton  ; //! Shower parameters of MC photons vs transerse momentum 
-  TH3F * fhPhotTrueE;  // MC truth E vs Recons E vs. Lambda of the cluster for MC photons
-  TH1F * fhRatioEPhoton  ; //! Reco/MC E distribution for photons      
-  
-  TH1F * fhEMCPi0;   //! Number of identified  (single shower) Pi0 
-  TH2F * fhPhiMCPi0;  //! Phi of identified Pi0
-  TH2F * fhEtaMCPi0;  //! eta of identified Pi0        
-  TH3F * fhLambdaMCPi0  ; //! Shower parameters of MC Pi0 vs transerse momentum   
-  TH3F * fhPi0TrueE;  // MC truth E vs Recons E vs. Lambda of the cluster for MC Pi0s
-  TH2F * fhProductionDistance; //! Distance from beam to production of the Pi0
-  TH2F * fhProductionRadius; //! R from beam to Pi0
-  TH2F * fhDecayAngle;//! Decay angle of the Pi0
-  TH1F * fhRatioEPi0  ; //! Reco/MC E distribution for Pi0s
-
-  TH1F * fhEMCPion;   //! Number of identified pions 
-  TH2F * fhPhiMCPion;  //! Phi of identified pions
-  TH2F * fhEtaMCPion;  //! eta of identified pions     
-  TH3F * fhLambdaMCPion; //! Shower parameters of MC pions vs transerse momentum   
-  TH3F * fhPionTrueE;  // MC truth E vs Recons E vs. Lambda of the cluster for MC pions
-  TH1F * fhRatioEPion  ; //! Reco/MC E distribution for charged pions      
-
-  TH1F * fhEMCProton;   //! Number of identified protons
-  TH2F * fhPhiMCProton;  //! Phi of identified protons
-  TH2F * fhEtaMCProton;  //! eta of identified protons 
-  TH3F * fhLambdaMCProton  ; //! Shower parameters of MC protons vs transerse momentum  
-  TH3F * fhProtonTrueE;  // MC truth E vs Recons E vs. Lambda of the cluster for MC protons
-  TH1F * fhRatioEProton  ; //! Reco/MC E distribution for protons       
-
-  TH1F * fhEMCAntiProton;   //! Number of identified antiprotons
-  TH2F * fhPhiMCAntiProton;  //! Phi of identified antiprotons
-  TH2F * fhEtaMCAntiProton;  //! eta of identified antiprotons 
-  TH3F * fhLambdaMCAntiProton  ; //! Shower parameters of MC antiprotons vs transerse momentum 
-  TH3F * fhAntiProtonTrueE;  // MC truth E vs Recons E vs. Lambda of the cluster for MC antiprotons
-  TH1F * fhRatioEAntiProton  ; //! Reco/MC E distribution for antiprotons        
-
-  TH1F * fhEMCNeutron;   //! Number of identified neutrons
-  TH2F * fhPhiMCNeutron;  //! Phi of identified neutrons
-  TH2F * fhEtaMCNeutron;  //! eta of identified neutrons       
-  TH3F * fhLambdaMCNeutron  ; //! Shower parameters of MC neutrons vs transerse momentum
-  TH3F * fhNeutronTrueE;  // MC truth E vs Recons E vs. Lambda of the cluster for MC Neutrons
-  TH1F * fhRatioENeutron  ; //! Reco/MC E distribution for photons         
-
-  TH1F * fhEMCEta;   //! Number of identified etas
-  TH2F * fhPhiMCEta;  //! Phi of identified etas
-  TH2F * fhEtaMCEta;  //! eta of identified etas       
-  TH3F * fhLambdaMCEta  ; //! Shower parameters of MC etas vs transerse momentum 
-
-  TH1D *  fhPrimPt; //! Pi0 pT spectrum truth.
-
-   ClassDef(AliAnaShowerParameter,1)
-
-} ;
-
-#endif//AliAnaShowerParameter_H
-
-
-
+#ifndef ALIANASHOWERPARAMETER_H\r
+#define ALIANASHOWERPARAMETER_H\r
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice     */\r
+/* $Id: AliAnaShowerParameter.h 27413 2008-07-18 13:28:12Z gconesab $ */\r
+\r
+//_________________________________________________________________________\r
+//\r
+// Class cloned from AliAnaPhoton, main aim is shower shape studies\r
+// \r
+// \r
+//\r
+//-- Author: Jocelyn Mlynarz (WSU) and Gustavo Conesa (LPSC)\r
+\r
+// --- ROOT system ---\r
+class TH3F;\r
+class TH2F ;\r
+class TH1F;\r
+class TString ;\r
+\r
+// --- ANALYSIS system ---\r
+#include "AliAnaPartCorrBaseClass.h"\r
+//#include "AliStack.h"\r
+//#include "TParticle.h"\r
+class AliStack;\r
+class TParticle;\r
+\r
+class TList ;\r
+class AliEMCALGeoUtils;\r
+class AliAnaShowerParameter : public AliAnaPartCorrBaseClass {\r
+\r
+public: \r
+\r
+  AliAnaShowerParameter() ; // default ctor\r
+  AliAnaShowerParameter(const AliAnaShowerParameter & g) ; // cpy ctor\r
+  AliAnaShowerParameter & operator = (const AliAnaShowerParameter & g) ;//cpy assignment\r
+  virtual ~AliAnaShowerParameter() ; //virtual dtor\r
+  \r
+  TList *  GetCreateOutputObjects();\r
+\r
+  void Init();\r
+\r
+  TObjString* GetAnalysisCuts();\r
+\r
+  void MakeAnalysisFillAOD()  ;\r
+    \r
+  void MakeAnalysisFillHistograms() ; \r
+  \r
+  void Print(const Option_t * opt)const;\r
+  \r
+  TString GetCalorimeter()   const {return fCalorimeter ; }\r
+  void SetCalorimeter(TString det)    {fCalorimeter = det ; }\r
+       \r
+  void InitParameters();\r
+\r
+  void SetTimeCut(Double_t min, Double_t max) {fTimeCutMin = min; fTimeCutMax = max;}\r
+  Double_t GetTimeCutMin() const {return fTimeCutMin;}\r
+  Double_t GetTimeCutMax() const {return fTimeCutMax;} \r
+  \r
+  void SetNCellsCut(Double_t min, Double_t max) {fNCellsCutMin = min; fNCellsCutMax = max;}\r
+  Double_t GetNCellsCutMin() const {return fNCellsCutMin;}\r
+  Double_t GetNCellsCutMax() const {return fNCellsCutMax;}     \r
+       \r
+  private:\r
\r
+  TString fCalorimeter ; // Calorimeter where the gamma is searched;\r
+  Float_t fNCellsCutMin ;\r
+  Float_t fNCellsCutMax ;\r
+  Float_t fLambdaCut ;\r
+  Double_t fTimeCutMin  ;    // Remove clusters/cells with time smaller than this value, in ns\r
+  Double_t fTimeCutMax  ;    // Remove clusters/cells with time larger than this value, in ns\r
+\r
+  //Histograms   \r
+  TH1F * fhNClusters  ; \r
+  TH2F * fhNCellCluster;\r
+  TH3F * fhEtaPhiPtCluster   ; \r
+  TH2F * fhLambdaPtCluster  ; \r
+\r
+  //MC\r
+  TH2F * fhLambdaPtPhoton ;    \r
+  TH2F * fhLambdaPtPi0 ;    \r
+  TH2F * fhLambdaPtPion ; \r
+  TH1D * fhPtTruthPi0 ;\r
+\r
+   ClassDef(AliAnaShowerParameter,1)\r
+\r
+} ;\r
\r
+\r
+#endif//AliAnaShowerParameter_H\r
+\r
+\r
+\r