///////////////////////////////////////////////////////////////////////////
#include <TChain.h>
+#include <TH1D.h>
#include <AliLog.h>
#include <AliAODHandler.h>
//_________________________________________________________________________________
AliAnalysisTaskDielectronFilter::AliAnalysisTaskDielectronFilter() :
- AliAnalysisTaskSE(),
- fDielectron(0),
- fSelectPhysics(kTRUE)
+ AliAnalysisTaskSE(),
+ fDielectron(0),
+ fSelectPhysics(kTRUE),
+ fTriggerMask(AliVEvent::kMB),
+ fEventStat(0x0)
{
//
// Constructor
//_________________________________________________________________________________
AliAnalysisTaskDielectronFilter::AliAnalysisTaskDielectronFilter(const char *name) :
- AliAnalysisTaskSE(name),
- fDielectron(0),
- fSelectPhysics(kTRUE)
+ AliAnalysisTaskSE(name),
+ fDielectron(0),
+ fSelectPhysics(kTRUE),
+ fTriggerMask(AliVEvent::kMB),
+ fEventStat(0x0)
{
//
// Constructor
//
DefineInput(0,TChain::Class());
DefineOutput(1, THashList::Class());
+ DefineOutput(2, TH1D::Class());
}
//_________________________________________________________________________________
// AddAODBranch("AliDielectronCandidates",fDielectron->GetPairArraysPointer(),"deltaAOD.Dielectron.root");
}
+//_________________________________________________________________________________
+void AliAnalysisTaskDielectronFilter::UserCreateOutputObjects()
+{
+ //
+ // Initilise histograms
+ //
+ if (!fEventStat){
+ fEventStat=new TH1D("hEventStat","Event statistics",5,0,5);
+ fEventStat->GetXaxis()->SetBinLabel(1,"Before Phys. Sel.");
+ fEventStat->GetXaxis()->SetBinLabel(2,"After Phys. Sel.");
+ fEventStat->GetXaxis()->SetBinLabel(3,"After Cand. Sel.");
+ }
+
+ PostData(2,fEventStat);
+}
+
//_________________________________________________________________________________
void AliAnalysisTaskDielectronFilter::UserExec(Option_t *)
{
}
// Was event selected ?
AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
- Bool_t isSelected = kTRUE;
+ UInt_t isSelected = AliVEvent::kAny;
if( fSelectPhysics && inputHandler && inputHandler->GetEventSelection() ) {
isSelected = inputHandler->IsEventSelected();
+ isSelected&=fTriggerMask;
}
-
- if (!isSelected) return;
+
+ //Before physics selection
+ fEventStat->Fill(0.);
+ if (isSelected==0) {
+ PostData(2,fEventStat);
+ return;
+ }
+ //after physics selection
+ fEventStat->Fill(1.);
//bz for AliKF
Double_t bz = InputEvent()->GetMagneticField();
AliAODExtension *extDielectron = dynamic_cast<AliAODHandler*>
((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler())->GetFilteredAOD("AliAOD.Dielectron.root");
extDielectron->SelectEvent();
+ //after candidate selection
+ fEventStat->Fill(2.);
+
//see if dielectron candidate branch exists, if not create is
TTree *t=extDielectron->GetTree();
if (!t->GetBranch("dielectrons")){
}
PostData(1, const_cast<THashList*>(fDielectron->GetHistogramList()));
+ PostData(2,fEventStat);
}
#include "AliAnalysisTaskSE.h"
class AliDielectron;
+class TH1D;
class AliAnalysisTaskDielectronFilter : public AliAnalysisTaskSE {
virtual void UserExec(Option_t *option);
virtual void Init();
+ virtual void UserCreateOutputObjects();
virtual void LocalInit() {Init();}
void UsePhysicsSelection(Bool_t phy=kTRUE) {fSelectPhysics=phy;}
+ void SetTriggerMask(UInt_t mask) {fTriggerMask=mask;}
+ UInt_t GetTriggerMask() const { return fTriggerMask; }
void SetDielectron(AliDielectron * const die) { fDielectron = die; }
AliDielectron *fDielectron; // J/psi framework object
Bool_t fSelectPhysics; // Whether to use physics selection
+ UInt_t fTriggerMask; // Event trigger mask
+
+ TH1D *fEventStat; //! Histogram with event statistics
AliAnalysisTaskDielectronFilter(const AliAnalysisTaskDielectronFilter &c);
AliAnalysisTaskDielectronFilter& operator= (const AliAnalysisTaskDielectronFilter &c);
///////////////////////////////////////////////////////////////////////////
#include <TChain.h>
+#include <TH1D.h>
#include <AliCFContainer.h>
#include <AliVEvent.h>
+#include <AliInputEventHandler.h>
+#include <AliESDInputHandler.h>
+#include <AliAnalysisManager.h>
#include "AliDielectron.h"
#include "AliDielectronHistos.h"
//_________________________________________________________________________________
AliAnalysisTaskDielectronSE::AliAnalysisTaskDielectronSE() :
AliAnalysisTaskSE(),
- fDielectron(0)
+ fDielectron(0),
+ fSelectPhysics(kFALSE),
+ fTriggerMask(AliVEvent::kMB),
+ fEventStat(0x0)
{
//
// Constructor
//_________________________________________________________________________________
AliAnalysisTaskDielectronSE::AliAnalysisTaskDielectronSE(const char *name) :
AliAnalysisTaskSE(name),
- fDielectron(0)
+ fDielectron(0),
+ fSelectPhysics(kFALSE),
+ fTriggerMask(AliVEvent::kMB),
+ fEventStat(0x0)
{
//
// Constructor
DefineInput(0,TChain::Class());
DefineOutput(1, THashList::Class());
DefineOutput(2, AliCFContainer::Class());
+ DefineOutput(3, TH1D::Class());
}
//_________________________________________________________________________________
return;
}
fDielectron->Init();
+ if (fDielectron->GetHistogramList()){
+ PostData(1, const_cast<THashList*>(fDielectron->GetHistogramList()));
+ }
+ if (fDielectron->GetCFManagerPair()){
+ PostData(2, const_cast<AliCFContainer*>(fDielectron->GetCFManagerPair()->GetContainer()));
+ }
+
+ if (!fEventStat){
+ fEventStat=new TH1D("hEventStat","Event statistics",5,0,5);
+ fEventStat->GetXaxis()->SetBinLabel(1,"Before Phys. Sel.");
+ fEventStat->GetXaxis()->SetBinLabel(2,"After Phys. Sel.");
+ }
+
+ PostData(3,fEventStat);
+
}
//_________________________________________________________________________________
if (!fDielectron) return;
+ AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+ AliESDInputHandler *esdHandler=0x0;
+ if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){
+ AliDielectronVarManager::SetESDpid(esdHandler->GetESDpid());
+ } else {
+ //load esd pid bethe bloch parameters depending on the existance of the MC handler
+ // yes: MC parameters
+ // no: data parameters
+ if (!AliDielectronVarManager::GetESDpid()){
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronVarManager::InitESDpid();
+ } else {
+ AliDielectronVarManager::InitESDpid(1);
+ }
+ }
+ }
+ // Was event selected ?
+ AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
+ UInt_t isSelected = AliVEvent::kAny;
+ if( fSelectPhysics && inputHandler && inputHandler->GetEventSelection() ) {
+ isSelected = inputHandler->IsEventSelected();
+ isSelected&=fTriggerMask;
+ }
+
+ //Before physics selection
+ fEventStat->Fill(0.);
+ if (isSelected==0) {
+ PostData(3,fEventStat);
+ return;
+ }
+ //after physics selection
+ fEventStat->Fill(1.);
+
//bz for AliKF
Double_t bz = InputEvent()->GetMagneticField();
AliKFParticle::SetField( bz );
if (fDielectron->GetCFManagerPair()){
PostData(2, const_cast<AliCFContainer*>(fDielectron->GetCFManagerPair()->GetContainer()));
}
+ PostData(3,fEventStat);
}
#include "AliAnalysisTaskSE.h"
class AliDielectron;
+class TH1D;
class AliAnalysisTaskDielectronSE : public AliAnalysisTaskSE {
virtual void UserExec(Option_t *option);
virtual void UserCreateOutputObjects();
+
+ void UsePhysicsSelection(Bool_t phy=kTRUE) {fSelectPhysics=phy;}
+ void SetTriggerMask(UInt_t mask) {fTriggerMask=mask;}
+ UInt_t GetTriggerMask() const { return fTriggerMask; }
void SetDielectron(AliDielectron * const die) { fDielectron = die; }
AliDielectron *fDielectron; // Dielectron framework object
+ Bool_t fSelectPhysics; // Whether to use physics selection
+ UInt_t fTriggerMask; // Event trigger mask
+
+ TH1D *fEventStat; //! Histogram with event statistics
+
AliAnalysisTaskDielectronSE(const AliAnalysisTaskDielectronSE &c);
AliAnalysisTaskDielectronSE& operator= (const AliAnalysisTaskDielectronSE &c);
///////////////////////////////////////////////////////////////////////////
#include <TChain.h>
+#include <TH1D.h>
#include <AliCFContainer.h>
#include <AliInputEventHandler.h>
fListDielectron(),
fListHistos(),
fListCF(),
- fSelectPhysics(kFALSE)
+ fSelectPhysics(kFALSE),
+ fTriggerMask(AliVEvent::kMB),
+ fEventStat(0x0)
{
//
// Constructor
fListDielectron(),
fListHistos(),
fListCF(),
- fSelectPhysics(kFALSE)
+ fSelectPhysics(kFALSE),
+ fTriggerMask(AliVEvent::kMB),
+ fEventStat(0x0)
{
//
// Constructor
DefineInput(0,TChain::Class());
DefineOutput(1, TList::Class());
DefineOutput(2, TList::Class());
+ DefineOutput(3, TH1D::Class());
fListHistos.SetName("Dielectron_Histos_Multi");
fListCF.SetName("Dielectron_CF_Multi");
}
if (die->GetHistogramList()) fListHistos.Add(const_cast<THashList*>(die->GetHistogramList()));
if (die->GetCFManagerPair()) fListCF.Add(const_cast<AliCFContainer*>(die->GetCFManagerPair()->GetContainer()));
}
+
+ if (!fEventStat){
+ fEventStat=new TH1D("hEventStat","Event statistics",5,0,5);
+ fEventStat->GetXaxis()->SetBinLabel(1,"Before Phys. Sel.");
+ fEventStat->GetXaxis()->SetBinLabel(2,"After Phys. Sel.");
+ }
+
+ PostData(1, &fListHistos);
+ PostData(2, &fListCF);
+ PostData(3, fEventStat);
}
//_________________________________________________________________________________
}
// Was event selected ?
AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
- Bool_t isSelected = kTRUE;
+ UInt_t isSelected = AliVEvent::kAny;
if( fSelectPhysics && inputHandler && inputHandler->GetEventSelection() ) {
isSelected = inputHandler->IsEventSelected();
+ isSelected&=fTriggerMask;
}
- if (!isSelected) return;
+ //Before physics selection
+ fEventStat->Fill(0.);
+ if (isSelected==0) {
+ PostData(3,fEventStat);
+ return;
+ }
+ //after physics selection
+ fEventStat->Fill(1.);
//bz for AliKF
Double_t bz = InputEvent()->GetMagneticField();
PostData(1, &fListHistos);
PostData(2, &fListCF);
+ PostData(3,fEventStat);
}
//_________________________________________________________________________________
#include "AliAnalysisTaskSE.h"
class AliDielectron;
+class TH1D;
class AliAnalysisTaskMultiDielectron : public AliAnalysisTaskSE {
virtual void FinishTaskOutput();
void UsePhysicsSelection(Bool_t phy=kTRUE) {fSelectPhysics=phy;}
+ void SetTriggerMask(UInt_t mask) {fTriggerMask=mask;}
+ UInt_t GetTriggerMask() const { return fTriggerMask; }
void AddDielectron(AliDielectron * const die) { fListDielectron.Add(die); }
TList fListCF; //! List with CF Managers
Bool_t fSelectPhysics; // Whether to use physics selection
+ UInt_t fTriggerMask; // Event trigger mask
+
+ TH1D *fEventStat; //! Histogram with event statistics
AliAnalysisTaskMultiDielectron(const AliAnalysisTaskMultiDielectron &c);
AliAnalysisTaskMultiDielectron& operator= (const AliAnalysisTaskMultiDielectron &c);
TNamed("AliDielectron","AliDielectron"),
fEventFilter("EventFilter"),
fTrackFilter("TrackFilter"),
+ fPairPreFilter("PairPreFilter"),
+ fPairPreFilterLegs("PairPreFilterLegs"),
fPairFilter("PairFilter"),
fPdgMother(443),
fPdgLeg1(11),
TNamed(name,title),
fEventFilter("EventFilter"),
fTrackFilter("TrackFilter"),
+ fPairPreFilter("PairPreFilter"),
+ fPairPreFilterLegs("PairPreFilterLegs"),
fPairFilter("PairFilter"),
fPdgMother(443),
fPdgLeg1(11),
//in case we have MC load the MC event and process the MC particles
if (AliDielectronMC::Instance()->HasMC()) {
- AliDielectronMC::Instance()->ConnectMCEvent();
+ if (!AliDielectronMC::Instance()->ConnectMCEvent()){
+ AliError("Could not properly connect the MC event, skipping this event!");
+ return;
+ }
ProcessMC();
}
// Fill Histogram information for tracks and pairs
//
- TString className;
+ TString className,className2;
Double_t values[AliDielectronVarManager::kNMaxValues];
//Fill event information
AliDielectronVarManager::Fill(ev, values);
- fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, values);
+ if (fHistos->GetHistogramList()->FindObject("Event"))
+ fHistos->FillClass("Event", AliDielectronVarManager::kNMaxValues, values);
//Fill track information, separately for the track array candidates
for (Int_t i=0; i<4; ++i){
className.Form("Track_%s",fgkTrackClassNames[i]);
+ if (!fHistos->GetHistogramList()->FindObject(className.Data())) continue;
Int_t ntracks=fTracks[i].GetEntriesFast();
for (Int_t itrack=0; itrack<ntracks; ++itrack){
AliDielectronVarManager::Fill(fTracks[i].UncheckedAt(itrack), values);
}
}
- //Fill Pair information, separately for all pair candidate arrays
+ //Fill Pair information, separately for all pair candidate arrays and the legs
+ TObjArray arrLegs(100);
for (Int_t i=0; i<10; ++i){
className.Form("Pair_%s",fgkPairClassNames[i]);
+ className2.Form("Track_Legs_%s",fgkPairClassNames[i]);
+ Bool_t pairClass=fHistos->GetHistogramList()->FindObject(className.Data())!=0x0;
+ Bool_t legClass=fHistos->GetHistogramList()->FindObject(className2.Data())!=0x0;
+ if (!pairClass&&!legClass) continue;
Int_t ntracks=PairArray(i)->GetEntriesFast();
for (Int_t ipair=0; ipair<ntracks; ++ipair){
- AliDielectronVarManager::Fill(PairArray(i)->UncheckedAt(ipair), values);
- fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
+ AliDielectronPair *pair=static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair));
+
+ //fill pair information
+ if (pairClass){
+ AliDielectronVarManager::Fill(pair, values);
+ fHistos->FillClass(className, AliDielectronVarManager::kNMaxValues, values);
+ }
+
+ //fill leg information, don't fill the information twice
+ if (legClass){
+ AliVParticle *d1=pair->GetFirstDaughter();
+ AliVParticle *d2=pair->GetSecondDaughter();
+ if (!arrLegs.FindObject(d1)){
+ AliDielectronVarManager::Fill(d1, values);
+ fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
+ arrLegs.Add(d1);
+ }
+ if (!arrLegs.FindObject(d2)){
+ AliDielectronVarManager::Fill(d2, values);
+ fHistos->FillClass(className2, AliDielectronVarManager::kNMaxValues, values);
+ arrLegs.Add(d2);
+ }
+ }
}
+ if (legClass) arrLegs.Clear();
}
}
}
}
+//________________________________________________________________
+void AliDielectron::PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2) {
+ //
+ // Prefilter tracks from pairs
+ // Needed for datlitz rejections
+ // remove all tracks from the Single track arrays that pass the cuts in this filter
+ //
+ Int_t pairIndex=GetPairIndex(arr1,arr2);
+
+ Int_t ntrack1=arrTracks1.GetEntriesFast();
+ Int_t ntrack2=arrTracks2.GetEntriesFast();
+
+ AliDielectronPair candidate;
+
+ UInt_t selectedMask=(1<<fPairPreFilter.GetCuts()->GetEntries())-1;
+
+ for (Int_t itrack1=0; itrack1<ntrack1; ++itrack1){
+ Int_t end=ntrack2;
+ if (arr1==arr2) end=itrack1;
+ Bool_t accepted=kFALSE;
+ for (Int_t itrack2=0; itrack2<end; ++itrack2){
+ AliVTrack *track1=static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1));
+ AliVTrack *track2=static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2));
+ if (!track1 || !track2) continue;
+ //create the pair
+ candidate.SetTracks(track1, fPdgLeg1,
+ track2, fPdgLeg2);
+ candidate.SetType(pairIndex);
+ candidate.SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(&candidate,fPdgMother));
+
+ //pair cuts
+ UInt_t cutMask=fPairPreFilter.IsSelected(&candidate);
+
+ //apply cut
+ if (cutMask!=selectedMask) continue;
+ accepted=kTRUE;
+ //remove the tracks from the Track arrays
+ arrTracks2.AddAt(0x0,itrack2);
+ //in case of like sign remove the track from both arrays!
+ if (arr1==arr2) arrTracks1.AddAt(0x0, itrack2);
+ }
+ if ( accepted ) arrTracks1.AddAt(0x0,itrack1);
+ }
+ //compress the track arrays
+ arrTracks1.Compress();
+ arrTracks2.Compress();
+
+ //apply leg cuts after the pre filter
+ if ( fPairPreFilterLegs.GetCuts()->GetEntries()>0 ) {
+ selectedMask=(1<<fPairPreFilterLegs.GetCuts()->GetEntries())-1;
+ //loop over tracks from array 1
+ for (Int_t itrack=0; itrack<arrTracks1.GetEntriesFast();++itrack){
+ //test cuts
+ UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks1.UncheckedAt(itrack));
+ //apply cut
+ if (cutMask!=selectedMask) arrTracks1.AddAt(0x0,itrack);;
+ }
+ arrTracks1.Compress();
+
+ //in case of like sign don't loop over second array
+ if (arr1==arr2) {
+ arrTracks2=arrTracks1;
+ } else {
+
+ //loop over tracks from array 2
+ for (Int_t itrack=0; itrack<arrTracks2.GetEntriesFast();++itrack){
+ //test cuts
+ UInt_t cutMask=fPairPreFilterLegs.IsSelected(arrTracks2.UncheckedAt(itrack));
+ //apply cut
+ if (cutMask!=selectedMask) arrTracks2.AddAt(0x0,itrack);
+ }
+ arrTracks2.Compress();
+
+ }
+ }
+}
+
//________________________________________________________________
void AliDielectron::FillPairArrays(Int_t arr1, Int_t arr2) {
//
// select pairs and fill pair candidate arrays
//
+
+ TObjArray arrTracks1=fTracks[arr1];
+ TObjArray arrTracks2=fTracks[arr2];
+
+ //process pre filter if set
+ if ( fPairPreFilter.GetCuts()->GetEntries()>0 ) PairPreFilter(arr1, arr2, arrTracks1, arrTracks2);
+
Int_t pairIndex=GetPairIndex(arr1,arr2);
- Int_t ntrack1=fTracks[arr1].GetEntriesFast();
- Int_t ntrack2=fTracks[arr2].GetEntriesFast();
+ Int_t ntrack1=arrTracks1.GetEntriesFast();
+ Int_t ntrack2=arrTracks2.GetEntriesFast();
AliDielectronPair *candidate=new AliDielectronPair;
if (arr1==arr2) end=itrack1;
for (Int_t itrack2=0; itrack2<end; ++itrack2){
//create the pair
- candidate->SetTracks(static_cast<AliVTrack*>(fTracks[arr1].UncheckedAt(itrack1)), fPdgLeg1,
- static_cast<AliVTrack*>(fTracks[arr2].UncheckedAt(itrack2)), fPdgLeg2);
+ candidate->SetTracks(static_cast<AliVTrack*>(arrTracks1.UncheckedAt(itrack1)), fPdgLeg1,
+ static_cast<AliVTrack*>(arrTracks2.UncheckedAt(itrack2)), fPdgLeg2);
candidate->SetType(pairIndex);
candidate->SetLabel(AliDielectronMC::Instance()->GetLabelMotherWithPdg(candidate,fPdgMother));
const AliAnalysisFilter& GetTrackFilter() const { return fTrackFilter; }
const AliAnalysisFilter& GetPairFilter() const { return fPairFilter; }
- AliAnalysisFilter& GetEventFilter() { return fEventFilter; }
- AliAnalysisFilter& GetTrackFilter() { return fTrackFilter; }
- AliAnalysisFilter& GetPairFilter() { return fPairFilter; }
-
+ AliAnalysisFilter& GetEventFilter() { return fEventFilter; }
+ AliAnalysisFilter& GetTrackFilter() { return fTrackFilter; }
+ AliAnalysisFilter& GetPairFilter() { return fPairFilter; }
+ AliAnalysisFilter& GetPairPreFilter() { return fPairPreFilter; }
+ AliAnalysisFilter& GetPairPreFilterLegs() { return fPairPreFilterLegs; }
+
void SetMotherPdg( Int_t pdgMother ) { fPdgMother=pdgMother; }
void SetLegPdg(Int_t pdgLeg1, Int_t pdgLeg2) { fPdgLeg1=pdgLeg1; fPdgLeg2=pdgLeg2; }
Int_t GetMotherPdg() const { return fPdgMother; }
private:
- AliAnalysisFilter fEventFilter; // Event cuts
- AliAnalysisFilter fTrackFilter; // leg cuts
- AliAnalysisFilter fPairFilter; // pair cuts
+ AliAnalysisFilter fEventFilter; // Event cuts
+ AliAnalysisFilter fTrackFilter; // leg cuts
+ AliAnalysisFilter fPairPreFilter; // pair prefilter cuts
+ AliAnalysisFilter fPairPreFilterLegs; // Leg filter after the pair prefilter cuts
+ AliAnalysisFilter fPairFilter; // pair cuts
Int_t fPdgMother; // pdg code of mother tracks
Int_t fPdgLeg1; // pdg code leg1
AliDielectronDebugTree *fDebugTree; // Debug tree output
void FillTrackArrays(AliVEvent * const ev, Int_t eventNr=0);
+ void PairPreFilter(Int_t arr1, Int_t arr2, TObjArray &arrTracks1, TObjArray &arrTracks2);
void FillPairArrays(Int_t arr1, Int_t arr2);
Int_t GetPairIndex(Int_t arr1, Int_t arr2) const {return arr1>=arr2?arr1*(arr1+1)/2+arr2:arr2*(arr2+1)/2+arr1;}
AliDielectron(const AliDielectron &c);
AliDielectron &operator=(const AliDielectron &c);
- ClassDef(AliDielectron,3);
+ ClassDef(AliDielectron,4);
};
inline void AliDielectron::InitPairCandidateArrays()
///////////////////////////////////////////////////////////////////////////
#include <TList.h>
+#include <TObjArray.h>
+#include <TVectorD.h>
#include <AliCFContainer.h>
#include <AliAnalysisFilter.h>
TNamed("DielectronCF","DielectronCF"),
fNSteps(0),
fNVars(0),
+ fVarBinLimits(0x0),
fNVarsLeg(0),
fNCuts(0),
fValues(0x0),
//
for (Int_t i=0; i<AliDielectronVarManager::kNMaxValues; ++i){
fVariables[i]=0;
+ fVariablesLeg[i]=0;
}
for (Int_t i=0; i<kNmaxAddSteps; ++i){
fNBinsLeg[i]=0;
fVarLoLimitLeg[i]=0.;
fVarUpLimitLeg[i]=0.;
+ fVarBinLimits[i]=0x0;
fStepMasks[i]=0xFFFFFF;
}
}
TNamed(name, title),
fNSteps(0),
fNVars(0),
+ fVarBinLimits(0x0),
fNVarsLeg(0),
fNCuts(0),
fValues(0x0),
//
for (Int_t i=0; i<AliDielectronVarManager::kNMaxValues; ++i){
fVariables[i]=0;
+ fVariablesLeg[i]=0;
}
for (Int_t i=0; i<kNmaxAddSteps; ++i){
//
if (!leg){
+ if (max<min){
+ Double_t tmp=min;
+ min=max;
+ max=tmp;
+ }
fVariables[fNVars] = (UInt_t)type;
fVarLoLimit[fNVars] = min;
fVarUpLimit[fNVars] = max;
}
}
+//________________________________________________________________
+void AliDielectronCF::AddVariable(AliDielectronVarManager::ValueTypes type, TVectorD * const binLimits)
+{
+ //
+ // Add a variable to the CF configuration
+ //
+ if (!binLimits) return;
+ if (!fVarBinLimits) fVarBinLimits=new TObjArray(kNmaxAddSteps);
+ fVarBinLimits->AddAt(binLimits, fNVars);
+ fVariables[fNVars] = (UInt_t)type;
+ fVarLoLimit[fNVars] = 2.;
+ fVarUpLimit[fNVars] = 1.;
+ fNBins[fNVars] = binLimits->GetNrows()-1;
+ ++fNVars;
+}
+
//________________________________________________________________
void AliDielectronCF::InitialiseContainer(const AliAnalysisFilter& filter)
{
Int_t nBins = fNBins[iVar];
Double_t loLim = fVarLoLimit[iVar];
Double_t upLim = fVarUpLimit[iVar];
- Double_t *binLim = new Double_t[nBins+1];
-
- // set the bin limits
- for(Int_t iBin=0; iBin<=nBins; iBin++) binLim[iBin] = loLim + (upLim-loLim) / nBins*(Double_t)iBin;
-
+ Double_t *binLim = 0x0;
+
+ if (upLim<loLim) {
+ if (!fVarBinLimits) {
+ AliError(Form("Expected variable bin limits for variable %d (%s)",
+ iVar,AliDielectronVarManager::GetValueName(type)));
+ continue;
+ }
+ if (!fVarBinLimits->At(iVar)) {
+ AliError(Form("Expected variable bin limits for variable %d (%s)",
+ iVar,AliDielectronVarManager::GetValueName(type)));
+ continue;
+ }
+ binLim=(static_cast<TVectorD*>(fVarBinLimits->At(iVar)))->GetMatrixArray();
+ } else {
+ binLim=new Double_t[nBins+1];
+ // set the bin limits
+ for(Int_t iBin=0; iBin<=nBins; iBin++) binLim[iBin] = loLim + (upLim-loLim) / nBins*(Double_t)iBin;
+ }
fCfContainer->SetBinLimits(iVar, binLim);
fCfContainer->SetVarTitle(iVar, AliDielectronVarManager::GetValueName(type));
- delete [] binLim;
+ if ( upLim>loLim ) delete [] binLim;
}
// initialize the variables and their bin limits for the Legs
#include <TNamed.h>
+#include <TVectorDfwd.h>
#include "AliDielectronVarManager.h"
class AliAnalysisCuts;
class AliAnalysisFilter;
class AliCFContainer;
class AliDielectronPair;
+class TObjArray;
class AliDielectronCF : public TNamed {
public:
void AddStepMask(UInt_t mask) { fStepMasks[fNStepMasks++]=mask; }
void AddVariable(AliDielectronVarManager::ValueTypes type, Int_t nbins, Double_t min, Double_t max, Bool_t leg=kFALSE);
-
+ void AddVariable(AliDielectronVarManager::ValueTypes type, TVectorD * const binLimits);
+
void InitialiseContainer(const AliAnalysisFilter& filter);
// void Fill(UInt_t mask, const TObject *particle);
Int_t fNBins[kNmaxAddSteps]; // array of numbers ob bins of the vars
Double_t fVarLoLimit[kNmaxAddSteps]; // array of the lower limits of the vars
Double_t fVarUpLimit[kNmaxAddSteps]; // array of the upper limits of the vars
+ TObjArray *fVarBinLimits; // alternative array of bin limits
Int_t fNVarsLeg; // number of variables for the legs
Int_t fNBinsLeg[kNmaxAddSteps]; // array of numbers ob bins of the vars for the legs
Int_t nstep=0;
while ( (o=next()) ){
AliCFContainer *cf=dynamic_cast<AliCFContainer*>(o);
- if (!o) continue;
+ if (!cf) continue;
nstep+=cf->GetNStep();
}
+ if (nstep==0) return;
Int_t nbins[1]={1};
fCfContainer=new AliCFContainer(GetName(), GetTitle(), nstep, 1, nbins);
TFile f(filename);
TList *l=f.GetListOfKeys();
- Int_t entries=l->GetEntries();
- if (entries==0) return;
-
- TKey *k=(TKey*)l->At(0);
- if (!k) return;
- TObject *o=k->ReadObj();
- if (o->IsA()->InheritsFrom(TSeqCollection::Class())){
- TSeqCollection *arr=static_cast<TSeqCollection*>(o);
- SetCFContainers(arr);
- } else if (o->IsA()==AliCFContainer::Class()){
- fCfContainer=static_cast<AliCFContainer*>(o);
- if (fEffGrid) delete fEffGrid;
- fEffGrid=new AliCFEffGrid("eff","eff",*fCfContainer);
+ TIter nextKey(l);
+ TKey *k=0x0;
+ while ( (k=static_cast<TKey*>(nextKey())) ){
+ TObject *o=k->ReadObj();
+ if (o->IsA()->InheritsFrom(TSeqCollection::Class())){
+ TSeqCollection *arr=static_cast<TSeqCollection*>(o);
+ SetCFContainers(arr);
+ } else if (o->IsA()==AliCFContainer::Class()){
+ fCfContainer=static_cast<AliCFContainer*>(o);
+ if (fEffGrid) delete fEffGrid;
+ fEffGrid=new AliCFEffGrid("eff","eff",*fCfContainer);
+ }
}
}
//only in local mode!!!
AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
- if (man && man->GetAnalysisType()!=AliAnalysisManager::kLocalAnalysis) return;
+ if (!man) return;
+ if (man->GetAnalysisType()!=AliAnalysisManager::kLocalAnalysis) return;
//Get File and event information
TObjString fileName;
#include <TCollection.h>
#include <THashList.h>
#include <TString.h>
+#include <TObjString.h>
#include <TObjArray.h>
#include <TFile.h>
#include <TError.h>
//
// Destructor
//
- fHistoList.Delete();
+ fHistoList.Clear();
delete fReservedWords;
}
//
TString drawStr(option);
+ TObjArray *arr=drawStr.Tokenize(";");
+ arr->SetOwner();
+ TIter nextOpt(arr);
+
+ TString drawClasses;
+ TObjString *ostr=0x0;
+
+ TString currentOpt;
+ TString testOpt;
+ while ( (ostr=(TObjString*)nextOpt()) ){
+ currentOpt=ostr->GetString();
+ currentOpt.Remove(TString::kBoth,'\t');
+ currentOpt.Remove(TString::kBoth,' ');
+
+ testOpt="classes=";
+ if ( currentOpt.Contains(testOpt.Data()) ){
+ drawClasses=currentOpt(testOpt.Length(),currentOpt.Length());
+ }
+ }
+
+ delete arr;
drawStr.ToLower();
//options
// Bool_t same=drawOpt.Contains("same"); //FIXME not yet implemented
return;
}
c=gPad->GetCanvas();
+ c->cd();
// c=new TCanvas;
}
TIter nextClass(&fHistoList);
THashList *classTable=0;
- Bool_t first=kTRUE;
+// Bool_t first=kTRUE;
while ( (classTable=(THashList*)nextClass()) ){
+ //test classes option
+ if (!drawClasses.IsNull() && !drawClasses.Contains(classTable->GetName())) continue;
//optimised division
Int_t nPads = classTable->GetEntries();
Int_t nCols = (Int_t)TMath::Ceil( TMath::Sqrt(nPads) );
if (!c) c=new TCanvas(canvasName.Data(),Form("%s: %s",GetName(),classTable->GetName()));
c->Clear();
} else {
- if (first){
- first=kFALSE;
- } else {
- c->Clear();
- }
+// if (first){
+// first=kFALSE;
+// if (nPads>1) gVirtualPS->NewPage();
+// } else {
+ if (nPads>1) c->Clear();
+// }
}
if (nCols>1||nRows>1) c->Divide(nCols,nRows);
histOpt.ReplaceAll("logz","");
h->Draw(drawOpt.Data());
}
- if (gVirtualPS) c->Update();
+ if (gVirtualPS) {
+ c->Update();
+ }
+
}
// if (gVirtualPS) delete c;
}
// set histogram classes and histograms to this instance. It will take onwnership!
//
ResetHistogramList();
- SetName(list.GetName());
+ TString name(GetName());
+ if (name == "AliDielectronHistos") SetName(list.GetName());
TIter next(&list);
TObject *o;
while ( (o=next()) ){
AliDielectronPID::AliDielectronPID() :
AliAnalysisCuts(),
fNcuts(0),
- fRequirePIDbit(kTRUE),
fESDpid(0x0)
{
//
fExclude[icut]=kFALSE;
fFunUpperCut[icut]=0x0;
fFunLowerCut[icut]=0x0;
+ fRequirePIDbit[icut]=0;
}
}
AliDielectronPID::AliDielectronPID(const char* name, const char* title) :
AliAnalysisCuts(name, title),
fNcuts(0),
- fRequirePIDbit(kTRUE),
fESDpid(0x0)
{
//
fExclude[icut]=kFALSE;
fFunUpperCut[icut]=0x0;
fFunLowerCut[icut]=0x0;
+ fRequirePIDbit[icut]=0;
}
}
//______________________________________________
void AliDielectronPID::AddCut(DetType det, AliPID::EParticleType type, Double_t nSigmaLow, Double_t nSigmaUp/*=-99999.*/,
- Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/)
+ Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/,
+ UInt_t pidBitType/*=AliDielectronPID::kRequire*/)
{
//
// Add a pid nsigma cut
fPmin[fNcuts]=pMin;
fPmax[fNcuts]=pMax;
fExclude[fNcuts]=exclude;
+ fRequirePIDbit[fNcuts]=pidBitType;
++fNcuts;
}
//______________________________________________
void AliDielectronPID::AddCut(DetType det, AliPID::EParticleType type, Double_t nSigmaLow, TF1 * const funUp,
- Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/)
+ Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/,
+ UInt_t pidBitType/*=AliDielectronPID::kRequire*/)
{
//
// cut using a TF1 as upper cut
return;
}
fFunUpperCut[fNcuts]=funUp;
- AddCut(det,type,nSigmaLow,0.,pMin,pMax,exclude);
+ AddCut(det,type,nSigmaLow,0.,pMin,pMax,exclude,pidBitType);
}
//______________________________________________
void AliDielectronPID::AddCut(DetType det, AliPID::EParticleType type, TF1 * const funLow, Double_t nSigmaUp,
- Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/)
+ Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/,
+ UInt_t pidBitType/*=AliDielectronPID::kRequire*/)
{
//
// cut using a TF1 as lower cut
return;
}
fFunLowerCut[fNcuts]=funLow;
- AddCut(det,type,0.,nSigmaUp,pMin,pMax,exclude);
+ AddCut(det,type,0.,nSigmaUp,pMin,pMax,exclude,pidBitType);
}
//______________________________________________
void AliDielectronPID::AddCut(DetType det, AliPID::EParticleType type, TF1 * const funLow, TF1 * const funUp,
- Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/)
+ Double_t pMin/*=0*/, Double_t pMax/*=0*/, Bool_t exclude/*=kFALSE*/,
+ UInt_t pidBitType/*=AliDielectronPID::kRequire*/)
{
//
// cut using a TF1 as lower and upper cut
}
fFunUpperCut[fNcuts]=funUp;
fFunLowerCut[fNcuts]=funLow;
- AddCut(det,type,0.,0.,pMin,pMax,exclude);
+ AddCut(det,type,0.,0.,pMin,pMax,exclude,pidBitType);
}
//______________________________________________
if (part->IsA()==AliESDtrack::Class()){
// ESD case in case the PID bit is not set, don't use this track!
AliESDtrack *track=static_cast<AliESDtrack*>(part);
- if (fRequirePIDbit&&!(track->GetStatus()&AliESDtrack::kITSpid)) return kFALSE;
+ if (fRequirePIDbit[icut]==AliDielectronPID::kRequire&&!(track->GetStatus()&AliESDtrack::kITSpid)) return kFALSE;
+ if (fRequirePIDbit[icut]==AliDielectronPID::kIfAvailable&&!(track->GetStatus()&AliESDtrack::kITSpid)) return kTRUE;
numberOfSigmas=fESDpid->NumberOfSigmasITS(track, fPartType[icut]);
}else{
if (part->IsA()==AliESDtrack::Class()){
// ESD case in case the PID bit is not set, don't use this track!
AliESDtrack *track=static_cast<AliESDtrack*>(part);
- if (fRequirePIDbit&&!(track->GetStatus()&AliESDtrack::kTPCpid)) return kFALSE;
+ if (fRequirePIDbit[icut]==AliDielectronPID::kRequire&&!(track->GetStatus()&AliESDtrack::kTPCpid)) return kFALSE;
+ if (fRequirePIDbit[icut]==AliDielectronPID::kIfAvailable&&!(track->GetStatus()&AliESDtrack::kTPCpid)) return kTRUE;
numberOfSigmas=fESDpid->NumberOfSigmasTPC(track, fPartType[icut]);
}else{
if (part->IsA()==AliESDtrack::Class()){
// ESD case in case the PID bit is not set, don't use this track!
AliESDtrack *track=static_cast<AliESDtrack*>(part);
- if (fRequirePIDbit&&!(track->GetStatus()&AliESDtrack::kTOFpid)) return kFALSE;
+ if (fRequirePIDbit[icut]==AliDielectronPID::kRequire&&!(track->GetStatus()&AliESDtrack::kTOFpid)) return kFALSE;
+ if (fRequirePIDbit[icut]==AliDielectronPID::kIfAvailable&&!(track->GetStatus()&AliESDtrack::kTOFpid)) return kTRUE;
numberOfSigmas=fESDpid->NumberOfSigmasTOF(track, fPartType[icut], fESDpid->GetTOFResponse().GetTimeZero());
}else{
AddCut(kTOF,AliPID::kKaon,-3.,3.,0.,1.,kTRUE);
AddCut(kTOF,AliPID::kProton,-6.,6.,0.,1.,kTRUE);
AddCut(kTOF,AliPID::kProton,-3.,3.,1.,2.,kTRUE);
- fRequirePIDbit=kFALSE;
} else if (def==5) {
AddCut(kTPC,AliPID::kElectron,-0.5,3);
AddCut(kTOF,AliPID::kElectron,-3,3,0,1.5);
AddCut(kTPC,AliPID::kElectron,lowerCut,3.);
AddCut(kTOF,AliPID::kElectron,-3,3,0,1.5);
} else if (def==7) {
- // lower cut TPC: parametrisation by HFE
- // upper cut TPC: 3 sigma
+ // wide TPC cut
// TOF ele band 3sigma 0<p<1.5GeV
AddCut(kTPC,AliPID::kElectron,10.);
AddCut(kTOF,AliPID::kElectron,-3,3,0,1.5);
+ } else if (def==8) {
+ // TOF 5 sigma inclusion if TOFpid available
+ // this should reduce K,p,Pi to a large extent
+ AddCut(kTOF,AliPID::kElectron,-5,5,0,200,kFALSE,AliDielectronPID::kIfAvailable);
+ } else if (def==9) {
+ // lower cut TPC: parametrisation by HFE
+ // upper cut TPC: 3 sigma
+ // TOF 5 sigma inclusion if TOFpid available
+ // this should reduce K,p,Pi to a large extent
+ TF1 *lowerCut=new TF1("lowerCut", "[0] * TMath::Exp([1]*x)", 0, 100);
+ lowerCut->SetParameters(-2.7,-0.4357);
+ AddCut(kTPC,AliPID::kElectron,lowerCut,3.);
+ AddCut(kTOF,AliPID::kElectron,-5,5,0,200,kFALSE,AliDielectronPID::kIfAvailable);
+ } else if (def==10) {
+ AddCut(kTOF,AliPID::kElectron,-5,5,0,200,kFALSE,AliDielectronPID::kIfAvailable);
+ AddCut(kTPC,AliPID::kElectron,3.);
+ AddCut(kTPC,AliPID::kPion,-3.,3.,0.,0.,kTRUE);
+ AddCut(kTPC,AliPID::kProton,-3.,3.,0.,0.,kTRUE);
+
}
}
class AliDielectronPID : public AliAnalysisCuts {
public:
enum DetType {kITS, kTPC, kTRD, kTOF};
+ enum PIDbitTupe {kIgnore=0, kRequire, kIfAvailable};
AliDielectronPID();
AliDielectronPID(const char*name, const char* title);
virtual ~AliDielectronPID();
void AddCut(DetType det, AliPID::EParticleType type, Double_t nSigmaLow, Double_t nSigmaUp=-99999.,
- Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE);
+ Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE, UInt_t pidBitType=AliDielectronPID::kRequire);
void AddCut(DetType det, AliPID::EParticleType type, Double_t nSigmaLow, TF1 * const funUp,
- Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE);
+ Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE, UInt_t pidBitType=AliDielectronPID::kRequire);
void AddCut(DetType det, AliPID::EParticleType type, TF1 * const funLow, Double_t nSigmaUp,
- Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE);
+ Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE, UInt_t pidBitType=AliDielectronPID::kRequire);
void AddCut(DetType det, AliPID::EParticleType type, TF1 * const funLow, TF1 * const funUp,
- Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE);
+ Double_t pMin=0, Double_t pMax=0, Bool_t exclude=kFALSE, UInt_t pidBitType=AliDielectronPID::kRequire);
void SetDefaults(Int_t def);
- void SetRequirePIDbit(Bool_t require) {fRequirePIDbit=require;}
//
//Analysis cuts interface
//const
TF1 *fFunUpperCut[kNmaxPID];//use function as upper cut
TF1 *fFunLowerCut[kNmaxPID];//use function as lower cut
UChar_t fNcuts; //number of cuts
-
- Bool_t fRequirePIDbit; //If to require the PID bits to be set.
+ UChar_t fRequirePIDbit[kNmaxPID]; //How to make use of the pid bit (see)
AliESDpid *fESDpid; //! esd pid object
AliDielectronPID(const AliDielectronPID &c);
AliDielectronPID &operator=(const AliDielectronPID &c);
- ClassDef(AliDielectronPID,1) // Dielectron PID
+ ClassDef(AliDielectronPID,2) // Dielectron PID
};
}
}
+//______________________________________________
+void AliDielectronPair::GetThetaPhiCM(Double_t &thetaHE, Double_t &phiHE, Double_t &thetaCS, Double_t &phiCS) const
+{
+ //
+ // Calculate theta and phi in helicity and Collins-Soper coordinate frame
+ //
+ const Double_t kBeamEnergy = 3500.;
+ Double_t pxyz1[3]={0,0,0};
+ Double_t pxyz2[3]={0,0,0};
+ Double_t eleMass=AliPID::ParticleMass(AliPID::kElectron);
+ Double_t proMass=AliPID::ParticleMass(AliPID::kProton);
+
+ AliVParticle *d1 = static_cast<AliVParticle*>(fRefD1.GetObject());
+ AliVParticle *d2 = static_cast<AliVParticle*>(fRefD2.GetObject());
+
+ d1->PxPyPz(pxyz1);
+ d2->PxPyPz(pxyz2);
+
+ TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+ TLorentzVector targMom(0.,0., kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+
+ // first & second daughter 4-mom
+ TLorentzVector p1Mom(pxyz1[0],pxyz1[1],pxyz1[2],
+ TMath::Sqrt(pxyz1[0]*pxyz1[0]+pxyz1[1]*pxyz1[1]+pxyz1[2]*pxyz1[2]+eleMass*eleMass));
+ TLorentzVector p2Mom(pxyz2[0],pxyz2[1],pxyz2[2],
+ TMath::Sqrt(pxyz2[0]*pxyz2[0]+pxyz2[1]*pxyz2[1]+pxyz2[2]*pxyz2[2]+eleMass*eleMass));
+ // J/Psi 4-momentum vector
+ TLorentzVector motherMom=p1Mom+p2Mom;
+
+ // boost all the 4-mom vectors to the mother rest frame
+ TVector3 beta = (-1.0/motherMom.E())*motherMom.Vect();
+ p1Mom.Boost(beta);
+ p2Mom.Boost(beta);
+ projMom.Boost(beta);
+ targMom.Boost(beta);
+
+ // x,y,z axes
+ TVector3 zAxisHE = (motherMom.Vect()).Unit();
+ TVector3 zAxisCS = ((projMom.Vect()).Unit()-(targMom.Vect()).Unit()).Unit();
+ TVector3 yAxis = ((projMom.Vect()).Cross(targMom.Vect())).Unit();
+ TVector3 xAxisHE = (yAxis.Cross(zAxisHE)).Unit();
+ TVector3 xAxisCS = (yAxis.Cross(zAxisCS)).Unit();
+
+ // fill theta and phi
+ if(d1->Charge()>0){
+ thetaHE = zAxisHE.Dot((p1Mom.Vect()).Unit());
+ thetaCS = zAxisCS.Dot((p1Mom.Vect()).Unit());
+ phiHE = TMath::ATan2((p1Mom.Vect()).Dot(yAxis), (p1Mom.Vect()).Dot(xAxisHE));
+ phiCS = TMath::ATan2((p1Mom.Vect()).Dot(yAxis), (p1Mom.Vect()).Dot(xAxisCS));
+ } else {
+ thetaHE = zAxisHE.Dot((p2Mom.Vect()).Unit());
+ thetaCS = zAxisCS.Dot((p2Mom.Vect()).Unit());
+ phiHE = TMath::ATan2((p2Mom.Vect()).Dot(yAxis), (p2Mom.Vect()).Dot(xAxisHE));
+ phiCS = TMath::ATan2((p2Mom.Vect()).Dot(yAxis), (p2Mom.Vect()).Dot(xAxisCS));
+ }
+}
+
//______________________________________________
Double_t AliDielectronPair::ThetaPhiCM(const AliVParticle* d1, const AliVParticle* d2,
- const Bool_t isHE, const Bool_t isTheta) {
- // The function calculates theta and phi in the mother rest frame with
+ const Bool_t isHE, const Bool_t isTheta)
+{
+ // The function calculates theta and phi in the mother rest frame with
// respect to the helicity coordinate system and Collins-Soper coordinate system
// TO DO: generalize for different decays (only J/Psi->e+e- now)
// Laboratory frame 4-vectors:
// projectile beam & target beam 4-mom
- const Double_t kBeamEnergy = 3500.; //TODO: need to retrieve the beam energy from somewhere
- TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton)));
- TLorentzVector targMom(0.,0.,kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton)));
+ // TODO: need to retrieve the beam energy from somewhere
+ const Double_t kBeamEnergy = 3500.;
+ Double_t px1=d1->Px();
+ Double_t py1=d1->Py();
+ Double_t pz1=d1->Pz();
+ Double_t px2=d2->Px();
+ Double_t py2=d2->Py();
+ Double_t pz2=d2->Pz();
+ Double_t eleMass=AliPID::ParticleMass(AliPID::kElectron);
+ Double_t proMass=AliPID::ParticleMass(AliPID::kProton);
+
+ TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+ TLorentzVector targMom(0.,0., kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
// first & second daughter 4-mom
- TLorentzVector p1Mom(d1->Px(),d1->Py(),d1->Pz(),TMath::Sqrt(d1->Px()*d1->Px()+d1->Py()*d1->Py()+d1->Pz()*d1->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
- TLorentzVector p2Mom(d2->Px(),d2->Py(),d2->Pz(),TMath::Sqrt(d2->Px()*d2->Px()+d2->Py()*d2->Py()+d2->Pz()*d2->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
+ TLorentzVector p1Mom(px1,py1,pz1,TMath::Sqrt(px1*px1+py1*py1+pz1*pz1+eleMass*eleMass));
+ TLorentzVector p2Mom(px2,py2,pz2,TMath::Sqrt(px2*px2+py2*py2+pz2*pz2+eleMass*eleMass));
// J/Psi 4-momentum vector
TLorentzVector motherMom=p1Mom+p2Mom;
// Laboratory frame 4-vectors:
// projectile beam & target beam 4-mom
- const Double_t kBeamEnergy = 3500.; //TODO: need to retrieve the beam energy from somewhere
- TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton)));
- TLorentzVector targMom(0.,0.,kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+AliPID::ParticleMass(AliPID::kProton)*AliPID::ParticleMass(AliPID::kProton)));
-
- // first & second daughter 4-mom
AliVParticle *d1 = dynamic_cast<AliVParticle*>(fRefD1.GetObject());
AliVParticle *d2 = dynamic_cast<AliVParticle*>(fRefD2.GetObject());
- TLorentzVector p1Mom(d1->Px(),d1->Py(),d1->Pz(),TMath::Sqrt(d1->Px()*d1->Px()+d1->Py()*d1->Py()+d1->Pz()*d1->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
- TLorentzVector p2Mom(d2->Px(),d2->Py(),d2->Pz(),TMath::Sqrt(d2->Px()*d2->Px()+d2->Py()*d2->Py()+d2->Pz()*d2->Pz()+AliPID::ParticleMass(AliPID::kElectron)*AliPID::ParticleMass(AliPID::kElectron)));
+
+ const Double_t kBeamEnergy = 3500.;
+ Double_t px1=d1->Px();
+ Double_t py1=d1->Py();
+ Double_t pz1=d1->Pz();
+ Double_t px2=d2->Px();
+ Double_t py2=d2->Py();
+ Double_t pz2=d2->Pz();
+ Double_t eleMass=AliPID::ParticleMass(AliPID::kElectron);
+ Double_t proMass=AliPID::ParticleMass(AliPID::kProton);
+
+ TLorentzVector projMom(0.,0.,-kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+ TLorentzVector targMom(0.,0., kBeamEnergy,TMath::Sqrt(kBeamEnergy*kBeamEnergy+proMass*proMass));
+
+ // first & second daughter 4-mom
+ // first & second daughter 4-mom
+ TLorentzVector p1Mom(px1,py1,pz1,TMath::Sqrt(px1*px1+py1*py1+pz1*pz1+eleMass*eleMass));
+ TLorentzVector p2Mom(px2,py2,pz2,TMath::Sqrt(px2*px2+py2*py2+pz2*pz2+eleMass*eleMass));
// J/Psi 4-momentum vector
TLorentzVector motherMom=p1Mom+p2Mom;
virtual Double_t Xv() const { return fPair.GetX(); }
virtual Double_t Yv() const { return fPair.GetY(); }
virtual Double_t Zv() const { return fPair.GetZ(); }
- virtual Bool_t XvYvZv(Double_t x[3]) const { x[0]=Xv(); x[1]=Xv(); x[2]=Zv(); return kTRUE; }
+ virtual Bool_t XvYvZv(Double_t x[3]) const { x[0]=Xv(); x[1]=Yv(); x[2]=Zv(); return kTRUE; }
virtual Double_t OneOverPt() const { return Pt()>0.?1./Pt():0.; } //TODO: check
virtual Double_t Phi() const { return fPair.GetPhi();}
void SetLabel(Int_t label) {fLabel=label;}
//inter leg information
- Double_t OpeningAngle() const { return fD1.GetAngle(fD2); }
- Double_t DistanceDaughters() const { return fD1.GetDistanceFromParticle(fD2); }
- Double_t DistanceDaughtersXY() const { return fD1.GetDistanceFromParticleXY(fD2); }
- Double_t DeviationDaughters() const { return fD1.GetDeviationFromParticle(fD2); }
- Double_t DeviationDaughtersXY() const { return fD1.GetDeviationFromParticleXY(fD2); }
+ Double_t OpeningAngle() const { return fD1.GetAngle(fD2); }
+ Double_t DistanceDaughters() const { return fD1.GetDistanceFromParticle(fD2); }
+ Double_t DistanceDaughtersXY() const { return fD1.GetDistanceFromParticleXY(fD2); }
+ Double_t DeviationDaughters() const { return fD1.GetDeviationFromParticle(fD2); }
+ Double_t DeviationDaughtersXY() const { return fD1.GetDeviationFromParticleXY(fD2); }
+ Double_t DeltaEta() const { return TMath::Abs(fD1.GetEta()-fD2.GetEta()); }
+ Double_t DeltaPhi() const { Double_t dphi=TMath::Abs(fD1.GetPhi()-fD2.GetPhi());
+ return (dphi>TMath::Pi())?dphi-TMath::Pi():dphi; }
// calculate cos(theta*) and phi* in HE and CS pictures
+ void GetThetaPhiCM(Double_t &thetaHE, Double_t &phiHE, Double_t &thetaCS, Double_t &phiCS) const;
+
Double_t ThetaPhiCM(Bool_t isHE, Bool_t isTheta) const;
static Double_t ThetaPhiCM(const AliVParticle* d1, const AliVParticle* d2,
- const Bool_t isHE, const Bool_t isTheta);
-
+ const Bool_t isHE, const Bool_t isTheta);
// internal KF particle
const AliKFParticle& GetKFParticle() const { return fPair; }
const AliKFParticle& GetKFFirstDaughter() const { return fD1; }
if (!sig||!back||!combined) {
AliError("Both, signal and background function need to be set!");
+ return;
}
fSignal=sig;
fBackground=back;
"NclsTPC",
"NFclsTPC",
"TPCsignalN",
+ "TPCchi2PerCluster",
+ "TrackStatus",
+
"NclsTRD",
"TRDntracklets",
"TRDpidQuality",
"PhiCS",
"LegDistance",
"LegDistanceXY",
+ "DeltaEta",
+ "DeltaPhi",
"Merr",
"DCA",
"PairType",
kNclsTPC, // number of clusters assigned in the TPC
kNFclsTPC, // number of findable clusters in the TPC
kTPCsignalN, // number of points used for dEdx
+ kTPCchi2Cl, // chi2/cl in TPC
+ kTrackStatus, // track status bits
kNclsTRD, // number of clusters assigned in the TRD
kTRDntracklets, // number of TRD tracklets used for tracking/PID TODO: correct getter
kPhiCS, // phi in mother's rest frame in Collins-Soper picture
kLegDist, // distance of the legs
kLegDistXY, // distance of the legs in XY
+ kDeltaEta, // Absolute value of Delta Eta for the legs
+ kDeltaPhi, // Absolute value of Delta Phi for the legs
kMerr, // error of mass calculation
kDCA, // distance of closest approach TODO: not implemented yet
kPairType, // type of the pair, like like sign ++ unlikesign ...
Double_t pidProbs[AliPID::kSPECIES];
// Fill AliESDtrack interface specific information
+ Double_t tpcNcls=particle->GetNcls(1);
values[AliDielectronVarManager::kNclsITS] = particle->GetNcls(0); // TODO: get rid of the plain numbers
- values[AliDielectronVarManager::kNclsTPC] = particle->GetNcls(1); // TODO: get rid of the plain numbers
+ values[AliDielectronVarManager::kNclsTPC] = tpcNcls; // TODO: get rid of the plain numbers
values[AliDielectronVarManager::kNFclsTPC] = particle->GetTPCNclsF();
values[AliDielectronVarManager::kTPCsignalN] = particle->GetTPCsignalN();
values[AliDielectronVarManager::kNclsTRD] = particle->GetNcls(2); // TODO: get rid of the plain numbers
values[AliDielectronVarManager::kTRDntracklets] = particle->GetTRDntracklets(); // TODO: GetTRDtracklets/GetTRDntracklets?
values[AliDielectronVarManager::kTRDpidQuality] = particle->GetTRDpidQuality();
-
+ values[AliDielectronVarManager::kTrackStatus] = (Double_t)particle->GetStatus();
+
+
+ values[AliDielectronVarManager::kTPCchi2Cl] = -1;
+ if (tpcNcls>0) values[AliDielectronVarManager::kTPCchi2Cl] = particle->GetTPCchi2() / tpcNcls;
//TRD pidProbs
particle->GetTRDpid(pidProbs);
values[AliDielectronVarManager::kTRDprobEle] = pidProbs[AliPID::kElectron];
values[AliDielectronVarManager::kNclsTRD] = 0;
values[AliDielectronVarManager::kTRDntracklets] = 0;
values[AliDielectronVarManager::kTRDpidQuality] = 0;
-
+
+ values[AliDielectronVarManager::kTPCchi2Cl] = -1;
+ values[AliDielectronVarManager::kTrackStatus] = (Double_t)particle->GetStatus();
+
//TRD pidProbs
//TODO: set correctly
values[AliDielectronVarManager::kTRDprobEle] = 0;
// Fill AliDielectronPair specific information
const AliKFParticle &kfPair = pair->GetKFParticle();
-
+
+ Double_t thetaHE=0;
+ Double_t phiHE=0;
+ Double_t thetaCS=0;
+ Double_t phiCS=0;
+
+ pair->GetThetaPhiCM(thetaHE,phiHE,thetaCS,phiCS);
+
values[AliDielectronVarManager::kChi2NDF] = kfPair.GetChi2()/kfPair.GetNDF();
values[AliDielectronVarManager::kDecayLength] = kfPair.GetDecayLength();
values[AliDielectronVarManager::kR] = kfPair.GetR();
values[AliDielectronVarManager::kOpeningAngle] = pair->OpeningAngle();
- values[AliDielectronVarManager::kThetaHE] = pair->ThetaPhiCM(kTRUE, kTRUE);
- values[AliDielectronVarManager::kPhiHE] = pair->ThetaPhiCM(kTRUE, kFALSE);
- values[AliDielectronVarManager::kThetaCS] = pair->ThetaPhiCM(kFALSE, kTRUE);
- values[AliDielectronVarManager::kPhiCS] = pair->ThetaPhiCM(kFALSE, kFALSE);
+ values[AliDielectronVarManager::kThetaHE] = thetaHE;
+ values[AliDielectronVarManager::kPhiHE] = phiHE;
+ values[AliDielectronVarManager::kThetaCS] = thetaCS;
+ values[AliDielectronVarManager::kPhiCS] = phiCS;
values[AliDielectronVarManager::kLegDist] = pair->DistanceDaughters();
values[AliDielectronVarManager::kLegDistXY] = pair->DistanceDaughtersXY();
+ values[AliDielectronVarManager::kDeltaEta] = pair->DeltaEta();
+ values[AliDielectronVarManager::kDeltaPhi] = pair->DeltaPhi();
values[AliDielectronVarManager::kMerr] = kfPair.GetErrMass()>1e-30&&kfPair.GetMass()>1e-30?kfPair.GetErrMass()/kfPair.GetMass():1000000;
values[AliDielectronVarManager::kPairType] = pair->GetType();
}
TString configFile("$ALICE_ROOT/PWG3/dielectron/macros/ConfigJpsi2eeData.C");
- if (mgr->GetInputEventHandler()->IsA()==AliAODInputHandler::Class()){
- ::Info("AddTaskJPSI", "Using AOD configuration");
- configFile="$ALICE_ROOT/PWG3/dielectron/macros/ConfigJpsi2eeDataAOD.C";
- }
+// if (mgr->GetInputEventHandler()->IsA()==AliAODInputHandler::Class()){
+// ::Info("AddTaskJPSI", "Using AOD configuration");
+// configFile="$ALICE_ROOT/PWG3/dielectron/macros/ConfigJpsi2eeDataAOD.C";
+// }
//create task and add it to the manager
AliAnalysisTaskMultiDielectron *task=new AliAnalysisTaskMultiDielectron("MultiDie");
//create output container
AliAnalysisDataContainer *cOutputHist1 =
- mgr->CreateContainer("QA", TList::Class(), AliAnalysisManager::kOutputContainer,
+ mgr->CreateContainer("jpsi_QA", TList::Class(), AliAnalysisManager::kOutputContainer,
containerName.Data());
AliAnalysisDataContainer *cOutputHist2 =
- mgr->CreateContainer("CF", TList::Class(), AliAnalysisManager::kOutputContainer,
+ mgr->CreateContainer("jpsi_CF", TList::Class(), AliAnalysisManager::kOutputContainer,
+ containerName.Data());
+
+ AliAnalysisDataContainer *cOutputHist3 =
+ mgr->CreateContainer("jpsi_EventStat", TH1D::Class(), AliAnalysisManager::kOutputContainer,
containerName.Data());
mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
mgr->ConnectOutput(task, 1, cOutputHist1);
mgr->ConnectOutput(task, 2, cOutputHist2);
+ mgr->ConnectOutput(task, 3, cOutputHist3);
return task;
}
//create output container
AliAnalysisDataContainer *cOutputHist1 =
- mgr->CreateContainer("QA", TList::Class(), AliAnalysisManager::kOutputContainer,
+ mgr->CreateContainer("jpsi_FilterQA",
+ THashList::Class(),
+ AliAnalysisManager::kOutputContainer,
containerName.Data());
AliAnalysisDataContainer *cOutputHist2 =
- mgr->CreateContainer("CF", TList::Class(), AliAnalysisManager::kOutputContainer,
+ mgr->CreateContainer("jpsi_FilterEventStat",
+ TH1D::Class(),
+ AliAnalysisManager::kOutputContainer,
containerName.Data());
mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
void SetupPairCuts(AliDielectron *die, Int_t cutDefinition);
AliESDtrackCuts *SetupESDtrackCuts(Int_t cutDefinition);
-/*
-trackQ+Pt>1GeV+|TPCnSigE|<3
-trackQ+PID1
-trackQ+PID6
-trackQ+PID7
-trackQ+Pt>.5GeV
-*/
-TString names=("trackQ+Pt>1GeV+|TPCnSigE|<3;trackQ+PID1;trackQ+PID6;trackQ+PID7;trackQ+Pt>.5GeV");
+
+TString names=("basicQ+SPDfirst+pt>.6+PID");
+
TObjArray *arrNames=names.Tokenize(";");
const Int_t nDie=arrNames->GetEntries();
AliDielectron *die =
new AliDielectron(Form("%s",name.Data()),
Form("Track cuts: %s",name.Data()));
+
// cut setup
SetupTrackCuts(die,cutDefinition);
//
// histogram setup
// only if an AliDielectronHistos object is attached to the
- // dielectron framework the QA histograms will be filled
+ // dielectron framework histograms will be filled
//
InitHistograms(die,cutDefinition);
// the last definition uses no cuts and only the QA histograms should be filled!
-// if (cutDefinition<nDie-1) InitCF(die,cutDefinition);
-
+// if (cutDefinition<nDie-1)
+ InitCF(die,cutDefinition);
+
return die;
}
//ESD quality cuts
die->GetTrackFilter().AddCuts(SetupESDtrackCuts(cutDefinition));
-
- //QA no CF
- if (cutDefinition==nDie-1) {
- //Pt cut
- AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>.5","Pt>.5");
- pt->AddCut(AliDielectronVarManager::kPt,.5,1e30);
+ if(cutDefinition==0) {
+ //Pt cut ----------------------------------------------------------
+ AliDielectronVarCuts *pt = new AliDielectronVarCuts("ptCut","pt cut");
+ pt->AddCut(AliDielectronVarManager::kPt,0.6,1e30);
die->GetTrackFilter().AddCuts(pt);
- return;
- }
-
- // pt + 3 sigma ele TPC
- if (cutDefinition==0){
- AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>1+|TPCnSigE|<3","Pt>1+|TPCnSigE|<3");
- pt->AddCut(AliDielectronVarManager::kPt,1.,1e30);
- pt->AddCut(AliDielectronVarManager::kTPCnSigmaEle, -3., 3.);
- die->GetTrackFilter().AddCuts(pt);
- }
-
- //PID 1
- if (cutDefinition==1){
- AliDielectronPID *pid=new AliDielectronPID("pid1","pid1");
- pid->SetDefaults(1);
+ // PID cuts --------------------------------------------------------
+ AliDielectronPID *pid = new AliDielectronPID("PID10","TPC nSigma |e|<3 + |Pi|>3 + |P|>3 + TOF nSigma |e|<3");
+ pid->SetDefaults(10);
die->GetTrackFilter().AddCuts(pid);
}
-
- //PID 6
- if (cutDefinition==2){
- AliDielectronPID *pid=new AliDielectronPID("pid6","pid6");
- pid->SetDefaults(6);
- die->GetTrackFilter().AddCuts(pid);
- }
-
- //PID 7
- if (cutDefinition==3){
- AliDielectronPID *pid=new AliDielectronPID("pid7","pid7");
- pid->SetDefaults(7);
- die->GetTrackFilter().AddCuts(pid);
- }
-
}
//______________________________________________________________________________________
// Setup the pair cuts
//
-
- // reject conversions
- // and select mass region
- AliDielectronVarCuts *openingAngleCut=new AliDielectronVarCuts("OpeningAngle","Opening angle > 35mrad");
- openingAngleCut->AddCut(AliDielectronVarManager::kOpeningAngle,.035,4.);
- openingAngleCut->AddCut(AliDielectronVarManager::kM,2.,4.);
- die->GetPairFilter().AddCuts(openingAngleCut);
+ //Invariant mass selection
+ AliDielectronVarCuts *invMassCut=new AliDielectronVarCuts("2<M<4+|Y|<.8","2<M<4 + |Y|<.8");
+ invMassCut->AddCut(AliDielectronVarManager::kM,2.,4.);
+ invMassCut->AddCut(AliDielectronVarManager::kY,-0.8,0.8);
+ die->GetPairFilter().AddCuts(invMassCut);
}
//______________________________________________________________________________________
//
AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts;
+ // basic track quality cuts (basicQ)
esdTrackCuts->SetMaxDCAToVertexZ(3.0);
- esdTrackCuts->SetMaxDCAToVertexXY(.07);
- esdTrackCuts->SetRequireTPCRefit(kTRUE);
- esdTrackCuts->SetRequireITSRefit(kTRUE);
+ esdTrackCuts->SetMaxDCAToVertexXY(1.0);
+
+ esdTrackCuts->SetEtaRange( -0.8 , 0.8 );
+
esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
- esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
+ esdTrackCuts->SetRequireITSRefit(kTRUE);
+ esdTrackCuts->SetRequireTPCRefit(kTRUE);
- esdTrackCuts->SetMinNClustersTPC(110);
+ esdTrackCuts->SetMinNClustersTPC(90);
esdTrackCuts->SetMaxChi2PerClusterTPC(4);
+ esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
+
return esdTrackCuts;
}
+
//______________________________________________________________________________________
void InitHistograms(AliDielectron *die, Int_t cutDefinition)
{
//
- // Initialise the QA histograms
+ // Initialise the histograms
//
-
- //Setup QA histograms
- AliDielectronHistos *histos=
- new AliDielectronHistos(die->GetName(),
- die->GetTitle());
+
+ //Setup histogram Manager
+ AliDielectronHistos *histos=new AliDielectronHistos(die->GetName(),die->GetTitle());
//Initialise histogram classes
histos->SetReservedWords("Track;Pair");
- //Event class (only for last QA)
- if (cutDefinition==nDie-1) histos->AddClass("Event");
-
- //Track classes, only first event
+ //Track classes
+ //to fill also track info from 2nd event loop until 2
for (Int_t i=0; i<2; ++i){
histos->AddClass(Form("Track_%s",AliDielectron::TrackClassName(i)));
}
- //Pair classes, only first event
+ //Pair classes
+ // to fill also mixed event histograms loop until 10
for (Int_t i=0; i<3; ++i){
histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(i)));
}
-
- //Event histograms
- if (cutDefinition==nDie-1){
- //add histograms to event class
- histos->UserHistogram("Event","nEvents","Number of processed events after cuts;Number events",
- 1,0.,1.,AliDielectronVarManager::kNevents);
+
+ //legs from pair
+ for (Int_t i=0; i<3; ++i){
+ histos->AddClass(Form("Track_Legs_%s",AliDielectron::PairClassName(i)));
}
//add histograms to Track classes
- histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",200,0,20.,AliDielectronVarManager::kPt);
+ histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",400,0,20.,AliDielectronVarManager::kPt);
+ histos->UserHistogram("Track","TPCnCls","Number of Clusters TPC;TPC number clusteres;#tracks",160,-0.5,159.5,AliDielectronVarManager::kNclsTPC);
+
+ histos->UserHistogram("Track","dXY","dXY;dXY [cm];#tracks",500,-1.,1.,AliDielectronVarManager::kImpactParXY);
+ histos->UserHistogram("Track","dZ","dZ;dZ [cm];#tracks",600,-3.,3.,AliDielectronVarManager::kImpactParZ);
+ histos->UserHistogram("Track","Eta_Phi","Eta Phi Map; Eta; Phi;#tracks",
+ 200,-1,1,200,0,6.285,AliDielectronVarManager::kEta,AliDielectronVarManager::kPhi);
+
histos->UserHistogram("Track","dEdx_P","dEdx;P [GeV];TPC signal (arb units);#tracks",
400,0.2,20.,200,0.,200.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCsignal,kTRUE);
histos->UserHistogram("Track","TPCnSigmaEle_P","TPC number of sigmas Electrons;P [GeV];TPC number of sigmas Electrons;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
- histos->UserHistogram("Track","TPCnSigmaPio_P","TPC number of sigmas Pions;P [GeV];TPC number of sigmas Pions;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaPio,kTRUE);
- histos->UserHistogram("Track","dXY","dXY;dXY [cm];#tracks",200,-2.,2.,AliDielectronVarManager::kImpactParXY);
-
+ 400,0.2,20.,200,-10.,10.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
+
//add histograms to Pair classes
histos->UserHistogram("Pair","InvMass","Inv.Mass;Inv. Mass [GeV];#pairs",
- 500,0.,4.,AliDielectronVarManager::kM);
+ 201,-.01,4.01,AliDielectronVarManager::kM);
histos->UserHistogram("Pair","Rapidity","Rapidity;Rapidity;#pairs",
- 100,-2.,2.,AliDielectronVarManager::kY);
+ 100,-1.,1.,AliDielectronVarManager::kY);
histos->UserHistogram("Pair","OpeningAngle","Opening angle;angle",
100,0.,3.15,AliDielectronVarManager::kOpeningAngle);
- histos->UserHistogram("Pair","Chi2/NDF","#Chi^{2}/NDF;#Chi^{2}/NDF",
- 100, 0., 20., AliDielectronVarManager::kChi2NDF);
die->SetHistogramManager(histos);
}
//
// Setupd the CF Manager if needed
//
+
AliDielectronCF *cf=new AliDielectronCF(die->GetName(),die->GetTitle());
//pair variables
- cf->AddVariable(AliDielectronVarManager::kPt,20,0,10);
+ TVectorD *binLimPt=new TVectorD(6);
+ (*binLimPt)[0]=0.0; (*binLimPt)[1]=0.8; (*binLimPt)[2]=1.4; (*binLimPt)[3]=2.8; (*binLimPt)[4]=4.2; (*binLimPt)[5]=9.9;
+ cf->AddVariable(AliDielectronVarManager::kPt,binLimPt);
cf->AddVariable(AliDielectronVarManager::kY,40,-2,2);
- cf->AddVariable(AliDielectronVarManager::kM,200,-0.01,3.99);
+ cf->AddVariable(AliDielectronVarManager::kM,150,0.,150*.027); //27Mev Steps
cf->AddVariable(AliDielectronVarManager::kPairType,10,0,10);
//leg variables
- cf->AddVariable(AliDielectronVarManager::kPt,20,0,10,kTRUE);
-
- //only in this case write MC truth info
- cf->SetStepsForCutsIncreasing();
- if (cutDefinition==0){
- cf->SetStepForMCtruth();
- }
+ cf->AddVariable(AliDielectronVarManager::kPt,2,0.8,1.2,kTRUE);
+ cf->AddVariable(AliDielectronVarManager::kNclsTPC,3,90,120,kTRUE);
die->SetCFManagerPair(cf);
+
}
+
die->GetTrackFilter().AddCuts(SetupESDtrackCuts());
//Pt cut
- AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>.5+60<dEdx<100","Pt>.5 && 60<dEdx<100");
- pt->AddCut(AliDielectronVarManager::kPt,.5,1e30);
+ AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>.5+60<dEdx<100","Pt>.6 && 60<dEdx<100");
+ pt->AddCut(AliDielectronVarManager::kPt,.6,1e30);
pt->AddCut(AliDielectronVarManager::kTPCsignal,60.,100.);
die->GetTrackFilter().AddCuts(pt);
//Invarian mass selection
AliDielectronVarCuts *invMassCut=new AliDielectronVarCuts("InvMass","2<M<4");
- invMassCut->AddCut(AliDielectronVarManager::kM,2.,4.);
+ invMassCut->AddCut(AliDielectronVarManager::kM,2.,1e30);
// invMassCut->AddCut(AliDielectronVarManager::kPairType,1.);
die->GetPairFilter().AddCuts(invMassCut);
esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
- esdTrackCuts->SetMinNClustersTPC(100);
+ esdTrackCuts->SetMinNClustersTPC(90);
esdTrackCuts->SetMaxChi2PerClusterTPC(4);
return esdTrackCuts;
//Initialise histogram classes
histos->SetReservedWords("Track;Pair");
- //Event class
- histos->AddClass("Event");
//Track classes
//to fill also track info from 2nd event loop until 2
histos->AddClass(Form("Pair_%s",AliDielectron::PairClassName(i)));
}
- //add histograms to event class
- histos->UserHistogram("Event","nEvents","Number of processed events after cuts;Number events",
- 1,0.,1.,AliDielectronVarManager::kNevents);
-
//add histograms to Track classes
- histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",200,0,20.,AliDielectronVarManager::kPt);
+ histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",400,0,20.,AliDielectronVarManager::kPt);
+ histos->UserHistogram("Track","TPCnCls","Number of Clusters TPC;TPC number clusteres;#tracks",160,-0.5,159.5,AliDielectronVarManager::kNclsTPC);
+
+ histos->UserHistogram("Track","dXY","dXY;dXY [cm];#tracks",500,-1.,1.,AliDielectronVarManager::kImpactParXY);
+ histos->UserHistogram("Track","dZ","dZ;dZ [cm];#tracks",600,-3.,3.,AliDielectronVarManager::kImpactParZ);
+ histos->UserHistogram("Track","Eta_Phi","Eta Phi Map; Eta; Phi;#tracks",
+ 200,-1,1,200,0,6.285,AliDielectronVarManager::kEta,AliDielectronVarManager::kPhi);
+
histos->UserHistogram("Track","dEdx_P","dEdx;P [GeV];TPC signal (arb units);#tracks",
400,0.2,20.,200,0.,200.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCsignal,kTRUE);
histos->UserHistogram("Track","TPCnSigmaEle_P","TPC number of sigmas Electrons;P [GeV];TPC number of sigmas Electrons;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
- histos->UserHistogram("Track","TPCnSigmaKao_P","TPC number of sigmas Kaons;P [GeV];TPC number of sigmas Kaons;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaKao,kTRUE);
- histos->UserHistogram("Track","TPCnSigmaPio_P","TPC number of sigmas Pions;P [GeV];TPC number of sigmas Pions;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaPio,kTRUE);
-
-
- histos->UserHistogram("Track","TOFnSigmaKao_P","TOF number of sigmas Kaons;P [GeV];TOF number of sigmas Kaons;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTOFnSigmaKao,kTRUE);
- histos->UserHistogram("Track","TOFnSigmaPro_P","TOF number of sigmas Protons;P [GeV];TOF number of sigmas Protons;#tracks",
- 400,0.2,20.,100,-5.,5.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTOFnSigmaPro,kTRUE);
-
- histos->UserHistogram("Track","dXY","dXY;dXY [cm];#tracks",200,-2.,2.,AliDielectronVarManager::kImpactParXY);
-
- histos->UserHistogram("Track","TPCnCls","Number of Clusters TPC;TPC number clusteres;#tracks",159,0.,159.,AliDielectronVarManager::kNclsTPC);
+ 400,0.2,20.,200,-10.,10.,AliDielectronVarManager::kPIn,AliDielectronVarManager::kTPCnSigmaEle,kTRUE);
//add histograms to Pair classes
histos->UserHistogram("Pair","InvMass","Inv.Mass;Inv. Mass [GeV];#pairs",
201,-.01,4.01,AliDielectronVarManager::kM);
histos->UserHistogram("Pair","Rapidity","Rapidity;Rapidity;#pairs",
- 100,-2.,2.,AliDielectronVarManager::kY);
+ 100,-1.,1.,AliDielectronVarManager::kY);
histos->UserHistogram("Pair","OpeningAngle","Opening angle;angle",
100,0.,3.15,AliDielectronVarManager::kOpeningAngle);
- histos->UserHistogram("Pair","Chi2/NDF","#Chi^{2}/NDF;#Chi^{2}/NDF",
- 100, 0., 20., AliDielectronVarManager::kChi2NDF);
die->SetHistogramManager(histos);
}
+c->Draw("Leg1_TOF_nSigma_Electrons:Leg1_P_InnerParam>>nSigE","cut","colz")
+c->Draw("Leg2_TOF_nSigma_Electrons:Leg2_P_InnerParam>>+nSigE","cut","colz")
+
c->Draw("Leg1_TOF_nSigma_Protons:Leg1_P_InnerParam>>nSigP","cut","goff")
c->Draw("Leg2_TOF_nSigma_Protons:Leg2_P_InnerParam>>+nSigP","cut","colz")
c->SetAlias("cut","Leg2_P_InnerParam>1.5&&cutE")
c->SetMarkerStyle(20);
c->SetMarkerColor(kRed);
+c->Draw("M>>hM3(50,1.99,3.99)","cut","esame");
+
+
+
+
+//=================== pass 1 pass 2 comparison =======================
+
+
+//------ cuts ------------
+
+
+c->SetAlias("cutPro","Leg1_TPC_signal>50*exp(-.5*(Leg1_P_InnerParam-2))&&Leg2_TPC_signal>50*exp(-.5*(Leg2_P_InnerParam-2))&&Leg1_TPC_signal<85&&Leg2_TPC_signal<85");
+
+c->SetAlias("Pcut","Leg1_P_InnerParam>1.3&&Leg2_P_InnerParam>1.3")
+c->SetAlias("Ptcut","Leg1_Pt>1&&Leg2_Pt>1")
+c->SetAlias("TOFcut","abs(Leg1_TOF_nSigma_Electrons)<3&&abs(Leg2_TOF_nSigma_Electrons)<3");
+c->SetAlias("TOFcut2","(Leg1_P_InnerParam<1.3&&abs(Leg1_TOF_nSigma_Electrons)<3||Leg1_P_InnerParam>=1.3)&&(Leg2_P_InnerParam<1.3&&abs(Leg2_TOF_nSigma_Electrons)<3||Leg2_P_InnerParam>=1.3)");
+c->SetAlias("TPCcut","(Leg1_TPC_signal>70*(1-exp(-1*(Leg1_P_InnerParam+2))))&&(Leg2_TPC_signal>70*(1-exp(-1*(Leg2_P_InnerParam+2))))")
+c->SetAlias("NClcut","Leg1_NclsTPC>120&&Leg2_NclsTPC>120");
+
+c->SetAlias("eleParam","Leg1_TPC_nSigma_Electrons<5&&Leg2_TPC_nSigma_Electrons<5&&Leg1_TPC_nSigma_Electrons>-2.65*exp(-0.6757*Leg1_P_InnerParam)&&Leg2_TPC_nSigma_Electrons>-2.65*exp(-0.6757*Leg2_P_InnerParam)")
+c->SetAlias("cut","PairType==1&&eleParam")
+c->SetAlias("cut","1==1")
+c->SetAlias("cut","NClcut")
+
+c->SetAlias("cut","TOFcut&&TPCcut&&NClcut")
+
+c->SetAlias("cut","TOFcut2&&TPCcut&&NClcut")
+
+c->SetAlias("cut","cutPro&&TPCcut&&NClcut")
+
+c->SetAlias("cut","Pcut&&TPCcut&&NClcut")
+
+c->SetAlias("cut","Ptcut&&TPCcut&&NClcut")
+
+
+
+//------------ plots --------------
+
+//no cut
+c->SetAlias("cut","1==1")
+c1->SetLogx(0)
+c1->SetLogz(0)
+c->SetLineColor(kBlack);
+c->Draw("M>>hM(50,1.99,3.99)","cut","e");
+hM->SetTitle(";Inv. Mass [GeV]; Pair (e+e-) / 40GeV")
+c1->Modified()
+c1->Update();
+c1->Print("pics/M_noCut.png");
+
+//=============
+//ncl: No cut
+//=============
+c->SetAlias("cut","1==1")
+c1->SetLogx(1)
+c1->SetLogz(0)
+gStyle->SetOptStat(0);
+c->Draw("Leg1_NclsTPC:Leg1_TPC_signal:Leg1_P_InnerParam>>test(1000,0,40,200,60,100)","cut","profcolz")
+c->Draw("Leg2_NclsTPC:Leg2_TPC_signal:Leg2_P_InnerParam>>+test","cut","profcolz")
+test->SetMinimum(80)
+test->SetTitle("mean TPC number of clusters;P_{TPC} [GeV]; TPC signal [arb. units]")
+c1->Modified()
+c1->Update();
+c1->Print("pics/TPCnCl_P.png");
+
+//=============
+//TPCsignal: ncl cut
+//=============
+c->SetAlias("cut","NClcut")
+c->Draw("Leg1_NclsTPC:Leg1_TPC_signal:Leg1_P_InnerParam>>test(1000,0,40,200,60,100)","cut","profcolz")
+c->Draw("Leg2_NclsTPC:Leg2_TPC_signal:Leg2_P_InnerParam>>+test","cut","profcolz")
+test->SetMinimum(80)
+test->SetTitle("mean TPC number of clusters;P_{TPC} [GeV]; TPC signal [arb. units]")
+c1->Modified()
+c1->Update();
+c1->Print("pics/TPCnCl_P_cutNcl.png");
+
+
+//=============
+//tpc signal + signal cut
+//=============
+c1->SetLogx(1)
+c1->SetLogy(0)
+c1->SetLogz(1)
+h.GetHistogram("TPCsignal","sigTPC")->GetYaxis()->SetRangeUser(60,100);
+c->SetAlias("cut","NClcut")
+c->Draw("Leg1_TPC_signal:Leg1_P_InnerParam>>sigTPC","cut","colz")
+c->Draw("Leg2_TPC_signal:Leg2_P_InnerParam>>+sigTPC","cut","colz")
+TF1 f("f1","[0]*(1-exp(-[1]*(x-[2])))",0.3,40);
+f.SetParameters(70,1,-2);
+f.Draw("same");
+c1->Modified();
+c1->Update();
+c1->Print("pics/TPCsignal_P_cutNcl.png");
+
+//------- Mass
+
+c1->SetLogx(0)
+c1->SetLogy(1)
+c1->SetLogz(0)
+c->SetAlias("cut","1==1")
+c->SetLineColor(kBlack);
+c->SetMarkerColor(kBlack);
+c->Draw("M>>hM(50,1.99,3.99)","cut","e");
+hM->SetTitle(";Inv. Mass [GeV]; Pair (e+e-) / 40GeV")
+hM->SetMinimum(5e2);
+c->SetAlias("cut","NClcut")
+c->SetLineColor(kRed);
+c->SetMarkerColor(kRed);
c->Draw("M>>hM2(50,1.99,3.99)","cut","esame");
+c1->Modified();
+c1->Update();
+c1->Print("pics/M_nClCut.png");
+
+
+//==========
+//tpc signal: ncl + tpc cut
+//==========
+c1->SetLogx(1)
+c1->SetLogy(0)
+c1->SetLogz(1)
+c->SetAlias("cut","TPCcut&&NClcut")
+c->Draw("Leg1_TPC_signal:Leg1_P_InnerParam>>sigTPC","cut","colz")
+c->Draw("Leg2_TPC_signal:Leg2_P_InnerParam>>+sigTPC","cut","colz")
+c1->Modified();
+c1->Update();
+c1->Print("pics/TPCsignal_P_cutNcl_tpc.png");
+
+/--- Mass
+
+c1->SetLogx(0)
+c1->SetLogy(1)
+c1->SetLogz(0)
+c->SetAlias("cut","1==1")
+c->SetLineColor(kBlack);
+c->SetMarkerColor(kBlack);
+c->Draw("M>>hM(50,1.99,3.99)","cut","e");
+hM->SetTitle(";Inv. Mass [GeV]; Pair (e+e-) / 40GeV")
+hM->SetMinimum(5);
+c->SetAlias("cut","NClcut")
+c->SetLineColor(kRed);
+c->SetMarkerColor(kRed);
+c->Draw("M>>hM2(50,1.99,3.99)","cut","esame");
+c->SetAlias("cut","TPCcut&&NClcut")
+c->SetLineColor(kGreen);
+c->SetMarkerColor(kGreen);
+c->Draw("M>>hM3(50,1.99,3.99)","cut","esame");
+c1->Modified();
+c1->Update();
+c1->Print("pics/M_nClCut_tpc.png");
+
+
+//========
+//TPC signal: ncl + tpc + tof cut
+//=======
+c1->SetLogx(1)
+c1->SetLogy(0)
+c1->SetLogz(1)
+c->SetAlias("cut","TOFcut2&&TPCcut&&NClcut")
+c->Draw("Leg1_TPC_signal:Leg1_P_InnerParam>>sigTPC","cut","colz")
+c->Draw("Leg2_TPC_signal:Leg2_P_InnerParam>>+sigTPC","cut","colz")
+c1->Modified();
+c1->Update();
+c1->Print("pics/TPCsignal_P_cutNcl_tpc.png");
+
+//--- Mass
+
+c1->SetLogx(0)
+c1->SetLogy(0)
+c1->SetLogz(0)
+c->SetAlias("cut","1==1")
+c->SetAlias("cut","TPCcut&&NClcut")
+c->SetLineColor(kGreen);
+c->Draw("M>>hM(50,1.99,3.99)","cut","e");
+hM->SetTitle(";Inv. Mass [GeV]; Pair (e+e-) / 40GeV")
+hM->SetMinimum(.1);
+c->SetAlias("cut","TOFcut2&&TPCcut&&NClcut")
+c->SetLineColor(kBlue);
+c->Draw("M>>hM2(50,1.99,3.99)","cut","esame");
+c1->Modified();
+c1->Update();
+c1->Print("pics/M_nClCut_tpc.png");
+
+//========
+//Inv Mass: different cuts
+//=======
+
+c1->SetLogx(0)
+c1->SetLogy(0)
+c1->SetLogz(0)
+c->SetAlias("cut","Ptcut&&TPCcut&&NClcut")
+c->SetLineColor(kMagenta);
+c->SetMarkerColor(kMagenta);
+c->SetMarkerStyle(22);
+c->Draw("M>>hM(50,1.99,3.99)","cut","e");
+
+c->SetAlias("cut","Pcut&&TPCcut&&NClcut")
+c->SetLineColor(kCyan+1);
+c->SetMarkerColor(kCyan+1);
+c->SetMarkerStyle(21);
+c->Draw("M>>hM2(50,1.99,3.99)","cut","esame");
+
+c->SetAlias("cut","TOFcut2&&TPCcut&&NClcut")
+c->SetMarkerStyle(20);
+c->SetLineColor(kBlue);
+c->SetMarkerColor(kBlue);
+c->Draw("M>>hM3(50,1.99,3.99)","cut","esame");
+
+c1->Modified();
+c1->Update();
+c1->Print("pics/M_nClCut_tpc_tof.png");
+
+
+
*/
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff