#include "AliESDpid.h"
#include "AliCentrality.h"
#include "AliTracker.h"
+#include "AliAODInputHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODTrack.h"
#include "AliAnalysisNetParticleDistribution.h"
fESDHandler(NULL),
fPIDResponse(NULL),
fESD(NULL),
+ fAODHandler(NULL),
+ fAOD(NULL),
fIsMC(kFALSE),
fMCEvent(NULL),
fStack(NULL),
fESDTrackCuts(NULL),
fEtaMax(0.9),
fPtRange(NULL),
+ fAODtrackCutBit(1024),
fNp(NULL),
fNCorrNp(2),
fCorrNp(NULL),
*/
//________________________________________________________________________
-Int_t AliAnalysisNetParticleDistribution::Initialize(AliAnalysisNetParticleHelper* helper, AliESDtrackCuts* cuts, Bool_t isMC, Float_t *ptRange, Float_t etaMax) {
+Int_t AliAnalysisNetParticleDistribution::Initialize(AliAnalysisNetParticleHelper* helper, AliESDtrackCuts* cuts, Bool_t isMC, Float_t *ptRange, Float_t etaMax, Int_t trackCutBit) {
// -- Initialize
fHelper = helper;
fIsMC = isMC;
fPtRange = ptRange;
fEtaMax = etaMax;
+ fAODtrackCutBit = trackCutBit;
// ------------------------------------------------------------------
// -- N particles / N anti-particles
}
//________________________________________________________________________
-Int_t AliAnalysisNetParticleDistribution::SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent) {
+Int_t AliAnalysisNetParticleDistribution::SetupEvent(AliESDInputHandler *esdHandler, AliAODInputHandler *aodHandler, AliMCEvent *mcEvent) {
// -- Setup Event
ResetEvent();
// -- Get ESD objects
- fESDHandler = esdHandler;
- fPIDResponse = esdHandler->GetPIDResponse();
- fESD = fESDHandler->GetEvent();
+ if(esdHandler){
+ fESDHandler = esdHandler;
+ fPIDResponse = esdHandler->GetPIDResponse();
+ fESD = fESDHandler->GetEvent();
+ }
+
+ // -- Get AOD objects
+ else if(aodHandler){
+ fAODHandler = aodHandler;
+ fPIDResponse = aodHandler->GetPIDResponse();
+ fAOD = fAODHandler->GetEvent();
+ }
// -- Get MC objects
fMCEvent = mcEvent;
// -- Reset ESD Event
fESD = NULL;
+ // -- Reset AOD Event
+ fAOD = NULL;
+
// -- Reset MC Event
if (fIsMC)
fMCEvent = NULL;
//________________________________________________________________________
Int_t AliAnalysisNetParticleDistribution::Process() {
// -- Process NetParticle Distributions
-
+
// -- Fill ESD tracks
- ProcessESDTracks();
+ if (fESD) ProcessESDTracks();
+
+ // -- Fill AOD tracks
+ else if (fAOD) ProcessAODTracks();
- // -- Fill MC truth particles
+ // -- Fill MC truth particles (missing for AOD XXX)
if (fIsMC) {
ProcessStackParticles();
ProcessStackControlParticles();
return 0;
}
+//________________________________________________________________________
+Int_t AliAnalysisNetParticleDistribution::ProcessAODTracks() {
+ // -- Process ESD tracks and fill histograms
+
+ for (Int_t idxTrack = 0; idxTrack < fAOD->GetNumberOfTracks(); ++idxTrack) {
+ AliAODTrack *track = (AliAODTrack*)fAOD->GetTrack(idxTrack);
+
+ // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
+ // -- Check track cuts
+ // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
+
+ // -- Check if charged track is accepted for basic parameters
+ if (!fHelper->IsTrackAcceptedBasicCharged(track))
+ continue;
+
+ // -- Check if accepted
+ if(!track->TestFilterBit(fAODtrackCutBit))
+ continue;
+
+ // -- Check if accepted in rapidity window
+ Double_t yP;
+ if (!fHelper->IsTrackAcceptedRapidity(track, yP))
+ continue;
+
+ // -- Check if accepted bt PID from TPC or TPC+TOF
+ Double_t pid[2];
+ if (!fHelper->IsTrackAcceptedPID(track, pid))
+ continue;
+
+ // -- Check if accepted with thighter DCA cuts XXX
+ // if (fHelper->IsTrackAcceptedDCA(track))
+ // continue;
+
+ // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
+ // -- Fill Probe Particle
+ // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
+
+ Double_t aTrack[6] = {
+ Double_t(fHelper->GetCentralityBin()), // 0 centrality
+ track->Pt(), // 1 pt
+ track->Charge(), // 2 sign
+ track->Eta(), // 3 eta
+ track->Phi(), // 4 phi
+ yP // 5 rapidity
+ };
+
+ fHnTrackUnCorr->Fill(aTrack);
+
+ // -- Count particle / anti-particle
+ // ------------------------------------------------------------------
+ // idxPart = 0 -> anti particle
+ // idxPart = 1 -> particle
+
+ Int_t idxPart = (track->Charge() < 0) ? 0 : 1;
+ fNp[idxPart] += 1.;
+
+ for (Int_t ii = 0; ii < fNCorrNp; ++ii)
+ fCorrNp[ii][idxPart] += fHelper->GetTrackbyTrackCorrectionFactor(aTrack, ii);
+
+ } // for (Int_t idxTrack = 0; idxTrack < fESD->GetNumberOfTracks(); ++idxTrack) {
+
+ // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
+ // -- Fill Particle Fluctuation Histograms
+ // -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
+
+ // -- Uncorrected
+ FillHistSet("fHDist", fNp);
+
+ // -- Corrected
+ for (Int_t ii = 0 ; ii < fNCorrNp; ++ii)
+ FillHistSet(Form("fHDistCorr%d", ii), fCorrNp[ii]);
+
+ return 0;
+}
+
//________________________________________________________________________
Int_t AliAnalysisNetParticleDistribution::ProcessStackParticles() {
// -- Process primary particles from the stack and fill histograms
class AliPIDResponse;
class AliESDInputHandler;
class AliESDtrackCuts;
+class AliAODInputHandler;
class AliAnalysisNetParticleDistribution : public TNamed {
*/
/** Initialize */
- Int_t Initialize(AliAnalysisNetParticleHelper* helper, AliESDtrackCuts* cuts, Bool_t isMC, Float_t *ptRange, Float_t etaMax);
+ Int_t Initialize(AliAnalysisNetParticleHelper* helper, AliESDtrackCuts* cuts, Bool_t isMC, Float_t *ptRange, Float_t etaMax, Int_t trackCutBit);
/** Add histograms to outlist */
void CreateHistograms(TList *outList);
/** Setup Event */
- Int_t SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent);
+ Int_t SetupEvent(AliESDInputHandler *esdHandler, AliAODInputHandler *aodHandler, AliMCEvent *mcEvent);
/** Resre Event */
void ResetEvent();
/** Process ESD tracks and fill histograms */
Int_t ProcessESDTracks();
+ /** Process AOD tracks and fill histograms */
+ Int_t ProcessAODTracks();
+
/** Process primary particles from the stack and fill histograms */
Int_t ProcessStackParticles();
AliPIDResponse *fPIDResponse; //! Ptr to PID response Object
AliESDEvent *fESD; //! Ptr to ESD event
+ AliAODInputHandler *fAODHandler; //! Ptr to AOD handler
+ AliAODEvent *fAOD; //! Ptr to AOD event
+
Bool_t fIsMC; // Is MC event
AliMCEvent *fMCEvent; //! Ptr to MC event
// -----------------------------------------------------------------------
Float_t fEtaMax; // Max, absolut eta
Float_t *fPtRange; // Array of pt [min,max]
+
+ Int_t fAODtrackCutBit; // Track filter bit for AOD tracks
// -----------------------------------------------------------------------
Float_t *fNp; // Array of particle/anti-particle counts
#include "TMath.h"
#include "TAxis.h"
+#include "TDatabasePDG.h"
#include "AliESDEvent.h"
#include "AliESDInputHandler.h"
#include "AliStack.h"
#include "AliMCEvent.h"
-
#include "AliESDtrackCuts.h"
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
+#include "AliAODMCParticle.h"
#include "AliAnalysisNetParticleEffCont.h"
fESD(NULL),
fESDTrackCuts(NULL),
+ fAOD(NULL),
+ fArrayMC(NULL),
+
fCentralityBin(-1.),
fNTracks(0),
+ fAODtrackCutBit(1024),
+
fStack(NULL),
fMCEvent(NULL),
*/
//________________________________________________________________________
-void AliAnalysisNetParticleEffCont::Initialize(AliESDtrackCuts *cuts , AliAnalysisNetParticleHelper* helper) {
+void AliAnalysisNetParticleEffCont::Initialize(AliESDtrackCuts *cuts , AliAnalysisNetParticleHelper* helper, Int_t trackCutBit) {
// -- ESD track cuts
// -------------------
for (Int_t ii = 0; ii < 2 ; ++ii)
fLabelsRec[ii] = NULL;
+ // -- AOD track filter bit
+ // -------------------------
+ fAODtrackCutBit = trackCutBit;
+
+
// -- Create THnSparse Histograms
// --------------------------------
CreateHistograms();
return 0;
}
+//________________________________________________________________________
+Int_t AliAnalysisNetParticleEffCont::SetupEvent(AliAODInputHandler *aodHandler) {
+ // -- Setup event
+
+ // -- Reset of event
+ // -------------------
+ ResetEvent();
+
+ // -- Setup of event
+ // -------------------
+ fAOD = aodHandler->GetEvent();
+ fNTracks = fAOD->GetNumberOfTracks();
+
+ fArrayMC = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
+ if (!fArrayMC)
+ AliFatal("No array of MC particles found !!!");
+
+ fCentralityBin = fHelper->GetCentralityBin();
+
+ // -- Create label arrays
+ // ------------------------
+ fLabelsRec[0] = new Int_t[fNTracks];
+ if(!fLabelsRec[0]) {
+ AliError("Cannot create fLabelsRec[0]");
+ return -1;
+ }
+
+ fLabelsRec[1] = new Int_t[fNTracks];
+ if(!fLabelsRec[1]) {
+ AliError("Cannot create fLabelsRec[1] for TPC");
+ return -1;
+ }
+
+ for(Int_t ii = 0; ii < fNTracks; ++ii) {
+ fLabelsRec[0][ii] = 0;
+ fLabelsRec[1][ii] = 0;
+ }
+
+ return 0;
+}
+
//________________________________________________________________________
void AliAnalysisNetParticleEffCont::ResetEvent() {
// -- Reset event
// -- Setup (clean, create and fill) MC labels
// ---------------------------------------------
- FillMCLabels();
-
+ if(fESD)
+ FillMCLabels();
+ else if(fAOD)
+ FillMCLabelsAOD();
+
// -- Fill MC histograms for efficiency studies
// -----------------------------------------------
- FillMCEffHist();
+ if(fESD)
+ FillMCEffHist();
+ else if(fAOD)
+ FillMCEffHistAOD();
return;
}
return;
}
+//________________________________________________________________________
+void AliAnalysisNetParticleEffCont::FillMCLabelsAOD() {
+ // Fill MC labels for AOD analysis
+ // Loop over ESD tracks and fill arrays with MC lables
+ // fLabelsRec[0] : all Tracks
+ // fLabelsRec[1] : all Tracks accepted by PID of TPC
+ // Check every accepted track if correctly identified
+ // otherwise check for contamination
+
+ for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
+ AliAODTrack *track = fAOD->GetTrack(idxTrack);
+
+ // -- Check if track is accepted for basic parameters
+ if (!fHelper->IsTrackAcceptedBasicCharged(track))
+ continue;
+
+ // -- Check if accepted
+ if (!track->TestFilterBit(fAODtrackCutBit))
+ continue;
+
+ // -- Check if accepted in rapidity window
+ Double_t yP;
+ if (!fHelper->IsTrackAcceptedRapidity(track, yP))
+ continue;
+
+ // -- Check if accepted with thighter DCA cuts (not yet done for AODs XXX)
+ //if (!fHelper->IsTrackAcceptedDCA(track))
+ // continue;
+
+ Int_t label = TMath::Abs(track->GetLabel());
+
+ // -- Fill Label of all reconstructed
+ fLabelsRec[0][idxTrack] = label;
+
+ // -- Check if accepted by PID from TPC or TPC+TOF
+ Double_t pid[2];
+ if (!fHelper->IsTrackAcceptedPID(track, pid))
+ continue;
+
+ // -- Fill Label of all reconstructed && recPid_TPC+TOF
+ fLabelsRec[1][idxTrack] = label;
+
+ // -- Check for contamination and fill contamination THnSparse
+ CheckContTrackAOD(label, track->Charge(), idxTrack);
+
+ } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
+
+ return;
+}
+
//________________________________________________________________________
void AliAnalysisNetParticleEffCont::CheckContTrack(Int_t label, Float_t sign, Int_t idxTrack) {
- // Check if particle is contamination or correctly identified
+ // Check if particle is contamination or correctly identified for ESDs
// Fill contamination THnSparse
-
+
TParticle* particle = fStack->Particle(label);
if (!particle)
return;
fHnCont->Fill(hnCont);
}
+//________________________________________________________________________
+void AliAnalysisNetParticleEffCont::CheckContTrackAOD(Int_t label, Float_t sign, Int_t idxTrack) {
+ // Check if particle is contamination or correctly identified for AODs
+ // Fill contamination THnSparse
+
+ AliAODMCParticle* particle = (AliAODMCParticle*)fArrayMC->At(label);
+
+ if (!particle)
+ return;
+
+ Bool_t isSecondaryFromWeakDecay = kFALSE;
+ Bool_t isSecondaryFromMaterial = kFALSE;
+
+ // -- Check if correctly identified
+ if (particle->GetPdgCode() == (sign*fPdgCode)) {
+
+ // Check if is physical primary -> all ok
+ if(particle->IsPhysicalPrimary())
+ return;
+
+ // -- Check if secondaries from material or weak decay
+ isSecondaryFromWeakDecay = particle->IsSecondaryFromWeakDecay();
+ isSecondaryFromMaterial = particle->IsSecondaryFromMaterial();
+
+ }
+
+ // -- Get MC pdg
+ Float_t contSign = 0.;
+ if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() == 0.) contSign = 0.;
+ else if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() < 0.) contSign = -1.;
+ else if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() > 0.) contSign = 1.;
+
+ // -- Get contaminating particle
+ Float_t contPart = 0;
+ if (isSecondaryFromWeakDecay) contPart = 7; // probeParticle from WeakDecay
+ else if (isSecondaryFromMaterial) contPart = 8; // probeParticle from Material
+ else {
+ if (TMath::Abs(particle->GetPdgCode()) == 211) contPart = 1; // pion
+ else if (TMath::Abs(particle->GetPdgCode()) == 321) contPart = 2; // kaon
+ else if (TMath::Abs(particle->GetPdgCode()) == 2212) contPart = 3; // proton
+ else if (TMath::Abs(particle->GetPdgCode()) == 11) contPart = 4; // electron
+ else if (TMath::Abs(particle->GetPdgCode()) == 13) contPart = 5; // muon
+ else contPart = 6; // other
+ }
+
+ // -- Get Reconstructed values
+ Float_t etaRec = 0.;
+ Float_t phiRec = 0.;
+ Float_t ptRec = 0.;
+
+ AliAODTrack *track = fAOD->GetTrack(idxTrack);
+
+ if (track) {
+ // if no track present (which should not happen)
+ // -> pt = 0. , which is not in the looked at range
+
+ // -- Get Reconstructed eta/phi/pt
+ etaRec = track->Eta();
+ phiRec = track->Phi();
+ ptRec = track->Pt();
+ }
+
+ // -- Fill THnSparse
+ Double_t hnCont[14] = {fCentralityBin, particle->Eta(), particle->Y(), particle->Phi(), particle->Pt(), sign,
+ contPart, contSign, etaRec, phiRec, ptRec,
+ particle->Eta()-etaRec, particle->Phi()-phiRec, particle->Pt()-ptRec};
+ fHnCont->Fill(hnCont);
+}
+
//________________________________________________________________________
void AliAnalysisNetParticleEffCont::FillMCEffHist() {
- // Fill efficiency THnSparse
+ // Fill efficiency THnSparse for ESDs
Int_t nPart = fStack->GetNprimary();
return;
}
+
+//________________________________________________________________________
+void AliAnalysisNetParticleEffCont::FillMCEffHistAOD() {
+ // Fill efficiency THnSparse for AODs
+
+ Int_t nPart = fArrayMC->GetEntriesFast();
+
+ for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
+
+ AliAODMCParticle* particle = (AliAODMCParticle*)fArrayMC->At(idxMC);
+
+ // -- Check basic MC properties -> charged physical primary
+ if (!fHelper->IsParticleAcceptedBasicCharged(particle))
+ continue;
+
+ // -- Check rapidity window
+ Double_t yMC;
+ if (!fHelper->IsParticleAcceptedRapidity((TParticle*)particle, yMC))
+ continue;
+
+ // -- Check if probeParticle / anti-probeParticle
+ if (TMath::Abs(particle->GetPdgCode()) != fPdgCode)
+ continue;
+
+ // -- Get sign of particle
+ Float_t sign = (particle->GetPdgCode() < 0) ? -1. : 1.;
+
+ // -- Get if particle is findable
+ Float_t findable = 1.;// Float_t(fHelper->IsParticleFindable(idxMC)); XXX
+
+ // -- Get recStatus and pidStatus
+ Float_t recStatus = 0.;
+ Float_t recPid = 0.;
+
+ // -- Get Reconstructed values
+ Float_t etaRec = 0.;
+ Float_t phiRec = 0.;
+ Float_t ptRec = 0.;
+
+ // -- Loop over all labels
+ for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
+ if (idxMC == fLabelsRec[0][idxRec]) {
+ recStatus = 1.;
+
+ if (idxMC == fLabelsRec[1][idxRec])
+ recPid = 1.;
+
+ AliVTrack *track = NULL;
+ if(fESD)
+ track = fESD->GetTrack(idxRec);
+ else if(fAOD)
+ track = fAOD->GetTrack(idxRec);
+
+ if (track) {
+ // if no track present (which should not happen)
+ // -> pt = 0. , which is not in the looked at range
+
+ // -- Get Reconstructed eta/phi/pt
+ etaRec = track->Eta();
+ phiRec = track->Phi();
+ ptRec = track->Pt();
+ }
+ break;
+ }
+ } // for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
+
+ // -- Fill THnSparse
+ Double_t hnEff[15] = {fCentralityBin, particle->Eta(), particle->Y(), particle->Phi(), particle->Pt(), sign,
+ findable, recStatus, recPid, etaRec, phiRec, ptRec,
+ particle->Eta()-etaRec, particle->Phi()-phiRec, particle->Pt()-ptRec};
+ fHnEff->Fill(hnEff);
+
+ } // for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
+
+ return;
+}
class AliESDEvent;
class AliESDInputHandler;
+class AliAODInputHandler;
class AliMCEvent;
class AliAnalysisNetParticleEffCont : public TNamed {
*/
/** Initialize */
- void Initialize(AliESDtrackCuts *cuts, AliAnalysisNetParticleHelper* helper);
+ void Initialize(AliESDtrackCuts *cuts, AliAnalysisNetParticleHelper* helper, Int_t trackCutBit);
/** Setup Event */
- Int_t SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent);
+ Int_t SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent);
+ Int_t SetupEvent(AliAODInputHandler *esdHandler); // MC particles are stored in AOD
/** Reset Event */
void ResetEvent();
/** Fill MC labels */
void FillMCLabels();
+ void FillMCLabelsAOD();
/** Fill efficiency THnSparse */
void FillMCEffHist();
+ void FillMCEffHistAOD();
/** Check if particle is contamination */
void CheckContTrack(Int_t label, Float_t sign, Int_t idxTrack);
+ void CheckContTrackAOD(Int_t label, Float_t sign, Int_t idxTrack);
/*
* ---------------------------------------------------------------------------------
AliESDEvent *fESD; //! ESD object
AliESDtrackCuts *fESDTrackCuts; //! ESD cuts
+ // --- AOD only ----------------------------------------------------------
+
+ AliAODEvent *fAOD; //! AOD object
+ TClonesArray *fArrayMC; //! array of MC particles
+
// -----------------------------------------------------------------------
Float_t fCentralityBin; // Centrality of current event
Int_t fNTracks; // N Tracks in the current event
+ Int_t fAODtrackCutBit; // Track filter bit for AOD tracks
+
// --- MC only -----------------------------------------------------------
AliStack *fStack; //! Ptr to stack
#include "TSystem.h"
#include "TFile.h"
#include "TPRegexp.h"
+#include "TDatabasePDG.h"
#include "AliStack.h"
#include "AliMCEvent.h"
#include "AliESDtrackCuts.h"
#include "AliESDInputHandler.h"
#include "AliESDpid.h"
+#include "AliAODInputHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODMCParticle.h"
#include "AliCentrality.h"
#include "AliTracker.h"
fESDHandler(NULL),
fPIDResponse(NULL),
fESD(NULL),
+ fAODHandler(NULL),
+ fAOD(NULL),
fMCEvent(NULL),
fStack(NULL),
}
//________________________________________________________________________
-Int_t AliAnalysisNetParticleHelper::SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent) {
+Int_t AliAnalysisNetParticleHelper::SetupEvent(AliESDInputHandler *esdHandler, AliAODInputHandler *aodHandler, AliMCEvent *mcEvent) {
// -- Setup Event
// -- Get ESD objects
- fESDHandler = esdHandler;
- fPIDResponse = esdHandler->GetPIDResponse();
- fESD = fESDHandler->GetEvent();
+ if(esdHandler){
+ fESDHandler = esdHandler;
+ fPIDResponse = esdHandler->GetPIDResponse();
+ fESD = fESDHandler->GetEvent();
+ }
+
+ // -- Get AOD objects
+ else if(aodHandler){
+ fAODHandler = aodHandler;
+ fPIDResponse = aodHandler->GetPIDResponse();
+ fAOD = fAODHandler->GetEvent();
+ }
// -- Get MC objects
fMCEvent = mcEvent;
// > 0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
// > 0 1 2 3 4 5 6 7 8 9
- AliCentrality *centrality = fESD->GetCentrality();
-
+ AliCentrality *centrality = NULL;
+
+ if(esdHandler)
+ centrality = fESD->GetCentrality();
+ else if(aodHandler)
+ centrality = fAOD->GetHeader()->GetCentralityP();
+
Int_t centBin = centrality->GetCentralityClass10("V0M");
if (centBin == 0)
fCentralityBin = centrality->GetCentralityClass5("V0M");
fCentralityBin = centBin + 1;
else
fCentralityBin = -2;
-
+
// -- Stay within the max centrality bin
if (fCentralityBin >= fCentralityBinMax)
fCentralityBin = -2;
-
+
fCentralityPercentile = centrality->GetCentralityPercentile("V0M");
- // -- Update TPC pid with eta correction
- UpdateEtaCorrectedTPCPid();
+ if(esdHandler){
+ // -- Update TPC pid with eta correction (only for ESDs?) XXX
+ UpdateEtaCorrectedTPCPid();
+ }
return 0;
}
for (Int_t ii = 0; ii < fNTriggers; ++ii)
aTriggerFired[ii] = kFALSE;
- if ((fESDHandler->IsEventSelected() & AliVEvent::kMB)) aTriggerFired[0] = kTRUE;
- if ((fESDHandler->IsEventSelected() & AliVEvent::kCentral)) aTriggerFired[1] = kTRUE;
- if ((fESDHandler->IsEventSelected() & AliVEvent::kSemiCentral)) aTriggerFired[2] = kTRUE;
- if ((fESDHandler->IsEventSelected() & AliVEvent::kEMCEJE)) aTriggerFired[3] = kTRUE;
- if ((fESDHandler->IsEventSelected() & AliVEvent::kEMCEGA)) aTriggerFired[4] = kTRUE;
-
+ if(fESDHandler){
+ if ((fESDHandler->IsEventSelected() & AliVEvent::kMB)) aTriggerFired[0] = kTRUE;
+ if ((fESDHandler->IsEventSelected() & AliVEvent::kCentral)) aTriggerFired[1] = kTRUE;
+ if ((fESDHandler->IsEventSelected() & AliVEvent::kSemiCentral)) aTriggerFired[2] = kTRUE;
+ if ((fESDHandler->IsEventSelected() & AliVEvent::kEMCEJE)) aTriggerFired[3] = kTRUE;
+ if ((fESDHandler->IsEventSelected() & AliVEvent::kEMCEGA)) aTriggerFired[4] = kTRUE;
+ }
+ else if(fAODHandler){
+ if ((fAODHandler->IsEventSelected() & AliVEvent::kMB)) aTriggerFired[0] = kTRUE;
+ if ((fAODHandler->IsEventSelected() & AliVEvent::kCentral)) aTriggerFired[1] = kTRUE;
+ if ((fAODHandler->IsEventSelected() & AliVEvent::kSemiCentral)) aTriggerFired[2] = kTRUE;
+ if ((fAODHandler->IsEventSelected() & AliVEvent::kEMCEJE)) aTriggerFired[3] = kTRUE;
+ if ((fAODHandler->IsEventSelected() & AliVEvent::kEMCEGA)) aTriggerFired[4] = kTRUE;
+ }
+
Bool_t isTriggered = kFALSE;
for (Int_t ii=0; ii<fNTriggers; ++ii) {
// -- 2 - No Vertex
++iCut;
- const AliESDVertex* vtxESD = fESD->GetPrimaryVertexTracks();
- if (!vtxESD)
- aEventCuts[iCut] = 1;
+ const AliESDVertex* vtxESD = NULL;
+ const AliAODVertex* vtxAOD = NULL;
+ if (fESD){
+ vtxESD = fESD->GetPrimaryVertexTracks();
+ if (!vtxESD)
+ aEventCuts[iCut] = 1;
+ }
+ else if (fAOD){
+ vtxAOD = fAOD->GetPrimaryVertex();
+ if (!vtxAOD)
+ aEventCuts[iCut] = 1;
+ }
// -- 3 - Vertex z outside cut window
++iCut;
if(TMath::Abs(vtxESD->GetZv()) > fVertexZMax)
aEventCuts[iCut] = 1;
}
+ else if(vtxAOD){
+ if(TMath::Abs(vtxAOD->GetZ()) > fVertexZMax)
+ aEventCuts[iCut] = 1;
+ }
else
aEventCuts[iCut] = 1;
-
+
// -- 4 - Centrality = -1 (no centrality or not hadronic)
++iCut;
if(fCentralityBin == -1.)
return kTRUE;
}
//________________________________________________________________________
+Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedBasicCharged(AliAODMCParticle *particle) {
+ // -- Check if MC particle is accepted for basic parameters
+
+ if (!particle)
+ return kFALSE;
+
+ // -- check if PDF code exists
+ if (!(TDatabasePDG::Instance()->GetParticle(particle->PdgCode())))
+ return kFALSE;
+
+ // -- check if charged
+ if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() == 0.0)
+ return kFALSE;
+
+ // -- check if physical primary
+ if(!particle->IsPhysicalPrimary())
+ return kFALSE;
+
+ return kTRUE;
+}
+//________________________________________________________________________
Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedBasicNeutral(TParticle *particle, Int_t idxMC) {
// -- Check if MC particle is accepted for basic parameters
return kTRUE;
}
+//________________________________________________________________________
+Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedBasicNeutral(AliAODMCParticle *particle) {
+ // -- Check if MC particle is accepted for basic parameters
+
+ if (!particle)
+ return kFALSE;
+ // -- check if PDF code exists
+ if (!(TDatabasePDG::Instance()->GetParticle(particle->PdgCode())))
+ return kFALSE;
+
+ // -- check if charged
+ if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() != 0.0)
+ return kFALSE;
+
+ // -- check if physical primary
+ if(!particle->IsPhysicalPrimary())
+ return kFALSE;
+
+ return kTRUE;
+}
//________________________________________________________________________
Bool_t AliAnalysisNetParticleHelper::IsParticleAcceptedRapidity(TParticle *particle, Double_t &yP) {
// -- Check if particle is accepted
* Accept Track Methods - public
* ---------------------------------------------------------------------------------
*/
-
+
//________________________________________________________________________
Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedBasicCharged(AliESDtrack *track) {
// -- Check if track is accepted
// > for basic parameters
-
+
if (!track)
return kFALSE;
if (!track->GetInnerParam())
return kFALSE;
+ return kTRUE;
+}
+
+//________________________________________________________________________
+Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedBasicCharged(AliAODTrack *track) {
+ // -- Check if track is accepted
+ // > for basic parameters
+
+ if (!track)
+ return kFALSE;
+
+ if (track->Charge() == 0)
+ return kFALSE;
+
+ //// -- Get momentum for dEdx --> returns always ZERO XXX
+ //if (!track->GetInnerParam())
+ // return kFALSE;
+
return kTRUE;
}
//________________________________________________________________________
-Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedRapidity(AliESDtrack *track, Double_t &yP) {
+Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedRapidity(AliVTrack *track, Double_t &yP) {
// -- Check if track is accepted
// > in rapidity
// > return 0 if not accepted
Double_t pvec[3];
track->GetPxPyPz(pvec);
- Double_t p = track->GetP();
+ Double_t p = track->P();
Double_t eP = TMath::Sqrt(p*p + mP*mP);
yP = 0.5 * TMath::Log((eP + pvec[2]) / (eP - pvec[2]));
}
//________________________________________________________________________
-Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedDCA(AliESDtrack *track) {
+Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedDCA(AliESDtrack *track) { //ONLY FOR ESDs so far XXX
// -- Check if track is accepted
// > for DCA, if both SPD layers have hits
}
//________________________________________________________________________
-Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedPID(AliESDtrack *track, Double_t* pid) {
+Bool_t AliAnalysisNetParticleHelper::IsTrackAcceptedPID(AliVTrack *track, Double_t* pid) {
// -- Check if track is accepted
// > provides TPC and TOF nSigmas to the argument
// -- Fill event statistics
for (Int_t idx = 0; idx < fHEventStatMax ; ++idx) {
+
if (aEventCuts[idx])
isRejected = kTRUE;
else
class AliStack;
class AliPIDResponse;
class AliESDInputHandler;
+class AliAODInputHandler;
+class AliAODEvent;
+class AliAODTrack;
+class AliAODMCParticle;
class AliAnalysisNetParticleHelper : public TNamed {
Int_t Initialize(Bool_t isMC);
/** Setup Event */
- Int_t SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent);
+ Int_t SetupEvent(AliESDInputHandler *esdHandler, AliAODInputHandler *aodHandler, AliMCEvent *mcEvent);
/*
* ---------------------------------------------------------------------------------
*/
/** Check if charged MC particle is accepted for basic parameters */
+ /** NOT possible for AODs (AliAODMCParticle NOT from TParticle)*/
Bool_t IsParticleAcceptedBasicCharged(TParticle *particle, Int_t idxMC);
+ Bool_t IsParticleAcceptedBasicCharged(AliAODMCParticle *particle);
/** Check if neutral MC particle is accepted for basic parameters */
+ /** NOT possible for AODs (AliAODMCParticle NOT from TParticle)*/
Bool_t IsParticleAcceptedBasicNeutral(TParticle *particle, Int_t idxMC);
+ Bool_t IsParticleAcceptedBasicNeutral(AliAODMCParticle *particle);
/** Check if MC particle is accepted for Rapidity */
Bool_t IsParticleAcceptedRapidity(TParticle *particle, Double_t &yP);
*/
/** Check if track is accepted for basic parameters */
+ /** NOT possible with AliVTrack (GetInnerParam returns NULL) */
Bool_t IsTrackAcceptedBasicCharged(AliESDtrack *track);
+ Bool_t IsTrackAcceptedBasicCharged(AliAODTrack *track);
/** Check if track is accepted for Rapidity */
- Bool_t IsTrackAcceptedRapidity(AliESDtrack *track, Double_t &yP);
+ Bool_t IsTrackAcceptedRapidity(AliVTrack *track, Double_t &yP);
/** Check if track is accepted for DCA */
Bool_t IsTrackAcceptedDCA(AliESDtrack *track);
/** Check if track is accepted for PID */
- Bool_t IsTrackAcceptedPID(AliESDtrack *track, Double_t *pid);
+ Bool_t IsTrackAcceptedPID(AliVTrack *track, Double_t *pid);
/*
* ---------------------------------------------------------------------------------
AliESDInputHandler *fESDHandler; //! Ptr to ESD handler
AliPIDResponse *fPIDResponse; //! Ptr to PID response Object
AliESDEvent *fESD; //! Ptr to ESD event
+ AliAODInputHandler *fAODHandler; //! Ptr to AOD handler
+ AliAODEvent *fAOD; //! Ptr to AOD event
AliMCEvent *fMCEvent; //! Ptr to MC event
AliStack *fStack; //! Ptr to stack
#include "AliAnalysisTaskNetParticle.h"
#include "AliGenEventHeader.h"
#include "AliCentrality.h"
+#include "AliAODEvent.h"
+#include "AliAODInputHandler.h"
using namespace std;
fESDTrackCutsBkg(NULL),
fESDTrackCutsEff(NULL),
+ fAOD(NULL),
+ fAODHandler(NULL),
+
fIsMC(kFALSE),
+ fIsAOD(kFALSE),
fESDTrackCutMode(0),
fModeEffCreation(0),
fModeDCACreation(0),
fEtaMax(0.9),
fPtRange(new Float_t[2]),
- fPtRangeEff(new Float_t[2]) {
+ fPtRangeEff(new Float_t[2]),
+
+ fAODtrackCutBit(1024){
// Constructor
AliLog::SetClassDebugLevel("AliAnalysisTaskNetParticle",10);
// ------------------------------------------------------------------
// -- Add histograms from efficiency/contamination class
// ------------------------------------------------------------------
- if (fIsMC && fModeEffCreation == 1) {
+ if ((fIsAOD||fIsMC) && fModeEffCreation == 1) {
fOutListEff->Add(fEffCont->GetHnEff());
fOutListCont->Add(fEffCont->GetHnCont());
}
// -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
// -- Process Efficiency / Contamination Determination
// -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
- if (fIsMC && fModeEffCreation == 1)
+ if ((fIsMC||fIsAOD) && fModeEffCreation == 1)
fEffCont->Process();
// -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
// ------------------------------------------------------------------
// -- Create / Initialize Efficiency/Contamination
// ------------------------------------------------------------------
- if (fIsMC && fModeEffCreation == 1) {
+ if ((fIsMC||fIsAOD) && fModeEffCreation == 1) {
fEffCont = new AliAnalysisNetParticleEffCont;
- fEffCont->Initialize(fESDTrackCutsEff, fHelper);
+ fEffCont->Initialize(fESDTrackCutsEff, fHelper,fAODtrackCutBit);
}
// ------------------------------------------------------------------
// ------------------------------------------------------------------
if (fModeDistCreation == 1) {
fDist = new AliAnalysisNetParticleDistribution;
- fDist->Initialize(fHelper, fESDTrackCuts, fIsMC, fPtRange, fEtaMax);
+ fDist->Initialize(fHelper, fESDTrackCuts, fIsMC, fPtRange, fEtaMax,fAODtrackCutBit);
}
// ------------------------------------------------------------------
// -- ESD Event
// ------------------------------------------------------------------
- if (SetupESDEvent() < 0) {
- AliError("Setup ESD Event failed");
- return -1;
+ if(!fIsAOD){
+ if (SetupESDEvent() < 0) {
+ AliError("Setup ESD Event failed");
+ return -1;
+ }
}
+ // -- AOD Event
+ // ------------------------------------------------------------------
+ else {
+ if (SetupAODEvent() < 0) {
+ AliError("Setup AOD Event failed");
+ return -1;
+ }
+ }
+
// -- Setup MC Event
// ------------------------------------------------------------------
if (fIsMC && SetupMCEvent() < 0) {
// -- Setup Event for Helper / EffCont / DCA / Dist classes
// ------------------------------------------------------------------
- fHelper->SetupEvent(fESDHandler, fMCEvent);
+ fHelper->SetupEvent(fESDHandler, fAODHandler, fMCEvent);
if (fIsMC && fModeEffCreation)
fEffCont->SetupEvent(fESDHandler, fMCEvent);
+ if (fIsAOD && fModeEffCreation)
+ fEffCont->SetupEvent(fAODHandler);
+
if (fModeDCACreation == 1)
fDCA->SetupEvent(fESDHandler, fMCEvent);
if (fModeDistCreation == 1)
- fDist->SetupEvent(fESDHandler, fMCEvent);
+ fDist->SetupEvent(fESDHandler, fAODHandler, fMCEvent);
// -- Evaluate Event cuts
// ------------------------------------------------------------------
return 0;
}
+//________________________________________________________________________
+Int_t AliAnalysisTaskNetParticle::SetupAODEvent() {
+ // -- Setup AOD Event
+ // > return 0 for success
+ // > return -1 for failed setup
+
+ fAODHandler= dynamic_cast<AliAODInputHandler*>
+ (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+ if (!fAODHandler) {
+ AliError("Could not get AOD input handler");
+ return -1;
+ }
+
+ fAOD = fAODHandler->GetEvent();
+ if (!fAOD) {
+ AliError("Could not get AOD event");
+ return -1;
+ }
+
+ // -- Check PID response
+ // ------------------------------------------------------------------
+ if (!fAODHandler->GetPIDResponse()) {
+ AliError("Could not get PID response");
+ return -1;
+ }
+
+ // -- Check Vertex
+ // ------------------------------------------------------------------
+ if (!fAOD->GetPrimaryVertex()) {
+ AliError("Could not get primary vertex");
+ return -1;
+ }
+
+ // -- Check Centrality
+ // ------------------------------------------------------------------
+ if (!fAOD->GetHeader()->GetCentralityP()) {
+ AliError("Could not get centrality");
+ return -1;
+ }
+
+ return 0;
+}
+
//________________________________________________________________________
Int_t AliAnalysisTaskNetParticle::SetupMCEvent() {
// -- Setup MC Event
// -- Reset ESD Event
fESD = NULL;
+ // -- Reset ESD Event
+ fAOD = NULL;
+
// -- Reset MC Event
if (fIsMC)
fMCEvent = NULL;
class AliESDEvent;
class AliESDInputHandler;
+class AliAODEvent;
+class AliAODInputHandler;
class AliMCEvent;
class AliStack;
class AliPIDResponse;
*/
void SetIsMC() {fIsMC = kTRUE;}
+ void SetIsAOD(Bool_t b) {fIsAOD = b;}
void SetUseQATHnSparse(Bool_t b) {fUseQATHnSparse = b;}
void SetESDTrackCutMode(Int_t i) {fESDTrackCutMode = i;}
void SetModeEffCreation(Int_t i) {fModeEffCreation = i;}
void SetModeDCACreation(Int_t i) {fModeDCACreation = i;}
void SetModeDistCreation(Int_t i) {fModeDistCreation = i;}
+
void SetEtaMax(Float_t f) {fEtaMax = f;}
void SetPtRange(Float_t f1, Float_t f2) {fPtRange[0] = f1; fPtRange[1] = f2;}
void SetPtRangeEff(Float_t f1, Float_t f2) {fPtRangeEff[0] = f1; fPtRangeEff[1] = f2;}
+ void SetTrackFilterBit(Int_t i) {fAODtrackCutBit = i;}
+
/*
* ---------------------------------------------------------------------------------
* Setter - Pass-Through
/** Setup ESD Event */
Int_t SetupESDEvent();
+ /** Setup AOD Event */
+ Int_t SetupAODEvent();
+
/** Setup MC Event */
Int_t SetupMCEvent();
AliESDtrackCuts *fESDTrackCuts; //! ESD cuts
AliESDtrackCuts *fESDTrackCutsBkg; //! ESD cuts for Bkg
AliESDtrackCuts *fESDTrackCutsEff; //! ESD cuts for efficiency determination -> larger pt Range
+
+ // --- AOD only ----------------------------------------------------------
+
+ AliAODEvent *fAOD; //! Ptr to AOD event
+ AliAODInputHandler *fAODHandler; //! Ptr to AOD Handler
// --- Flags -------------------------------------------------------------
Bool_t fIsMC; // Is MC event
+ Bool_t fIsAOD; // analysis mode : 0 = ESDs | 1 = AODs
Int_t fESDTrackCutMode; // ESD track cut mode : 0 = clean | 1 dirty
Int_t fModeEffCreation ; // Correction creation mode : 1 = on | 0 = off
Int_t fModeDCACreation; // DCA creation mode : 1 = on | 0 = off
Float_t *fPtRange; // Array of pt [min,max]
Float_t *fPtRangeEff; // Array of pt [min,max] for efficiency
+ Int_t fAODtrackCutBit; // Track filter bit for AOD tracks
+
// -----------------------------------------------------------------------
ClassDef(AliAnalysisTaskNetParticle, 1);
AliAnalysisTask *AddTaskNetParticle(const Char_t * name = "jthaeder_NetProton",
Bool_t isModeDist, Bool_t isModeEff, Bool_t isModeDCA, Bool_t useQAThnSparse = kFALSE,
- Bool_t isCreateCSC = kFALSE) {
+ Bool_t isCreateCSC = kFALSE, Bool_t isModeAOD = kFALSE) {
TString sName(name);
task->SetParticleSpecies(AliPID::kKaon);
task->SetControlParticleSpecies(3122, kTRUE, "Lambda"); /// maybe something else ...
minPt = 0.2; maxPt = 0.4;
- minPtEff = 0.1; maxPtEff = 0.8;
+ minPtEff = 0.1; maxPtEff = 2.5;
minPtForTOF = 0.8;
}
else {
task->SetModeDCACreation(1); // => 1 = on | 0 = off (default)
if (isModeDist)
task->SetModeDistCreation(1); // => 1 = on | 0 = off (default)
-
+ if(isModeAOD){
+ task->SetIsAOD(1); // => 1 = AODs | 0 = ESDs
+ task->SetTrackFilterBit(1024); // 1024 = RAA cuts
+ }
// -- Enable QA plots
if (useQAThnSparse)
task->SetUseQATHnSparse(kTRUE);