/************************************************************************** * 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 iGetEntriesFast(s 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$ */ //_________________________________________________________________________ // Class for the prompt gamma analysis, isolation cut // // Class created from old AliPHOSGammaJet // (see AliRoot versions previous Release 4-09) // // -- Author: Gustavo Conesa (LNF-INFN) ////////////////////////////////////////////////////////////////////////////// // --- ROOT system --- #include #include #include #include "Riostream.h" // --- Analysis system --- #include "AliAnaGammaDirect.h" #include "AliLog.h" #include "AliCaloTrackReader.h" #include "AliIsolationCut.h" #include "AliNeutralMesonSelection.h" ClassImp(AliAnaGammaDirect) //____________________________________________________________________________ AliAnaGammaDirect::AliAnaGammaDirect() : AliAnaPartCorrBaseClass(), fDetector(""), fMakeIC(0), fReMakeIC(0), fMakeSeveralIC(0), fMakeInvMass(0), fhPtGamma(0),fhPhiGamma(0),fhEtaGamma(0), fhConeSumPt(0), //Several IC fNCones(0),fNPtThresFrac(0), fConeSizes(), fPtThresholds(), fPtFractions(), fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(), //MC fhPtPrompt(0),fhPhiPrompt(0),fhEtaPrompt(0), fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(), fhPtFragmentation(0),fhPhiFragmentation(0),fhEtaFragmentation(0), fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(), fhPtPi0Decay(0),fhPhiPi0Decay(0),fhEtaPi0Decay(0), fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(), fhPtOtherDecay(0),fhPhiOtherDecay(0),fhEtaOtherDecay(0), fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(), fhPtConversion(0),fhPhiConversion(0),fhEtaConversion(0), fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(), fhPtUnknown(0),fhPhiUnknown(0),fhEtaUnknown(0), fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown() { //default ctor //Initialize parameters InitParameters(); for(Int_t i = 0; i < 5 ; i++){ fConeSizes[i] = 0 ; fhPtSumIsolated[i] = 0 ; fhPtSumIsolatedPrompt[i] = 0 ; fhPtSumIsolatedFragmentation[i] = 0 ; fhPtSumIsolatedPi0Decay[i] = 0 ; fhPtSumIsolatedOtherDecay[i] = 0 ; fhPtSumIsolatedConversion[i] = 0 ; fhPtSumIsolatedUnknown[i] = 0 ; for(Int_t j = 0; j < 5 ; j++){ fhPtThresIsolated[i][j] = 0 ; fhPtFracIsolated[i][j] = 0 ; fhPtThresIsolatedPrompt[i][j] = 0 ; fhPtThresIsolatedFragmentation[i][j] = 0 ; fhPtThresIsolatedPi0Decay[i][j] = 0 ; fhPtThresIsolatedOtherDecay[i][j] = 0 ; fhPtThresIsolatedConversion[i][j] = 0 ; fhPtThresIsolatedUnknown[i][j] = 0 ; fhPtFracIsolatedPrompt[i][j] = 0 ; fhPtFracIsolatedFragmentation[i][j] = 0 ; fhPtFracIsolatedPi0Decay[i][j] = 0 ; fhPtFracIsolatedOtherDecay[i][j] = 0 ; fhPtFracIsolatedConversion[i][j] = 0 ; fhPtFracIsolatedUnknown[i][j] = 0 ; } } for(Int_t i = 0; i < 5 ; i++){ fPtFractions[i]= 0 ; fPtThresholds[i]= 0 ; } } //____________________________________________________________________________ AliAnaGammaDirect::AliAnaGammaDirect(const AliAnaGammaDirect & g) : AliAnaPartCorrBaseClass(g), fDetector(g.fDetector), fMakeIC(g.fMakeIC), fReMakeIC(g.fReMakeIC), fMakeSeveralIC(g.fMakeSeveralIC), fMakeInvMass(g.fMakeInvMass), fhPtGamma(g.fhPtGamma),fhPhiGamma(g.fhPhiGamma), fhEtaGamma(g.fhEtaGamma), fhConeSumPt(g.fhConeSumPt), //Several IC fNCones(g.fNCones),fNPtThresFrac(g.fNPtThresFrac), fConeSizes(), fPtThresholds(), fPtFractions(), fhPtThresIsolated(), fhPtFracIsolated(), fhPtSumIsolated(), //MC fhPtPrompt(g.fhPtPrompt),fhPhiPrompt(g.fhPhiPrompt),fhEtaPrompt(g.fhEtaPrompt), fhPtThresIsolatedPrompt(), fhPtFracIsolatedPrompt(), fhPtSumIsolatedPrompt(), fhPtFragmentation(g.fhPtFragmentation),fhPhiFragmentation(g.fhPhiFragmentation),fhEtaFragmentation(g.fhEtaFragmentation), fhPtThresIsolatedFragmentation(), fhPtFracIsolatedFragmentation(), fhPtSumIsolatedFragmentation(), fhPtPi0Decay(g.fhPtPi0Decay),fhPhiPi0Decay(g.fhPhiPi0Decay),fhEtaPi0Decay(g.fhEtaPi0Decay), fhPtThresIsolatedPi0Decay(), fhPtFracIsolatedPi0Decay(), fhPtSumIsolatedPi0Decay(), fhPtOtherDecay(g.fhPtOtherDecay),fhPhiOtherDecay(g.fhPhiOtherDecay),fhEtaOtherDecay(g.fhEtaOtherDecay), fhPtThresIsolatedOtherDecay(), fhPtFracIsolatedOtherDecay(), fhPtSumIsolatedOtherDecay(), fhPtConversion(g. fhPtConversion),fhPhiConversion(g.fhPhiConversion),fhEtaConversion(g.fhEtaConversion), fhPtThresIsolatedConversion(), fhPtFracIsolatedConversion(), fhPtSumIsolatedConversion(), fhPtUnknown(g.fhPtUnknown),fhPhiUnknown(g.fhPhiUnknown),fhEtaUnknown(g.fhEtaUnknown), fhPtThresIsolatedUnknown(), fhPtFracIsolatedUnknown(), fhPtSumIsolatedUnknown() { // cpy ctor //Several IC for(Int_t i = 0; i < fNCones ; i++){ fConeSizes[i] = g.fConeSizes[i]; fhPtSumIsolated[i] = g.fhPtSumIsolated[i]; fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i]; fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i]; fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i]; fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i]; fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i]; fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i]; for(Int_t j = 0; j < fNPtThresFrac ; j++){ fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j]; fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j]; fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j]; fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j]; fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j]; fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j]; fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j]; fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j]; fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j]; fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j]; fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j]; fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j]; fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j]; fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j]; } } for(Int_t i = 0; i < fNPtThresFrac ; i++){ fPtFractions[i]= g.fPtFractions[i]; fPtThresholds[i]= g.fPtThresholds[i]; } } //_________________________________________________________________________ AliAnaGammaDirect & AliAnaGammaDirect::operator = (const AliAnaGammaDirect & g) { // assignment operator if(&g == this) return *this; fMakeIC = g.fMakeIC ; fReMakeIC = g.fReMakeIC ; fMakeSeveralIC = g.fMakeSeveralIC ; fMakeInvMass = g.fMakeInvMass ; fDetector = g.fDetector ; fhPtGamma = g.fhPtGamma ; fhPhiGamma = g.fhPhiGamma ; fhEtaGamma = g.fhEtaGamma ; fhConeSumPt = g.fhConeSumPt ; fhPtPrompt = g.fhPtPrompt; fhPhiPrompt = g.fhPhiPrompt; fhEtaPrompt = g.fhEtaPrompt; fhPtFragmentation = g.fhPtFragmentation; fhPhiFragmentation = g.fhPhiFragmentation; fhEtaFragmentation = g.fhEtaFragmentation; fhPtPi0Decay = g.fhPtPi0Decay; fhPhiPi0Decay = g.fhPhiPi0Decay; fhEtaPi0Decay = g.fhEtaPi0Decay; fhPtOtherDecay = g.fhPtOtherDecay; fhPhiOtherDecay = g.fhPhiOtherDecay; fhEtaOtherDecay = g.fhEtaOtherDecay; fhPtConversion = g. fhPtConversion; fhPhiConversion = g.fhPhiConversion; fhEtaConversion = g.fhEtaConversion; fhPtUnknown = g.fhPtUnknown; fhPhiUnknown = g.fhPhiUnknown; fhEtaUnknown = g.fhEtaUnknown; //Several IC fNCones = g.fNCones ; fNPtThresFrac = g.fNPtThresFrac ; for(Int_t i = 0; i < fNCones ; i++){ fConeSizes[i] = g.fConeSizes[i]; fhPtSumIsolated[i] = g.fhPtSumIsolated[i] ; fhPtSumIsolatedPrompt[i] = g.fhPtSumIsolatedPrompt[i]; fhPtSumIsolatedFragmentation[i] = g.fhPtSumIsolatedFragmentation[i]; fhPtSumIsolatedPi0Decay[i] = g.fhPtSumIsolatedPi0Decay[i]; fhPtSumIsolatedOtherDecay[i] = g.fhPtSumIsolatedOtherDecay[i]; fhPtSumIsolatedConversion[i] = g.fhPtSumIsolatedConversion[i]; fhPtSumIsolatedUnknown[i] = g.fhPtSumIsolatedUnknown[i]; for(Int_t j = 0; j < fNPtThresFrac ; j++){ fhPtThresIsolated[i][j] = g.fhPtThresIsolated[i][j] ; fhPtFracIsolated[i][j] = g.fhPtFracIsolated[i][j] ; fhPtThresIsolatedPrompt[i][j] = g.fhPtThresIsolatedPrompt[i][j]; fhPtThresIsolatedFragmentation[i][j] = g.fhPtThresIsolatedFragmentation[i][j]; fhPtThresIsolatedPi0Decay[i][j] = g.fhPtThresIsolatedPi0Decay[i][j]; fhPtThresIsolatedOtherDecay[i][j] = g.fhPtThresIsolatedOtherDecay[i][j]; fhPtThresIsolatedConversion[i][j] = g.fhPtThresIsolatedConversion[i][j]; fhPtThresIsolatedUnknown[i][j] = g.fhPtThresIsolatedUnknown[i][j]; fhPtFracIsolatedPrompt[i][j] = g.fhPtFracIsolatedPrompt[i][j]; fhPtFracIsolatedFragmentation[i][j] = g.fhPtFracIsolatedFragmentation[i][j]; fhPtFracIsolatedPi0Decay[i][j] = g.fhPtFracIsolatedPi0Decay[i][j]; fhPtFracIsolatedOtherDecay[i][j] = g.fhPtFracIsolatedOtherDecay[i][j]; fhPtFracIsolatedConversion[i][j] = g.fhPtFracIsolatedConversion[i][j]; fhPtFracIsolatedUnknown[i][j] = g.fhPtFracIsolatedUnknown[i][j]; } } for(Int_t i = 0; i < fNPtThresFrac ; i++){ fPtThresholds[i]= g.fPtThresholds[i]; fPtFractions[i]= g.fPtFractions[i]; } return *this; } //____________________________________________________________________________ AliAnaGammaDirect::~AliAnaGammaDirect() { //dtor //do not delete histograms delete [] fConeSizes ; delete [] fPtThresholds ; delete [] fPtFractions ; } //_________________________________________________________________________ Bool_t AliAnaGammaDirect::CheckInvMass(const Int_t icalo,const TLorentzVector mom, Double_t *vertex, TClonesArray * pl){ //Search if there is a companion decay photon to the candidate // and discard it in such case TLorentzVector mom2 ; for(Int_t jcalo = 0; jcalo < pl->GetEntriesFast(); jcalo++){ if(icalo==jcalo) continue ; AliAODCaloCluster * calo = dynamic_cast (pl->At(jcalo)); //Cluster selection, not charged, with photon id and in fidutial cut, fill TLorentz if(!SelectCluster(calo, vertex, mom2)) continue ; //Select good pair (good phit, pt cuts, aperture and invariant mass) if(GetNeutralMesonSelection()->SelectPair(mom, mom2)){ if(GetDebug()>1)printf("Selected gamma pair: pt %f, phi %f, eta%f",(mom+mom2).Pt(), (mom+mom2).Phi(), (mom+mom2).Eta()); return kTRUE ; } }//loop return kFALSE; } //_________________________________________________________________________ Int_t AliAnaGammaDirect::CheckOrigin(const Int_t label){ //Play with the MC stack if available //Check origin of the candidates, good for PYTHIA AliStack * stack = GetMCStack() ; 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(label < 8 ) printf("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 printf("Parent with label %d\n",iParent); //return tag if(mPdg == 22){ if(mStatus == 1){ if(iParent < 8) { if(pPdg == 22) return kPrompt; else return kFragmentation; } else if(pStatus == 11){ if(pPdg == 111) return kPi0Decay ; else if (pPdg == 321) return kEtaDecay ; else return kOtherDecay ; } }//Status 1 : Pythia generated else if(mStatus == 0){ if(pPdg ==22 || pPdg ==11) return kConversion ; if(pPdg == 111) return kPi0Decay ; else if (pPdg == 221) return kEtaDecay ; else return kOtherDecay ; }//status 0 : geant generated }//Mother gamma else if(mPdg == 111) return kPi0 ; else if(mPdg == 221) return kEta ; else if(mPdg ==11){ if(mStatus == 0) return kConversion ; else return kElectron ; } else return kUnknown; }//Good label value else{ if(label < 0 ) printf("*** bad label or no stack ***: label %d \n", label); if(label >= stack->GetNtrack()) printf("*** large label ***: label %d, n tracks %d \n", label, stack->GetNtrack()); return kUnknown; }//Bad label return kUnknown; } //________________________________________________________________________ TList * AliAnaGammaDirect::GetCreateOutputObjects() { // Create histograms to be saved in output file and // store them in outputContainer TList * outputContainer = new TList() ; outputContainer->SetName("DirectGammaHistos") ; //Histograms of highest gamma identified in Event fhPtGamma = new TH1F("hPtGamma","Number of #gamma over calorimeter",240,0,120); fhPtGamma->SetYTitle("N"); fhPtGamma->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtGamma) ; fhPhiGamma = new TH2F ("hPhiGamma","#phi_{#gamma}",200,0,120,200,0,7); fhPhiGamma->SetYTitle("#phi"); fhPhiGamma->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiGamma) ; fhEtaGamma = new TH2F ("hEtaGamma","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaGamma->SetYTitle("#eta"); fhEtaGamma->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaGamma) ; if(!fMakeSeveralIC){ fhConeSumPt = new TH2F ("hConePtSum","#Sigma p_{T} in cone ",200,0,120,100,0,100); fhConeSumPt->SetYTitle("#Sigma p_{T}"); fhConeSumPt->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhConeSumPt) ; } if(IsDataMC()){ fhPtPrompt = new TH1F("hPtPrompt","Number of #gamma over calorimeter",240,0,120); fhPtPrompt->SetYTitle("N"); fhPtPrompt->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtPrompt) ; fhPhiPrompt = new TH2F ("hPhiPrompt","#phi_{#gamma}",200,0,120,200,0,7); fhPhiPrompt->SetYTitle("#phi"); fhPhiPrompt->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiPrompt) ; fhEtaPrompt = new TH2F ("hEtaPrompt","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaPrompt->SetYTitle("#eta"); fhEtaPrompt->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaPrompt) ; fhPtFragmentation = new TH1F("hPtFragmentation","Number of #gamma over calorimeter",240,0,120); fhPtFragmentation->SetYTitle("N"); fhPtFragmentation->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtFragmentation) ; fhPhiFragmentation = new TH2F ("hPhiFragmentation","#phi_{#gamma}",200,0,120,200,0,7); fhPhiFragmentation->SetYTitle("#phi"); fhPhiFragmentation->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiFragmentation) ; fhEtaFragmentation = new TH2F ("hEtaFragmentation","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaFragmentation->SetYTitle("#eta"); fhEtaFragmentation->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaFragmentation) ; fhPtPi0Decay = new TH1F("hPtPi0Decay","Number of #gamma over calorimeter",240,0,120); fhPtPi0Decay->SetYTitle("N"); fhPtPi0Decay->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtPi0Decay) ; fhPhiPi0Decay = new TH2F ("hPhiPi0Decay","#phi_{#gamma}",200,0,120,200,0,7); fhPhiPi0Decay->SetYTitle("#phi"); fhPhiPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiPi0Decay) ; fhEtaPi0Decay = new TH2F ("hEtaPi0Decay","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaPi0Decay->SetYTitle("#eta"); fhEtaPi0Decay->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaPi0Decay) ; fhPtOtherDecay = new TH1F("hPtOtherDecay","Number of #gamma over calorimeter",240,0,120); fhPtOtherDecay->SetYTitle("N"); fhPtOtherDecay->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtOtherDecay) ; fhPhiOtherDecay = new TH2F ("hPhiOtherDecay","#phi_{#gamma}",200,0,120,200,0,7); fhPhiOtherDecay->SetYTitle("#phi"); fhPhiOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiOtherDecay) ; fhEtaOtherDecay = new TH2F ("hEtaOtherDecay","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaOtherDecay->SetYTitle("#eta"); fhEtaOtherDecay->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaOtherDecay) ; fhPtConversion = new TH1F("hPtConversion","Number of #gamma over calorimeter",240,0,120); fhPtConversion->SetYTitle("N"); fhPtConversion->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtConversion) ; fhPhiConversion = new TH2F ("hPhiConversion","#phi_{#gamma}",200,0,120,200,0,7); fhPhiConversion->SetYTitle("#phi"); fhPhiConversion->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiConversion) ; fhEtaConversion = new TH2F ("hEtaConversion","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaConversion->SetYTitle("#eta"); fhEtaConversion->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaConversion) ; fhPtUnknown = new TH1F("hPtUnknown","Number of #gamma over calorimeter",240,0,120); fhPtUnknown->SetYTitle("N"); fhPtUnknown->SetXTitle("p_{T #gamma}(GeV/c)"); outputContainer->Add(fhPtUnknown) ; fhPhiUnknown = new TH2F ("hPhiUnknown","#phi_{#gamma}",200,0,120,200,0,7); fhPhiUnknown->SetYTitle("#phi"); fhPhiUnknown->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhPhiUnknown) ; fhEtaUnknown = new TH2F ("hEtaUnknown","#phi_{#gamma}",200,0,120,200,-0.8,0.8); fhEtaUnknown->SetYTitle("#eta"); fhEtaUnknown->SetXTitle("p_{T #gamma} (GeV/c)"); outputContainer->Add(fhEtaUnknown) ; }//Histos with MC if(fMakeSeveralIC){ char name[128]; char title[128]; for(Int_t icone = 0; iconeSetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolated[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolated[icone]) ; if(IsDataMC()){ sprintf(name,"hPtSumIsolatedPrompt_Cone_%d",icone); sprintf(title,"Candidate Prompt cone sum p_{T} for cone size %d vs candidate p_{T}",icone); fhPtSumIsolatedPrompt[icone] = new TH2F(name, title,240,0,120,120,0,10); fhPtSumIsolatedPrompt[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolatedPrompt[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolatedPrompt[icone]) ; sprintf(name,"hPtSumIsolatedFragmentation_Cone_%d",icone); sprintf(title,"Candidate Fragmentation cone sum p_{T} for cone size %d vs candidate p_{T}",icone); fhPtSumIsolatedFragmentation[icone] = new TH2F(name, title,240,0,120,120,0,10); fhPtSumIsolatedFragmentation[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolatedFragmentation[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolatedFragmentation[icone]) ; sprintf(name,"hPtSumIsolatedPi0Decay_Cone_%d",icone); sprintf(title,"Candidate Pi0Decay cone sum p_{T} for cone size %d vs candidate p_{T}",icone); fhPtSumIsolatedPi0Decay[icone] = new TH2F(name, title,240,0,120,120,0,10); fhPtSumIsolatedPi0Decay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolatedPi0Decay[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolatedPi0Decay[icone]) ; sprintf(name,"hPtSumIsolatedOtherDecay_Cone_%d",icone); sprintf(title,"Candidate OtherDecay cone sum p_{T} for cone size %d vs candidate p_{T}",icone); fhPtSumIsolatedOtherDecay[icone] = new TH2F(name, title,240,0,120,120,0,10); fhPtSumIsolatedOtherDecay[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolatedOtherDecay[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolatedOtherDecay[icone]) ; sprintf(name,"hPtSumIsolatedConversion_Cone_%d",icone); sprintf(title,"Candidate Conversion cone sum p_{T} for cone size %d vs candidate p_{T}",icone); fhPtSumIsolatedConversion[icone] = new TH2F(name, title,240,0,120,120,0,10); fhPtSumIsolatedConversion[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolatedConversion[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolatedConversion[icone]) ; sprintf(name,"hPtSumIsolatedUnknown_Cone_%d",icone); sprintf(title,"Candidate Unknown cone sum p_{T} for cone size %d vs candidate p_{T}",icone); fhPtSumIsolatedUnknown[icone] = new TH2F(name, title,240,0,120,120,0,10); fhPtSumIsolatedUnknown[icone]->SetYTitle("#Sigma p_{T} (GeV/c)"); fhPtSumIsolatedUnknown[icone]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtSumIsolatedUnknown[icone]) ; }//Histos with MC for(Int_t ipt = 0; iptSetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolated[icone][ipt]) ; sprintf(name,"hPtFracIsol_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolated[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolated[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolated[icone][ipt]) ; if(IsDataMC()){ sprintf(name,"hPtThresIsolPrompt_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtThresIsolatedPrompt[icone][ipt] = new TH1F(name, title,240,0,120); fhPtThresIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolatedPrompt[icone][ipt]) ; sprintf(name,"hPtFracIsolPrompt_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Prompt p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolatedPrompt[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolatedPrompt[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolatedPrompt[icone][ipt]) ; sprintf(name,"hPtThresIsolFragmentation_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtThresIsolatedFragmentation[icone][ipt] = new TH1F(name, title,240,0,120); fhPtThresIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolatedFragmentation[icone][ipt]) ; sprintf(name,"hPtFracIsolFragmentation_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Fragmentation p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolatedFragmentation[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolatedFragmentation[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolatedFragmentation[icone][ipt]) ; sprintf(name,"hPtThresIsolPi0Decay_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtThresIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,240,0,120); fhPtThresIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolatedPi0Decay[icone][ipt]) ; sprintf(name,"hPtFracIsolPi0Decay_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Pi0Decay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolatedPi0Decay[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolatedPi0Decay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolatedPi0Decay[icone][ipt]) ; sprintf(name,"hPtThresIsolOtherDecay_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtThresIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,240,0,120); fhPtThresIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolatedOtherDecay[icone][ipt]) ; sprintf(name,"hPtFracIsolOtherDecay_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate OtherDecay p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolatedOtherDecay[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolatedOtherDecay[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolatedOtherDecay[icone][ipt]) ; sprintf(name,"hPtThresIsolConversion_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtThresIsolatedConversion[icone][ipt] = new TH1F(name, title,240,0,120); fhPtThresIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolatedConversion[icone][ipt]) ; sprintf(name,"hPtFracIsolConversion_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Conversion p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolatedConversion[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolatedConversion[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolatedConversion[icone][ipt]) ; sprintf(name,"hPtThresIsolUnknown_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtThresIsolatedUnknown[icone][ipt] = new TH1F(name, title,240,0,120); fhPtThresIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtThresIsolatedUnknown[icone][ipt]) ; sprintf(name,"hPtFracIsolUnknown_Cone_%d_Pt%d",icone,ipt); sprintf(title,"Isolated candidate Unknown p_{T} distribution for cone size %d and p_{T}^{th} %d",icone,ipt); fhPtFracIsolatedUnknown[icone][ipt] = new TH1F(name, title,240,0,120); fhPtFracIsolatedUnknown[icone][ipt]->SetXTitle("p_{T} (GeV/c)"); outputContainer->Add(fhPtFracIsolatedUnknown[icone][ipt]) ; }//Histos with MC }//icone loop }//ipt loop } //Keep neutral meson selection histograms if requiered //Setting done in AliNeutralMesonSelection if(fMakeInvMass && GetNeutralMesonSelection()){ TList * nmsHistos = GetNeutralMesonSelection()->GetCreateOutputObjects() ; cout<<"NMSHistos "<< nmsHistos<AreNeutralMesonSelectionHistosKept()) for(Int_t i = 0; i < nmsHistos->GetEntries(); i++) outputContainer->Add(nmsHistos->At(i)) ; } return outputContainer ; } //__________________________________________________________________ void AliAnaGammaDirect::MakeAnalysisFillAOD() { //Do analysis and fill aods //Search for the isolated photon in fDetector with pt > GetMinPt() //Fill CaloClusters if working with ESDs //if(GetReader()->GetDataType() == AliCaloTrackReader::kESD) ConnectAODCaloClusters(); Int_t n = 0, nfrac = 0; Bool_t isolated = kFALSE ; Float_t coneptsum = 0 ; TClonesArray * pl = new TClonesArray; //Get vertex for photon momentum calculation Double_t vertex[]={0,0,0} ; //vertex ; if(!GetReader()->GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(vertex); //Select the detector of the photon if(fDetector == "PHOS") pl = GetAODPHOS(); else if (fDetector == "EMCAL") pl = GetAODEMCAL(); //cout<<"Number of entries "<GetEntriesFast()<GetEntriesFast(); icalo++){ AliAODCaloCluster * calo = dynamic_cast (pl->At(icalo)); //Cluster selection, not charged, with photon id and in fidutial cut if(!SelectCluster(calo,vertex,mom)) continue ; //If too small pt, skip if(mom.Pt() < GetMinPt() || mom.Pt() > GetMaxPt() ) continue ; //Play with the MC stack if available Int_t tag =-1; if(IsDataMC()){ //Check origin of the candidates tag = CheckOrigin(calo->GetLabel(0)); if(GetDebug() > 0) printf("Origin of candidate %d\n",tag); }//Work with stack also //Check invariant mass Bool_t decay = kFALSE ; if(fMakeInvMass) decay = CheckInvMass(icalo,mom,vertex,pl); if(decay) continue ; //Create AOD for analysis AliAODParticleCorrelation ph = AliAODParticleCorrelation(mom); ph.SetLabel(calo->GetLabel(0)); ph.SetPdg(AliCaloPID::kPhoton); ph.SetDetector(fDetector); ph.SetTag(tag); if(fMakeIC){ n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0; GetIsolationCut()->MakeIsolationCut((TSeqCollection*)GetAODCTS(), (TSeqCollection*)pl, vertex, kTRUE, &ph,icalo,-1, n,nfrac,coneptsum, isolated); if(isolated){ //cout<<"Isolated : E "< 1) printf("End fill AODs "); } //__________________________________________________________________ void AliAnaGammaDirect::MakeAnalysisFillHistograms() { //Do analysis and fill histograms Int_t n = 0, nfrac = 0; Bool_t isolated = kFALSE ; Float_t coneptsum = 0 ; //cout<<"MakeAnalysisFillHistograms"<GetDataType()== AliCaloTrackReader::kMC) GetReader()->GetVertex(v); //Loop on stored AOD photons Int_t naod = GetAODBranch()->GetEntriesFast(); if(GetDebug() > 0) printf("histo aod branch entries %d\n", naod); for(Int_t iaod = 0; iaod < naod ; iaod++){ AliAODParticleCorrelation* ph = dynamic_cast (GetAODBranch()->At(iaod)); Int_t pdg = ph->GetPdg(); //Only isolate photons in detector fDetector if(pdg != AliCaloPID::kPhoton || ph->GetDetector() != fDetector) continue; if(fMakeSeveralIC) { //Analysis of multiple IC at same time MakeSeveralICAnalysis(ph,v); continue; } else if(fReMakeIC){ //In case a more strict IC is needed in the produced AOD n=0; nfrac = 0; isolated = kFALSE; coneptsum = 0; GetIsolationCut()->MakeIsolationCut((TSeqCollection*)ph->GetRefIsolationConeTracks(), (TSeqCollection*)ph->GetRefIsolationConeClusters(), v, kFALSE, ph,-1,-1, n,nfrac,coneptsum, isolated); } //Fill histograms (normal case) if(!fReMakeIC || (fReMakeIC && isolated)){ //printf("Isolated Gamma: pt %f, phi %f, eta %f", ph->Pt(),ph->Phi(),ph->Eta()) ; //Fill prompt gamma histograms Float_t ptcluster = ph->Pt(); Float_t phicluster = ph->Phi(); Float_t etacluster = ph->Eta(); fhPtGamma ->Fill(ptcluster); fhPhiGamma ->Fill(ptcluster,phicluster); fhEtaGamma ->Fill(ptcluster,etacluster); fhConeSumPt->Fill(ptcluster,coneptsum); if(IsDataMC()){ Int_t tag =ph->GetTag(); if(tag == kPrompt){ fhPtPrompt ->Fill(ptcluster); fhPhiPrompt ->Fill(ptcluster,phicluster); fhEtaPrompt ->Fill(ptcluster,etacluster); } else if(tag==kFragmentation) { fhPtFragmentation ->Fill(ptcluster); fhPhiFragmentation ->Fill(ptcluster,phicluster); fhEtaFragmentation ->Fill(ptcluster,etacluster); } else if(tag==kPi0Decay) { fhPtPi0Decay ->Fill(ptcluster); fhPhiPi0Decay ->Fill(ptcluster,phicluster); fhEtaPi0Decay ->Fill(ptcluster,etacluster); } else if(tag==kEtaDecay || tag==kOtherDecay) { fhPtOtherDecay ->Fill(ptcluster); fhPhiOtherDecay ->Fill(ptcluster,phicluster); fhEtaOtherDecay ->Fill(ptcluster,etacluster); } else if(tag==kConversion) { fhPtConversion ->Fill(ptcluster); fhPhiConversion ->Fill(ptcluster,phicluster); fhEtaConversion ->Fill(ptcluster,etacluster); } else{ fhPtUnknown ->Fill(ptcluster); fhPhiUnknown ->Fill(ptcluster,phicluster); fhEtaUnknown ->Fill(ptcluster,etacluster); } }//Histograms with MC } }// aod loop } //____________________________________________________________________________ void AliAnaGammaDirect::InitParameters() { //Initialize the parameters of the analysis. fDetector = "PHOS" ; fMakeIC = kTRUE; fReMakeIC = kFALSE ; fMakeSeveralIC = kFALSE ; fMakeInvMass = kFALSE ; //----------- Several IC----------------- fNCones = 5 ; fNPtThresFrac = 6 ; fConeSizes[0] = 0.1; fConeSizes[1] = 0.2; fConeSizes[2] = 0.3; fConeSizes[3] = 0.4; fConeSizes[4] = 0.5; fPtThresholds[0]=0.; fPtThresholds[1]=1.; fPtThresholds[2]=2.; fPtThresholds[3]=3.; fPtThresholds[4]=4.;fPtThresholds[5]=5.; fPtFractions[0]=0.05; fPtFractions[1]=0.075; fPtFractions[2]=0.1; fPtFractions[3]=1.25; fPtFractions[4]=1.5;fPtFractions[5]=2.; } //__________________________________________________________________ void AliAnaGammaDirect::MakeSeveralICAnalysis(AliAODParticleCorrelation* ph, Double_t v[3]) { //Isolation Cut Analysis for both methods and different pt cuts and cones Float_t ptC = ph->Pt(); Float_t phiC = ph->Phi(); Float_t etaC = ph->Eta(); fhPtGamma->Fill(ptC); fhPhiGamma->Fill(ptC,ph->Phi()); fhEtaGamma->Fill(ptC,ph->Eta()); Int_t tag =ph->GetTag(); if(IsDataMC()){ if(tag == kPrompt){ fhPtPrompt ->Fill(ptC); fhPhiPrompt ->Fill(ptC,phiC); fhEtaPrompt ->Fill(ptC,etaC); } else if(tag==kFragmentation) { fhPtFragmentation ->Fill(ptC); fhPhiFragmentation ->Fill(ptC,phiC); fhEtaFragmentation ->Fill(ptC,etaC); } else if(tag==kPi0Decay) { fhPtPi0Decay ->Fill(ptC); fhPhiPi0Decay ->Fill(ptC,phiC); fhEtaPi0Decay ->Fill(ptC,etaC); } else if(tag==kEtaDecay || tag==kOtherDecay) { fhPtOtherDecay ->Fill(ptC); fhPhiOtherDecay ->Fill(ptC,phiC); fhEtaOtherDecay ->Fill(ptC,etaC); } else if(tag==kConversion) { fhPtConversion ->Fill(ptC); fhPhiConversion ->Fill(ptC,phiC); fhEtaConversion ->Fill(ptC,etaC); } else{ fhPtUnknown ->Fill(ptC); fhPhiUnknown ->Fill(ptC,phiC); fhEtaUnknown ->Fill(ptC,etaC); } }//Histograms with MC //Keep original setting used when filling AODs, reset at end of analysis Float_t ptthresorg = GetIsolationCut()->GetPtThreshold(); Float_t ptfracorg = GetIsolationCut()->GetPtFraction(); Float_t rorg = GetIsolationCut()->GetConeSize(); Float_t coneptsum = 0 ; Int_t n[10][10];//[fNCones][fNPtThresFrac]; Int_t nfrac[10][10];//[fNCones][fNPtThresFrac]; Bool_t isolated = kFALSE; for(Int_t icone = 0; iconeSetConeSize(fConeSizes[icone]); coneptsum = 0 ; for(Int_t ipt = 0; iptSetPtThreshold(fPtThresholds[ipt]); GetIsolationCut()->MakeIsolationCut((TSeqCollection*)ph->GetRefIsolationConeTracks(), (TSeqCollection*)ph->GetRefIsolationConeClusters(), v, kFALSE, ph,-1,-1, n[icone][ipt],nfrac[icone][ipt],coneptsum, isolated); if(n[icone][ipt] == 0) { fhPtThresIsolated[icone][ipt]->Fill(ptC); if(IsDataMC()){ if(tag==kPrompt) fhPtThresIsolatedPrompt[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kConversion) fhPtThresIsolatedConversion[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kFragmentation) fhPtThresIsolatedFragmentation[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kPi0Decay) fhPtThresIsolatedPi0Decay[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kOtherDecay || tag==kEtaDecay) fhPtThresIsolatedOtherDecay[icone][ipt]->Fill(ptC,coneptsum) ; else fhPtThresIsolatedUnknown[icone][ipt]->Fill(ptC,coneptsum) ; } } if(nfrac[icone][ipt] == 0) { fhPtFracIsolated[icone][ipt]->Fill(ptC); if(IsDataMC()){ if(tag==kPrompt) fhPtFracIsolatedPrompt[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kConversion) fhPtFracIsolatedConversion[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kFragmentation) fhPtFracIsolatedFragmentation[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kPi0Decay) fhPtFracIsolatedPi0Decay[icone][ipt]->Fill(ptC,coneptsum) ; else if(tag==kOtherDecay || tag==kEtaDecay) fhPtFracIsolatedOtherDecay[icone][ipt]->Fill(ptC,coneptsum) ; else fhPtFracIsolatedUnknown[icone][ipt]->Fill(ptC,coneptsum) ; } } }//pt thresh loop fhPtSumIsolated[icone]->Fill(ptC,coneptsum) ; if(IsDataMC()){ if(tag==kPrompt) fhPtSumIsolatedPrompt[icone]->Fill(ptC,coneptsum) ; else if(tag==kConversion) fhPtSumIsolatedConversion[icone]->Fill(ptC,coneptsum) ; else if(tag==kFragmentation) fhPtSumIsolatedFragmentation[icone]->Fill(ptC,coneptsum) ; else if(tag==kPi0Decay) fhPtSumIsolatedPi0Decay[icone]->Fill(ptC,coneptsum) ; else if(tag==kOtherDecay || tag==kEtaDecay) fhPtSumIsolatedOtherDecay[icone]->Fill(ptC,coneptsum) ; else fhPtSumIsolatedUnknown[icone]->Fill(ptC,coneptsum) ; } }//cone size loop //Reset original parameters for AOD analysis GetIsolationCut()->SetPtThreshold(ptthresorg); GetIsolationCut()->SetPtFraction(ptfracorg); GetIsolationCut()->SetConeSize(rorg); } //__________________________________________________________________ void AliAnaGammaDirect::Print(const Option_t * opt) const { //Print some relevant parameters set for the analysis if(! opt) return; printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ; printf("Min Gamma pT = %2.2f\n", GetMinPt()) ; printf("Max Gamma pT = %3.2f\n", GetMaxPt()) ; // if(IsCaloPIDOn())printf("Check PID \n") ; // if(IsCaloPIDRecalculationOn())printf("Recalculate PID \n") ; // if(IsFidutialCutOn())printf("Check Fidutial cut \n") ; printf("Make Isolation = %d \n", fMakeIC) ; printf("ReMake Isolation = %d \n", fReMakeIC) ; printf("Make Several Isolation = %d \n", fMakeSeveralIC) ; if(fMakeSeveralIC){ printf("N Cone Sizes = %d\n", fNCones) ; printf("Cone Sizes = \n") ; for(Int_t i = 0; i < fNCones; i++) printf(" %1.2f;", fConeSizes[i]) ; printf(" \n") ; printf("N pT thresholds/fractions = %d\n", fNPtThresFrac) ; printf(" pT thresholds = \n") ; for(Int_t i = 0; i < fNPtThresFrac; i++) printf(" %2.2f;", fPtThresholds[i]) ; printf(" \n") ; printf(" pT fractions = \n") ; for(Int_t i = 0; i < fNPtThresFrac; i++) printf(" %2.2f;", fPtFractions[i]) ; } printf(" \n") ; } //____________________________________________________________________________ Bool_t AliAnaGammaDirect::SelectCluster(AliAODCaloCluster * calo, Double_t vertex[3], TLorentzVector & mom){ //Select cluster depending on its pid and acceptance selections //Skip matched clusters with tracks if(calo->GetNTracksMatched() > 0) return kFALSE ; //Check PID calo->GetMomentum(mom,vertex);//Assume that come from vertex in straight line Int_t pdg = AliCaloPID::kPhoton; if(IsCaloPIDOn()){ //Get most probable PID, 2 options check PID weights (in MC this option is mandatory) //or redo PID, recommended option for EMCal. if(!IsCaloPIDRecalculationOn() || GetReader()->GetDataType() == AliCaloTrackReader::kMC ) pdg = GetCaloPID()->GetPdg(fDetector,calo->PID(),mom.E());//PID with weights else pdg = GetCaloPID()->GetPdg(fDetector,mom,calo->GetM02(),0,0,0,0);//PID recalculated if(GetDebug() > 1) printf("PDG of identified particle %d\n",pdg); //If it does not pass pid, skip if(pdg != AliCaloPID::kPhoton) return kFALSE ; } //Check acceptance selection if(IsFidutialCutOn()){ Bool_t in = GetFidutialCut()->IsInFidutialCut(mom,fDetector) ; if(! in ) return kFALSE ; } if(GetDebug() > 1) printf("Correlation photon selection cuts passed: pT %3.2f, pdg %d\n",mom.Pt(), pdg); return kTRUE; }