ana.C: Included possibility to read MC data directly from galice.root

AliAnalysisTaskParticleCorrelation.cxx: During initialization add generators particles missing in PDG
AliCaloTrackMCReader.cxx: Remove skip of particles with pdg not in Database
AliCaloPID.cxx,.h: CheckOrigin depends on MC generator, HERWIG added. Electrons where given a photon identification, corrected
AliAnaScale.cxx: Scale only keys which are histograms, begin with "TH"

diff --git a/PWG4/PartCorrBase/AliAnaScale.cxx b/PWG4/PartCorrBase/AliAnaScale.cxx
index b2cbc20515a8a0ec77525c82be4db3f54b42f929..6899facbd04f6ea87482f48da09dfa1eb767fc24 100755 (executable)
--- a/PWG4/PartCorrBase/AliAnaScale.cxx
+++ b/PWG4/PartCorrBase/AliAnaScale.cxx
@@ -29,6 +29,7 @@
 
 #include "AliAnaScale.h" 
 #include "AliAnalysisManager.h"
 
 #include "AliAnaScale.h" 
 #include "AliAnalysisManager.h"

+#include "Riostream.h"

 
 //______________________________________________________________________________
 AliAnaScale::AliAnaScale() : 
 
 //______________________________________________________________________________
 AliAnaScale::AliAnaScale() : 


