[u/mrichter/AliRoot.git] / PWGHF / hfe / AliAnalysisTaskElecHadronCorrel.cxx

/**************************************************************************
 * 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 is 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.                  *               
 **************************************************************************/


////////////////////////////////////////////////////////////////////////
//                                                                    //
//  Task for Heavy Flavour Electron-Hadron DeltaPhi Correlation       //
//  Non-Photonic Electron identified with Invariant mass              //
//  analysis methos in function  SelectPhotonicElectron               //
//  DeltaPhi calculated in function  ElectronHadCorrel                // 
//                                                                    //
//  Author: Deepa Thomas (Utrecht University)                         //
//                                                                    //
////////////////////////////////////////////////////////////////////////

#include "TChain.h"
#include "TTree.h"
#include "TH2F.h"
#include "TMath.h"
#include "TCanvas.h"
#include "THnSparse.h"
#include "TLorentzVector.h"
#include "TString.h"
#include "TFile.h"

#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"

#include "AliESDEvent.h"
#include "AliESDHandler.h"
#include "AliAODEvent.h"
#include "AliAODHandler.h"

#include "AliAnalysisTaskElecHadronCorrel.h"
#include "TGeoGlobalMagField.h"
#include "AliLog.h"
#include "AliAnalysisTaskSE.h"
#include "TRefArray.h"
#include "TVector.h"
#include "AliESDInputHandler.h"
#include "AliESDpid.h"
#include "AliESDtrackCuts.h"
#include "AliPhysicsSelection.h"
#include "AliESDCaloCluster.h"
#include "AliAODCaloCluster.h"
#include "AliEMCALRecoUtils.h"
#include "AliEMCALGeometry.h"
#include "AliGeomManager.h"
#include "stdio.h"
#include "TGeoManager.h"
#include "iostream"
#include "fstream"

#include "AliEMCALTrack.h"
#include "AliMagF.h"

#include "AliKFParticle.h"
#include "AliKFVertex.h"

#include "AliPID.h"
#include "AliPIDResponse.h"
#include "AliHFEcontainer.h"
#include "AliHFEcuts.h"
#include "AliHFEpid.h"
#include "AliHFEpidBase.h"
#include "AliHFEpidQAmanager.h"
#include "AliHFEtools.h"
#include "AliCFContainer.h"
#include "AliCFManager.h"

ClassImp(AliAnalysisTaskElecHadronCorrel)
//________________________________________________________________________
AliAnalysisTaskElecHadronCorrel::AliAnalysisTaskElecHadronCorrel(const char *name) 
  : AliAnalysisTaskSE(name)
  ,fESD(0)
  ,fGeom(0)
  ,fOutputList(0)
  ,fTrackCuts1(0)
  ,fTrackCuts2(0)
  ,fCuts(0)
  ,fIdentifiedAsOutInz(kFALSE)
  ,fPassTheEventCut(kFALSE)
  ,fRejectKinkMother(kFALSE)
  ,fVz(0.0)
  ,fCFM(0)	
  ,fPID(0)
  ,fPIDqa(0)	       
  ,fOpeningAngleCut(0.1)
  ,fInvmassCut(0.01)	
  ,fNoEvents(0)
  ,fTrkpt(0)
  ,fTrkEovPBef(0)	 
  ,fTrkEovPAft(0)	
  ,fdEdxBef(0)	 
  ,fdEdxAft(0)	 
  ,fInvmassLS(0)		
  ,fInvmassULS(0)		
  ,fOpeningAngleLS(0)	
  ,fOpeningAngleULS(0)	
  ,fSemiIncElecDphi(0) 	
  ,fPhotElecDphi(0)  	
  ,fInclusiveElecDphi(0)  	
  ,fDphiMassHigh(0)		
  ,fDphiULSMassLow(0)	
  ,fDphiLSMassLow(0)
  ,fPhotoElecPt(0)
  ,fSemiInclElecPt(0)
  ,fTrackPtBefTrkCuts(0)	 
  ,fTrackPtAftTrkCuts(0)
  ,fTPCnsigma(0)	 
{
  //Named constructor
  
  fPID = new AliHFEpid("hfePid");
  fTrackCuts1 = new AliESDtrackCuts();
  fTrackCuts2 = new AliESDtrackCuts();  

  // Define input and output slots here
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 id reserved by the base class for AOD
  // Output slot #1 writes into a TH1 container
  // DefineOutput(1, TH1I::Class());
  DefineOutput(1, TList::Class());
  //  DefineOutput(3, TTree::Class());
}

//________________________________________________________________________
AliAnalysisTaskElecHadronCorrel::AliAnalysisTaskElecHadronCorrel() 
  : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisElecHadCorrel")
  ,fESD(0)
  ,fGeom(0)  
  ,fOutputList(0)
  ,fTrackCuts1(0)
  ,fTrackCuts2(0)
  ,fCuts(0)
  ,fIdentifiedAsOutInz(kFALSE)
  ,fPassTheEventCut(kFALSE)
  ,fRejectKinkMother(kFALSE)
  ,fVz(0.0)
  ,fCFM(0)	
  ,fPID(0)       
  ,fPIDqa(0)	       
  ,fOpeningAngleCut(0.1)
  ,fInvmassCut(0.01)	
  ,fNoEvents(0)
  ,fTrkpt(0)
  ,fTrkEovPBef(0)	 
  ,fTrkEovPAft(0)	 
  ,fdEdxBef(0)	 
  ,fdEdxAft(0)	 
  ,fInvmassLS(0)		
  ,fInvmassULS(0)		
  ,fOpeningAngleLS(0)	
  ,fOpeningAngleULS(0)	
  ,fSemiIncElecDphi(0) 	
  ,fPhotElecDphi(0)  	
  ,fInclusiveElecDphi(0)  	
  ,fDphiMassHigh(0)		
  ,fDphiULSMassLow(0)	
  ,fDphiLSMassLow(0)
  ,fPhotoElecPt(0)
  ,fSemiInclElecPt(0)
  ,fTrackPtBefTrkCuts(0)	 
  ,fTrackPtAftTrkCuts(0)	 	  
  ,fTPCnsigma(0)	 
{
	//Default constructor
	fPID = new AliHFEpid("hfePid");

	fTrackCuts1 = new AliESDtrackCuts();
	fTrackCuts2 = new AliESDtrackCuts();

	// Constructor
	// Define input and output slots here
	// Input slot #0 works with a TChain
	DefineInput(0, TChain::Class());
	// Output slot #0 id reserved by the base class for AOD
	// Output slot #1 writes into a TH1 container
	// DefineOutput(1, TH1I::Class());
	DefineOutput(1, TList::Class());
	//DefineOutput(3, TTree::Class());
}
//_________________________________________

