-// $Id$
+// $Id$
//**************************************************************************
//* This file is property of and copyright by the ALICE Project *
/// @brief D0 selection for D0-HFE correlation
///
#include "AliDxHFEParticleSelectionEl.h"
+#include "AliSelectNonHFE.h"
+#include "AliReducedParticle.h"
+#include "AliESDtrackCuts.h"
#include "AliVParticle.h"
#include "AliVEvent.h"
#include "AliPID.h"
#include "AliHFEtools.h"
#include "AliHFEcuts.h"
#include "AliAODTrack.h"
+#include "AliAODEvent.h"
#include "AliAnalysisDataSlot.h"
#include "AliAnalysisDataContainer.h"
#include "AliAnalysisManager.h"
+#include "AliExternalTrackParam.h"
+#include "AliAODVertex.h"
+#include "AliHFEextraCuts.h"
#include "AliCFManager.h"
#include "THnSparse.h"
+#include "AliLog.h"
#include "TH1F.h"
+#include "TH2F.h"
#include "TAxis.h"
+#include "TList.h"
#include "TObjArray.h"
+#include "AliKFParticle.h"
+#include "AliKFVertex.h"
#include <iostream>
#include <cerrno>
#include <memory>
AliDxHFEParticleSelectionEl::AliDxHFEParticleSelectionEl(const char* opt)
: AliDxHFEParticleSelection("electron", opt)
- , fPID(NULL)
+ , fPIDTOFTPC(NULL)
+ , fPIDTOF(NULL)
+ , fPIDTPC(NULL)
, fElectronProperties(NULL)
- , fWhichCut(NULL)
+ , fHistoList(NULL)
+ , fCutPidList(NULL)
+ , fPIDResponse(NULL)
, fCuts(NULL)
+ , fSelNHFE(NULL)
+ , fTrackCuts(NULL)
, fCFM(NULL)
+ , fFinalCutStep(kPIDTOFTPC)
+ , fInvMassLow(0.1)
+ , fUseInvMassCut(kNoInvMass)
+ , fSystem(0)
+ , fTrackNum(0)
+ , fImpactParamCutRadial(3)
+ , fEtaCut(0.8)
+ , fSurvivedCutStep(kNotSelected)
+ , fStoreCutStepInfo(kFALSE)
+
{
// constructor
//
//
//
//
+ ParseArguments(opt);
}
AliDxHFEParticleSelectionEl::~AliDxHFEParticleSelectionEl()
delete fElectronProperties;
fElectronProperties=NULL;
}
+ if(fHistoList){
+ delete fHistoList;
+ fHistoList=NULL;
+ }
if(fCFM){
delete fCFM;
fCFM=NULL;
}
- if(fWhichCut){
- delete fWhichCut;
- fWhichCut=NULL;
+ if(fPIDResponse){
+ delete fPIDResponse;
+ fPIDResponse=NULL;
}
-
// NOTE: external objects fPID and fCuts are not deleted here
- fPID=NULL;
+ fPIDTOFTPC=NULL;
+ fPIDTOF=NULL;
+ fPIDTPC=NULL;
fCuts=NULL;
+ if(fSelNHFE){
+ delete fSelNHFE;
+ fSelNHFE=NULL;
+ }
+ if(fTrackCuts){
+ delete fTrackCuts;
+ fTrackCuts=NULL;
+ }
}
const char* AliDxHFEParticleSelectionEl::fgkCutBinNames[]={
"kHFEcuts",
"kHFEcutsTOFPID",
"kHFEcutsTPCPID",
- "kPID",
+ "kPIDTOF",
+ "kPIDTPC",
+ "kPIDTPCTOF",
+ "kINVMASS",
"Selected e"
};
// init function
//
int iResult=0;
+ // TODO: think about initializing fSelNHFE and fTrackCuts from the macro
+ //Initialization of invariant mass cut function (AliSelectNonHFE)
+ fSelNHFE= new AliSelectNonHFE("IM","IM");
+ // Cuts for associated track in Inv Mass cut
+ fTrackCuts = new AliESDtrackCuts();
+ fTrackCuts->SetAcceptKinkDaughters(kFALSE);
+ fTrackCuts->SetRequireTPCRefit(kTRUE);
+ fTrackCuts->SetEtaRange(-0.9,0.9);
+ fTrackCuts->SetRequireSigmaToVertex(kTRUE);
+ fTrackCuts->SetMaxChi2PerClusterTPC(4.0);
+ fTrackCuts->SetMinNClustersTPC(80);
+ fTrackCuts->SetPtRange(0.3,1e10);
+
+ fSelNHFE->SetTrackCuts(-3, 3, fTrackCuts);
+ fSelNHFE->SetInvariantMassCut(fInvMassLow);
+ fSelNHFE->SetAlgorithm("KF");
+ fSelNHFE->SetAODanalysis(kTRUE);
+
// Implicit call to InitControlObjects() before setting up CFM and fCuts
// (if not there)
// TODO: error handling?
fCuts->Initialize(fCFM);
+ //Setting up TPC PID
+ //Add settings for asymmetric cut on nSigma TPC
+ //TODO: have this completely set up from addtask
+ const int paramSize=4;
+ Double_t params[paramSize];
+ memset(params, 0, sizeof(Double_t)*paramSize);
+ params[0]=-1.;
+ fPIDTPC = new AliHFEpid("hfePidTPC");
+ if(!fPIDTPC->GetNumberOfPIDdetectors()) {
+ fPIDTPC->AddDetector("TPC",1);
+ }
+ fPIDTPC->ConfigureTPCdefaultCut(NULL, params, 3.);
+ fPIDTPC->InitializePID();
+
+ if(fUseInvMassCut>kNoInvMass) AliDebug(2,Form("Setting up with invariant mass cut of %f",fInvMassLow));
+ if(fStoreCutStepInfo) AliDebug(2,"will store cut step info");
+ AliDebug(2,Form("Stopping selection after step: %d", fFinalCutStep));
+
return 0;
}
{
//
// Defines the THnSparse. For now, only calls CreatControlTHnSparse
+ // TODO: Add also invariant mass? here and in correlation, to do cut afterwards..
// here is the only place to change the dimension
const int thnSize = 3;
fElectronProperties=DefineTHnSparse();
AddControlObject(fElectronProperties);
- //
- fWhichCut= new TH1F("fWhichCut","effective cut for a rejected particle",kNCutLabels,-0.5,kNCutLabels-0.5);
- for (int iLabel=0; iLabel<kNCutLabels; iLabel++)
- fWhichCut->GetXaxis()->SetBinLabel(iLabel+1, fgkCutBinNames[iLabel]);
- AddControlObject(fWhichCut);
+ fHistoList=new TList;
+ fHistoList->SetName("HFE Histograms");
+ fHistoList->SetOwner();
+
+ // Histogram storing which cuts have been applied to the tracks
+ // TODO: revise the names of the histograms, starting with 'f'
+ // confuses the reader to think of class members
+ fHistoList->Add(CreateControlHistogram("fWhichCut","effective cut for a rejected particle", kNCutLabels, fgkCutBinNames));
+ fHistoList->Add(CreateControlHistogram("fTPCnClAOD","Nr TPC clusters/track for all tracks", 160, 0., 159.));
+ fHistoList->Add(CreateControlHistogram("fTPCnClSingleTrackCuts","Nr TPC clusters/track After SingleTrackCuts", 160, 0., 159.));
+ fHistoList->Add(CreateControlHistogram("fTPCnClTPCTOFPID","Nr TPC clusters/track After TPCTOF PID", 160, 0., 159.));
+
+ double dEdxBins[6]={1000,0.,10.,200,0.,200.};
+ double nSigBins[6]={1000,0.,10.,200,-10.,10.};
+
+ // dEdx plots, TPC signal vs momentum
+ fHistoList->Add(CreateControl2DHistogram("fdEdx", "dEdx before cuts", dEdxBins,"momentum (GeV/c)","dE/dx in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fdEdxCut", "dEdx after cuts",dEdxBins,"momentum (GeV/c)","dE/dx in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fdEdxPidTOF", "dEdx after TOF pid",dEdxBins,"momentum (GeV/c)","dE/dx in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fdEdxPidTPC", "dEdx after TPC pid",dEdxBins,"momentum (GeV/c)","dE/dx in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fdEdxPid", "dEdx after pid",dEdxBins,"momentum (GeV/c)","dE/dx in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fdEdxIM", "dEdx after Inv Mass",dEdxBins,"momentum (GeV/c)","dE/dx in TPC (a.u.)"));
+
+ // nSigmaTPC vs momentum
+ fHistoList->Add(CreateControl2DHistogram("fnSigTPC", "nSigTPC before cuts",nSigBins,"momentum (GeV/c)","nSigma in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTPCCut", "nSigmaTPC after cuts",nSigBins,"momentum (GeV/c)","nSigma in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTPCPidTPC", "nSigmaTPC after TPC PID", nSigBins,"momentum (GeV/c)","nSigma in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTPCPidTOF", "nSigmaTPC after TOF PID", nSigBins,"momentum (GeV/c)","nSigma in TPC (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTPCPid", "nSigmaTPC after PID", nSigBins,"momentum (GeV/c)","nSigma in TPC (a.u.)"));
+
+ // nSigmaTOF vs momentum
+ fHistoList->Add(CreateControl2DHistogram("fnSigTOF", "nSigmaTOF before cuts",nSigBins,"momentum (GeV/c)","nSigma in TOF (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTOFCut", "nSigmaTOF after cuts",nSigBins,"momentum (GeV/c)","nSigma in TOF (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTOFPidTOF", "nSigmaTOF after TOF PID", nSigBins,"momentum (GeV/c)","nSigma in TOF (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTOFPidTPC", "nSigmaTOF after TPC PID", nSigBins,"momentum (GeV/c)","nSigma in TOF (a.u.)"));
+ fHistoList->Add(CreateControl2DHistogram("fnSigTOFPid", "nSigmaTOF after PID", nSigBins,"momentum (GeV/c)","nSigma in TOF (a.u.)"));
+
+ // Invariant mass LS and ULS without cut
+ fHistoList->Add(CreateControlHistogram("fInvMassLS", "Invariant mass LS", 1000, 0., 0.5));
+ fHistoList->Add(CreateControlHistogram("fInvMassULS", "Invariant mass ULS", 1000, 0., 0.5));
+
+ // Invariant mass LS and ULS without cut, extended x-axis
+ fHistoList->Add(CreateControlHistogram("fInvMass2SelLS", "Invariant mass two selected particles LS", 1000, 0., 5.));
+ fHistoList->Add(CreateControlHistogram("fInvMass2SelULS", "Invariant mass two selected particles ULS", 1000, 0., 5.));
+
+ // Invariant mass LS and ULS with cut
+ // TODO: Remove if remove first try at invariant mass
+ fHistoList->Add(CreateControlHistogram("fInvMass2SelLScut", "Invariant mass two selected particles LS (cut)", 1000, 0., 0.5));
+ fHistoList->Add(CreateControlHistogram("fInvMass2SelULScut", "Invariant mass two selected particles ULS (cut)", 1000, 0., 0.5));
+
+ fSelNHFE->SetHistMass((TH1F*)fHistoList->FindObject("fInvMassULS"));
+ fSelNHFE->SetHistMassBack((TH1F*)fHistoList->FindObject("fInvMassLS"));
+
+ AddControlObject(fHistoList);
return AliDxHFEParticleSelection::InitControlObjects();
}
}
+TObjArray* AliDxHFEParticleSelectionEl::Select(const AliVEvent* pEvent)
+{
+ /// create selection from 'Tracks' member of the event,
+ /// array contains only pointers but does not own the objects
+ /// object array needs to be deleted by caller
+ if (!pEvent) return NULL;
+
+ TList* selectedTracks=new TList;
+ fSelNHFE->SetPIDresponse(fPIDResponse);
+ fPIDTPC->SetPIDResponse(fPIDResponse);
+ TObjArray* finalTracks=new TObjArray;
+ if (!finalTracks) return NULL;
+ finalTracks->SetOwner(kFALSE); // creating new track objects below
+ TObjArray* finalTracksIM=new TObjArray;
+ if (!finalTracksIM) return NULL;
+ finalTracksIM->SetOwner(kFALSE); // creating new track objects below
+
+ int nofTracks=pEvent->GetNumberOfTracks();
+ fTrackNum=0;
+ int selectionCode=0;
+ // int origin
+ for (int itrack=0; itrack<nofTracks; itrack++) {
+ selectionCode=0;
+ AliVParticle* track=pEvent->GetTrack(itrack);
+ selectionCode=IsSelected(track,pEvent);
+ fTrackNum++;
+ // If fUseInvMassCut is set to true, only call HistogramParticleProperties
+ // (here) for those who didn't pass single track cuts and pid.
+ // if fUseInvMassCut is not set, then call the function for all tracks
+ // TODO: Change this when have added invariant mass i thnsparse
+ //CHANGED FOR OLD INV MASS
+ if (selectionCode==0) continue;
+ AliReducedParticle *trackReduced2 =(AliReducedParticle*)CreateParticle(track);
+ finalTracks->Add(trackReduced2);
+ if(fUseInvMassCut!=kInvMassTwoSelected || selectionCode==0)
+ HistogramParticleProperties(trackReduced2, selectionCode);
+ selectedTracks->Add(track);
+ }
+
+ // Calculating invariant mass for electron pairs with same selection criteria
+ // TODO: Could be removed if it is verified that looser cuts is better, but keep
+ // for now
+ if(selectedTracks->GetSize()>0)
+ {
+ Bool_t* selTrackIndex=new Bool_t[selectedTracks->GetSize()];
+ if (selTrackIndex) {
+ InvMassFilter(selectedTracks, selTrackIndex);
+
+ //Only remove electrons if fUseInvMassCut is set
+ if(fUseInvMassCut==kInvMassTwoSelected){
+ for(Int_t k=0; k<selectedTracks->GetSize(); k++)
+ {
+ selectionCode=0; //Reset selectionCode here to be based on invariant mass cut
+ //On the basis that selectedTracks and finalTracks contain same particles
+ AliAODTrack *trackIM=(AliAODTrack*)selectedTracks->At(k);
+ AliReducedParticle *trackReduced=(AliReducedParticle*)finalTracks->At(k);
+ if(selTrackIndex[k]==kTRUE)
+ {
+ ((TH2D*)fHistoList->FindObject("fdEdxIM"))->Fill(trackIM->GetTPCmomentum(), trackIM->GetTPCsignal());
+ finalTracksIM->Add(trackReduced);
+ selectionCode=1;
+ }
+ HistogramParticleProperties(trackReduced, selectionCode);
+ }
+ }
+ delete [] selTrackIndex;
+ }
+
+ }
+
+ HistogramEventProperties(AliDxHFEParticleSelection::kHistoNrTracksPrEvent,finalTracks->GetEntries());
+ if(fUseInvMassCut!=kInvMassTwoSelected)
+ return finalTracks;
+ else
+ return finalTracksIM;
+
+}
+
int AliDxHFEParticleSelectionEl::IsSelected(AliVParticle* pEl, const AliVEvent* pEvent)
{
/// select El candidates
- // TODO: How to handle MC? would be too much duplicate code if copy entire IsSelected.
+ if(!pEvent){
+ AliError("No event information");
+ return 0;
+ }
+ fSurvivedCutStep=kNotSelected;
AliAODTrack *track=(AliAODTrack*)pEl;
fCFM->SetRecEventInfo(pEvent);
+
+ ((TH2D*)fHistoList->FindObject("fdEdx"))->Fill(track->GetTPCmomentum(), track->GetTPCsignal());
+ ((TH2D*)fHistoList->FindObject("fnSigTPC"))->Fill(track->GetTPCmomentum(), fPIDResponse->NumberOfSigmasTPC(track,AliPID::kElectron));
+ ((TH2D*)fHistoList->FindObject("fnSigTOF"))->Fill(track->P(), fPIDResponse->NumberOfSigmasTOF(track,AliPID::kElectron));
+ ((TH1D*)fHistoList->FindObject("fTPCnClAOD"))->Fill(track->GetTPCNcls());
+
+
+ if(fFinalCutStep==kNoCuts){
+ Float_t radial=999;
+ //Float_t z=999;
+ Bool_t acceptTrack=kFALSE;
+
+ //Cut on Eta
+ if(TMath::Abs(track->Eta()) < fEtaCut) {
+ acceptTrack=kTRUE;
+ }
+
+ const Double_t kBeampiperadius=fImpactParamCutRadial;
+ Double_t dcaD[2]={-999.,-999.}, covD[3]={-999.,-999.,-999.};
+
+ AliAODEvent *aodevent = (AliAODEvent*)(pEvent);
+ if(!aodevent) {
+ AliDebug(1, "No aod event available\n");
+ return 0;
+ }
+
+ AliAODVertex *vtxAODSkip = aodevent->GetPrimaryVertex();
+ if(!vtxAODSkip) return 0;
+ AliExternalTrackParam etp; etp.CopyFromVTrack(track);
+ if(etp.PropagateToDCA(vtxAODSkip, aodevent->GetMagneticField(), kBeampiperadius, dcaD, covD)) {
+ radial = dcaD[0];
+ //z = dcaD[1];
+ }
+
+ // Cut on z
+ if(acceptTrack && TMath::Abs(radial) < fImpactParamCutRadial){
+ acceptTrack=kTRUE;
+ }
+
+ if(acceptTrack)
+ fSurvivedCutStep=kNoCuts;
+
+ if(!acceptTrack) return 0;
+
+ // Should only return at this point if you only want to store all tracks
+ AliDebug(2,"Returns after kNoCuts");
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected);
+ return 1;
+ }
//--------track cut selection-----------------------
//Using AliHFECuts:
// RecKine: ITSTPC cuts
if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)){
- fWhichCut->Fill(kRecKineITSTPC);
+ AliDebug(4,"Cut: kStepRecKineITSTPC");
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kRecKineITSTPC);
return 0;
}
-
+ if(fStoreCutStepInfo) fSurvivedCutStep=kRecKineITSTPC;
+ if(fFinalCutStep==kRecKineITSTPC) {AliDebug(2,"Returns after kStepRecKineITSTPC"); ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
+
// RecPrim
if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) {
- fWhichCut->Fill(kRecPrim);
- return 0;
+ AliDebug(4,"Cut: kStepRecPrim");
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kRecPrim);
+ if(!fStoreCutStepInfo) return 0;
+ else return 1; //return 1 because it passed cuts above, but not this (no need to go further)
}
-
+ if(fStoreCutStepInfo) fSurvivedCutStep=kRecPrim;
+ if(fFinalCutStep==kRecPrim) {AliDebug(2,"Returns after kRecPrim"); ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
+
// HFEcuts: ITS layers cuts
if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) {
- fWhichCut->Fill(kHFEcutsITS);
- return 0;
+ AliDebug(4,"Cut: kStepHFEcutsITS");
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kHFEcutsITS);
+ if(!fStoreCutStepInfo) return 0;
+ else return 1; //return 1 because it passed cuts above, but not this (no need to go further)
+
}
+ if(fStoreCutStepInfo) fSurvivedCutStep=kHFEcutsITS;
+ if(fFinalCutStep==kHFEcutsITS) {AliDebug(2,"Returns after kHFEcutsITS "); ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
// HFE cuts: TOF PID and mismatch flag
if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTOF, track)) {
- fWhichCut->Fill(kHFEcutsTOF);
- return 0;
+ AliDebug(4,"Cut: kStepHFEcutsTOF");
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kHFEcutsTOF);
+ if(!fStoreCutStepInfo) return 0;
+ else return 1; //return 1 because it passed cuts above, but not this (no need to go further)
}
+ if(fStoreCutStepInfo) fSurvivedCutStep=kHFEcutsTOF;
+ if(fFinalCutStep==kHFEcutsTOF) {AliDebug(2,"Returns after kHFEcutsTOF"); ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
// HFE cuts: TPC PID cleanup
if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)){
- fWhichCut->Fill(kHFEcutsTPC);
- return 0;
+ AliDebug(4,"Cut: kStepHFEcutsTPC");
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kHFEcutsTPC);
+ if(!fStoreCutStepInfo) return 0;
+ else return 1; //return 1 because it passed cuts above, but not this (no need to go further)
}
-
- // HFEcuts: Nb of tracklets TRD0
- //if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTRD, track)) continue;
+ if(fStoreCutStepInfo) fSurvivedCutStep=kHFEcutsTPC;
+
+ ((TH2D*)fHistoList->FindObject("fdEdxCut"))->Fill(track->GetTPCmomentum(), track->GetTPCsignal());
+ ((TH2D*)fHistoList->FindObject("fnSigTPCCut"))->Fill(track->GetTPCmomentum(), fPIDResponse->NumberOfSigmasTPC(track,AliPID::kElectron));
+ ((TH2D*)fHistoList->FindObject("fnSigTOFCut"))->Fill(track->P(), fPIDResponse->NumberOfSigmasTOF(track,AliPID::kElectron));
+ ((TH1D*)fHistoList->FindObject("fTPCnClSingleTrackCuts"))->Fill(track->GetTPCNcls());
+
+ if(fFinalCutStep==kHFEcutsTPC) {AliDebug(2,"Returns after track cuts"); ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
//--------PID selection-----------------------
hfetrack.SetRecTrack(track);
// TODO: configurable colliding system
- //if(IsPbPb()) hfetrack.SetPbPb();
- hfetrack.SetPP();
+ // TODO: Check problem with PbPb, for now workaround with setting multiplicity to -1
+ if(GetSystem()==1) {
+ hfetrack.SetPbPb();
+ hfetrack.SetCentrality(-1);
+ }
+ else hfetrack.SetPP();
+ // hfetrack.SetMulitplicity(ncontribVtx);
+
+ // TODO: Put this into a while-loop instead, looping over the number of pid objects in the cut-list?
+ // This needs a bit of thinking and finetuning (wrt histogramming)
+ if(fPIDTOF && fPIDTOF->IsSelected(&hfetrack)) {
+ ((TH2D*)fHistoList->FindObject("fdEdxPidTOF"))->Fill(track->GetTPCmomentum(), track->GetTPCsignal());
+ ((TH2D*)fHistoList->FindObject("fnSigTPCPidTOF"))->Fill(track->GetTPCmomentum(), fPIDResponse->NumberOfSigmasTPC(track,AliPID::kElectron));
+ ((TH2D*)fHistoList->FindObject("fnSigTOFPidTOF"))->Fill(track->P(), fPIDResponse->NumberOfSigmasTOF(track,AliPID::kElectron));
+ if(fStoreCutStepInfo){fSurvivedCutStep=kPIDTOF; }
+ if(fFinalCutStep==kPIDTOF) {AliDebug(2,"Returns at PIDTOF");((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
+ }
+ else{
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kPIDTOF);
+ }
+
+ if(fFinalCutStep==kPIDTOF) {AliDebug(2,"Returns at PIDTOF"); return 0;}
+
- if(fPID && fPID->IsSelected(&hfetrack)) {
+ if(fPIDTPC && fPIDTPC->IsSelected(&hfetrack)) {
+ ((TH2D*)fHistoList->FindObject("fdEdxPidTPC"))->Fill(track->GetTPCmomentum(), track->GetTPCsignal());
+ ((TH2D*)fHistoList->FindObject("fnSigTPCPidTPC"))->Fill(track->GetTPCmomentum(), fPIDResponse->NumberOfSigmasTPC(track,AliPID::kElectron));
+ ((TH2D*)fHistoList->FindObject("fnSigTOFPidTPC"))->Fill(track->P(), fPIDResponse->NumberOfSigmasTOF(track,AliPID::kElectron));
+ if(fStoreCutStepInfo){fSurvivedCutStep=kPIDTPC; }
+ if(fFinalCutStep==kPIDTPC) {AliDebug(2,"Returns at PIDTPC"); ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected); return 1;}
+ }
+ else{
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kPIDTPC);
+ //return 0; //if rejected by TOF, will also be rejected by TPCTOF
+ }
+ if(fFinalCutStep==kPIDTPC) {AliDebug(2,"Returns at PIDTTPC"); return 0;}
+
+ //Combined tof & tpc pid
+ if(fPIDTOFTPC && fPIDTOFTPC->IsSelected(&hfetrack)) {
AliDebug(3,"Inside FilldPhi, electron is selected");
- fWhichCut->Fill(kSelected);
- return 1;
+ ((TH2D*)fHistoList->FindObject("fdEdxPid"))->Fill(track->GetTPCmomentum(), track->GetTPCsignal());
+ ((TH2D*)fHistoList->FindObject("fnSigTPCPid"))->Fill(track->GetTPCmomentum(), fPIDResponse->NumberOfSigmasTPC(track,AliPID::kElectron));
+ ((TH2D*)fHistoList->FindObject("fnSigTOFPid"))->Fill(track->P(), fPIDResponse->NumberOfSigmasTOF(track,AliPID::kElectron));
+ ((TH1D*)fHistoList->FindObject("fTPCnClTPCTOFPID"))->Fill(track->GetTPCNcls());
+ // Only do this if
+ if(fStoreCutStepInfo){
+ fSurvivedCutStep=kPIDTOFTPC;
+
+ }
}
else{
- fWhichCut->Fill(kPID);
- return 0;
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kPIDTOFTPC);
+ AliDebug(4,"Cut: kTPCTOFPID");
+ if(!fStoreCutStepInfo) return 0;
+ else { return 1; }//return 1 because it passed cuts above, but not this (no need to go further)
}
+
+ if(fUseInvMassCut==kInvMassSingleSelected)
+ {
+ AliDebug(4,"invmass check");
+ fSelNHFE->FindNonHFE(fTrackNum, track, const_cast<AliVEvent*>(pEvent));
+ if(fSelNHFE->IsLS() || fSelNHFE->IsULS())
+ {
+ //Not selected
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kINVMASS);
+ AliDebug(4,"Cut: Invmass");
+ if(!fStoreCutStepInfo) return 0;
+ }
+ else
+ if(fStoreCutStepInfo) fSurvivedCutStep=kINVMASS;
+
+ }
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected);
+ return 1;
+
}
void AliDxHFEParticleSelectionEl::SetCuts(TObject* cuts, int level)
return;
}
- if (level==kCutPID) {
- fPID=dynamic_cast<AliHFEpid*>(cuts);
- if (!fPID && cuts) {
+ if (level==kCutPIDTOFTPC) {
+ fPIDTOFTPC=dynamic_cast<AliHFEpid*>(cuts);
+ if (!fPIDTOFTPC && cuts) {
+ AliError(Form("cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
+ }
+ return;
+ }
+
+ if (level==kCutPIDTOF) {
+ fPIDTOF=dynamic_cast<AliHFEpid*>(cuts);
+ if (!fPIDTOF && cuts) {
AliError(Form("cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
}
return;
}
+ if(level==kCutList){
+ fCutPidList=dynamic_cast<TList*>(cuts);
+ if (!fCutPidList && cuts) {
+ AliError(Form("cuts object is not of required type TList but %s", cuts->ClassName()));
+ }
+ else{
+ // TODO: Could be done more elegantly, at the moment requires that the cut and pid objects are in
+ // a specific order..
+ TObject *obj=NULL;
+ int iii=0;
+ TIter next(fCutPidList);
+ while((obj = next())){
+ iii++;
+ if(iii==1) {
+ fCuts=dynamic_cast<AliHFEcuts*>(obj);
+ if (!fCuts)
+ AliError(Form("Cut object is not of required type AliHFEcuts but %s", cuts->ClassName()));
+ }
+ if(iii==2){
+ fPIDTOFTPC=dynamic_cast<AliHFEpid*>(obj);
+ if (!fPIDTOFTPC)
+ AliError(Form("(TOFTPC) cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
+ }
+ if(iii==3){
+ fPIDTOF=dynamic_cast<AliHFEpid*>(obj);
+ if (!fPIDTOF)
+ AliError(Form("(TOF) cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
+ }
+ /*if(iii=4){
+ fPIDTPC=dynamic_cast<AliHFEpid*>(obj);
+ if (!fPIDTPC)
+ AliError(Form("(TPC) cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
+ }*/
+ }
+
+ }
+ return;
+ }
}
//________________________________________________________________________
//if(!fCuts->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
return kTRUE;
}
+
+int AliDxHFEParticleSelectionEl::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("elreco=")){
+ argument.ReplaceAll("elreco=", "");
+ int selectionStep=kPIDTOFTPC; //Default
+ if(argument.CompareTo("alltracks")==0) selectionStep=kNoCuts;
+ if(argument.CompareTo("afterreckineitstpc")==0) selectionStep=kRecKineITSTPC;
+ if(argument.CompareTo("afterrecprim")==0) selectionStep=kRecPrim;
+ if(argument.CompareTo("afterhfeits")==0) selectionStep=kHFEcutsITS;
+ if(argument.CompareTo("afterhfetof")==0) selectionStep=kHFEcutsTOF;
+ if(argument.CompareTo("afterhfetpc")==0) selectionStep=kHFEcutsTPC;
+ if(argument.CompareTo("aftertrackcuts")==0) selectionStep=kHFEcutsTPC;
+ if(argument.CompareTo("aftertofpid")==0) selectionStep=kPIDTOF;
+ if(argument.CompareTo("aftertpcpid")==0) selectionStep=kPIDTPC;
+ if(argument.CompareTo("afterfullpid")==0) selectionStep=kPIDTOFTPC;
+ SetFinalCutStep(selectionStep);
+ continue;
+ }
+ if(argument.BeginsWith("twoselectedinvmasscut")){
+ fUseInvMassCut=kInvMassTwoSelected;
+ AliInfo("Using Invariant mass cut for two selected particles");
+ continue;
+ }
+ if(argument.BeginsWith("useinvmasscut")){
+ fUseInvMassCut=kInvMassSingleSelected;
+ AliInfo("Using Invariant mass cut for single selected particle and looser cuts on partner");
+ continue;
+ }
+ if(argument.BeginsWith("invmasscut=")){
+ argument.ReplaceAll("invmasscut=", "");
+ fInvMassLow=argument.Atof();
+ AliInfo(Form("Using invariant mass cut: %f",fInvMassLow));
+ // fUseInvMassCut=kInvMassSingleSelected;
+ continue;
+ }
+ if(argument.BeginsWith("impactparamcut=")){
+ argument.ReplaceAll("impactparamcut=", "");
+ fImpactParamCutRadial=argument.Atof();
+ AliInfo(Form("Using impact parameter cut: %f",fImpactParamCutRadial));
+ continue;
+ }
+ if(argument.BeginsWith("etacut=")){
+ argument.ReplaceAll("etacut=", "");
+ fEtaCut=argument.Atof();
+ AliInfo(Form("Using Eta cut: %f",fEtaCut));
+ continue;
+ }
+ if(argument.BeginsWith("storelastcutstep")){
+ AliInfo("Stores the last cut step");
+ SetStoreLastCutStep(kTRUE);
+ continue;
+ }
+ // forwarding of single argument works, unless key-option pairs separated
+ // by blanks are introduced
+ AliDxHFEParticleSelection::ParseArguments(argument);
+ }
+
+ return 0;
+}
+
+
+void AliDxHFEParticleSelectionEl::InvMassFilter(TList *elList, Bool_t *selIndx)
+{
+
+ //Function for getting invariant mass of electron pairs
+
+ for(int i=0; i<((elList->GetSize())-1); i++)
+ {
+ AliAODTrack *trackAsso=(AliAODTrack*)elList->At(i);
+ for(int j=i+1; j<elList->GetSize(); j++)
+ {
+
+ AliAODTrack *track=(AliAODTrack*)elList->At(j);
+ if(trackAsso && track)
+ {
+
+ Double_t mass=-999., width = -999;
+ Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
+
+ Int_t chargeAsso = trackAsso->Charge();
+ Int_t charge = track->Charge();
+
+ Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
+ if(charge>0) fPDGe1 = -11;
+ if(chargeAsso>0) fPDGe2 = -11;
+
+ if(charge == chargeAsso) fFlagLS = kTRUE;
+ if(charge != chargeAsso) fFlagULS = kTRUE;
+
+ AliKFParticle ge1 = AliKFParticle(*track, fPDGe1);
+ AliKFParticle ge2 = AliKFParticle(*trackAsso, fPDGe2);
+ AliKFParticle recg(ge1, ge2);
+
+ if(recg.GetNDF()<1) continue;
+ Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
+ if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
+
+ recg.GetMass(mass,width);
+
+ if(fFlagLS) {
+ // if(mass<0.5)
+ // ((TH1D*)fHistoList->FindObject("fInvMass2SelLScut"))->Fill(mass);
+ ((TH1D*)fHistoList->FindObject("fInvMass2SelLS"))->Fill(mass);
+ if(mass>fInvMassLow)
+ {
+ ((TH1D*)fHistoList->FindObject("fInvMass2SelLScut"))->Fill(mass);
+ selIndx[i]=kTRUE;
+ selIndx[j]=kTRUE;
+ }
+ }
+ if(fFlagULS) {
+ //if(mass<0.5)
+ // ((TH1D*)fHistoList->FindObject("fInvMass2SelULScut"))->Fill(mass);
+ ((TH1D*)fHistoList->FindObject("fInvMass2SelULS"))->Fill(mass);
+ if(mass>fInvMassLow)
+ {
+ ((TH1D*)fHistoList->FindObject("fInvMass2SelULScut"))->Fill(mass);
+ selIndx[i]=kTRUE;
+ selIndx[j]=kTRUE;
+ }
+
+ }
+
+ }
+ }
+ }
+}