#include "TH1F.h"
#include "TAxis.h"
#include "AliAODTrack.h"
+#include "AliReducedParticle.h"
#include <iostream>
#include <cerrno>
#include <memory>
AliDxHFEParticleSelectionMCEl::AliDxHFEParticleSelectionMCEl(const char* opt)
: AliDxHFEParticleSelectionEl(opt)
, fMCTools()
- , fElectronProperties(NULL)
- , fWhichCut(NULL)
+ , fPDGnotMCElectron(NULL)
+ , fPDGNotHFMother(NULL)
, fOriginMother(0)
, fResultMC(0)
+ , fUseKine(kFALSE)
+ , fMotherPDGs(0)
+ , fUseMCReco(kFALSE)
+ , fSelectionStep(AliDxHFEParticleSelectionEl::kNoCuts)
+ , fStoreCutStepInfo(kFALSE)
{
// constructor
//
fMCTools.~AliDxHFEToolsMC();
+ ParseArguments(opt);
+
+ // This is also checked in AddTasks, but just for safety!
+ if(fUseMCReco && fUseKine) AliFatal("CAN'T SET BOTH usekine AND elmcreco AT THE SAME TIME");
+
// TODO: argument scan, build tool options accordingly
// e.g. set mc mode first/last, skip control histograms
TString toolopt("pdg=11 mc-last");
+ if(fUseKine) toolopt+=" usekine";
new (&fMCTools) AliDxHFEToolsMC(toolopt);
}
"others"
};
+
+const char* AliDxHFEParticleSelectionMCEl::fgkPDGBinLabels[]={
+ "positron",
+ "electron",
+ "#mu+",
+ "#mu-",
+ "#pi+",
+ "#pi-",
+ "K+",
+ "K-",
+ "proton",
+ "antiproton",
+ "others"
+};
+
AliDxHFEParticleSelectionMCEl::~AliDxHFEParticleSelectionMCEl()
{
// destructor
+
+ if(fPDGnotMCElectron){
+ delete fPDGnotMCElectron;
+ fPDGnotMCElectron=NULL;
+ }
+
+}
+
+int AliDxHFEParticleSelectionMCEl::Init()
+{
+ //
+ // init function
+ //
+
+ // Particles considered HFE background. can be expanded
+ fMotherPDGs.push_back(AliDxHFEToolsMC::kPDGpi0);
+ fMotherPDGs.push_back(AliDxHFEToolsMC::kPDGeta);
+ fMotherPDGs.push_back(AliDxHFEToolsMC::kPDGgamma);
+ fMotherPDGs.push_back(AliDxHFEToolsMC::kPDGJpsi);
+
+ int iResult=0;
+ iResult=AliDxHFEParticleSelectionEl::Init();
+ if (iResult<0) return iResult;
+
+ // Histo containing PDG of track which was not MC truth electron
+ // TODO: Add to list in baseclass
+ fPDGnotMCElectron=CreateControlHistogram("fPDGnotMCElectron","PDG of track not MC truth electron",AliDxHFEToolsMC::kNofPDGLabels,fgkPDGBinLabels);
+ fPDGNotHFMother=CreateControlHistogram("fPDGNotHFMother","PDG of mother not HF",5000);
+ AddControlObject(fPDGnotMCElectron);
+ AddControlObject(fPDGNotHFMother);
+ return 0;
}
THnSparse* AliDxHFEParticleSelectionMCEl::DefineTHnSparse()
//
// Defines the THnSparse.
- const int thnSize = 6;
- InitTHnSparseArray(thnSize);
const double Pi=TMath::Pi();
TString name;
+ THnSparse* thn=NULL;
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 };
- const char* thnNames[thnSize]={
- "Pt",
- "Phi",
- "Eta",
- "MC statistics electron",
- "Mother",
- "PDG of selected electron"
- };
-
- return CreateControlTHnSparse(name,thnSize,thnBins,thnMin,thnMax,thnNames);
+ if(fStoreCutStepInfo){
+ const int thnSizeExt =5;
+
+ InitTHnSparseArray(thnSizeExt);
+
+ // TODO: Redo binning of distributions more?
+ // 0 1 2 3
+ // Pt Phi Eta mother
+ int thnBinsExt[thnSizeExt] = { 100, 100, 100, 15, kNCutLabels };
+ double thnMinExt [thnSizeExt] = { 0, 0, -1., -1.5, kRecKineITSTPC };
+ double thnMaxExt [thnSizeExt] = { 10, 2*Pi, 1., 13.5, kSelected };
+ const char* thnNamesExt[thnSizeExt]={
+ "Pt",
+ "Phi",
+ "Eta",
+ "Mother", //bin==-1: Not MC truth electron
+ "Last survived cut step"
+ };
+ thn=(THnSparse*)CreateControlTHnSparse(name,thnSizeExt,thnBinsExt,thnMinExt,thnMaxExt,thnNamesExt);
+ }
+ else{
+
+ const int thnSize =4;
+ InitTHnSparseArray(thnSize);
+
+ // TODO: Redo binning of distributions more?
+ // 0 1 2 3
+ // Pt Phi Eta mother
+ int thnBins[thnSize] = { 100, 100, 100, 15 };
+ double thnMin [thnSize] = { 0, 0, -1., -1.5 };
+ double thnMax [thnSize] = { 10, 2*Pi, 1., 13.5 };
+ const char* thnNames[thnSize]={
+ "Pt",
+ "Phi",
+ "Eta",
+ "Mother", //bin==-1: Not MC truth electron
+ };
+ thn=(THnSparse*)CreateControlTHnSparse(name,thnSize,thnBins,thnMin,thnMax,thnNames);
+
+ }
+ return thn;
}
int AliDxHFEParticleSelectionMCEl::FillParticleProperties(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;
+ if (!p) return -ENODATA;
+ // handle different types of tracks, can be extended
+ AliReducedParticle *trRP=dynamic_cast<AliReducedParticle*>(p);
int i=0;
if (dimension!=GetDimTHnSparse()) {
// 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?
-
+ memset(data, 0, dimension*sizeof(data[0]));
+ data[i++]=p->Pt();
+ data[i++]=p->Phi();
+ data[i++]=p->Eta();
+ if (trRP) data[i]=trRP->GetOriginMother();
+ else data[i]=fOriginMother;
+ i++; // take out of conditionals to be save
+ if (i<dimension) {
+ if(fStoreCutStepInfo) data[i]=GetLastSurvivedCutsStep();
+ i++; // take out of conditionals to be save
+ }
return i;
}
if (!p || !pEvent){
return -EINVAL;
}
+ fOriginMother=-1;
+
+ if(!fUseKine && !fUseMCReco){
// step 1:
// optional MC selection before the particle selection
if (fMCTools.MCFirst() && (iResult=CheckMC(p, pEvent))==0) {
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
// result stored to be filled into THnSparse
// TODO: strictly speaken the particles should be rejected
// if not mc selected, however skip this for the moment, because of
// the logic outside
+ // This line will run always when running directly over the stack or over MC reconstructed tracks
+ // For MC reconstructed tracks, should consider adding more constrictions on the tracks,
+ // maybe do the track selection first (which means could call AliDxHFEParticleSelectionEl::IsSelected()
+ // and don't do PID
+
+ if(fUseMCReco) {// && ( fSelectionStep > AliDxHFEParticleSelectionEl::kNoCuts )){
+ // always check the base class method in mode fUseMCReco
+ iResult=AliDxHFEParticleSelectionEl::IsSelected(p, pEvent);
+ if(iResult == 0) return iResult;
+ }
fResultMC=CheckMC(p, pEvent);
+
+ if(fUseKine || fUseMCReco)
+ return fResultMC;
return iResult;
+
}
int AliDxHFEParticleSelectionMCEl::CheckMC(AliVParticle* p, const AliVEvent* pEvent)
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)) {
+
+ if(fUseKine){
+ // If run on kinematical level, check if particle is electron (only ones who are relevant)
+ // and if pass IsPhysicalPrimary()
+ Int_t test = fMCTools.CheckMCParticle(p);
+ if(test==1) return 0;
+ // Add additional contraints on the electrons?
+ // remove delta e? also remove E<300MeV?
+ // Should also mark dalitz decay and gamma conversion..
+ AliAODMCParticle *mcp=dynamic_cast<AliAODMCParticle*>(p);
+ if(!mcp || !mcp->IsPhysicalPrimary()) return 0;
+
+ }
+
+ int pdgParticle=-1;
+ if (!fUseKine && fMCTools.RejectByPDG(p,false, &pdgParticle)) {
// rejected by pdg
- // would also like to use this info? or not meaningful?
- // NEED TO CHANGE THIS!
+ // TODO: Move this to fMCTools???? Can this be part of the statistics in the MC class?
+ fPDGnotMCElectron->Fill(fMCTools.MapPDGLabel(pdgParticle));
return 0;
}
int pdgMother=0;
+ // Find PDG of first mother
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<int> 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?
+ // Check if first mother is counted as background
+ Bool_t isNotBackground=fMCTools.RejectByPDG(pdgMother,fMotherPDGs);
+
+ if(!isNotBackground){
+ // Set fOriginMother if mother is counted as background
+ // TODO: 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::kPDGJpsi): fOriginMother=AliDxHFEToolsMC::kNrOrginMother+3;break;
}
}
+ else{
+ // If loop over Stack, also checks if First mother is a HF meson
+ Bool_t isHFmeson =fMCTools.TestMotherHFMeson(TMath::Abs(pdgMother));
- /*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;
- }*/
+ if(isHFmeson){
+ // If first mother is a HF meson, loops back to find
+ // original quark + if there was a gluon. Result is
+ // stored in fOriginMother
+ pdgMother=fMCTools.FindMotherPDG(p,AliDxHFEToolsMC::kGetOriginMother);
+ fOriginMother=fMCTools.GetOriginMother();
+ }
+ else{
+ //NotHFmother
+ fPDGNotHFMother->Fill(pdgMother);
+ fOriginMother=AliDxHFEToolsMC::kNrOrginMother+4;
+ }
+
+ }
return 1;
}
/// clear internal memory
fMCTools.Clear(option);
}
+
+AliVParticle *AliDxHFEParticleSelectionMCEl::CreateParticle(AliVParticle* track)
+{
+
+ AliReducedParticle *part = new AliReducedParticle(track->Eta(), track->Phi(), track->Pt(),track->Charge(),fOriginMother);
+
+ return part;
+
+}
+
+int AliDxHFEParticleSelectionMCEl::ParseArguments(const char* arguments)
+{
+ // parse arguments and set internal flags
+ TString strArguments(arguments);
+ auto_ptr<TObjArray> tokens(strArguments.Tokenize(" "));
+ if (!tokens.get()) return 0;
+
+ AliInfo(strArguments);
+ TIter next(tokens.get());
+ TObject* token;
+ while ((token=next())) {
+ TString argument=token->GetName();
+ if (argument.BeginsWith("usekine") ){
+ fUseKine=kTRUE;
+ continue;
+ }
+ if (argument.BeginsWith("elmcreco")){
+ fUseMCReco=kTRUE;
+ if(argument.BeginsWith("elmcreco=")){
+ argument.ReplaceAll("elmcreco=", "");
+ if(argument.CompareTo("alltracks")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kNoCuts;
+ if(argument.CompareTo("afterreckineitstpc")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kRecKineITSTPC;
+ if(argument.CompareTo("afterrecprim")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kRecPrim;
+ if(argument.CompareTo("afterhfeits")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kHFEcutsITS;
+ if(argument.CompareTo("afterhfetof")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kHFEcutsTOF;
+ if(argument.CompareTo("afterhfetpc")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kHFEcutsTPC;
+ if(argument.CompareTo("aftertrackcuts")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kHFEcutsTPC;
+ if(argument.CompareTo("aftertofpid")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kPIDTOF;
+ if(argument.CompareTo("aftertpcpid")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kPIDTPC;
+ if(argument.CompareTo("afterfullpid")==0) fSelectionStep=AliDxHFEParticleSelectionEl::kPIDTOFTPC;
+
+ AliDxHFEParticleSelectionEl::SetFinalCutStep(fSelectionStep);
+ }
+ continue;
+ }
+ if(argument.BeginsWith("storelastcutstep")){
+ AliInfo("Stores the last cut step");
+ fUseMCReco=kTRUE;
+ fStoreCutStepInfo=kTRUE;
+ AliDxHFEParticleSelectionEl::SetStoreLastCutStep(kTRUE);
+ continue;
+ }
+ // forwarding of single argument works, unless key-option pairs separated
+ // by blanks are introduced
+ AliDxHFEParticleSelection::ParseArguments(argument);
+ }
+
+ return 0;
+}