AliAnalysisTaskElecHadronCorrel::~AliAnalysisTaskElecHadronCorrel()
{
  //Destructor 
  
  delete fOutputList;
  delete fGeom;
  delete fPID;
  delete fCFM;
  delete fPIDqa;
  delete fTrackCuts1;
  delete fTrackCuts2;
}
//_________________________________________

void AliAnalysisTaskElecHadronCorrel::UserExec(Option_t*)
{
  //Main loop
  //Called for each event
  
  // create pointer to event
  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
  if (!fESD) {
    printf("ERROR: fESD not available\n");
    return;
  }
  
  if(!fCuts){
    AliError("HFE cuts not available");
    return;
  }
  
  if(!fPID->IsInitialized()){ 
    // Initialize PID with the given run number
    AliWarning("PID not initialised, get from Run no");
    fPID->InitializePID(fESD->GetRunNumber());
  }
  
  //-------trigger selection
  UInt_t res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
  if (res==0)
    return;
  
  if( (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kFastOnly) )
    return;
  
  if(!(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kEMC1) ) return;
  //if(!(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB) ) return;
  
  Int_t fNOtrks =  fESD->GetNumberOfTracks();
  const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
  
  Double_t pVtxZ = -999;
  pVtxZ = pVtx->GetZ();
  
  // Event cut
  //	if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
  
  //Check if its an LED event 
  if(IsLEDEvent()) return;
  
  if(TMath::Abs(pVtxZ)>10) return;
  fNoEvents->Fill(0);
  
  if(fNOtrks<2) return;
  
  AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
  if(!pidResponse){
    AliDebug(1, "Using default PID Response");
    pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class()); 
  }
  
  fPID->SetPIDResponse(pidResponse);
  
  fCFM->SetRecEventInfo(fESD);
  
  // Track loop 
  for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
    AliESDtrack* track = fESD->GetTrack(iTracks);
    if (!track) {
      printf("ERROR: Could not receive track %d\n", iTracks);
      continue;
    }
    
    fTrackPtBefTrkCuts->Fill(track->Pt());		
    // RecKine: ITSTPC cuts  
    if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
    
    //RecKink
    if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
      if(track->GetKinkIndex(0) != 0) continue;
    } 
    
    // RecPrim
    if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
    
    // HFEcuts: ITS layers cuts
    if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
    
    // HFE cuts: TPC PID cleanup
    if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
    
    fTrackPtAftTrkCuts->Fill(track->Pt());		
    
    Double_t fClsE = -999, p = -999, fEovP=-999, pt = -999, dEdx=-999, fTPCnSigma=0;
    // Track extrapolation
    Int_t fClsId = track->GetEMCALcluster();
    if(fClsId == -1) continue;
    AliESDCaloCluster *cluster = fESD->GetCaloCluster(fClsId);
    
    pt = track->Pt();
    if(pt<1) continue;
    fTrkpt->Fill(pt);
    fClsE = cluster->E();
    p = track->P();
    dEdx = track->GetTPCsignal();
    fEovP = fClsE/p;
    fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
    fdEdxBef->Fill(p,dEdx);
    fTPCnsigma->Fill(p,fTPCnSigma);
    
    if(fTPCnSigma >= 1.5 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,fEovP);
    Int_t pidpassed = 1;
    //--- track accepted, do PID
    AliHFEpidObject hfetrack;
    hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
    hfetrack.SetRecTrack(track);
    hfetrack.SetPP();
    if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
    
    if(pidpassed==0) continue;
    fTrkEovPAft->Fill(pt,fEovP);
    fdEdxAft->Fill(p,dEdx);
    
    Bool_t fFlagPhotonicElec = kFALSE;
    // select photonic electron
    SelectPhotonicElectron(iTracks,track,fFlagPhotonicElec);
    
    //Inclusive electron-hadron correlation
    ElectronHadCorrel(iTracks, track, fInclusiveElecDphi);
    
    // photonic electron
    if(fFlagPhotonicElec){
      //Electron hadron correlation
      ElectronHadCorrel(iTracks, track, fPhotElecDphi);
      fPhotoElecPt->Fill(pt);
    }
    
    // Semi inclusive electron 
    if(!fFlagPhotonicElec){
      //Electron hadron correlation
      ElectronHadCorrel(iTracks, track, fSemiIncElecDphi);
      fSemiInclElecPt->Fill(pt);
    }
  }
  
  
  PostData(1, fOutputList);
}
//_________________________________________
void AliAnalysisTaskElecHadronCorrel::UserCreateOutputObjects()
{
  //Create histograms
  TGeoManager::Import("geometry.root");
  fGeom =  AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
  
  //--------Initialize PID
  fPID->SetHasMCData(kFALSE);
  if(!fPID->GetNumberOfPIDdetectors()) 
    {
      fPID->AddDetector("TPC", 0);
      fPID->AddDetector("EMCAL", 1);
    }
  
  fPID->SortDetectors(); 
  fPIDqa = new AliHFEpidQAmanager();
  fPIDqa->Initialize(fPID);
  
  //--------Initialize correction Framework and Cuts
  fCFM = new AliCFManager;
  const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
  fCFM->SetNStepParticle(kNcutSteps);
  for(Int_t istep = 0; istep < kNcutSteps; istep++)
    fCFM->SetParticleCutsList(istep, NULL);
  
  if(!fCuts){
    AliWarning("Cuts not available. Default cuts will be used");
    fCuts = new AliHFEcuts;
    fCuts->CreateStandardCuts();
  }
  fCuts->Initialize(fCFM);
  
  //---------Output Tlist
  fOutputList = new TList();
  fOutputList->SetOwner();
  fOutputList->Add(fPIDqa->MakeList("PIDQA"));
  
  fNoEvents = new TH1F("fNoEvents","",1,0,1) ;
  fOutputList->Add(fNoEvents);
  
  fTrkpt = new TH1F("fTrkpt","track pt",1000,0,50);
  fOutputList->Add(fTrkpt);
  
  fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",1000,0,50);
  fOutputList->Add(fTrackPtBefTrkCuts);
  
  fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",1000,0,50);
  fOutputList->Add(fTrackPtAftTrkCuts);
  
  fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",1000,0,50,200,-10,10);
  fOutputList->Add(fTPCnsigma);
  
  fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",1000,0,50,100,0,2);
  fOutputList->Add(fTrkEovPBef);
  
  fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",1000,0,50,100,0,2);
  fOutputList->Add(fTrkEovPAft);
  
  fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",1000,0,50,150,0,150);
  fOutputList->Add(fdEdxBef);
  
  fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",1000,0,50,150,0,150);
  fOutputList->Add(fdEdxAft);
  
  fInvmassLS = new TH1F("fInvmassLS", "Inv mass of LS (e,e); mass(GeV/c^2); counts;", 1000,0,0.5);
  fOutputList->Add(fInvmassLS);
  
  fInvmassULS = new TH1F("fInvmassULS", "Inv mass of ULS (e,e); mass(GeV/c^2); counts;", 1000,0,0.5);
  fOutputList->Add(fInvmassULS);
  
  fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
  fOutputList->Add(fOpeningAngleLS);
  
  fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
  fOutputList->Add(fOpeningAngleULS);
  
  fSemiIncElecDphi = new TH2F("fSemiIncElecDphi", "Semi Inclusive elec-had Dphi correlation",200,0,20,100,-3.14,3.14);
  fOutputList->Add(fSemiIncElecDphi);
  
  fPhotElecDphi = new TH2F("fPhotElecDphi", "Photon elec-had Dphi correlation",200,0,20,100,-3.14,3.14);
  fOutputList->Add(fPhotElecDphi);
  
  fInclusiveElecDphi = new TH2F("fInclusiveElecDphi", "Inclusive elec-had Dphi correlation",200,0,20,100,-3.14,3.14);
  fOutputList->Add(fInclusiveElecDphi);
  
  fDphiMassHigh = new TH2F("fDphiMassHigh", "e-h Dphi LS+ULS, mass>0.01",200,0,20,100,-3.14,3.14);
  fOutputList->Add(fDphiMassHigh);
  
  fDphiULSMassLow = new TH2F("fDphiULSMassLow", "e-h Dphi ULS, mass<0.01",200,0,20,100,-3.14,3.14);
  fOutputList->Add(fDphiULSMassLow);
  
  fDphiLSMassLow = new TH2F("fDphiLSMassLow", "e-h Dphi LS, mass<0.01",200,0,20,100,-3.14,3.14);
  fOutputList->Add(fDphiLSMassLow);
  
  fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",1000,0,100);
  fOutputList->Add(fPhotoElecPt);
  
  fSemiInclElecPt = new TH1F("fSemiInclElecPt", "Semi-inclusive electron pt",1000,0,100);
  fOutputList->Add(fSemiInclElecPt);
  PostData(1,fOutputList);
}

