-// $Id$
+// $Id$
//**************************************************************************
//* This file is property of and copyright by the ALICE Project *
#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)
, fElectronProperties(NULL)
, fHistoList(NULL)
, fPIDResponse(NULL)
, fCuts(NULL)
, fCFM(NULL)
+ , fFinalCutStep(kPIDTOFTPC)
+ , fInvMassLow(0.1)
+ , fUseInvMassCut(kFALSE)
{
// constructor
//
//
//
//
+ ParseArguments(opt);
}
AliDxHFEParticleSelectionEl::~AliDxHFEParticleSelectionEl()
fPIDResponse=NULL;
}
// NOTE: external objects fPID and fCuts are not deleted here
- fPID=NULL;
+ fPIDTOFTPC=NULL;
fPIDTOF=NULL;
fCuts=NULL;
}
{
//
// 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);
- double dEdxBins[6]={1000,0.,10.,200,0.,200.};
- double nSigBins[6]={1000,0.,10.,200,-10.,10.};
-
fHistoList=new TList;
fHistoList->SetName("HFE Histograms");
fHistoList->SetOwner();
// Histogram storing which cuts have been applied to the tracks
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("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("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("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("fInvMassLSExt", "Invariant mass LS", 1000, 0., 5.));
+ fHistoList->Add(CreateControlHistogram("fInvMassULSExt", "Invariant mass ULS", 1000, 0., 5.));
+
+ // Invariant mass LS and ULS with cut
+ fHistoList->Add(CreateControlHistogram("fInvMassLScut", "Invariant mass LS (cut)", 1000, 0., 0.5));
+ fHistoList->Add(CreateControlHistogram("fInvMassULScut", "Invariant mass ULS (cut)", 1000, 0., 0.5));
+
AddControlObject(fHistoList);
}
+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;
+
+ TObjArray* finalTracks=new TObjArray;
+ if (!finalTracks) return NULL;
+ finalTracks->SetOwner(kFALSE); // creating new track objects below
+
+ int nofTracks=pEvent->GetNumberOfTracks();
+ int selectionCode=0;
+ for (int itrack=0; itrack<nofTracks; itrack++) {
+ selectionCode=0;
+ AliVParticle* track=pEvent->GetTrack(itrack);
+ selectionCode=IsSelected(track,pEvent);
+ // 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
+ if(!fUseInvMassCut || selectionCode==0)
+ HistogramParticleProperties(track, selectionCode);
+ if (selectionCode==0) continue;
+ selectedTracks->Add(track);
+ //if fUseInvMassCut, fill the array below
+ if(!fUseInvMassCut) finalTracks->Add(CreateParticle(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()];
+ InvMassFilter(selectedTracks, selTrackIndex);
+
+ //Only remove electrons if fUseInvMassCut is set
+ if(fUseInvMassCut){
+ for(Int_t k=0; k<selectedTracks->GetSize(); k++)
+ {
+ selectionCode=0; //Reset selectionCode here to be based on invariant mass cut
+ AliAODTrack *trackIM=(AliAODTrack*)selectedTracks->At(k);
+ if(selTrackIndex[k]==kTRUE)
+ {
+ ((TH2D*)fHistoList->FindObject("fdEdxIM"))->Fill(trackIM->GetTPCmomentum(), trackIM->GetTPCsignal());
+ finalTracks->Add(CreateParticle(trackIM));
+ selectionCode=1;
+ }
+ HistogramParticleProperties(trackIM, selectionCode);
+ }
+ }
+
+ }
+ return finalTracks;
+}
+
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;
+ }
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());
+
//--------track cut selection-----------------------
//Using AliHFECuts:
// RecKine: ITSTPC cuts
((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kHFEcutsTPC);
return 0;
}
-
- // HFEcuts: Nb of tracklets TRD0
- //if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTRD, track)) continue;
+
((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"); return 1;}
//--------PID selection-----------------------
hfetrack.SetRecTrack(track);
// TODO: configurable colliding system
- //if(IsPbPb()) hfetrack.SetPbPb();
- hfetrack.SetPP();
+ if(GetSystem()==1) hfetrack.SetPbPb();
+ else hfetrack.SetPP();
// 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)
((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(fFinalCutStep==kPIDTOF) {AliDebug(2,"Returns at PIDTOF"); return 1;}
}
else{
((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kPIDTOF);
- return 0;
}
//Combined tof & tpc pid
- if(fPID && fPID->IsSelected(&hfetrack)) {
+ if(fPIDTOFTPC && fPIDTOFTPC->IsSelected(&hfetrack)) {
AliDebug(3,"Inside FilldPhi, electron is selected");
((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kSelected);
((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());
+
return 1;
}
else{
- ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kPID);
+ ((TH1D*)fHistoList->FindObject("fWhichCut"))->Fill(kPIDTOFTPC);
return 0;
}
return;
}
- if (level==kCutPID) {
- fPID=dynamic_cast<AliHFEpid*>(cuts);
- if (!fPID && cuts) {
+ if (level==kCutPIDTOFTPC) {
+ fPID
TOFTPC=dynamic_cast<AliHFEpid*>(cuts);
+ if (!fPIDTOFTPC && cuts) {
AliError(Form("cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
}
return;
AliError(Form("Cut object is not of required type AliHFEcuts but %s", cuts->ClassName()));
}
if(iii==2){
- fPID=dynamic_cast<AliHFEpid*>(obj);
- if (!fPID)
+ fPID
TOFTPC=dynamic_cast<AliHFEpid*>(obj);
+ if (!fPIDTOFTPC)
AliError(Form("cuts object is not of required type AliHFEpid but %s", cuts->ClassName()));
}
if(iii==3){
//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("aftertrackcuts")==0) selectionStep=kHFEcutsTPC;
+ if(argument.CompareTo("aftertofpid")==0) selectionStep=kPIDTOF;
+ if(argument.CompareTo("afterfullpid")==0) selectionStep=kPIDTOFTPC;
+ SetFinalCutStep(selectionStep);
+ continue;
+ }
+ if(argument.BeginsWith("useinvmasscut")){
+ fUseInvMassCut=kTRUE;
+ AliInfo("Using Invariant mass cut");
+ continue;
+ }
+ if(argument.BeginsWith("invmasscut=")){
+ argument.ReplaceAll("invmasscut=", "");
+ fInvMassLow=argument.Atoi();
+ fUseInvMassCut=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("fInvMassLS"))->Fill(mass);
+ ((TH1D*)fHistoList->FindObject("fInvMassLSExt"))->Fill(mass);
+ if(mass>fInvMassLow)
+ {
+ ((TH1D*)fHistoList->FindObject("fInvMassLScut"))->Fill(mass);
+ selIndx[i]=kTRUE;
+ selIndx[j]=kTRUE;
+ }
+ }
+ if(fFlagULS) {
+ if(mass<0.5)
+ ((TH1D*)fHistoList->FindObject("fInvMassULS"))->Fill(mass);
+ ((TH1D*)fHistoList->FindObject("fInvMassULSExt"))->Fill(mass);
+ if(mass>fInvMassLow)
+ {
+ ((TH1D*)fHistoList->FindObject("fInvMassULScut"))->Fill(mass);
+ selIndx[i]=kTRUE;
+ selIndx[j]=kTRUE;
+ }
+
+ }
+
+ }
+ }
+ }
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff