From d731501a24d5c7e2102e9cdff83d35eb584e7d45 Mon Sep 17 00:00:00 2001 From: prino Date: Tue, 11 Sep 2012 13:23:15 +0000 Subject: [PATCH] MC implementations for D0-HFE correlation (Matthias) --- PWGHF/CMakelibPWGHFcorrelationHF.pkg | 3 + .../AliAnalysisTaskDxHFECorrelation.cxx | 14 +- .../AliAnalysisTaskDxHFEParticleSelection.cxx | 32 +- .../AliAnalysisTaskDxHFEParticleSelection.h | 3 + PWGHF/correlationHF/AliDxHFECorrelation.cxx | 4 +- .../AliDxHFEParticleSelection.cxx | 105 ++++- .../correlationHF/AliDxHFEParticleSelection.h | 32 +- .../AliDxHFEParticleSelectionD0.cxx | 131 ++++-- .../AliDxHFEParticleSelectionD0.h | 10 + .../AliDxHFEParticleSelectionEl.cxx | 137 ++++-- .../AliDxHFEParticleSelectionEl.h | 21 +- .../AliDxHFEParticleSelectionMCD0.cxx | 197 +++++++++ .../AliDxHFEParticleSelectionMCD0.h | 62 +++ .../AliDxHFEParticleSelectionMCEl.cxx | 235 +++++++++++ .../AliDxHFEParticleSelectionMCEl.h | 71 ++++ PWGHF/correlationHF/AliDxHFEToolsMC.cxx | 399 ++++++++++++++++++ PWGHF/correlationHF/AliDxHFEToolsMC.h | 196 +++++++++ 17 files changed, 1530 insertions(+), 122 deletions(-) create mode 100644 PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.cxx create mode 100644 PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.h create mode 100644 PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.cxx create mode 100644 PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.h create mode 100644 PWGHF/correlationHF/AliDxHFEToolsMC.cxx create mode 100644 PWGHF/correlationHF/AliDxHFEToolsMC.h diff --git a/PWGHF/CMakelibPWGHFcorrelationHF.pkg b/PWGHF/CMakelibPWGHFcorrelationHF.pkg index 1751c7249d6..ba0893b3222 100644 --- a/PWGHF/CMakelibPWGHFcorrelationHF.pkg +++ b/PWGHF/CMakelibPWGHFcorrelationHF.pkg @@ -30,6 +30,9 @@ set ( CLASS_HDRS AliDxHFEParticleSelection.h AliDxHFEParticleSelectionD0.h AliDxHFEParticleSelectionEl.h + AliDxHFEParticleSelectionMCD0.h + AliDxHFEParticleSelectionMCEl.h + AliDxHFEToolsMC.h AliDxHFECorrelation.h AliAnalysisTaskDxHFEParticleSelection.h AliAnalysisTaskDxHFECorrelation.h diff --git a/PWGHF/correlationHF/AliAnalysisTaskDxHFECorrelation.cxx b/PWGHF/correlationHF/AliAnalysisTaskDxHFECorrelation.cxx index 260154f7046..7203f1ce256 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskDxHFECorrelation.cxx +++ b/PWGHF/correlationHF/AliAnalysisTaskDxHFECorrelation.cxx @@ -26,7 +26,9 @@ #include "AliAnalysisTaskDxHFECorrelation.h" #include "AliDxHFECorrelation.h" #include "AliDxHFEParticleSelectionD0.h" +#include "AliDxHFEParticleSelectionMCD0.h" #include "AliDxHFEParticleSelectionEl.h" +#include "AliDxHFEParticleSelectionMCEl.h" #include "AliAnalysisManager.h" #include "AliLog.h" #include "AliESDInputHandler.h" @@ -143,15 +145,17 @@ void AliAnalysisTaskDxHFECorrelation::UserCreateOutputObjects() default: selectionD0Options+="FillD0D0bar "; } - fD0s=new AliDxHFEParticleSelectionD0(selectionD0Options); + if(fUseMC) fD0s=new AliDxHFEParticleSelectionMCD0(selectionD0Options); + else fD0s=new AliDxHFEParticleSelectionD0(selectionD0Options); fD0s->SetCuts(fCutsD0); - fD0s->InitControlObjects(); + fD0s->Init(); //Electrons - fElectrons=new AliDxHFEParticleSelectionEl; + if(fUseMC) fElectrons=new AliDxHFEParticleSelectionMCEl; + else fElectrons=new AliDxHFEParticleSelectionEl; fElectrons->SetCuts(fPID, AliDxHFEParticleSelectionEl::kCutPID); fElectrons->SetCuts(fCutsHFE, AliDxHFEParticleSelectionEl::kCutHFE); - fElectrons->InitControlObjects(); + fElectrons->Init(); //Correlation fCorrelation=new AliDxHFECorrelation; @@ -160,7 +164,7 @@ void AliAnalysisTaskDxHFECorrelation::UserCreateOutputObjects() // Fix for merging: // Retrieving the individual objects created - // and storing them instead of fD0s,fElectrons etc.. + // and storing them instead of fD0s, fElectrons etc.. TList *list =(TList*)fD0s->GetControlObjects(); TObject *obj=NULL; diff --git a/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.cxx b/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.cxx index 71e4662c39e..7d39fd65a6b 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.cxx +++ b/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.cxx @@ -26,6 +26,7 @@ #include "AliAnalysisTaskDxHFEParticleSelection.h" #include "AliDxHFEParticleSelection.h" #include "AliDxHFEParticleSelectionD0.h" +#include "AliDxHFEParticleSelectionMCD0.h" #include "AliAnalysisManager.h" #include "AliAnalysisCuts.h" #include "AliLog.h" @@ -52,6 +53,7 @@ AliAnalysisTaskDxHFEParticleSelection::AliAnalysisTaskDxHFEParticleSelection(con , fCuts(NULL) , fSelector(NULL) , fUseMC(kFALSE) + , fFillOnlyD0D0bar(0) { // constructor // @@ -101,22 +103,18 @@ void AliAnalysisTaskDxHFEParticleSelection::UserCreateOutputObjects() fOutput = new TList; fOutput->SetOwner(); - std::auto_ptr selector(new AliDxHFEParticleSelectionD0); - if (!selector.get()) return; - if (fCuts) { - AliRDHFCutsD0toKpi* cuts=dynamic_cast(fCuts); - if (!cuts) { - AliFatal(Form("cut object %s is of incorrect type %s, expecting AliRDHFCutsD0toKpi", fCuts->GetName(), fCuts->ClassName())); - return; - } - selector->SetCuts(fCuts); - AliInfo(Form("Initializing particle selection %s, using %d pt bins", selector->GetName(), cuts->GetNPtBins())); - } else { - AliWarning("no cut object available for particle selection"); + // setting up for D0s + TString selectionD0Options; + switch (fFillOnlyD0D0bar) { + case 1: selectionD0Options+="FillOnlyD0 "; break; + case 2: selectionD0Options+="FillOnlyD0bar "; break; + default: selectionD0Options+="FillD0D0bar "; } - selector->InitControlObjects(); - fSelector=selector.release(); + if(fUseMC) fSelector=new AliDxHFEParticleSelectionMCD0(selectionD0Options); + else fSelector=new AliDxHFEParticleSelectionD0(selectionD0Options); + fSelector->SetCuts(fCuts); + fSelector->Init(); // Retrieving the list containing histos and THnSparse // and storing them instead of fSelector @@ -186,7 +184,7 @@ void AliAnalysisTaskDxHFEParticleSelection::UserExec(Option_t* /*option*/) return; } - Int_t nInD0toKpi = inputArray->GetEntriesFast(); + // Int_t nInD0toKpi = inputArray->GetEntriesFast(); fSelector->HistogramEventProperties(AliDxHFEParticleSelection::kEventsSel); @@ -206,14 +204,12 @@ void AliAnalysisTaskDxHFEParticleSelection::UserExec(Option_t* /*option*/) //Test to see if I have read in D0s and retrieved them after selection Int_t nD0Selected = selectedTracks->GetEntriesFast(); - printf("Number of D0->Kpi Start: %d , End: %d\n",nInD0toKpi,nD0Selected); - fSelector->HistogramEventProperties(AliDxHFEParticleSelection::kEventsD0); + fSelector->HistogramEventProperties(AliDxHFEParticleSelection::kEventsWithParticle); for(Int_t iD0toKpi = 0; iD0toKpi < nD0Selected; iD0toKpi++) { AliAODRecoDecayHF2Prong *particle = (AliAODRecoDecayHF2Prong*)selectedTracks->UncheckedAt(iD0toKpi); if (!particle) continue; - cout << "D0s inv mass: " << particle->InvMassD0() << endl; } PostData(1, fOutput); diff --git a/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.h b/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.h index 1bb1bdd8a43..11ee0a311a6 100644 --- a/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.h +++ b/PWGHF/correlationHF/AliAnalysisTaskDxHFEParticleSelection.h @@ -49,6 +49,7 @@ class AliAnalysisTaskDxHFEParticleSelection : public AliAnalysisTaskSE { virtual void Terminate(Option_t*); void SetOption(const char* opt) { fOption = opt; } + void SetFillOnlyD0D0bar(Int_t flagfill){fFillOnlyD0D0bar=flagfill;} virtual void SetUseMC(Bool_t useMC){fUseMC=useMC;} virtual void SetCuts(AliAnalysisCuts* cuts){fCuts=cuts;} Bool_t GetUseMC() const {return fUseMC;} @@ -68,6 +69,8 @@ class AliAnalysisTaskDxHFEParticleSelection : public AliAnalysisTaskSE { AliAnalysisCuts* fCuts; // Cuts AliDxHFEParticleSelection* fSelector; // selector instance bool fUseMC; // use MC info + Int_t fFillOnlyD0D0bar; // flag to set what to fill (0 = both, 1 = D0 only, 2 = D0bar only) + ClassDef(AliAnalysisTaskDxHFEParticleSelection, 2); }; diff --git a/PWGHF/correlationHF/AliDxHFECorrelation.cxx b/PWGHF/correlationHF/AliDxHFECorrelation.cxx index c0fd9ab5358..446c21648ac 100644 --- a/PWGHF/correlationHF/AliDxHFECorrelation.cxx +++ b/PWGHF/correlationHF/AliDxHFECorrelation.cxx @@ -67,7 +67,7 @@ AliDxHFECorrelation::AliDxHFECorrelation(const char* name) const char* AliDxHFECorrelation::fgkEventControlBinNames[]={ "nEventsAll", "nEventsSelected", - "nEventsD0" + "nEventsD0", "nEventsD0e" }; @@ -162,7 +162,7 @@ int AliDxHFECorrelation::Init() } for (iLabel=0; iLabelGetXaxis()->SetBinLabel(iLabel, fgkEventControlBinNames[iLabel]); + hEventControl->GetXaxis()->SetBinLabel(iLabel+1, fgkEventControlBinNames[iLabel]); fCorrProperties=CorrProperties.release(); AddControlObject(fCorrProperties); diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelection.cxx b/PWGHF/correlationHF/AliDxHFEParticleSelection.cxx index 2a5201ca1a2..4d72677c14f 100644 --- a/PWGHF/correlationHF/AliDxHFEParticleSelection.cxx +++ b/PWGHF/correlationHF/AliDxHFEParticleSelection.cxx @@ -49,7 +49,7 @@ AliDxHFEParticleSelection::AliDxHFEParticleSelection(const char* name, const cha , fControlObjects(NULL) , fhEventControl(NULL) , fhTrackControl(NULL) - , fUseMC(false) + , fUseMC(0) , fVerbosity(0) { // constructor @@ -62,7 +62,18 @@ AliDxHFEParticleSelection::AliDxHFEParticleSelection(const char* name, const cha const char* AliDxHFEParticleSelection::fgkEventControlBinNames[]={ "nEventsAll", "nEventsSelected", - "nEventsD0" + "nEventsWParticle" +}; + +const char* AliDxHFEParticleSelection::fgkTrackControlBinNames[]={ + "nTrackAll", + "nTrackSelected", +}; + +const char* AliDxHFEParticleSelection::fgkTrackControlDimNames[]={ + "Pt", + "Phi", + "Eta" }; AliDxHFEParticleSelection::~AliDxHFEParticleSelection() @@ -76,9 +87,23 @@ AliDxHFEParticleSelection::~AliDxHFEParticleSelection() fhTrackControl=NULL; } +int AliDxHFEParticleSelection::Init() +{ + // + // Init part sel. Calls InitControlObjects() + // + + InitControlObjects(); + return 0; +} + + int AliDxHFEParticleSelection::InitControlObjects() { + // // init control objects + // TODO: Change to private now that have Init()? + // if (fVerbosity>0) { AliInfo("Setting up control objects"); } @@ -90,14 +115,63 @@ int AliDxHFEParticleSelection::InitControlObjects() fhEventControl=hEventControl.release(); for (int iLabel=0; iLabelGetXaxis()->SetBinLabel(iLabel, fgkEventControlBinNames[iLabel]); + fhEventControl->GetXaxis()->SetBinLabel(iLabel+1, fgkEventControlBinNames[iLabel]); AddControlObject(fhEventControl); fhTrackControl=hTrackControl.release(); + for (int iLabel=0; iLabelGetXaxis()->SetBinLabel(iLabel+1, fgkTrackControlBinNames[iLabel]); AddControlObject(fhTrackControl); return 0; } +THnSparse* AliDxHFEParticleSelection::CreateControlTHnSparse(const char* name, + int thnSize, + int* thnBins, + double* thnMin, + double* thnMax, + const char** binLabels) const +{ + // + // Creates THnSparse. + // + + AliInfo("Setting up THnSparse"); + + std::auto_ptr th(new THnSparseF(name, name, thnSize, thnBins, thnMin, thnMax)); + if (th.get()==NULL) { + return NULL; + } + for (int iLabel=0; iLabelGetAxis(iLabel)->SetTitle(binLabels[iLabel]); + + } + return th.release(); + +} + +THnSparse* AliDxHFEParticleSelection::DefineTHnSparse() const +{ + // + // Defines the THnSparse. For now, only calls CreatControlTHnSparse + + //TODO: Should make it more general. Or maybe one can use this here, and skip in PartSelEl? + + const int thnSize = 3; + const double Pi=TMath::Pi(); + TString name; + // 0 1 2 + // Pt Phi Eta + int thnBins[thnSize] = { 1000, 200, 500}; + double thnMin [thnSize] = { 0, 0, -1.}; + double thnMax [thnSize] = { 100, 2*Pi, 1.}; + + name.Form("%s info", GetName()); + + return CreateControlTHnSparse(name,thnSize,thnBins,thnMin,thnMax,fgkTrackControlDimNames); + +} + int AliDxHFEParticleSelection::AddControlObject(TObject* pObj) { /// add control object to list, the base class becomes owner of the object @@ -132,15 +206,34 @@ int AliDxHFEParticleSelection::HistogramEventProperties(int bin) int AliDxHFEParticleSelection::HistogramParticleProperties(AliVParticle* p, int selected) { /// histogram particle properties + if (!p) return -EINVAL; if (!fControlObjects) return 0; // TODO: use enums for the bins of the control histogram - fhTrackControl->Fill(0); - if (selected) fhTrackControl->Fill(1); + fhTrackControl->Fill(kTrackAll); + if (selected) fhTrackControl->Fill(kTrackSel); return 0; } +int AliDxHFEParticleSelection::DefineParticleProperties(AliVParticle* p, Double_t* data, int dimension) const +{ + // fill the data array from the particle data + if (!data) return -EINVAL; + int i=0; + // TODO: this corresponds to the THnSparse dimensions which is available in the same class + // use this consistently + const int requiredDimension=3; + if (dimension!=requiredDimension) { + // TODO: think about filling only the available data and throwing a warning + return -ENOSPC; + } + data[i++]=p->Pt(); + data[i++]=p->Phi(); + data[i++]=p->Eta(); + return i; +} + TObjArray* AliDxHFEParticleSelection::Select(const AliVEvent* pEvent) { /// create selection from 'Tracks' member of the event, @@ -152,7 +245,7 @@ TObjArray* AliDxHFEParticleSelection::Select(const AliVEvent* pEvent) int nofTracks=pEvent->GetNumberOfTracks(); for (int itrack=0; itrackGetTrack(itrack); - int selectionCode=IsSelected(track); + int selectionCode=IsSelected(track,pEvent); HistogramParticleProperties(track, selectionCode); if (selectionCode==0) continue; selectedTracks->Add(track); diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelection.h b/PWGHF/correlationHF/AliDxHFEParticleSelection.h index c6712e7e11d..62f29accb2d 100644 --- a/PWGHF/correlationHF/AliDxHFEParticleSelection.h +++ b/PWGHF/correlationHF/AliDxHFEParticleSelection.h @@ -45,9 +45,15 @@ class AliDxHFEParticleSelection : public TNamed { enum { kEventsAll = 0, kEventsSel, - kEventsD0, + kEventsWithParticle, kNEventPropertyLabels }; + + enum { + kTrackAll = 0, + kTrackSel, + kNTrackPropertyLabels + }; /// set options void SetOption(const char* opt) { fOption = opt; } @@ -55,6 +61,7 @@ class AliDxHFEParticleSelection : public TNamed { virtual Option_t* GetOption() const { return fOption;} /// init the control objects + virtual int Init(); virtual int InitControlObjects(); /// create selection from 'Tracks' member of the event, @@ -72,6 +79,9 @@ class AliDxHFEParticleSelection : public TNamed { /// histogram event properties virtual int HistogramEventProperties(int bin); + // TODO: function can be renamed to better describe what it's doing + virtual int DefineParticleProperties(AliVParticle* p, Double_t* date, int dimension) const; + /// check and add track to internal array int CheckAndAdd(AliVParticle* p); @@ -79,6 +89,8 @@ class AliDxHFEParticleSelection : public TNamed { // different types; a type cast check is implemented in the method virtual void SetCuts(TObject* /*cuts*/, int /*level*/=0) {} + // TODO: check whether that is needed, should be covered by the specific + // child implementation Bool_t GetUseMC() const {return fUseMC;} /// get selected tracks @@ -104,6 +116,16 @@ class AliDxHFEParticleSelection : public TNamed { void SetVerbosity(int verbosity) {fVerbosity=verbosity;} /// get verbosity int GetVerbosity() const {return fVerbosity;} + + /// create control THnSparse + THnSparse* CreateControlTHnSparse(const char* name, + int thnSize, + int* thnBins, + double* thnMin, + double* thnMax, + const char** binLabels) const; + + virtual THnSparse* DefineTHnSparse() const; protected: /// add control object to list, the base class becomes owner of the object @@ -112,12 +134,6 @@ class AliDxHFEParticleSelection : public TNamed { /// histogram particle properties virtual int HistogramParticleProperties(AliVParticle* p, int selected=1); - /// set from the specific child class implementation - /// TODO: think about removing this function - /// in principle that should not be necessary, the analysis should - /// be blind after the particle selection - void SetUseMC(bool useMC=true) {fUseMC=useMC;} - private: /// copy contructor prohibited AliDxHFEParticleSelection(const AliDxHFEParticleSelection&); @@ -135,6 +151,8 @@ class AliDxHFEParticleSelection : public TNamed { int fVerbosity; //! verbosity static const char* fgkEventControlBinNames[]; //! bin labels for event control histogram + static const char* fgkTrackControlBinNames[]; //! bin labels for track control histogram + static const char* fgkTrackControlDimNames[]; //! bin labels fro track thnsparse ClassDef(AliDxHFEParticleSelection, 2); diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.cxx b/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.cxx index c6c36bd142f..648d7da3477 100644 --- a/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.cxx +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.cxx @@ -51,6 +51,8 @@ AliDxHFEParticleSelectionD0::AliDxHFEParticleSelectionD0(const char* opt) , fD0Daughter1(NULL) , fCuts(NULL) , fFillOnlyD0D0bar(0) + , fD0InvMass(0.0) + , fPtBin(-1) { // constructor // @@ -85,35 +87,28 @@ AliDxHFEParticleSelectionD0::~AliDxHFEParticleSelectionD0() fCuts=NULL; } +const char* AliDxHFEParticleSelectionD0::fgkTrackControlBinNames[]={ + "Pt", + "Phi", + "Ptbin", + "D0InvMass", + "Eta" +}; + +const char* AliDxHFEParticleSelectionD0::fgkDgTrackControlBinNames[]={ + "Pt", + "Phi", + "Ptbin", + "D0InvMass", + "Eta" +}; + int AliDxHFEParticleSelectionD0::InitControlObjects() { /// init the control objects, can be overloaded by childs which should /// call AliDxHFEParticleSelection::InitControlObjects() explicitly - AliInfo("Setting up THnSparse"); - TString name; - const int thnSize = 4; - const double pi=TMath::Pi(); - - // TODO: theta? - // 0 1 2 3 - // mass Pt Phi Ptbin - int thnBins[thnSize] = { 200, 1000, 100, 100 }; - double thnMin [thnSize] = { 1.5648, 0, 0, 0 }; - double thnMax [thnSize] = { 2.1648, 100, (2*pi), 100 }; - - name.Form("%s info", GetName()); - std::auto_ptr D0Properties(new THnSparseF(name, name, thnSize, thnBins, thnMin, thnMax)); - - if (D0Properties.get()==NULL) { - return -ENOMEM; - } - int axis=0; - D0Properties->GetAxis(axis++)->SetTitle("D0 inv mass"); - D0Properties->GetAxis(axis++)->SetTitle("Pt"); - D0Properties->GetAxis(axis++)->SetTitle("Phi"); - D0Properties->GetAxis(axis++)->SetTitle("Ptbin"); - fD0Properties=D0Properties.release(); + fD0Properties=DefineTHnSparse(); AddControlObject(fD0Properties); //Adding control objects for the daughters @@ -123,9 +118,55 @@ int AliDxHFEParticleSelectionD0::InitControlObjects() return AliDxHFEParticleSelection::InitControlObjects(); } +THnSparse* AliDxHFEParticleSelectionD0::DefineTHnSparse() const +{ + // + // Defines the THnSparse. For now, only calls CreateControlTHnSparse + // TODO: remove pt?? (Have ptbin) + + const int thnSize2 = 5; + const double Pi=TMath::Pi(); + TString name; + name.Form("%s info", GetName()); + + // 0 1 2 3 4 + // Pt Phi Ptbin D0InvMass Eta + int thnBins[thnSize2] = { 1000, 200, 21, 200, 500 }; + double thnMin [thnSize2] = { 0, 0, 0, 1.5648, -1. }; + double thnMax [thnSize2] = { 100, 2*Pi, 20, 2.1648, 1. }; + + return CreateControlTHnSparse(name,thnSize2,thnBins,thnMin,thnMax,fgkTrackControlBinNames); + +} + +int AliDxHFEParticleSelectionD0::DefineParticleProperties(AliVParticle* p, Double_t* data, int dimension) const +{ + // fill the data array from the particle data + if (!data) return -EINVAL; + // AliAODTrack *track=(AliAODTrack*)p; + AliAODRecoDecayHF2Prong* track=dynamic_cast(p); + if (!track) return -ENODATA; + int i=0; + // TODO: this corresponds to the THnSparse dimensions which is available in the same class + // use this consistently + const int requiredDimension=5; + if (dimension!=requiredDimension) { + // TODO: think about filling only the available data and throwing a warning + return -ENOSPC; + } + data[i++]=track->Pt(); + data[i++]=track->Phi(); + data[i++]=fPtBin; + data[i++]=fD0InvMass; + data[i++]=track->Eta(); + + return i; +} + int AliDxHFEParticleSelectionD0::InitControlObjectsDaughters(TString name, int daughter) { - //Setting up Control objects for the daughters. + // Setting up Control objects for the daughters. + // Move to ParticleSelection?? AliInfo("Setting up daughter THnSparse"); const int thnSize2 = 5; @@ -141,12 +182,9 @@ int AliDxHFEParticleSelectionD0::InitControlObjectsDaughters(TString name, int d if (DaughterProperties.get()==NULL) { return -ENOMEM; } - int axis=0; - DaughterProperties->GetAxis(axis++)->SetTitle("Pt"); - DaughterProperties->GetAxis(axis++)->SetTitle("Phi"); - DaughterProperties->GetAxis(axis++)->SetTitle("Ptbin"); - DaughterProperties->GetAxis(axis++)->SetTitle("D0InvMass"); - DaughterProperties->GetAxis(axis++)->SetTitle("Eta"); + + for(int iLabel=0; iLabel< 5;iLabel++) + DaughterProperties->GetAxis(iLabel)->SetTitle(fgkDgTrackControlBinNames[iLabel]); if(daughter==0){ fD0Daughter0=DaughterProperties.release(); @@ -186,14 +224,17 @@ int AliDxHFEParticleSelectionD0::HistogramParticleProperties(AliVParticle* p, in // Only D0s are filled // TODO: Also include D0bar if ((selectionCode==1 || selectionCode==3) && fFillOnlyD0D0bar<2) { - Double_t invmassD0 = part->InvMassD0(); - AliDebug(3,Form("pt %f, nr pt bins %d",part->Pt(),fCuts->GetNPtBins())); - Int_t ptbin=fCuts->PtBin(part->Pt()); - Double_t D0Properties[] = {invmassD0,part->Pt(),part->Phi(),ptbin}; - Double_t KProperties[]={prongneg->Pt(),prongneg->Phi(),ptbin, invmassD0,prongneg->Eta()}; - Double_t piProperties[]={prongpos->Pt(),prongpos->Phi(),ptbin,invmassD0,prongpos->Eta()}; - - if(fD0Properties) fD0Properties->Fill(D0Properties); + fD0InvMass= part->InvMassD0(); + fPtBin=fCuts->PtBin(part->Pt()); + + // TODO: avoid repeated allocation of the arrays + Double_t KProperties[]={prongneg->Pt(),prongneg->Phi(),fPtBin, fD0InvMass,prongneg->Eta()}; + Double_t piProperties[]={prongpos->Pt(),prongpos->Phi(),fPtBin,fD0InvMass,prongpos->Eta()}; + + // TODO: make array a member, consistent dimensions for THnSparse and array + Double_t d0Properties[5]={0.0, 0.0, 0.0, 0.0, 0.0}; + DefineParticleProperties(p, d0Properties, 5); + if(fD0Properties) fD0Properties->Fill(d0Properties); if(fD0Daughter0) fD0Daughter0->Fill(piProperties); if(fD0Daughter1) fD0Daughter1->Fill(KProperties); } @@ -215,8 +256,8 @@ TObjArray* AliDxHFEParticleSelectionD0::Select(TObjArray* pTracks, const AliVEve if (!track) continue; int selectionCode=IsSelected(track,pEvent); HistogramParticleProperties(track, selectionCode); - //TODO: Also add selection for D0bar + //TODO: Also add selection for D0bar // Add track if it is either defined as D0(selectionCode==1) or both // D0bar and a D0 (selectionCode==3) if (! ((selectionCode==1 || selectionCode==3) && fFillOnlyD0D0bar<2)) continue; @@ -228,7 +269,9 @@ TObjArray* AliDxHFEParticleSelectionD0::Select(TObjArray* pTracks, const AliVEve int AliDxHFEParticleSelectionD0::IsSelected(AliVParticle* p, const AliVEvent* pEvent) { /// TODO: implement specific selection of D0 candidates - /// Could also return values based on where where selection "failed" + /// Could also return values based on where where selection "failed + /// Selected. Return 0 (none), 1(D0), 2(D0bar) or 3 (both) + int selectionCode=0; AliAODRecoDecayHF2Prong *d0 = dynamic_cast(p); @@ -241,9 +284,9 @@ int AliDxHFEParticleSelectionD0::IsSelected(AliVParticle* p, const AliVEvent* pE // AliRDHFCuts::IsSelected does not allow this AliRDHFCuts* cuts=const_cast(fCuts); if (!cuts) { - selectionCode=1; + selectionCode=0; } else if(cuts->IsInFiducialAcceptance(d0->Pt(),d0->Y(421)) ) { - // if(cuts->IsSelected(d0,AliRDHFCuts::kTracks,pEvent))fNentries->Fill(6); + Int_t ptbin=cuts->PtBin(d0->Pt()); if(ptbin==-1) { AliDebug(1,"Pt out of bounds"); @@ -255,7 +298,7 @@ int AliDxHFEParticleSelectionD0::IsSelected(AliVParticle* p, const AliVEvent* pE AliAODEvent* aod=NULL; if (pEvent) aod=dynamic_cast(const_cast(pEvent)); - // Selected. Return 0 (none), 1(D0), 2(D0bar) or 3 (both) + // Selected. Return 0 (none), 1 (D0), 2 (D0bar) or 3 (both) selectionCode=cuts->IsSelected(d0,AliRDHFCuts::kAll,aod); AliDebug(1,Form("Candidate is %d \n", selectionCode)); diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.h b/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.h index 93c43f151a3..f0bacb203d6 100644 --- a/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.h +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionD0.h @@ -33,6 +33,9 @@ class AliDxHFEParticleSelectionD0 : public AliDxHFEParticleSelection { /// overloaded from AliDxHFEParticleSelection: init the control objects virtual int InitControlObjects(); + virtual THnSparse* DefineTHnSparse() const; + // TODO: function can be renamed to better describe what it's doing + virtual int DefineParticleProperties(AliVParticle* p, Double_t* date, int dimension) const; //Function for daughter control objects //TODO: move to AliDxHFEParticleSelection to be used for several particles? @@ -49,6 +52,8 @@ class AliDxHFEParticleSelectionD0 : public AliDxHFEParticleSelection { //AliRDHFCutsD0toKpi GetCuts()const {return fCuts;} Int_t GetFillOnlyD0D0bar() const {return fFillOnlyD0D0bar;} + Int_t GetPtBin() const {return fPtBin;} + Double_t GetInvMass() const {return fD0InvMass;} protected: /// overloaded from AliDxHFEParticleSelection: histogram particle properties @@ -65,6 +70,11 @@ class AliDxHFEParticleSelectionD0 : public AliDxHFEParticleSelection { THnSparse* fD0Daughter1; // the particle properties of selected particles AliRDHFCuts* fCuts; //! pointer to external cuts object Int_t fFillOnlyD0D0bar; // flag to set what to fill (0 = both, 1 = D0 only, 2 = D0bar only) + Double_t fD0InvMass; // D0InvMass + Int_t fPtBin; // Pt Bin + + static const char* fgkTrackControlBinNames[]; //! bin labels for track control histogram + static const char* fgkDgTrackControlBinNames[]; //! bin labels for track control histogram // TODO: at the moment the dimensions of the different THnSparse objects are different // needs to be consolidated diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.cxx b/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.cxx index 0835d9740f6..27677052b49 100644 --- a/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.cxx +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.cxx @@ -76,42 +76,43 @@ AliDxHFEParticleSelectionEl::~AliDxHFEParticleSelectionEl() delete fCFM; fCFM=NULL; } + if(fWhichCut){ + delete fWhichCut; + fWhichCut=NULL; + } // NOTE: external objects fPID and fCuts are not deleted here fPID=NULL; fCuts=NULL; } - -int AliDxHFEParticleSelectionEl::InitControlObjects() -{ - /// init control and monitoring objects - AliInfo("Electron THnSparse"); - const int thnSize = 3; - const double Pi=TMath::Pi(); - TString name;// ="e information"; - // 0 1 2 - // Pt Phi Eta - int thnBins[thnSize] = { 1000, 200, 500}; - double thnMin [thnSize] = { 0, 0, -1.}; - double thnMax [thnSize] = { 100, 2*Pi, 1.}; - name.Form("%s info", GetName()); - std::auto_ptr ElectronProperties(new THnSparseF(name, name, thnSize, thnBins, thnMin, thnMax)); +const char* AliDxHFEParticleSelectionEl::fgkTrackControlBinNames[]={ + "Pt", + "Phi", + "Eta" +}; - if (ElectronProperties.get()==NULL) { - return -ENOMEM; - } - int axis=0; - ElectronProperties->GetAxis(axis++)->SetTitle("Pt"); - ElectronProperties->GetAxis(axis++)->SetTitle("Phi"); - ElectronProperties->GetAxis(axis++)->SetTitle("Eta"); +const char* AliDxHFEParticleSelectionEl::fgkCutBinNames[]={ + "kRecKineITSTPC", + "kRecPrim", + "kHFEcuts", + "kHFEcutsTOFPID", + "kHFEcutsTPCPID", + "kPID", + "Selected e" +}; - fElectronProperties=ElectronProperties.release(); +int AliDxHFEParticleSelectionEl::Init() +{ + // + // init function + // + int iResult=0; - AddControlObject(fElectronProperties); - - fWhichCut= new TH1F("fWhichCut","effective cut for a rejected particle",6,-0.5,5.5); - AddControlObject(fWhichCut); + // Implicit call to InitControlObjects() before setting up CFM and fCuts + // (if not there) + iResult=AliDxHFEParticleSelection::Init(); + if (iResult<0) return iResult; //--------Initialize correction Framework and Cuts // Consider moving this, either to separate function or @@ -137,6 +138,43 @@ int AliDxHFEParticleSelectionEl::InitControlObjects() fCuts->Initialize(fCFM); return 0; + +} + +THnSparse* AliDxHFEParticleSelectionEl::DefineTHnSparse() const +{ + // + // Defines the THnSparse. For now, only calls CreatControlTHnSparse + const int thnSize = 3; + const double Pi=TMath::Pi(); + TString name; + name.Form("%s info", GetName()); + + // 0 1 2 + // Pt Phi Eta + int thnBins[thnSize] = { 1000, 200, 500}; + double thnMin [thnSize] = { 0, 0, -1.}; + double thnMax [thnSize] = { 100, 2*Pi, 1.}; + + return CreateControlTHnSparse(name,thnSize,thnBins,thnMin,thnMax,fgkTrackControlBinNames); + +} + +int AliDxHFEParticleSelectionEl::InitControlObjects() +{ + /// init control and monitoring objects + AliInfo("Electron THnSparse"); + + fElectronProperties=DefineTHnSparse(); + AddControlObject(fElectronProperties); + + // + fWhichCut= new TH1F("fWhichCut","effective cut for a rejected particle",kNCutLabels,-0.5,kNCutLabels-0.5); + for (int iLabel=0; iLabelGetXaxis()->SetBinLabel(iLabel+1, fgkCutBinNames[iLabel]); + AddControlObject(fWhichCut); + + return AliDxHFEParticleSelection::InitControlObjects(); } int AliDxHFEParticleSelectionEl::HistogramParticleProperties(AliVParticle* p, int selectionCode) @@ -145,17 +183,40 @@ int AliDxHFEParticleSelectionEl::HistogramParticleProperties(AliVParticle* p, in if (!p) return -EINVAL; //if (!fControlObjects) return 0; if(selectionCode==0) return 0; - - AliAODTrack *track=(AliAODTrack*)p; - Double_t eProperties[]={track->Pt(),track->Phi(),track->Eta()}; + // TODO: make array a member, consistent dimensions for THnSparse and array + Double_t eProperties[3]={0.0, 0.0, 0.0}; + DefineParticleProperties(p, eProperties, 3); if(fElectronProperties) fElectronProperties->Fill(eProperties); - + return 0; } +int AliDxHFEParticleSelectionEl::DefineParticleProperties(AliVParticle* p, Double_t* data, int dimension) const +{ + // fill the data array from the particle data + if (!data) return -EINVAL; + AliAODTrack *track=(AliAODTrack*)p; + if (!track) return -ENODATA; + int i=0; + // TODO: this corresponds to the THnSparse dimensions which is available in the same class + // use this consistently + const int requiredDimension=3; + if (dimension!=requiredDimension) { + // TODO: think about filling only the available data and throwing a warning + return -ENOSPC; + } + data[i++]=track->Pt(); + data[i++]=track->Phi(); + data[i++]=track->Eta(); + return i; +} + + int AliDxHFEParticleSelectionEl::IsSelected(AliVParticle* pEl, const AliVEvent*) { /// select El candidates + // TODO: How to handle MC? would be too much duplicate code if copy entire IsSelected. + AliAODTrack *track=(AliAODTrack*)pEl; @@ -163,31 +224,31 @@ int AliDxHFEParticleSelectionEl::IsSelected(AliVParticle* pEl, const AliVEvent*) //Using AliHFECuts: // RecKine: ITSTPC cuts if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)){ - fWhichCut->Fill(0); + fWhichCut->Fill(kRecKineITSTPC); return 0; } // RecPrim if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) { - fWhichCut->Fill(1); + fWhichCut->Fill(kRecPrim); return 0; } // HFEcuts: ITS layers cuts if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) { - fWhichCut->Fill(2); + fWhichCut->Fill(kHFEcutsITS); return 0; } // HFE cuts: TOF PID and mismatch flag if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTOF, track)) { - fWhichCut->Fill(3); + fWhichCut->Fill(kHFEcutsTOF); return 0; } // HFE cuts: TPC PID cleanup if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)){ - fWhichCut->Fill(4); + fWhichCut->Fill(kHFEcutsTPC); return 0; } @@ -206,11 +267,11 @@ int AliDxHFEParticleSelectionEl::IsSelected(AliVParticle* pEl, const AliVEvent*) if(fPID && fPID->IsSelected(&hfetrack)) { AliDebug(3,"Inside FilldPhi, electron is selected"); - + fWhichCut->Fill(kSelected); return 1; } else{ - fWhichCut->Fill(5); + fWhichCut->Fill(kPID); return 0; } } diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.h b/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.h index e3a91b665cb..914baf85d5e 100644 --- a/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.h +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionEl.h @@ -44,16 +44,31 @@ class AliDxHFEParticleSelectionEl : public AliDxHFEParticleSelection { kNCuts }; + enum { + kRecKineITSTPC=0, + kRecPrim, + kHFEcutsITS, + kHFEcutsTOF, + kHFEcutsTPC, + kPID, + kSelected, + kNCutLabels + }; + + ///overloaded from AliDxHFEParticleSelection: Init + virtual int Init(); + /// overloaded from AliDxHFEParticleSelection: init the control objects virtual int InitControlObjects(); + virtual THnSparse* DefineTHnSparse() const; /// overloaded from AliDxHFEParticleSelection: check particle virtual int IsSelected(AliVParticle* p, const AliVEvent*); virtual int HistogramParticleProperties(AliVParticle* p, int selected); - // overloaded from AliDxHFEParticleSelection: specific for electrons - //virtual TObjArray* Select(const AliVEvent *pEvent); + // TODO: function can be renamed to better describe what it's doing + virtual int DefineParticleProperties(AliVParticle* p, Double_t* date, int dimension) const; /// set cuts object: a type cast check is implemented in the method virtual void SetCuts(TObject* /*cuts*/, int /*level*/=0); @@ -75,6 +90,8 @@ class AliDxHFEParticleSelectionEl : public AliDxHFEParticleSelection { AliHFEcuts* fCuts; //! Cuts for HF electrons AliCFManager* fCFM; //! Correction Framework Manager + static const char* fgkTrackControlBinNames[]; //! bin labels for track control histogram + static const char* fgkCutBinNames[]; //! bin labels for cuts histogram ClassDef(AliDxHFEParticleSelectionEl, 2); }; diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.cxx b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.cxx new file mode 100644 index 00000000000..1f843fe93e2 --- /dev/null +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.cxx @@ -0,0 +1,197 @@ +// $Id$ + +//************************************************************************** +//* This file is property of and copyright by the ALICE Project * +//* ALICE Experiment at CERN, All rights reserved. * +//* * +//* Primary Authors: Matthias Richter * +//* Hege Erdal * +//* * +//* 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. * +//************************************************************************** + +/// @file AliDxHFEParticleSelectionMCD0.cxx +/// @author Hege Erdal, Matthias Richter +/// @date 2012-07-19 +/// @brief MC D0 selection for D0-HFE correlation +/// + +#include "AliDxHFEParticleSelectionMCD0.h" +#include "AliAODRecoDecayHF2Prong.h" +#include "AliAODTrack.h" +#include "AliAODMCParticle.h" +#include "AliVParticle.h" +#include +#include +#include + +using namespace std; + +/// ROOT macro for the implementation of ROOT specific class methods +ClassImp(AliDxHFEParticleSelectionMCD0) + +// TODO: can that be a common definition +const char* AliDxHFEParticleSelectionMCD0::fgkTrackControlBinNames[]={ +"Pt", +"Phi", +"Ptbin", +"D0InvMass", +"Eta", +"Statistics D0", +"Mother of D0" +}; + +AliDxHFEParticleSelectionMCD0::AliDxHFEParticleSelectionMCD0(const char* opt) + : AliDxHFEParticleSelectionD0(opt) + , fMCTools() + , fResultMC(0) + , fOriginMother(0) +{ + // constructor + // + // + // + // + + // TODO: argument scan, pass only relevant arguments to tools + fMCTools.~AliDxHFEToolsMC(); + TString toolopt("pdg=421 mc-last"); + new (&fMCTools) AliDxHFEToolsMC(toolopt); +} + +AliDxHFEParticleSelectionMCD0::~AliDxHFEParticleSelectionMCD0() +{ + // destructor +} +THnSparse* AliDxHFEParticleSelectionMCD0::DefineTHnSparse() const +{ + // + // Defines the THnSparse. + // could/should remove Pt and leave Ptbin + + const int thnSize2 = 7; + const double Pi=TMath::Pi(); + TString name; + name.Form("%s info", GetName()); + + // 0 1 2 3 4 5 6 + // Pt Phi Ptbin D0InvMass Eta 'stat D0' mother + int thnBins[thnSize2] = { 1000, 200, 21, 200, 500, 2, 10 }; + double thnMin [thnSize2] = { 0, 0, 0, 1.5648, -1., -0.5, -1.5 }; + double thnMax [thnSize2] = { 100, 2*Pi, 20, 2.1648, 1., 1.5, 8.5 }; + + return CreateControlTHnSparse(name,thnSize2,thnBins,thnMin,thnMax,fgkTrackControlBinNames); + +} + +int AliDxHFEParticleSelectionMCD0::DefineParticleProperties(AliVParticle* p, Double_t* data, int dimension) const +{ + // fill the data array from the particle data + if (!data) return -EINVAL; + AliAODTrack *track=(AliAODTrack*)p; + if (!track) return -ENODATA; + int i=0; + // TODO: this corresponds to the THnSparse dimensions which is available in the same class + // use this consistently + const int requiredDimension=7; + if (dimension!=requiredDimension) { + // TODO: think about filling only the available data and throwing a warning + return -ENOSPC; + } + data[i++]=track->Pt(); + data[i++]=track->Phi(); + data[i++]=AliDxHFEParticleSelectionMCD0::GetPtBin(); + data[i++]=AliDxHFEParticleSelectionMCD0::GetInvMass(); + data[i++]=track->Eta(); + data[i++]=fResultMC; // stat electron (MC electron or not) + data[i++]=fOriginMother; // at the moment not included background. Should expand + + return i; +} + +int AliDxHFEParticleSelectionMCD0::IsSelected(AliVParticle* p, const AliVEvent* pEvent) +{ + /// overloaded from AliDxHFEParticleSelection: check particle + /// H: Have changed function. Now doing particle selection first, then run MC over + /// selected tracks. Could configure it to be configurable, but not sure if it + /// is needed. + /// result from normal track selection is returned, result from MC is stored in + /// THnSparse. + + int iResult=0; + fOriginMother=-1; + + // step 1: + // MC selection + if (fMCTools.MCFirst() && (iResult=CheckMC(p, pEvent))==0) { + // histograming? + return iResult; + } + + // step 2 or 1, depending on sequence: + // normal particle selection + iResult=AliDxHFEParticleSelectionD0::IsSelected(p, pEvent); + if (fMCTools.MCFirst() || iResult==0) return iResult; + + // step 2, only executed if MC check is last + // MC selection - > Should maybe also distinguish between D0 and D0bar + iResult=CheckMC(p, pEvent); + // TODO: why do we need to store the result in a member? + fResultMC=iResult; + + return iResult; +} + +int AliDxHFEParticleSelectionMCD0::CheckMC(AliVParticle* p, const AliVEvent* pEvent) +{ + /// check if MC criteria are fulfilled + // Check both D0 and D0bar (for now only D0) + + if (!p || !pEvent){ + return -EINVAL; + } + int iResult=0; + + if (!fMCTools.IsInitialized() && (iResult=fMCTools.InitMCParticles(pEvent))<0) { + // TODO: message? but has to be filtered in order to avoid message flood + return 0; // no meaningful filtering on mc possible + } + + AliAODRecoDecayHF2Prong *particle = dynamic_cast(p); + + if(!particle) return 0; + + Int_t pdgDgD0toKpi[2]={AliDxHFEToolsMC::kPDGkaon,AliDxHFEToolsMC::kPDGpion}; + + TClonesArray* fMCArray = dynamic_cast(fMCTools.GetMCArray()); + if(!fMCArray) {cout << "no array" << endl; return -1;} + + // find associated MC particle for D0->Kpi + Int_t MClabel=-9999; + + //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx). Checks both D0s and daughters + MClabel=particle->MatchToMC(AliDxHFEToolsMC::kPDGD0,fMCArray,2,pdgDgD0toKpi); + + if(MClabel<0){ + fOriginMother=-1; + return 0; + } + + fMCTools.SetMClabel(MClabel); + fMCTools.FindMotherPDG(p,AliDxHFEToolsMC::kGetOriginMother); + fOriginMother=fMCTools.GetOriginMother(); + + return 1; +} + +void AliDxHFEParticleSelectionMCD0::Clear(const char* option) +{ + /// clear internal memory + fMCTools.Clear(option); +} diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.h b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.h new file mode 100644 index 00000000000..6a3a9cdcfec --- /dev/null +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCD0.h @@ -0,0 +1,62 @@ +//-*- Mode: C++ -*- +// $Id$ + +//* This file is property of and copyright by the ALICE Project * +//* ALICE Experiment at CERN, All rights reserved. * +//* See cxx source for full Copyright notice * + +/// @file AliDxHFEParticleSelectionMCD0.h +/// @author Hege Erdal, Matthias Richter +/// @date 2012-07-19 +/// @brief D0 MC selection for D-HFE correlations +/// + +#ifndef ALIDXHFEPARTICLESELECTIONMCD0_H +#define ALIDXHFEPARTICLESELECTIONMCD0_H + +#include "AliDxHFEParticleSelectionD0.h" +#include "AliDxHFEToolsMC.h" + +/** + * @class AliDxHFEParticleSelectionMCD0 + * Monte Carlo D0 selection for D-HFE correlations, implements the specific + * selection criteria. + */ +class AliDxHFEParticleSelectionMCD0 : public AliDxHFEParticleSelectionD0 { + public: + /// constructor + AliDxHFEParticleSelectionMCD0(const char* opt=""); + /// destructor + virtual ~AliDxHFEParticleSelectionMCD0(); + + /// overloaded from AliDxHFEParticleSelection: check particle + virtual int IsSelected(AliVParticle* p, const AliVEvent *pEvent=NULL); + + virtual THnSparse* DefineTHnSparse() const; + // TODO: function can be renamed to better describe what it's doing + virtual int DefineParticleProperties(AliVParticle* p, Double_t* date, int dimension) const; + + /// check MC criteria + int CheckMC(AliVParticle* p, const AliVEvent* pEvent); + + /// clear internal memory + virtual void Clear(const char* option=""); + + protected: + + private: + /// copy contructor prohibited + AliDxHFEParticleSelectionMCD0(const AliDxHFEParticleSelectionMCD0&); + /// assignment operator prohibited + AliDxHFEParticleSelectionMCD0& operator=(const AliDxHFEParticleSelectionMCD0&); + + AliDxHFEToolsMC fMCTools; // MC selction tools + int fResultMC; // Result on MC check + int fOriginMother; // Holds info on the original mother particle + static const char* fgkTrackControlBinNames[]; //! bin labels for track control histogram + + + ClassDef(AliDxHFEParticleSelectionMCD0, 1); +}; + +#endif diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.cxx b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.cxx new file mode 100644 index 00000000000..2c7418f5e8d --- /dev/null +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.cxx @@ -0,0 +1,235 @@ +// $Id$ + +//************************************************************************** +//* This file is property of and copyright by the ALICE Project * +//* ALICE Experiment at CERN, All rights reserved. * +//* * +//* Primary Authors: Matthias Richter * +//* Hege Erdal * +//* * +//* 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. * +//************************************************************************** + +/// @file AliDxHFEParticleSelectionMCEl.cxx +/// @author Hege Erdal, Matthias Richter +/// @date 2012-07-19 +/// @brief MC El selection for D0-HFE correlation +/// + +#include "AliDxHFEParticleSelectionMCEl.h" +#include "AliVParticle.h" +#include "AliLog.h" +#include "THnSparse.h" +#include "AliAODMCParticle.h" +#include "TH1F.h" +#include "TAxis.h" +#include "AliAODTrack.h" +#include +#include +#include + +using namespace std; +using std::vector; + + +/// ROOT macro for the implementation of ROOT specific class methods +ClassImp(AliDxHFEParticleSelectionMCEl) + +AliDxHFEParticleSelectionMCEl::AliDxHFEParticleSelectionMCEl(const char* opt) + : AliDxHFEParticleSelectionEl(opt) + , fMCTools() + , fElectronProperties(NULL) + , fWhichCut(NULL) + , fOriginMother(0) + , fResultMC(0) +{ + // constructor + // + // + // + // + + fMCTools.~AliDxHFEToolsMC(); + + // TODO: argument scan, build tool options accordingly + // e.g. set mc mode first/last, skip control histograms + TString toolopt("pdg=11 mc-last"); + new (&fMCTools) AliDxHFEToolsMC(toolopt); +} + +//at the moment not used, keep for now (to be used when plotting) +const char* AliDxHFEParticleSelectionMCEl::fgkPDGMotherBinLabels[]={ + "d", + "u", + "s", + "c", + "b", + "gluon", + "gamma", + "#pi^{0}", + "#eta", + "proton", + "others" +}; + +const char* AliDxHFEParticleSelectionMCEl::fgkTrackControlBinNames[]={ + "Pt", + "Phi", + "Eta", + "MC statistics electron" + "Mother", + "PDG of selected electron", +}; + +AliDxHFEParticleSelectionMCEl::~AliDxHFEParticleSelectionMCEl() +{ + // destructor +} + +THnSparse* AliDxHFEParticleSelectionMCEl::DefineTHnSparse() const +{ + // + // Defines the THnSparse. + + const int thnSize = 6; + const double Pi=TMath::Pi(); + TString name; + name.Form("%s info", GetName()); + + // 0 1 2 3 4 5 + // Pt Phi Eta 'stat e' mother "pdg e" + int thnBins[thnSize] = { 1000, 200, 500, 2, 14, 10 }; + double thnMin [thnSize] = { 0, 0, -1., -0.5, -1.5, -0.5 }; + double thnMax [thnSize] = { 100, 2*Pi, 1., 1.5, 12.5, 9.5 }; + + return CreateControlTHnSparse(name,thnSize,thnBins,thnMin,thnMax,fgkTrackControlBinNames); + +} + +int AliDxHFEParticleSelectionMCEl::DefineParticleProperties(AliVParticle* p, Double_t* data, int dimension) const +{ + // fill the data array from the particle data + if (!data) return -EINVAL; + AliAODTrack *track=(AliAODTrack*)p; + if (!track) return -ENODATA; + int i=0; + // TODO: this corresponds to the THnSparse dimensions which is available in the same class + // use this consistently + const int requiredDimension=6; + if (dimension!=requiredDimension) { + // TODO: think about filling only the available data and throwing a warning + return -ENOSPC; + } + data[i++]=track->Pt(); + data[i++]=track->Phi(); + data[i++]=track->Eta(); + data[i++]=fResultMC; // stat electron (MC electron or not) + data[i++]=fOriginMother; // at the moment not included background. Should expand + data[i++]=1; // PDG e - not sure if needed, maybe only needed as separate histo? + + return i; +} + +int AliDxHFEParticleSelectionMCEl::IsSelected(AliVParticle* p, const AliVEvent* pEvent) +{ + /// overloaded from AliDxHFEParticleSelection: check particle + /// H: Have changed function. Now doing particle selection first, then run MC over + /// selected tracks. Could configure it to be configurable, but not sure if it + /// is needed. + /// Result from normal track selection is returned, result from MC is stored in + /// THnSparse. + + int iResult=0; + if (!p || !pEvent){ + return -EINVAL; + } + // step 1: + // optional MC selection before the particle selection + if (fMCTools.MCFirst() && (iResult=CheckMC(p, pEvent))==0) { + // histograming? + return iResult; + } + + // step 2 or 1, depending on sequence: + // normal particle selection + iResult=AliDxHFEParticleSelectionEl::IsSelected(p, pEvent); + if (fMCTools.MCFirst() || iResult==0) return iResult; + + // step 2, only executed if MC check is last + // optional MC selection after the particle selection + iResult=CheckMC(p, pEvent); + // TODO: why do we need to store the result in a member? + fResultMC=iResult; + + return iResult; +} + +int AliDxHFEParticleSelectionMCEl::CheckMC(AliVParticle* p, const AliVEvent* pEvent) +{ + /// check if MC criteria are fulfilled + if (!p || !pEvent){ + return -EINVAL; + } + fOriginMother=-1; + int iResult=0; + + if (!fMCTools.IsInitialized() && (iResult=fMCTools.InitMCParticles(pEvent))<0) { + // TODO: message? but has to be filtered in order to avoid message flood + return 0; // no meaningful filtering on mc possible + } + + if (fMCTools.RejectByPDG(p,false)) { + // rejected by pdg + // would also like to use this info? or not meaningful? + // NEED TO CHANGE THIS! + return 0; + } + + int pdgMother=0; + pdgMother=fMCTools.FindMotherPDG(p,AliDxHFEToolsMC::kGetFirstMother); + + // Particles considered HFE background. can be expanded + // Should be created only once + // TODO: that needs to be configured once to avoid performance penalty + vector motherPDGs; + motherPDGs.push_back(AliDxHFEToolsMC::kPDGpi0); + motherPDGs.push_back(AliDxHFEToolsMC::kPDGeta); + motherPDGs.push_back(AliDxHFEToolsMC::kPDGgamma); + motherPDGs.push_back(AliDxHFEToolsMC::kPDGJpsi); + + if(fMCTools.RejectByPDG(pdgMother,motherPDGs)){ + pdgMother=fMCTools.FindMotherPDG(p,AliDxHFEToolsMC::kGetOriginMother); + fOriginMother=fMCTools.GetOriginMother(); + } + else{ + //Could this be done in a more elegant way? + switch(pdgMother){ + case(AliDxHFEToolsMC::kPDGpi0): fOriginMother=AliDxHFEToolsMC::kNrOrginMother; break; + case(AliDxHFEToolsMC::kPDGeta): fOriginMother=AliDxHFEToolsMC::kNrOrginMother+1; break; + case(AliDxHFEToolsMC::kPDGgamma): fOriginMother=AliDxHFEToolsMC::kNrOrginMother+2;break; + case(AliDxHFEToolsMC::kPDGJpsi): fOriginMother=AliDxHFEToolsMC::kNrOrginMother+3;break; + } + } + + /*if (fMCTools.RejectByMotherPDG(p)) { + // rejected by pdg of original mother + // H: want pdg of origin process to be stored in THnSparse + // Not sure if this is needed... Need to check, using AliDxHFEToolsMC, who + // first mother are, and also what origin is. Use this info here. + return 0; + }*/ + + return 1; +} + +void AliDxHFEParticleSelectionMCEl::Clear(const char* option) +{ + /// clear internal memory + fMCTools.Clear(option); +} diff --git a/PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.h b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.h new file mode 100644 index 00000000000..c9780c76349 --- /dev/null +++ b/PWGHF/correlationHF/AliDxHFEParticleSelectionMCEl.h @@ -0,0 +1,71 @@ +//-*- Mode: C++ -*- +// $Id$ + +//* This file is property of and copyright by the ALICE Project * +//* ALICE Experiment at CERN, All rights reserved. * +//* See cxx source for full Copyright notice * + +/// @file AliDxHFEParticleSelectionMCEl.h +/// @author Hege Erdal, Matthias Richter +/// @date 2012-07-19 +/// @brief El MC selection for D-HFE correlations +/// + +#ifndef ALIDXHFEPARTICLESELECTIONMCEL_H +#define ALIDXHFEPARTICLESELECTIONMCEL_H + +#include "AliDxHFEParticleSelectionEl.h" +#include "AliDxHFEToolsMC.h" + +class THnSparse; +class TH1; + +/** + * @class AliDxHFEParticleSelectionMCEl + * Monte Carlo electron selection for D-HFE correlations, implements the specific + * selection criteria. + */ +class AliDxHFEParticleSelectionMCEl : public AliDxHFEParticleSelectionEl { + public: + /// constructor + AliDxHFEParticleSelectionMCEl(const char* opt=""); + /// destructor + virtual ~AliDxHFEParticleSelectionMCEl(); + + // Setting up control objects: overloaded from AliDxHFEParticleSelection + virtual THnSparse* DefineTHnSparse() const; + + // TODO: function can be renamed to better describe what it's doing + virtual int DefineParticleProperties(AliVParticle* p, Double_t* date, int dimension) const; + + /// overloaded from AliDxHFEParticleSelection: check particle + virtual int IsSelected(AliVParticle* p, const AliVEvent *pEvent=NULL); + + /// check MC criteria + int CheckMC(AliVParticle* p, const AliVEvent* pEvent); + + /// clear internal memory + virtual void Clear(const char* option=""); + + protected: + + private: + /// copy contructor prohibited + AliDxHFEParticleSelectionMCEl(const AliDxHFEParticleSelectionMCEl&); + /// assignment operator prohibited + AliDxHFEParticleSelectionMCEl& operator=(const AliDxHFEParticleSelectionMCEl&); + + /// TODO: check if the label definitions can be used from the ToolsMC + static const char* fgkPDGMotherBinLabels[]; + static const char* fgkTrackControlBinNames[]; //! bin labels for track control histogram + + AliDxHFEToolsMC fMCTools; // MC selection tools + THnSparse* fElectronProperties; // the particle properties of selected particles + TH1* fWhichCut; // effective cut for a rejected particle + int fOriginMother; // Holds the origin motherquark (process) + int fResultMC; // Holds information on check on MC + + ClassDef(AliDxHFEParticleSelectionMCEl, 1); +}; + +#endif diff --git a/PWGHF/correlationHF/AliDxHFEToolsMC.cxx b/PWGHF/correlationHF/AliDxHFEToolsMC.cxx new file mode 100644 index 00000000000..e8ec3a311f8 --- /dev/null +++ b/PWGHF/correlationHF/AliDxHFEToolsMC.cxx @@ -0,0 +1,399 @@ +// $Id$ + +//************************************************************************** +//* This file is property of and copyright by the ALICE Project * +//* ALICE Experiment at CERN, All rights reserved. * +//* * +//* Primary Authors: Matthias Richter * +//* Hege Erdal * +//* * +//* 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. * +//************************************************************************** + +/// @file AliDxHFEToolsMC.cxx +/// @author Hege Erdal, Matthias Richter +/// @date 2012-07-19 +/// @brief Common Tools for MC particle selection +/// + +#include "AliDxHFEToolsMC.h" +#include "AliAODMCParticle.h" +#include "AliVEvent.h" +#include "AliVParticle.h" +#include "AliLog.h" +#include "TObjArray.h" +#include "TString.h" +#include "TH1D.h" +#include +#include +#include + +using namespace std; + +/// ROOT macro for the implementation of ROOT specific class methods +ClassImp(AliDxHFEToolsMC) + +AliDxHFEToolsMC::AliDxHFEToolsMC(const char* option) + : fSequence(kMCLast) + , fMCParticles(NULL) + , fPDGs() + , fMotherPDGs() + , fHistPDG(NULL) + , fHistPDGMother(NULL) + , fOriginMother(kOriginNone) + , fMClabel(-1) +{ + // constructor + // + // + // + // + Init(option); +} + +const char* AliDxHFEToolsMC::fgkPDGBinLabels[]={ + "positron", + "electron", + "#mu+", + "#mu-", + "#pi+", + "#pi-", + "K+", + "K-", + "proton", + "antiproton", + "others" +}; + +const char* AliDxHFEToolsMC::fgkPDGMotherBinLabels[]={ + "d", + "u", + "s", + "c", + "b", + "gluon", + "gamma", + "#pi^{0}", + "#eta", + "proton", + "others" +}; + +const char* AliDxHFEToolsMC::fgkStatisticsBinLabels[]={ + "all", // all MC particles + "found", // selected particles with correct MC + "fake" // selected particles without corresponding MC +}; + +AliDxHFEToolsMC::~AliDxHFEToolsMC() +{ + // destructor + if (fHistPDG) delete fHistPDG; + fHistPDG=NULL; + if (fHistPDGMother) delete fHistPDGMother; + fHistPDGMother=NULL; +} + +int AliDxHFEToolsMC::Init(const char* option) +{ + // initialize according to options + TString strOption(option); + bool bControlHist=true; + std::auto_ptr tokens(strOption.Tokenize(" ")); + if (tokens.get() && tokens->GetEntriesFast()>0) { + for (int itoken=0; itokenGetEntriesFast(); itoken++) { + if (tokens->At(itoken)==NULL) continue; + TString arg=tokens->At(itoken)->GetName(); + const char* key=""; + key="pdg="; + if (arg.BeginsWith(key)) { + arg.ReplaceAll(key, ""); + fPDGs.push_back(arg.Atoi()); + } + key="mother-pdg="; + if (arg.BeginsWith(key)) { + arg.ReplaceAll(key, ""); + fMotherPDGs.push_back(arg.Atoi()); + } + key="control-hist="; + if (arg.BeginsWith(key)) { + arg.ReplaceAll(key, ""); + bControlHist=arg.CompareTo("off")!=0; + } + key="mc-first"; + if (arg.BeginsWith(key)) { + fSequence=kMCFirst; + } + key="mc-last"; + if (arg.BeginsWith(key)) { + fSequence=kMCLast; + } + } + } + + if (bControlHist) { + fHistPDG=CreateControlHistogram("histPDG", + "pdg code of selected particle", + sizeof(fgkPDGBinLabels)/sizeof(const char*), + fgkPDGBinLabels); + fHistPDGMother=CreateControlHistogram("histPDGMother", + "pdg code of first mother of selected particle", + sizeof(fgkPDGMotherBinLabels)/sizeof(fgkPDGMotherBinLabels[0]), + fgkPDGMotherBinLabels); + } + return 0; +} + +int AliDxHFEToolsMC::InitMCParticles(const AliVEvent* pEvent) +{ + + // init MC info from event object + if (!pEvent) return -EINVAL; + + // TODO: choose branch name depending on VEvent type; configurable? + TString branchname(AliAODMCParticle::StdBranchName()); + TObject* o=pEvent->FindListObject(branchname); + if (!o) { + AliWarningClass(Form("can not find MC info '%s' in event of type '%s'", branchname.Data(), pEvent->ClassName())); + return -ENOENT; + } + fMCParticles = dynamic_cast(o); + if (!fMCParticles) { + AliWarningClass(Form("ignoring MC info '%s' of wrong type '%s', expecting TObjArray", branchname.Data(), o->ClassName())); + return -ENODATA; + } + + return 0; +} + +bool AliDxHFEToolsMC::RejectByPDG(int pdg, const vector &list) const +{ + // check if pdg should be rejected, particle is not rejected + // if it is in the list, returns always false if list is empty + if (list.size()==0) return false; + for (vector::const_iterator i=list.begin(); + i!=list.end(); i++) { + if (*i==pdg) return false; + } + return true; +} + +bool AliDxHFEToolsMC::RejectByPDG(AliVParticle* p, bool doStatistics) +{ + // check if pdg should be rejected + // always false if not pdg list is initialized + if (!p) return false; + Int_t MClabel = p->GetLabel(); + if(MClabel<0) return 0; + + int pdgPart=-1; + // TODO: there might be more types of particles to be checked + AliAODMCParticle* aodmcp=0; + aodmcp=dynamic_cast(fMCParticles->At(MClabel)); + if (aodmcp) + pdgPart=TMath::Abs(aodmcp->GetPdgCode()); + if (pdgPart<0) return 0; + + bool bReject=RejectByPDG(pdgPart, fPDGs); + if (doStatistics && fHistPDG) { + // TODO: think about histogramming mode, e.g. histogramming of rejected particles? + if (!bReject) { + fHistPDG->Fill(MapPDGLabel(p->PdgCode())); + } + } + return bReject; +} + +bool AliDxHFEToolsMC::RejectByMotherPDG(AliVParticle* p, bool doStatistics) +{ + // check if pdg should be rejected by mother + // always false if not mother pdg list is initialized + // H: think maybe this is particle specific, and should be moved to PartSelMCEl + if (!p) return false; + int motherpdg=FindPdgOriginMother(p); + // TODO: This should be tweaked. Want to + bool bReject=RejectByPDG(motherpdg, fMotherPDGs); + if (doStatistics && fHistPDGMother) { + // TODO: think about histogramming mode, e.g. histogramming of rejected particles? + if (!bReject) { + fHistPDGMother->Fill(MapPDGMotherLabel(p->PdgCode())); + } + } + return bReject; +} + +int AliDxHFEToolsMC::FindMotherPDG(AliVParticle* p, bool bReturnFirstMother) +{ + // Will find and return pdg of the mother, either first or loop down to the + // initial quark + + // To reset fOriginMother. + fOriginMother=kOriginNone; + + if (!p) return false; + int motherpdg=FindPdgOriginMother(p, bReturnFirstMother); + + return motherpdg; +} + + +TH1* AliDxHFEToolsMC::CreateControlHistogram(const char* name, + const char* title, + int nBins, + const char** binLabels) const +{ + /// create control histogram + std::auto_ptr h(new TH1D(name, title, nBins, -0.5, nBins-0.5)); + if (!h.get()) return NULL; + for (int iLabel=0; iLabelGetXaxis()->SetBinLabel(iLabel+1, binLabels[iLabel]); + } + + return h.release(); +} + +int AliDxHFEToolsMC::MapPDGLabel(int pdg) const +{ + /// mapping of pdg code to enum + switch (pdg) { + case kPDGelectron : return kPDGLabelElectron; + case -kPDGelectron : return kPDGLabelPositron; + case kPDGmuon : return kPDGLabelMuPlus; + case -kPDGmuon : return kPDGLabelMuMinus; + case kPDGpion : return kPDGLabelPiPlus; + case -kPDGpion : return kPDGLabelPiMinus; + case kPDGkaon : return kPDGLabelKPlus; + case -kPDGkaon : return kPDGLabelKMinus; + case kPDGproton : return kPDGLabelProton; + case -kPDGproton : return kPDGLabelAntiproton; + default: + return kPDGLabelOthers; + } +} + +int AliDxHFEToolsMC::MapPDGMotherLabel(int pdg) const +{ + /// mapping of pdg code to enum + switch (pdg) { + case kPDGd : return kPDGMotherLabelD; + case kPDGu : return kPDGMotherLabelU; + case kPDGs : return kPDGMotherLabelS; + case kPDGc : return kPDGMotherLabelC; + case kPDGb : return kPDGMotherLabelB; + case kPDGgluon : return kPDGMotherLabelGluon; + case kPDGgamma : return kPDGMotherLabelGamma; + case kPDGpi0 : return kPDGMotherLabelPi0; + case kPDGeta : return kPDGMotherLabelEta; + case kPDGproton: return kPDGMotherLabelProton; + default: + return kPDGLabelOthers; + } +} + +int AliDxHFEToolsMC::FindPdgOriginMother(AliVParticle* p, bool bReturnFirstMother) +{ + // Return the pgd of original mother particle + // TODO: need also to have specific for D0, electron etc + // for instance to mark when you have gluon, charm or beauty + // among the mothers. Or maybe this will be the same anyway? + // TODO: implement tests on origin, if charm/beauty quark and if + // they came from gluon. use booleans to set this which can be accessed from + // outside? Something like fSequence. + + if (!p) return kPDGnone; + + Int_t imother=-1; + Int_t MClabel=0; + + // Either using MClabel set from outside (needed for Dmesons), or find it + if(fMClabel<0){ + MClabel = p->GetLabel(); + if(MClabel<0){ + return kPDGnone; + } + } + else MClabel=fMClabel; + + // try different classes, unfortunately there is no common base class + AliAODMCParticle* aodmcp=0; + aodmcp=dynamic_cast(fMCParticles->At(MClabel)); + if (!aodmcp) { + return kPDGnone; + } + imother = aodmcp->GetMother(); + + // Absolute or +/- on pgd ??? + Int_t pdg=TMath::Abs(aodmcp->GetPdgCode()); + + if (imother<0){ + // also check this particle + CheckOriginMother(pdg); + return pdg; + } + + if (!fMCParticles->At(imother)) { + AliErrorClass(Form("no mc particle with label %d", imother)); + return pdg; + } + + AliAODMCParticle * mother=dynamic_cast(fMCParticles->At(imother)); + if (!mother) { + AliErrorClass(Form("mc mother particle of wrong class type %s", fMCParticles->At(imother)->ClassName())); + return pdg; + } + + if(bReturnFirstMother){ + return mother->GetPdgCode(); + } + + CheckOriginMother(pdg); + + if (mother->GetPdgCode()==kPDGproton){ + // why? - H: This is the proton, can't get further back + // To be discussed whether to do this a different way + return pdg; + } + + //Reset fMClabel to find pdg of mother if looping on D + fMClabel=-1; + pdg=FindPdgOriginMother(mother); + + return pdg; +} + +void AliDxHFEToolsMC::CheckOriginMother(int pdg) +{ + + // Checking if the particle is a quark or gluon and setting fOriginMother accordingly. + // Right now only check on quark. Need also check on whether it comes from D or B meson? + + switch(pdg){ + case(kPDGc): + fOriginMother = kOriginCharm; break; + case(kPDGb): + fOriginMother = kOriginBeauty; break; + case(kPDGgluon): + if(fOriginMother==kOriginCharm) fOriginMother=kOriginGluonCharm; + else if(fOriginMother==kOriginBeauty) fOriginMother=kOriginGluonBeauty; + else fOriginMother=kOriginGluon; + break; + case(kPDGd): + fOriginMother=kOriginDown; break; + case(kPDGu): + fOriginMother=kOriginUp; break; + case(kPDGs): + fOriginMother=kOriginStrange; break; + } +} + +void AliDxHFEToolsMC::Clear(const char* /*option*/) +{ + // clear internal memory + fMCParticles=NULL; +} diff --git a/PWGHF/correlationHF/AliDxHFEToolsMC.h b/PWGHF/correlationHF/AliDxHFEToolsMC.h new file mode 100644 index 00000000000..5337452af32 --- /dev/null +++ b/PWGHF/correlationHF/AliDxHFEToolsMC.h @@ -0,0 +1,196 @@ +//-*- Mode: C++ -*- +// $Id$ + +//* This file is property of and copyright by the ALICE Project * +//* ALICE Experiment at CERN, All rights reserved. * +//* See cxx source for full Copyright notice * + +/// @file AliDxHFEToolsMC.h +/// @author Hege Erdal, Matthias Richter +/// @date 2012-07-19 +/// @brief Common Tools for MC particle selection +/// + +#ifndef ALIDXHFETOOLSMC_H +#define ALIDXHFETOOLSMC_H + +#include "TObject.h" +#include + +class AliVEvent; +class AliVParticle; +class TH1; + +using std::vector; + +/** + * @class AliDxHFEToolsMC + * Common Tools for MC particle selection. + */ +class AliDxHFEToolsMC { + public: + /// constructor + AliDxHFEToolsMC(const char* options=""); + /// destructor + virtual ~AliDxHFEToolsMC(); + + // different pdgs + enum{ + kPDGnone=0, + kPDGd=1, + kPDGu=2, + kPDGs=3, + kPDGc=4, + kPDGb=5, + kPDGelectron=11, + kPDGmuon=13, + kPDGgluon=21, + kPDGgamma=22, + kPDGpi0=111, + kPDGpion=211, + kPDGeta=221, + kPDGkaon=321, + kPDGD0=421, + kPDGJpsi=443, + kPDGproton=2212 + }; + + enum { + kPDGLabelPositron, + kPDGLabelElectron, + kPDGLabelMuPlus, + kPDGLabelMuMinus, + kPDGLabelPiPlus, + kPDGLabelPiMinus, + kPDGLabelKPlus, + kPDGLabelKMinus, + kPDGLabelProton, + kPDGLabelAntiproton, + kPDGLabelOthers, + kNofPDGLabels + }; + + enum { + kPDGMotherLabelD, + kPDGMotherLabelU, + kPDGMotherLabelS, + kPDGMotherLabelC, + kPDGMotherLabelB, + kPDGMotherLabelGluon, + kPDGMotherLabelGamma, + kPDGMotherLabelPi0, + kPDGMotherLabelEta, + kPDGMotherLabelProton, + kPDGMotherLabelOthers, + kNofPDGMotherLabels + }; + + enum { + kMCFirst = 0, + kMCLast + }; + + enum { + kOriginNone=0, + kOriginDown, + kOriginUp, + kOriginStrange, + kOriginCharm, + kOriginBeauty, + kOriginGluon, + kOriginGluonCharm, + kOriginGluonBeauty, + kNrOrginMother + }; + + enum { + kGetOriginMother=0, + kGetFirstMother + }; + + /// initialize according to options + int Init(const char* /*option*/); + + /// init MC info from event object + int InitMCParticles(const AliVEvent* pEvent); + + /// Returning MC Array + TObjArray* GetMCArray() const {return fMCParticles;} + + /// check state + bool IsInitialized() const {return fMCParticles!=NULL;} + + /// clear internal memory + virtual void Clear(const char* option=""); + + /// flag indicating to do the selection on MC first + bool MCFirst() const {return fSequence==kMCFirst;} + /// flag indicating to do the selection on MC last + bool MCLast() const {return fSequence==kMCLast;} + + /// Return the result on check on initial quark + int GetOriginMother() const {return fOriginMother;} + + /// check if pdg should be rejected + /// always false if pdg list is not initialized + bool RejectByPDG(AliVParticle* p, bool doStatistics=true); + + /// check if pdg should be rejected by mother + /// always false if mother pdg list is not initialized + bool RejectByMotherPDG(AliVParticle* p, bool doStatistics=true); + + /// Finds pdg of first or origin mother and returns value + int FindMotherPDG(AliVParticle* p, bool bReturnFirstMother=false); + + /// step through tree and find the original mother particle + /// TODO: this can possibly be const, however, a member is set inside + /// check whether this is really necessary + int FindPdgOriginMother(AliVParticle* p,bool bReturnFirstMother=false); + + // Compare pdg to quark and gluon + void CheckOriginMother(int pdg); + + // Setting MC label from outside + void SetMClabel(int mclab){fMClabel=mclab;} + + /// TODO: want to have this function to be private again, currently + /// used with an external vector + /// check if pdg should be rejected, particle is not rejected + /// if it is in the list, returns always false if list is empty + bool RejectByPDG(int pdg, const vector &list) const; + + protected: + + private: + /// copy contructor prohibited + AliDxHFEToolsMC(const AliDxHFEToolsMC&); + /// assignment operator prohibited + AliDxHFEToolsMC& operator=(const AliDxHFEToolsMC&); + + /// create control histogram + TH1* CreateControlHistogram(const char* name, + const char* title, + int nBins, + const char** binLabels) const; + + /// mapping of pdg code to enum + int MapPDGLabel(int pdg) const; + /// mapping of pdg code to enum + int MapPDGMotherLabel(int pdg) const; + + static const char* fgkPDGBinLabels[]; + static const char* fgkPDGMotherBinLabels[]; + static const char* fgkStatisticsBinLabels[]; + + int fSequence; // sequence of checks + TObjArray* fMCParticles; //! pointer to external array of MC particles + vector fPDGs; // PDGs to be selected + vector fMotherPDGs; // mother PDGs to be selected + TH1* fHistPDG; // control histogram pdg of selected particle + TH1* fHistPDGMother; // control histogram pdg of selected particle + int fOriginMother; // Holds the origin motherquark (process) + int fMClabel; // MClabel passed from outside (default =-1) + + ClassDef(AliDxHFEToolsMC, 1); +}; +#endif -- 2.43.5