//________________________________________________________________________
void AliAnalysisTaskElecHadronCorrel::Terminate(Option_t *)
{
  // Info("Terminate");
	AliAnalysisTaskSE::Terminate();
}

//________________________________________________________________________
Bool_t AliAnalysisTaskElecHadronCorrel::ProcessCutStep(Int_t cutStep, AliVParticle *track)
{
  // Check single track cuts for a given cut step
  const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
  if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
  return kTRUE;
}
//_________________________________________
void AliAnalysisTaskElecHadronCorrel::SelectPhotonicElectron(Int_t itrack, AliESDtrack *track, Bool_t &fFlagPhotonicElec)
{
  //Identify non-heavy flavour electrons using Invariant mass method
  
  fTrackCuts1->SetAcceptKinkDaughters(kFALSE);
  fTrackCuts1->SetRequireTPCRefit(kTRUE);
  fTrackCuts1->SetEtaRange(-0.9,0.9);
  fTrackCuts1->SetRequireSigmaToVertex(kTRUE);
  fTrackCuts1->SetMaxChi2PerClusterTPC(3.5);
  fTrackCuts1->SetMinNClustersTPC(80);
  
  const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
  
  Bool_t flagPhotonicElec = kFALSE;
  
  for(Int_t jTracks = itrack+1; jTracks<fESD->GetNumberOfTracks(); jTracks++){
    AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
    if (!trackAsso) {
      printf("ERROR: Could not receive track %d\n", jTracks);
      continue;
    }
    
    Double_t dEdxAsso = -999., ptAsso=-999., openingAngle = -999.;
    Double_t mass=999., width = -999;
    Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
    
    dEdxAsso = trackAsso->GetTPCsignal();
    ptAsso = trackAsso->Pt();
    Int_t chargeAsso = trackAsso->Charge();
    Int_t charge = track->Charge();
    
    if(ptAsso <0.3) continue;
    if(!fTrackCuts1->AcceptTrack(trackAsso)) continue;
    if(dEdxAsso <70 || dEdxAsso>100) continue; //11a pass1
    
    Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
    if(charge>0) fPDGe1 = -11;
    if(chargeAsso>0) fPDGe2 = -11;
    
    if(charge == chargeAsso) fFlagLS = kTRUE;
    if(charge != chargeAsso) fFlagULS = kTRUE;
    
    AliKFParticle ge1(*track, fPDGe1);
    AliKFParticle ge2(*trackAsso, fPDGe2);
    AliKFParticle recg(ge1, ge2);
    
    if(recg.GetNDF()<1) continue;
    Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
    if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
    
    AliKFVertex primV(*pVtx);
    primV += recg;
    recg.SetProductionVertex(primV);
    
    recg.SetMassConstraint(0,0.0001);
    
    openingAngle = ge1.GetAngle(ge2);
    if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
    if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
    
    if(openingAngle > fOpeningAngleCut) continue;
    
    recg.GetMass(mass,width);
    
    if(fFlagLS) fInvmassLS->Fill(mass);
    if(fFlagULS) fInvmassULS->Fill(mass);
    
    if(mass>fInvmassCut){
      ElectronHadCorrel(itrack,track,fDphiMassHigh);
    }
    if(mass>fInvmassCut){
      if(fFlagULS) ElectronHadCorrel(itrack,track,fDphiULSMassLow);
      if(fFlagLS) ElectronHadCorrel(itrack,track,fDphiLSMassLow);
    }
    
    if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
      flagPhotonicElec = kTRUE;
    }
    
  }
  fFlagPhotonicElec = flagPhotonicElec;
  
}
//_________________________________________
void AliAnalysisTaskElecHadronCorrel::ElectronHadCorrel(Int_t itrack, AliESDtrack *track, TH2F *DphiPt)
{
  //Construct Delta Phi between electrons and hadrons
  
  fTrackCuts2->SetAcceptKinkDaughters(kFALSE);
  fTrackCuts2->SetRequireTPCRefit(kTRUE);
  fTrackCuts2->SetRequireITSRefit(kTRUE);
  fTrackCuts2->SetEtaRange(-0.9,0.9);
  fTrackCuts2->SetRequireSigmaToVertex(kTRUE);
  fTrackCuts2->SetMaxChi2PerClusterTPC(3.5);
  fTrackCuts2->SetMinNClustersTPC(80);
  
  for(Int_t ktracks = 0; ktracks<fESD->GetNumberOfTracks(); ktracks++){
    AliESDtrack* trackHad = fESD->GetTrack(ktracks);
    if (!trackHad) {
      printf("ERROR: Could not receive track %d\n", ktracks);
      continue;
    }
    if(ktracks == itrack) continue; //do not select the same electron
    
    Double_t ptHad= -999, pHad=-999., dEdxHad = -999;
    Double_t ptEle = -999;
    Double_t phiEle = -999, phiHad = -999, Dphi = -999;
    Double_t pi = 3.14;
    
    dEdxHad = trackHad->GetTPCsignal();
    ptHad = trackHad->Pt();
    pHad = trackHad->P();
    
    if(ptHad <0.3) continue;
    if(!fTrackCuts2->AcceptTrack(trackHad)) continue;
    
    phiEle = track->Phi();
    phiHad = trackHad->Phi();
    Dphi = phiEle - phiHad;
    if (Dphi > pi)
      Dphi = Dphi - 2*pi;
    if (Dphi < -pi)
      Dphi = Dphi + 2*pi;
    
    ptEle = track->Pt();
    
    DphiPt->Fill(ptEle,Dphi);
  }
}
//_________________________________________
Bool_t AliAnalysisTaskElecHadronCorrel::IsLEDEvent() const
{
 // Identify LED events which had to be discarded

 AliESDCaloCells *cells = fESD->GetEMCALCells();
 Short_t nCells = cells->GetNumberOfCells();
 Int_t nCellCount[2] = {0,0};
 for(Int_t iCell=0; iCell<nCells; iCell++)
   {
     Int_t cellId = cells->GetCellNumber(iCell);
     Double_t cellE = cells->GetCellAmplitude(cellId);
     Int_t sMod = fGeom->GetSuperModuleNumber(cellId);

     if(sMod==3 || sMod==4)
       {
         if(cellE>2)
           nCellCount[sMod-3]++;
       }
   }

 if(nCellCount[0]>7 || nCellCount[1]>7)
   {
     //cout<<"Bad event!"<<endl;
     return kTRUE;
   }
 return kFALSE;
}
Commit	Line	Data
57189f04	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16
	17	////////////////////////////////////////////////////////////////////////
	18	// //
	19	// Task for Heavy Flavour Electron-Hadron DeltaPhi Correlation //
	20	// Non-Photonic Electron identified with Invariant mass //
	21	// analysis methos in function SelectPhotonicElectron //
	22	// DeltaPhi calculated in function ElectronHadCorrel //
	23	// //
	24	// Author: Deepa Thomas (Utrecht University) //
	25	// //
	26	////////////////////////////////////////////////////////////////////////
	27
	28	#include "TChain.h"
	29	#include "TTree.h"
	30	#include "TH2F.h"
	31	#include "TMath.h"
	32	#include "TCanvas.h"
	33	#include "THnSparse.h"
	34	#include "TLorentzVector.h"
	35	#include "TString.h"
	36	#include "TFile.h"
	37
	38	#include "AliAnalysisTask.h"
	39	#include "AliAnalysisManager.h"
	40
	41	#include "AliESDEvent.h"
	42	#include "AliESDHandler.h"
	43	#include "AliAODEvent.h"
	44	#include "AliAODHandler.h"
	45
	46	#include "AliAnalysisTaskElecHadronCorrel.h"
	47	#include "TGeoGlobalMagField.h"
	48	#include "AliLog.h"
	49	#include "AliAnalysisTaskSE.h"
	50	#include "TRefArray.h"
	51	#include "TVector.h"
	52	#include "AliESDInputHandler.h"
	53	#include "AliESDpid.h"
	54	#include "AliESDtrackCuts.h"
	55	#include "AliPhysicsSelection.h"
	56	#include "AliESDCaloCluster.h"
	57	#include "AliAODCaloCluster.h"
	58	#include "AliEMCALRecoUtils.h"
	59	#include "AliEMCALGeometry.h"
	60	#include "AliGeomManager.h"
	61	#include "stdio.h"
	62	#include "TGeoManager.h"
	63	#include "iostream"
	64	#include "fstream"