@@ -100,11 +101,10 @@ void AliAnaScale::Exec(Option_t *)
   TObject * h ; 
   while ( (h = next()) ) { 
     if(h){
   TObject * h ; 
   while ( (h = next()) ) { 
     if(h){

-      if ( strcmp(h->ClassName(),"TNtuple") ) {
+      if ( !strncmp(h->ClassName(),"TH",2) ) {

       char name[128] ; 
       sprintf(name, "%sScaled", h->GetName()) ; 
       TH1 * hout = dynamic_cast<TH1*> (h->Clone(name)) ; 
       char name[128] ; 
       sprintf(name, "%sScaled", h->GetName()) ; 
       TH1 * hout = dynamic_cast<TH1*> (h->Clone(name)) ; 

-     

       if(fSumw2) hout->Sumw2();
       hout->Scale(fScale) ;  
       fOutputList->Add(hout) ; 
       if(fSumw2) hout->Sumw2();
       hout->Scale(fScale) ;  
       fOutputList->Add(hout) ; 

diff --git a/PWG4/PartCorrBase/AliAnalysisTaskParticleCorrelation.cxx b/PWG4/PartCorrBase/AliAnalysisTaskParticleCorrelation.cxx
index 744405ffe08f1e496e538b1a252a85e86e8214d1..c3a0ceacf0da58c7d37e52357ef34df53721fd49 100755 (executable)
--- a/PWG4/PartCorrBase/AliAnalysisTaskParticleCorrelation.cxx
+++ b/PWG4/PartCorrBase/AliAnalysisTaskParticleCorrelation.cxx
@@ -44,6 +44,7 @@
 #include "AliAODHandler.h"
 #include "AliStack.h"
 #include "AliLog.h"
 #include "AliAODHandler.h"
 #include "AliStack.h"
 #include "AliLog.h"

+#include "AliPDG.h"

 
 ClassImp(AliAnalysisTaskParticleCorrelation)
 
 
 ClassImp(AliAnalysisTaskParticleCorrelation)
 


@@ -130,6 +131,10 @@ void AliAnalysisTaskParticleCorrelation::Init()
   if(!fAna)
     AliFatal("Analysis pointer not initialized, abort analysis!");
   
   if(!fAna)
     AliFatal("Analysis pointer not initialized, abort analysis!");
   

+  // Add different generator particles to PDG Data Base 
+  // to avoid problems when reading MC generator particles
+  AliPDG::AddParticlesToPdgDataBase();
+

   // Initialise analysis
   fAna->Init();
   
   // Initialise analysis
   fAna->Init();
   

diff --git a/PWG4/PartCorrBase/AliCaloPID.cxx b/PWG4/PartCorrBase/AliCaloPID.cxx
index 040dccacbb8cd484044dfdb1bc5bfbfc46ec47e3..023f0399b45feeb7de3b9e069a4a035ec9f18561 100755 (executable)
--- a/PWG4/PartCorrBase/AliCaloPID.cxx
+++ b/PWG4/PartCorrBase/AliCaloPID.cxx
@@ -48,7 +48,7 @@ fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.),
 fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , 
 fPHOSNeutralWeight(0.), fPHOSWeightFormula(0), 
 fPHOSPhotonWeightFormula(0x0), fPHOSPi0WeightFormula(0x0),
 fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , 
 fPHOSNeutralWeight(0.), fPHOSWeightFormula(0), 
 fPHOSPhotonWeightFormula(0x0), fPHOSPi0WeightFormula(0x0),

-fDispCut(0.),fTOFCut(0.), fDebug(-1)
+fDispCut(0.),fTOFCut(0.), fDebug(-1), fMCGenerator("")

 {
        //Ctor
        
 {
        //Ctor
        


@@ -72,7 +72,7 @@ fPHOSWeightFormula(pid.fPHOSWeightFormula),
 fPHOSPhotonWeightFormula(pid.fPHOSPhotonWeightFormula), 
 fPHOSPi0WeightFormula(pid.fPHOSPi0WeightFormula), 
 fDispCut(pid.fDispCut),fTOFCut(pid.fTOFCut),
 fPHOSPhotonWeightFormula(pid.fPHOSPhotonWeightFormula), 
 fPHOSPi0WeightFormula(pid.fPHOSPi0WeightFormula), 
 fDispCut(pid.fDispCut),fTOFCut(pid.fTOFCut),

-fDebug(pid.fDebug)
+fDebug(pid.fDebug),fMCGenerator(pid.fMCGenerator)

 {
        // cpy ctor
        
 {
        // cpy ctor
        


@@ -104,7 +104,8 @@ AliCaloPID & AliCaloPID::operator = (const AliCaloPID & pid)
        fDispCut  = pid.fDispCut;
        fTOFCut   = pid.fTOFCut;
        fDebug    = pid.fDebug;
        fDispCut  = pid.fDispCut;
        fTOFCut   = pid.fTOFCut;
        fDebug    = pid.fDebug;

-       
+       fMCGenerator = pid.fMCGenerator;
+

        return *this;
        
 }
        return *this;
        
 }


@@ -121,66 +122,95 @@ AliCaloPID::~AliCaloPID() {
 
 //_________________________________________________________________________
 Int_t AliCaloPID::CheckOrigin(const Int_t label, AliStack * stack) const {
 
 //_________________________________________________________________________
 Int_t AliCaloPID::CheckOrigin(const Int_t label, AliStack * stack) const {

-       //Play with the MC stack if available
-       //Check origin of the candidates, good for PYTHIA
-       
-       if(!stack) AliFatal("Stack is not available, check analysis settings in configuration file, STOP!!");
-       
-       if(label >= 0 && label <  stack->GetNtrack()){
-               //Mother
-               TParticle * mom = stack->Particle(label);
-               Int_t mPdg = TMath::Abs(mom->GetPdgCode());
-               Int_t mStatus =  mom->GetStatusCode() ;
-               Int_t iParent =  mom->GetFirstMother() ;
-               if(fDebug > 0 && label < 8 ) printf("AliCaloPID::CheckOrigin: Mother is parton %d\n",iParent);
-               
-               //GrandParent
-               TParticle * parent = new TParticle ;
-               Int_t pPdg = -1;
-               Int_t pStatus =-1;
-               if(iParent > 0){
-                       parent = stack->Particle(iParent);
-                       pPdg = TMath::Abs(parent->GetPdgCode());
-                       pStatus = parent->GetStatusCode();  
-               }
-               else if(fDebug > 0 ) printf("AliCaloPID::CheckOrigin: Parent with label %d\n",iParent);
-               
-               //return tag
-               if(mPdg == 22){
-                       if(mStatus == 1){
-                               if(iParent < 8) {
-                                       if(pPdg == 22) return kMCPrompt;
-                                       else  return kMCFragmentation;
-                               }
-                               else if(pStatus == 11){
-                                       if(pPdg == 111) return kMCPi0Decay ;
-                                       else if (pPdg == 321)  return kMCEtaDecay ;
-                                       else  return kMCOtherDecay ;
-                               }
-                       }//Status 1 : Pythia generated
-                       else if(mStatus == 0){
-                               if(pPdg ==22 || pPdg ==11) return kMCConversion ;
-                               if(pPdg == 111) return kMCPi0Decay ;
-                               else if (pPdg == 221)  return kMCEtaDecay ;
-                               else  return kMCOtherDecay ;
-                       }//status 0 : geant generated
-               }//Mother Photon
-               else if(mPdg == 111)  return kMCPi0 ;
-               else if(mPdg == 221)  return kMCEta ;
-               else if(mPdg ==11){
-                       if(mStatus == 0) return kMCConversion ;
-                       else return kMCElectron ;
-               }
-               else return kMCUnknown;
-       }//Good label value
-       else{
-               if(label < 0 ) printf("AliCaloPID::CheckOrigin: *** bad label or no stack ***:  label %d \n", label);
-               if(label >=  stack->GetNtrack()) printf("AliCaloPID::CheckOrigin: *** large label ***:  label %d, n tracks %d \n", label, stack->GetNtrack());
-               return kMCUnknown;
-       }//Bad label
-       
-       return kMCUnknown;
-       
+  //Play with the MC stack if available
+  //Check origin of the candidates, good for PYTHIA
+  
+  if(!stack) AliFatal("Stack is not available, check analysis settings in configuration file, STOP!!");
+  
+  if(label >= 0 && label <  stack->GetNtrack()){
+    //Mother
+    TParticle * mom = stack->Particle(label);
+    Int_t mPdg = TMath::Abs(mom->GetPdgCode());
+    Int_t mStatus =  mom->GetStatusCode() ;
+    Int_t iParent =  mom->GetFirstMother() ;
+    if(fDebug > 0 && label < 8 ) printf("AliCaloPID::CheckOrigin: Mother is parton %d\n",iParent);
+    
+    //GrandParent
+    TParticle * parent = new TParticle ;
+    Int_t pPdg = -1;
+    Int_t pStatus =-1;
+    if(iParent > 0){
+      parent = stack->Particle(iParent);
+      pPdg = TMath::Abs(parent->GetPdgCode());
+      pStatus = parent->GetStatusCode();  
+    }
+    else if(fDebug > 0 ) printf("AliCaloPID::CheckOrigin: Parent with label %d\n",iParent);
+    
+    //return tag
+    if(mPdg == 22){
+      if(mStatus == 1){
+       if(fMCGenerator == "PYTHIA"){
+         if(iParent < 8 && iParent > 5) {//outgoing partons
+           if(pPdg == 22) return kMCPrompt;
+           else  return kMCFragmentation;
+         }//Outgoing partons
+         else if(pStatus == 11){//Decay
+           if(pPdg == 111) return kMCPi0Decay ;
+           else if (pPdg == 321)  return kMCEtaDecay ;
+           else  return kMCOtherDecay ;
+         }//Decay
+         else return kMCISR; //Initial state radiation
+       }//PYTHIA
+
+       else if(fMCGenerator == "HERWIG"){        
+         if(pStatus < 197){//Not decay
+           while(1){
+             if(parent->GetFirstMother()<=5) break;
+             iParent = parent->GetFirstMother();
+             parent=stack->Particle(iParent);
+             pStatus= parent->GetStatusCode();
+             pPdg = parent->GetPdgCode();
+           }//Look for the parton
+           
+           if(iParent < 8 && iParent > 5) {
+             if(pPdg == 22) return kMCPrompt;
+             else  return kMCFragmentation;
+           }
+           return kMCISR;//Initial state radiation
+         }//Not decay
+         else{//Decay
+           if(pPdg == 111) return kMCPi0Decay ;
+           else if (pPdg == 321)  return kMCEtaDecay ;
+           else  return kMCOtherDecay ;
+         }//Decay
+       }//HERWIG
+       else return  kMCUnknown;
+      }//Status 1 : Pythia generated
+      else if(mStatus == 0){
+       if(pPdg ==22 || pPdg ==11) return kMCConversion ;
+       if(pPdg == 111) return kMCPi0Decay ;
+       else if (pPdg == 221)  return kMCEtaDecay ;
+       else  return kMCOtherDecay ;
+      }//status 0 : geant generated
+    }//Mother Photon
+    else if(mPdg == 111)  return kMCPi0 ;
+    else if(mPdg == 221)  return kMCEta ;
+    else if(mPdg ==11){
+      printf("Origin electron, pT %f\n",mom->Pt());
+
+      if(mStatus == 0) return kMCConversion ;
+      else return kMCElectron ;
+    }
+    else return kMCUnknown;
+  }//Good label value
+  else{
+    if(label < 0 ) printf("AliCaloPID::CheckOrigin: *** bad label or no stack ***:  label %d \n", label);
+    if(label >=  stack->GetNtrack()) printf("AliCaloPID::CheckOrigin: *** large label ***:  label %d, n tracks %d \n", label, stack->GetNtrack());
+    return kMCUnknown;
+  }//Bad label
+       
+  return kMCUnknown;
+  

 }
 
 //_______________________________________________________________
 }
 
 //_______________________________________________________________


@@ -210,6 +240,7 @@ void AliCaloPID::InitParameters()
     fDispCut  = 1.5;
        fTOFCut   = 5.e-9;
        fDebug = -1;
     fDispCut  = 1.5;
        fTOFCut   = 5.e-9;
        fDebug = -1;

+       fMCGenerator = "PYTHIA";

 }
 
 //_______________________________________________________________
 }
 
 //_______________________________________________________________


@@ -267,9 +298,9 @@ Int_t AliCaloPID::GetPdg(const TString calo, const Double_t * pid, const Float_t
                        pdg = kNeutralUnknown ;
        }
        else{//EMCAL
                        pdg = kNeutralUnknown ;
        }
        else{//EMCAL

-               //Temporal solution, electrons and photons not differenciated
-               if(pid[AliAODCluster::kPhoton] + pid[AliAODCluster::kElectron]  > wPh) pdg = kPhoton ;
+               if(pid[AliAODCluster::kPhoton]  > wPh) pdg = kPhoton ;

                else if(pid[AliAODCluster::kPi0] > wPi0) pdg = kPi0 ; 
                else if(pid[AliAODCluster::kPi0] > wPi0) pdg = kPi0 ; 

+               else if(pid[AliAODCluster::kElectron]  > wE) pdg = kElectron ;

                else if(chargedHadronWeight + neutralHadronWeight > wCh) pdg = kChargedHadron ;  
                else if(neutralHadronWeight + chargedHadronWeight > wNe) pdg = kNeutralHadron ; 
                else pdg =  kNeutralUnknown ;
                else if(chargedHadronWeight + neutralHadronWeight > wCh) pdg = kChargedHadron ;  
                else if(neutralHadronWeight + chargedHadronWeight > wNe) pdg = kNeutralHadron ; 
                else pdg =  kNeutralUnknown ;


@@ -374,7 +405,8 @@ void AliCaloPID::Print(const Option_t * opt) const
        printf("TOF cut        = %e\n",fTOFCut);
        printf("Dispersion cut = %2.2f\n",fDispCut);
        printf("Debug level    = %d\n",fDebug);
        printf("TOF cut        = %e\n",fTOFCut);
        printf("Dispersion cut = %2.2f\n",fDispCut);
        printf("Debug level    = %d\n",fDebug);

-       
+       printf("MC Generator   = %s\n",fMCGenerator.Data());
+

        printf(" \n");
        
 } 
        printf(" \n");
        
 } 

diff --git a/PWG4/PartCorrBase/AliCaloPID.h b/PWG4/PartCorrBase/AliCaloPID.h
index b55c1a231c3505806b4aff492deeea012634d9b2..519069918e4798bac98a692d164d9291eed50eb3 100755 (executable)
--- a/PWG4/PartCorrBase/AliCaloPID.h
+++ b/PWG4/PartCorrBase/AliCaloPID.h
@@ -44,7 +44,7 @@ public:
        };
        
        
        };
        
        

-       enum mcTypes {kMCPrompt, kMCFragmentation, kMCPi0Decay, kMCEtaDecay, kMCOtherDecay, kMCPi0, kMCEta, kMCElectron, kMCConversion, kMCUnknown};
+       enum mcTypes {kMCPrompt, kMCFragmentation, kMCISR, kMCPi0Decay, kMCEtaDecay, kMCOtherDecay, kMCPi0, kMCEta, kMCElectron, kMCConversion, kMCUnknown};

        
        void InitParameters();
        Int_t CheckOrigin(const Int_t label, AliStack *  stack) const ;
        
        void InitParameters();
        Int_t CheckOrigin(const Int_t label, AliStack *  stack) const ;


@@ -101,6 +101,9 @@ public:
        void SetDebug(Int_t deb) {fDebug=deb;}
        Int_t GetDebug() const {return fDebug;} 
        
        void SetDebug(Int_t deb) {fDebug=deb;}
        Int_t GetDebug() const {return fDebug;} 
        

+       void SetMCGenerator(TString mcgen) {fMCGenerator=mcgen;}
+       TString GetMCGenerator() const {return fMCGenerator;}   
+

 private:
        
        Float_t      fEMCALPhotonWeight; //Bayesian PID weight for photons in EMCAL 
 private:
        
        Float_t      fEMCALPhotonWeight; //Bayesian PID weight for photons in EMCAL 


@@ -122,8 +125,9 @@ private:
        Float_t fTOFCut;     //Cut on TOF, used in PID evaluation
        
        Int_t    fDebug; //Debug level
        Float_t fTOFCut;     //Cut on TOF, used in PID evaluation
        
        Int_t    fDebug; //Debug level

-       
-       ClassDef(AliCaloPID,2)
+       TString fMCGenerator; // MC geneator used to generate data in simulation
+
+       ClassDef(AliCaloPID,3)

 } ;
 
 
 } ;
 
 

diff --git a/PWG4/PartCorrBase/AliCaloTrackMCReader.cxx b/PWG4/PartCorrBase/AliCaloTrackMCReader.cxx
index 600a52abff1293642874aea56f7904ed18803809..85d44282103dca51f5d5d626fb80599ff9a8bfa2 100755 (executable)
--- a/PWG4/PartCorrBase/AliCaloTrackMCReader.cxx
+++ b/PWG4/PartCorrBase/AliCaloTrackMCReader.cxx
@@ -138,7 +138,6 @@ void AliCaloTrackMCReader::InitParameters()
 void  AliCaloTrackMCReader::FillCalorimeters(const Int_t iParticle, TParticle* particle, TLorentzVector momentum,
                                             Int_t &indexPHOS, Int_t &indexEMCAL) {
        //Fill AODCaloClusters or TParticles lists of PHOS or EMCAL
 void  AliCaloTrackMCReader::FillCalorimeters(const Int_t iParticle, TParticle* particle, TLorentzVector momentum,
                                             Int_t &indexPHOS, Int_t &indexEMCAL) {
        //Fill AODCaloClusters or TParticles lists of PHOS or EMCAL

-       

        //In PHOS
        if(fFillPHOS && fFidutialCut->IsInFidutialCut(momentum,"PHOS") && momentum.Pt() > fPHOSPtMin){
                
        //In PHOS
        if(fFillPHOS && fFidutialCut->IsInFidutialCut(momentum,"PHOS") && momentum.Pt() > fPHOSPtMin){
                


@@ -209,7 +208,7 @@ void AliCaloTrackMCReader::FillInputEvent()
                if(KeepParticleWithStatus(particle->GetStatusCode()) && (particle->Pt() > 0) ){
                
                    //Skip bizarre particles, they crash when charge is calculated
                if(KeepParticleWithStatus(particle->GetStatusCode()) && (particle->Pt() > 0) ){
                
                    //Skip bizarre particles, they crash when charge is calculated

-                       if(TMath::Abs(pdg) == 3124 || TMath::Abs(pdg) > 10000000) continue ;
+                 //    if(TMath::Abs(pdg) == 3124 || TMath::Abs(pdg) > 10000000) continue ;

                        
                        charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
                        particle->Momentum(momentum);
                        
                        charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
                        particle->Momentum(momentum);

diff --git a/PWG4/PartCorrDep/AliAnaPhoton.cxx b/PWG4/PartCorrDep/AliAnaPhoton.cxx
index 50221bce24849bf3042e264371d16c2a1761e5bb..10327487b1ac75274b0d489a03bbe8a7cef85052 100755 (executable)
--- a/PWG4/PartCorrDep/AliAnaPhoton.cxx
+++ b/PWG4/PartCorrDep/AliAnaPhoton.cxx
@@ -37,9 +37,9 @@ ClassImp(AliAnaPhoton)
 //____________________________________________________________________________
   AliAnaPhoton::AliAnaPhoton() : 
     AliAnaPartCorrBaseClass(), fCalorimeter(""), 
 //____________________________________________________________________________
   AliAnaPhoton::AliAnaPhoton() : 
     AliAnaPartCorrBaseClass(), fCalorimeter(""), 

-       fMinDist(0.),fMinDist2(0.),fMinDist3(0.),
-       fhPtPhoton(0),fhPhiPhoton(0),fhEtaPhoton(0),
-       //MC
+    fMinDist(0.),fMinDist2(0.),fMinDist3(0.),
+    fhPtPhoton(0),fhPhiPhoton(0),fhEtaPhoton(0),
+    //MC

     fhPtPrompt(0),fhPhiPrompt(0),fhEtaPrompt(0), 
     fhPtFragmentation(0),fhPhiFragmentation(0),fhEtaFragmentation(0), 
     fhPtPi0Decay(0),fhPhiPi0Decay(0),fhEtaPi0Decay(0), 
     fhPtPrompt(0),fhPhiPrompt(0),fhEtaPrompt(0), 
     fhPtFragmentation(0),fhPhiFragmentation(0),fhEtaFragmentation(0), 
     fhPtPi0Decay(0),fhPhiPi0Decay(0),fhEtaPi0Decay(0), 


@@ -143,13 +143,13 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
        outputContainer->Add(fhPtPhoton) ; 
        
        fhPhiPhoton  = new TH2F
        outputContainer->Add(fhPtPhoton) ; 
        
        fhPhiPhoton  = new TH2F

-    ("hPhiPhoton","#phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 
+         ("hPhiPhoton","#phi_{#gamma}",nptbins,ptmin,ptmax,nphibins,phimin,phimax); 

        fhPhiPhoton->SetYTitle("#phi");
        fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
        outputContainer->Add(fhPhiPhoton) ; 
        
        fhEtaPhoton  = new TH2F
        fhPhiPhoton->SetYTitle("#phi");
        fhPhiPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
        outputContainer->Add(fhPhiPhoton) ; 
        
        fhEtaPhoton  = new TH2F

-    ("hEtaPhoton","#phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 
+         ("hEtaPhoton","#phi_{#gamma}",nptbins,ptmin,ptmax,netabins,etamin,etamax); 

        fhEtaPhoton->SetYTitle("#eta");
        fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
        outputContainer->Add(fhEtaPhoton) ;
        fhEtaPhoton->SetYTitle("#eta");
        fhEtaPhoton->SetXTitle("p_{T #gamma} (GeV/c)");
        outputContainer->Add(fhEtaPhoton) ;


@@ -281,7 +281,7 @@ TList *  AliAnaPhoton::GetCreateOutputObjects()
        parList += GetCaloPID()->GetPIDParametersList() ;
        
        //Get parameters set in FidutialCut class (not available yet)
        parList += GetCaloPID()->GetPIDParametersList() ;
        
        //Get parameters set in FidutialCut class (not available yet)

-       //parlist += GetCaloPID()->GetFidCutParametersList() 
+       //parlist += GetFidCut()->GetFidCutParametersList() 

        
        TObjString *oString= new TObjString(parList) ;
        outputContainer->Add(oString);
        
        TObjString *oString= new TObjString(parList) ;
        outputContainer->Add(oString);


@@ -344,7 +344,8 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
                //Set the indeces of the original caloclusters  
                aodph.SetCaloLabel(calo->GetID(),-1);
                aodph.SetDetector(fCalorimeter);
                //Set the indeces of the original caloclusters  
                aodph.SetCaloLabel(calo->GetID(),-1);
                aodph.SetDetector(fCalorimeter);

-               if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: Min pt cut and fidutial cut passed: pt %3.2f, phi %2.2f, eta %1.2f\n",aodph.Pt(),aodph.Phi(),aodph.Eta());    
+               if(GetDebug() > 1) 
+                 printf("AliAnaPhoton::FillAOD: Min pt cut and fidutial 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->GetDistToBadChannel() ; //Distance to bad channel
                        
                //Check Distance to Bad channel, set bit.
                Double_t distBad=calo->GetDistToBadChannel() ; //Distance to bad channel


@@ -363,7 +364,7 @@ void  AliAnaPhoton::MakeAnalysisFillAOD()
                if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){
                        //Get most probable PID, check PID weights (in MC this option is mandatory)
                        aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,calo->PID(),mom.E()));//PID with weights
                if(GetReader()->GetDataType() == AliCaloTrackReader::kMC){
                        //Get most probable PID, check PID weights (in MC this option is mandatory)
                        aodph.SetPdg(GetCaloPID()->GetPdg(fCalorimeter,calo->PID(),mom.E()));//PID with weights

-                       if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PDG of identified particle %d\n",aodph.GetPdg());
+                       if(GetDebug() > 1) printf("AliAnaPhoton::FillAOD: PDG of identified particle %d\n",aodph.GetPdg());      

                        //If primary is not photon, skip it.
                        if(aodph.GetPdg() != AliCaloPID::kPhoton) continue ;
                }                                       
                        //If primary is not photon, skip it.
                        if(aodph.GetPdg() != AliCaloPID::kPhoton) continue ;
                }                                       


@@ -415,9 +416,9 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
        //Do analysis and fill histograms
 
        //Get vertex for photon momentum calculation
        //Do analysis and fill histograms
 
        //Get vertex for photon momentum calculation

-       Double_t v[]={0,0,0} ; //vertex ;
+       //Double_t v[]={0,0,0} ; //vertex ;

        //if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) 
        //if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) 

-       GetReader()->GetVertex(v);
+                //GetReader()->GetVertex(v);

 
        //Loop on stored AOD photons
        Int_t naod = GetOutputAODBranch()->GetEntriesFast();
 
        //Loop on stored AOD photons
        Int_t naod = GetOutputAODBranch()->GetEntriesFast();


@@ -427,14 +428,16 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
                AliAODPWG4Particle* ph =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
                Int_t pdg = ph->GetPdg();
 
                AliAODPWG4Particle* ph =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));
                Int_t pdg = ph->GetPdg();
 

-               if(GetDebug() > 3) printf("AliAnaPhoton::FillHistos: PDG %d, MC TAG %d, Calorimeter %s\n", ph->GetPdg(),ph->GetTag(), (ph->GetDetector()).Data()) ;
+               if(GetDebug() > 3) 
+                 printf("AliAnaPhoton::FillHistos: PDG %d, MC TAG %d, Calorimeter %s\n", ph->GetPdg(),ph->GetTag(), (ph->GetDetector()).Data()) ;

                
                //If PID used, fill histos with photons in Calorimeter fCalorimeter
                if(IsCaloPIDOn())  
                        if(pdg != AliCaloPID::kPhoton) continue; 
                if(ph->GetDetector() != fCalorimeter) continue;
                
                
                //If PID used, fill histos with photons in Calorimeter fCalorimeter
                if(IsCaloPIDOn())  
                        if(pdg != AliCaloPID::kPhoton) continue; 
                if(ph->GetDetector() != fCalorimeter) continue;
                

-               if(GetDebug() > 1) printf("AliAnaPhoton::FillHistos: ID Photon: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;
+               if(GetDebug() > 1) 
+                 printf("AliAnaPhoton::FillHistos: ID Photon: pt %f, phi %f, eta %f\n", ph->Pt(),ph->Phi(),ph->Eta()) ;

                
                //Fill photon histograms 
                Float_t ptcluster = ph->Pt();
                
                //Fill photon histograms 
                Float_t ptcluster = ph->Pt();


@@ -477,7 +480,6 @@ void  AliAnaPhoton::MakeAnalysisFillHistograms()
                                fhEtaConversion ->Fill(ptcluster,etacluster);
                        }
                        else{
                                fhEtaConversion ->Fill(ptcluster,etacluster);
                        }
                        else{

-                               

                                fhPtUnknown   ->Fill(ptcluster);
                                fhPhiUnknown ->Fill(ptcluster,phicluster);
                                fhEtaUnknown ->Fill(ptcluster,etacluster);
                                fhPtUnknown   ->Fill(ptcluster);
                                fhPhiUnknown ->Fill(ptcluster,phicluster);
                                fhEtaUnknown ->Fill(ptcluster,etacluster);

diff --git a/PWG4/macros/ana.C b/PWG4/macros/ana.C
index 4aaed1b330bb5233a08ab0bb0cc62f9398cf0ce0..dadbfe2e0f4a84916257431157325b7746214eb9 100644 (file)
--- a/PWG4/macros/ana.C
+++ b/PWG4/macros/ana.C
@@ -36,7 +36,9 @@ const Int_t kNumberOfEventsPerFile = 100;
 //---------------------------------------------------------------------------
 
 const Bool_t kMC = kTRUE; //With real data kMC = kFALSE
 //---------------------------------------------------------------------------
 
 const Bool_t kMC = kTRUE; //With real data kMC = kFALSE

-const TString kInputData = "ESD";
+const TString kInputData = "ESD";//ESD, AOD, MC
+TString kTreeName = "esdTree";
+

 void ana(Int_t mode=mLocal, TString configName = "ConfigAnalysisPhoton")
 {
   // Main
 void ana(Int_t mode=mLocal, TString configName = "ConfigAnalysisPhoton")
 {
   // Main


@@ -51,7 +53,15 @@ void ana(Int_t mode=mLocal, TString configName = "ConfigAnalysisPhoton")
   //-------------------------------------------------------------------------------------------------
   //Create chain from ESD and from cross sections files, look below for options.
   //------------------------------------------------------------------------------------------------- 
   //-------------------------------------------------------------------------------------------------
   //Create chain from ESD and from cross sections files, look below for options.
   //------------------------------------------------------------------------------------------------- 

-  TChain *chain   = new TChain("esdTree") ;
+  if(kInputData == "ESD") kTreeName = "esdTree" ;
+  else if(kInputData == "AOD") kTreeName = "aodTree" ;
+  else if (kInputData == "MC") kTreeName = "TE" ;
+  else {
+    cout<<"Wrong  data type "<<kInputData<<endl;
+    break;
+  }
+
+  TChain *chain   = new TChain(kTreeName) ;

   TChain * chainxs = new TChain("Xsection") ;
   CreateChain(mode, chain, chainxs);  
 
   TChain * chainxs = new TChain("Xsection") ;
   CreateChain(mode, chain, chainxs);  
 


@@ -63,7 +73,7 @@ void ana(Int_t mode=mLocal, TString configName = "ConfigAnalysisPhoton")
     //-------------------------------------
     AliAnalysisManager *mgr  = new AliAnalysisManager("Manager", "Manager");
     // MC handler
     //-------------------------------------
     AliAnalysisManager *mgr  = new AliAnalysisManager("Manager", "Manager");
     // MC handler

-    if(kMC){
+    if(kMC || kInputData == "MC"){

       AliMCEventHandler* mcHandler = new AliMCEventHandler();
       mcHandler->SetReadTR(kFALSE);//Do not search TrackRef file
       mgr->SetMCtruthEventHandler(mcHandler);
       AliMCEventHandler* mcHandler = new AliMCEventHandler();
       mcHandler->SetReadTR(kFALSE);//Do not search TrackRef file
       mgr->SetMCtruthEventHandler(mcHandler);


@@ -336,6 +346,11 @@ void CreateChain(const anaModes mode, TChain * chain, TChain * chainxs){
       cout<<"NEVENT : "<<kEvent<<endl;
       cout<<"PATTERN: " <<kPattern<<endl;
       
       cout<<"NEVENT : "<<kEvent<<endl;
       cout<<"PATTERN: " <<kPattern<<endl;
       

+      TString datafile="";
+      if(kInputData == "ESD") datafile = "AliESDs.root" ;
+      else if(kInputData == "AOD") datafile = "aod.root" ;
+      else if(kInputData == "MC")  datafile = "galice.root" ;
+      

       //Loop on ESD files, add them to chain
       Int_t event =0;
       Int_t skipped=0 ; 
       //Loop on ESD files, add them to chain
       Int_t event =0;
       Int_t skipped=0 ; 


@@ -343,12 +358,12 @@ void CreateChain(const anaModes mode, TChain * chain, TChain * chainxs){
       char filexs[120] ;
       
       for (event = 0 ; event < kEvent ; event++) {
       char filexs[120] ;
       
       for (event = 0 ; event < kEvent ; event++) {

-       sprintf(file, "%s/%s%d/AliESDs.root", kInDir,kPattern,event) ; 
+       sprintf(file, "%s/%s%d/%s", kInDir,kPattern,event,datafile.Data()) ; 

        sprintf(filexs, "%s/%s%d/%s", kInDir,kPattern,event,kXSFileName) ; 
        TFile * fESD = 0 ; 
        //Check if file exists and add it, if not skip it
        if ( fESD = TFile::Open(file)) {
        sprintf(filexs, "%s/%s%d/%s", kInDir,kPattern,event,kXSFileName) ; 
        TFile * fESD = 0 ; 
        //Check if file exists and add it, if not skip it
        if ( fESD = TFile::Open(file)) {

-         if ( fESD->Get("esdTree") ) { 
+         if ( fESD->Get(kTreeName) ) { 

            printf("++++ Adding %s\n", file) ;
            chain->AddFile(file);
            chainxs->Add(filexs) ; 
            printf("++++ Adding %s\n", file) ;
            chain->AddFile(file);
            chainxs->Add(filexs) ;