65
66	#include "AliEMCALTrack.h"
67	#include "AliMagF.h"
68
69	#include "AliKFParticle.h"
70	#include "AliKFVertex.h"
71
72	#include "AliPID.h"
73	#include "AliPIDResponse.h"
74	#include "AliHFEcontainer.h"
75	#include "AliHFEcuts.h"
76	#include "AliHFEpid.h"
77	#include "AliHFEpidBase.h"
78	#include "AliHFEpidQAmanager.h"
79	#include "AliHFEtools.h"
80	#include "AliCFContainer.h"
81	#include "AliCFManager.h"
82
83	ClassImp(AliAnalysisTaskElecHadronCorrel)
84	//________________________________________________________________________
85	AliAnalysisTaskElecHadronCorrel::AliAnalysisTaskElecHadronCorrel(const char *name)
86	: AliAnalysisTaskSE(name)
87	,fESD(0)
88	,fGeom(0)
89	,fOutputList(0)
90	,fTrackCuts1(0)
91	,fTrackCuts2(0)
92	,fCuts(0)
93	,fIdentifiedAsOutInz(kFALSE)
94	,fPassTheEventCut(kFALSE)
95	,fRejectKinkMother(kFALSE)
96	,fVz(0.0)
97	,fCFM(0)
98	,fPID(0)
99	,fPIDqa(0)
100	,fOpeningAngleCut(0.1)
101	,fInvmassCut(0.01)
102	,fNoEvents(0)
103	,fTrkpt(0)
104	,fTrkEovPBef(0)
105	,fTrkEovPAft(0)
106	,fdEdxBef(0)
107	,fdEdxAft(0)
108	,fInvmassLS(0)
109	,fInvmassULS(0)
110	,fOpeningAngleLS(0)
111	,fOpeningAngleULS(0)
112	,fSemiIncElecDphi(0)
113	,fPhotElecDphi(0)
114	,fInclusiveElecDphi(0)
115	,fDphiMassHigh(0)
116	,fDphiULSMassLow(0)
117	,fDphiLSMassLow(0)
118	,fPhotoElecPt(0)
119	,fSemiInclElecPt(0)
120	,fTrackPtBefTrkCuts(0)
121	,fTrackPtAftTrkCuts(0)
122	,fTPCnsigma(0)
123	{
124	//Named constructor
125
126	fPID = new AliHFEpid("hfePid");
127	fTrackCuts1 = new AliESDtrackCuts();
128	fTrackCuts2 = new AliESDtrackCuts();
129
130	// Define input and output slots here
131	// Input slot #0 works with a TChain
132	DefineInput(0, TChain::Class());
133	// Output slot #0 id reserved by the base class for AOD
134	// Output slot #1 writes into a TH1 container
135	// DefineOutput(1, TH1I::Class());
136	DefineOutput(1, TList::Class());
137	// DefineOutput(3, TTree::Class());
138	}
139
140	//________________________________________________________________________
141	AliAnalysisTaskElecHadronCorrel::AliAnalysisTaskElecHadronCorrel()
142	: AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisElecHadCorrel")
143	,fESD(0)
144	,fGeom(0)
145	,fOutputList(0)
146	,fTrackCuts1(0)
147	,fTrackCuts2(0)
148	,fCuts(0)
149	,fIdentifiedAsOutInz(kFALSE)
150	,fPassTheEventCut(kFALSE)
151	,fRejectKinkMother(kFALSE)
152	,fVz(0.0)
153	,fCFM(0)
154	,fPID(0)
155	,fPIDqa(0)
156	,fOpeningAngleCut(0.1)
157	,fInvmassCut(0.01)
158	,fNoEvents(0)
159	,fTrkpt(0)
160	,fTrkEovPBef(0)
161	,fTrkEovPAft(0)
162	,fdEdxBef(0)
163	,fdEdxAft(0)
164	,fInvmassLS(0)
165	,fInvmassULS(0)
166	,fOpeningAngleLS(0)
167	,fOpeningAngleULS(0)
168	,fSemiIncElecDphi(0)
169	,fPhotElecDphi(0)
170	,fInclusiveElecDphi(0)
171	,fDphiMassHigh(0)
172	,fDphiULSMassLow(0)
173	,fDphiLSMassLow(0)
174	,fPhotoElecPt(0)
175	,fSemiInclElecPt(0)
176	,fTrackPtBefTrkCuts(0)
177	,fTrackPtAftTrkCuts(0)
178	,fTPCnsigma(0)
179	{
180	//Default constructor
181	fPID = new AliHFEpid("hfePid");
182
183	fTrackCuts1 = new AliESDtrackCuts();
184	fTrackCuts2 = new AliESDtrackCuts();
185
186	// Constructor
187	// Define input and output slots here
188	// Input slot #0 works with a TChain
189	DefineInput(0, TChain::Class());
190	// Output slot #0 id reserved by the base class for AOD
191	// Output slot #1 writes into a TH1 container
192	// DefineOutput(1, TH1I::Class());
193	DefineOutput(1, TList::Class());
194	//DefineOutput(3, TTree::Class());
195	}
196	//_________________________________________
197
198	AliAnalysisTaskElecHadronCorrel::~AliAnalysisTaskElecHadronCorrel()
199	{
200	//Destructor
201
202	delete fOutputList;
203	delete fGeom;
204	delete fPID;
205	delete fCFM;
206	delete fPIDqa;
207	delete fTrackCuts1;
208	delete fTrackCuts2;
209	}
210	//_________________________________________
211
212	void AliAnalysisTaskElecHadronCorrel::UserExec(Option_t*)
213	{
214	//Main loop
215	//Called for each event
216
217	// create pointer to event
218	fESD = dynamic_cast<AliESDEvent*>(InputEvent());
219	if (!fESD) {
220	printf("ERROR: fESD not available\n");
221	return;
222	}
223
224	if(!fCuts){
225	AliError("HFE cuts not available");
226	return;
227	}
228
229	if(!fPID->IsInitialized()){
230	// Initialize PID with the given run number
231	AliWarning("PID not initialised, get from Run no");
232	fPID->InitializePID(fESD->GetRunNumber());
233	}
234
235	//-------trigger selection
236	UInt_t res = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
237	if (res==0)
238	return;
239
240	if( (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kFastOnly) )
241	return;
242
243	if(!(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kEMC1) ) return;
244	//if(!(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB) ) return;
245
246	Int_t fNOtrks = fESD->GetNumberOfTracks();
247	const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
248
249	Double_t pVtxZ = -999;
250	pVtxZ = pVtx->GetZ();
251
252	// Event cut
253	// if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
254
255	//Check if its an LED event
256	if(IsLEDEvent()) return;
257
258	if(TMath::Abs(pVtxZ)>10) return;
259	fNoEvents->Fill(0);
260
261	if(fNOtrks<2) return;
262
263	AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
264	if(!pidResponse){
265	AliDebug(1, "Using default PID Response");
266	pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
267	}
268
269	fPID->SetPIDResponse(pidResponse);
270
271	fCFM->SetRecEventInfo(fESD);
272
273	// Track loop
274	for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
275	AliESDtrack* track = fESD->GetTrack(iTracks);
276	if (!track) {
277	printf("ERROR: Could not receive track %d\n", iTracks);
278	continue;
279	}
280
281	fTrackPtBefTrkCuts->Fill(track->Pt());
282	// RecKine: ITSTPC cuts
283	if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
284
285	//RecKink
286	if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
287	if(track->GetKinkIndex(0) != 0) continue;
288	}
289
290	// RecPrim
291	if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
292
293	// HFEcuts: ITS layers cuts
294	if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
295
296	// HFE cuts: TPC PID cleanup
297	if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
298
299	fTrackPtAftTrkCuts->Fill(track->Pt());
300
301	Double_t fClsE = -999, p = -999, fEovP=-999, pt = -999, dEdx=-999, fTPCnSigma=0;
302	// Track extrapolation
303	Int_t fClsId = track->GetEMCALcluster();
304	if(fClsId == -1) continue;
305	AliESDCaloCluster *cluster = fESD->GetCaloCluster(fClsId);
306
307	pt = track->Pt();
308	if(pt<1) continue;
309	fTrkpt->Fill(pt);
310	fClsE = cluster->E();
311	p = track->P();
312	dEdx = track->GetTPCsignal();
313	fEovP = fClsE/p;
314	fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
315	fdEdxBef->Fill(p,dEdx);
316	fTPCnsigma->Fill(p,fTPCnSigma);
317
318	if(fTPCnSigma >= 1.5 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,fEovP);
319	Int_t pidpassed = 1;
320	//--- track accepted, do PID
321	AliHFEpidObject hfetrack;
322	hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
323	hfetrack.SetRecTrack(track);
324	hfetrack.SetPP();
325	if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
326
327	if(pidpassed==0) continue;
328	fTrkEovPAft->Fill(pt,fEovP);
329	fdEdxAft->Fill(p,dEdx);
330
331	Bool_t fFlagPhotonicElec = kFALSE;
332	// select photonic electron
333	SelectPhotonicElectron(iTracks,track,fFlagPhotonicElec);
334
335	//Inclusive electron-hadron correlation
336	ElectronHadCorrel(iTracks, track, fInclusiveElecDphi);
337
338	// photonic electron
339	if(fFlagPhotonicElec){
340	//Electron hadron correlation
341	ElectronHadCorrel(iTracks, track, fPhotElecDphi);
342	fPhotoElecPt->Fill(pt);
343	}
344
345	// Semi inclusive electron
346	if(!fFlagPhotonicElec){
347	//Electron hadron correlation
348	ElectronHadCorrel(iTracks, track, fSemiIncElecDphi);
349	fSemiInclElecPt->Fill(pt);
350	}
351	}
352
353
354	PostData(1, fOutputList);
355	}
356	//_________________________________________
357	void AliAnalysisTaskElecHadronCorrel::UserCreateOutputObjects()
358	{
359	//Create histograms
360	TGeoManager::Import("geometry.root");
361	fGeom = AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
362
363	//--------Initialize PID
364	fPID->SetHasMCData(kFALSE);
365	if(!fPID->GetNumberOfPIDdetectors())
366	{
367	fPID->AddDetector("TPC", 0);
368	fPID->AddDetector("EMCAL", 1);
369	}
370
371	fPID->SortDetectors();
372	fPIDqa = new AliHFEpidQAmanager();
373	fPIDqa->Initialize(fPID);
374
375	//--------Initialize correction Framework and Cuts
376	fCFM = new AliCFManager;
377	const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
378	fCFM->SetNStepParticle(kNcutSteps);
379	for(Int_t istep = 0; istep < kNcutSteps; istep++)
380	fCFM->SetParticleCutsList(istep, NULL);
381
382	if(!fCuts){
383	AliWarning("Cuts not available. Default cuts will be used");
384	fCuts = new AliHFEcuts;
385	fCuts->CreateStandardCuts();
386	}
387	fCuts->Initialize(fCFM);
388
389	//---------Output Tlist
390	fOutputList = new TList();
391	fOutputList->SetOwner();
392	fOutputList->Add(fPIDqa->MakeList("PIDQA"));
393
394	fNoEvents = new TH1F("fNoEvents","",1,0,1) ;
395	fOutputList->Add(fNoEvents);
396
397	fTrkpt = new TH1F("fTrkpt","track pt",1000,0,50);
398	fOutputList->Add(fTrkpt);
399
400	fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",1000,0,50);
401	fOutputList->Add(fTrackPtBefTrkCuts);
402
403	fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",1000,0,50);
404	fOutputList->Add(fTrackPtAftTrkCuts);
405
406	fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",1000,0,50,200,-10,10);
407	fOutputList->Add(fTPCnsigma);
408
409	fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",1000,0,50,100,0,2);
410	fOutputList->Add(fTrkEovPBef);
411
412	fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",1000,0,50,100,0,2);
413	fOutputList->Add(fTrkEovPAft);
414
415	fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",1000,0,50,150,0,150);
416	fOutputList->Add(fdEdxBef);
417
418	fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",1000,0,50,150,0,150);
419	fOutputList->Add(fdEdxAft);
420
421	fInvmassLS = new TH1F("fInvmassLS", "Inv mass of LS (e,e); mass(GeV/c^2); counts;", 1000,0,0.5);
422	fOutputList->Add(fInvmassLS);
423
424	fInvmassULS = new TH1F("fInvmassULS", "Inv mass of ULS (e,e); mass(GeV/c^2); counts;", 1000,0,0.5);
425	fOutputList->Add(fInvmassULS);
426
427	fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
428	fOutputList->Add(fOpeningAngleLS);
429
430	fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
431	fOutputList->Add(fOpeningAngleULS);
432
433	fSemiIncElecDphi = new TH2F("fSemiIncElecDphi", "Semi Inclusive elec-had Dphi correlation",200,0,20,100,-3.14,3.14);
434	fOutputList->Add(fSemiIncElecDphi);
435
436	fPhotElecDphi = new TH2F("fPhotElecDphi", "Photon elec-had Dphi correlation",200,0,20,100,-3.14,3.14);
437	fOutputList->Add(fPhotElecDphi);
438
439	fInclusiveElecDphi = new TH2F("fInclusiveElecDphi", "Inclusive elec-had Dphi correlation",200,0,20,100,-3.14,3.14);
440	fOutputList->Add(fInclusiveElecDphi);
441
442	fDphiMassHigh = new TH2F("fDphiMassHigh", "e-h Dphi LS+ULS, mass>0.01",200,0,20,100,-3.14,3.14);
443	fOutputList->Add(fDphiMassHigh);
444
445	fDphiULSMassLow = new TH2F("fDphiULSMassLow", "e-h Dphi ULS, mass<0.01",200,0,20,100,-3.14,3.14);
446	fOutputList->Add(fDphiULSMassLow);
447
448	fDphiLSMassLow = new TH2F("fDphiLSMassLow", "e-h Dphi LS, mass<0.01",200,0,20,100,-3.14,3.14);
449	fOutputList->Add(fDphiLSMassLow);
450
451	fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",1000,0,100);
452	fOutputList->Add(fPhotoElecPt);
453
454	fSemiInclElecPt = new TH1F("fSemiInclElecPt", "Semi-inclusive electron pt",1000,0,100);
455	fOutputList->Add(fSemiInclElecPt);
456	PostData(1,fOutputList);
457	}
458
459	//________________________________________________________________________
460	void AliAnalysisTaskElecHadronCorrel::Terminate(Option_t *)
461	{
462	// Info("Terminate");
463	AliAnalysisTaskSE::Terminate();
464	}
465
466	//________________________________________________________________________
467	Bool_t AliAnalysisTaskElecHadronCorrel::ProcessCutStep(Int_t cutStep, AliVParticle *track)
468	{
469	// Check single track cuts for a given cut step
470	const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
471	if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
472	return kTRUE;
473	}
474	//_________________________________________
475	void AliAnalysisTaskElecHadronCorrel::SelectPhotonicElectron(Int_t itrack, AliESDtrack *track, Bool_t &fFlagPhotonicElec)
476	{
477	//Identify non-heavy flavour electrons using Invariant mass method
478
479	fTrackCuts1->SetAcceptKinkDaughters(kFALSE);
480	fTrackCuts1->SetRequireTPCRefit(kTRUE);
481	fTrackCuts1->SetEtaRange(-0.9,0.9);
482	fTrackCuts1->SetRequireSigmaToVertex(kTRUE);
483	fTrackCuts1->SetMaxChi2PerClusterTPC(3.5);
484	fTrackCuts1->SetMinNClustersTPC(80);
485
486	const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
487
488	Bool_t flagPhotonicElec = kFALSE;
489
490	for(Int_t jTracks = itrack+1; jTracks<fESD->GetNumberOfTracks(); jTracks++){
491	AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
492	if (!trackAsso) {
493	printf("ERROR: Could not receive track %d\n", jTracks);
494	continue;
495	}
496
497	Double_t dEdxAsso = -999., ptAsso=-999., openingAngle = -999.;
498	Double_t mass=999., width = -999;
499	Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
500
501	dEdxAsso = trackAsso->GetTPCsignal();
502	ptAsso = trackAsso->Pt();
503	Int_t chargeAsso = trackAsso->Charge();
504	Int_t charge = track->Charge();
505
506	if(ptAsso <0.3) continue;
507	if(!fTrackCuts1->AcceptTrack(trackAsso)) continue;
508	if(dEdxAsso <70 \|\| dEdxAsso>100) continue; //11a pass1
509
510	Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
511	if(charge>0) fPDGe1 = -11;
512	if(chargeAsso>0) fPDGe2 = -11;
513
514	if(charge == chargeAsso) fFlagLS = kTRUE;
515	if(charge != chargeAsso) fFlagULS = kTRUE;
516
517	AliKFParticle ge1(*track, fPDGe1);
518	AliKFParticle ge2(*trackAsso, fPDGe2);
519	AliKFParticle recg(ge1, ge2);
520
521	if(recg.GetNDF()<1) continue;
522	Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
523	if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
524
525	AliKFVertex primV(*pVtx);
526	primV += recg;
527	recg.SetProductionVertex(primV);
528
529	recg.SetMassConstraint(0,0.0001);
530
531	openingAngle = ge1.GetAngle(ge2);
532	if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
533	if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
534
535	if(openingAngle > fOpeningAngleCut) continue;
536
537	recg.GetMass(mass,width);
538
539	if(fFlagLS) fInvmassLS->Fill(mass);
540	if(fFlagULS) fInvmassULS->Fill(mass);
541
542	if(mass>fInvmassCut){
543	ElectronHadCorrel(itrack,track,fDphiMassHigh);
544	}
545	if(mass>fInvmassCut){
546	if(fFlagULS) ElectronHadCorrel(itrack,track,fDphiULSMassLow);
547	if(fFlagLS) ElectronHadCorrel(itrack,track,fDphiLSMassLow);
548	}
549
550	if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
551	flagPhotonicElec = kTRUE;
552	}
553
554	}
555	fFlagPhotonicElec = flagPhotonicElec;
556
557	}
558	//_________________________________________
559	void AliAnalysisTaskElecHadronCorrel::ElectronHadCorrel(Int_t itrack, AliESDtrack track, TH2F DphiPt)
560	{
561	//Construct Delta Phi between electrons and hadrons
562
563	fTrackCuts2->SetAcceptKinkDaughters(kFALSE);
564	fTrackCuts2->SetRequireTPCRefit(kTRUE);
565	fTrackCuts2->SetRequireITSRefit(kTRUE);
566	fTrackCuts2->SetEtaRange(-0.9,0.9);
567	fTrackCuts2->SetRequireSigmaToVertex(kTRUE);
568	fTrackCuts2->SetMaxChi2PerClusterTPC(3.5);
569	fTrackCuts2->SetMinNClustersTPC(80);
570
571	for(Int_t ktracks = 0; ktracks<fESD->GetNumberOfTracks(); ktracks++){
572	AliESDtrack* trackHad = fESD->GetTrack(ktracks);
573	if (!trackHad) {
574	printf("ERROR: Could not receive track %d\n", ktracks);
575	continue;
576	}
577	if(ktracks == itrack) continue; //do not select the same electron
578
579	Double_t ptHad= -999, pHad=-999., dEdxHad = -999;
580	Double_t ptEle = -999;
581	Double_t phiEle = -999, phiHad = -999, Dphi = -999;
582	Double_t pi = 3.14;
583
584	dEdxHad = trackHad->GetTPCsignal();
585	ptHad = trackHad->Pt();
586	pHad = trackHad->P();
587
588	if(ptHad <0.3) continue;
589	if(!fTrackCuts2->AcceptTrack(trackHad)) continue;
590
591	phiEle = track->Phi();
592	phiHad = trackHad->Phi();
593	Dphi = phiEle - phiHad;
594	if (Dphi > pi)
595	Dphi = Dphi - 2*pi;
596	if (Dphi < -pi)
597	Dphi = Dphi + 2*pi;
598
599	ptEle = track->Pt();
600
601	DphiPt->Fill(ptEle,Dphi);
602	}
603	}
604	//_________________________________________
605	Bool_t AliAnalysisTaskElecHadronCorrel::IsLEDEvent() const
606	{
607	// Identify LED events which had to be discarded
608
609	AliESDCaloCells *cells = fESD->GetEMCALCells();
610	Short_t nCells = cells->GetNumberOfCells();
611	Int_t nCellCount[2] = {0,0};
612	for(Int_t iCell=0; iCell<nCells; iCell++)
613	{
614	Int_t cellId = cells->GetCellNumber(iCell);
615	Double_t cellE = cells->GetCellAmplitude(cellId);
616	Int_t sMod = fGeom->GetSuperModuleNumber(cellId);
617
618	if(sMod==3 \|\| sMod==4)
619	{
620	if(cellE>2)
621	nCellCount[sMod-3]++;
622	}
623	}
624
625	if(nCellCount[0]>7 \|\| nCellCount[1]>7)
626	{
627	//cout<<"Bad event!"<<endl;
628	return kTRUE;
629	}
630	return kFALSE;
631	}
632