--- /dev/null
+AliAnalysisTask *AddTaskJPSI(const char* config=""){
+ //get the current analysis manager
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+ if (!mgr) {
+ ::Error("AddTasJPSI", "No analysis manager found.");
+ return NULL;
+ }
+ if (!mgr->GetInputEventHandler()) {
+ ::Error("AddTaskJPSI", "This task requires an input event handler");
+ return NULL;
+ }
+
+ 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";
+ }
+
+ //create task and add it to the manager
+ AliAnalysisTaskMultiDielectron *task=new AliAnalysisTaskMultiDielectron("MultiDie");
+ mgr->AddTask(task);
+
+ //load dielectron configuration file
+ gROOT->LoadMacro(configFile.Data());
+
+ //add dielectron analysis with different cuts to the task
+ for (Int_t i=0; i<nDie; ++i){ //nDie defined in config file
+ AliDielectron *jpsi=ConfigJpsi2ee(i);
+ task->AddDielectron(jpsi);
+ }
+
+ //----------------------
+ //create data containers
+ //----------------------
+
+ //find input container
+ AliAnalysisDataContainer *cinput = mgr->GetCommonInputContainer();
+
+ TString containerName = mgr->GetCommonFileName();
+ containerName += ":PWG3_dielectron";
+
+ //create output container
+
+ AliAnalysisDataContainer *cOutputHist1 =
+ mgr->CreateContainer("QA", TList::Class(), AliAnalysisManager::kOutputContainer,
+ containerName.Data());
+
+ AliAnalysisDataContainer *cOutputHist2 =
+ mgr->CreateContainer("CF", TList::Class(), AliAnalysisManager::kOutputContainer,
+ containerName.Data());
+
+ mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer());
+ mgr->ConnectOutput(task, 1, cOutputHist1);
+ mgr->ConnectOutput(task, 2, cOutputHist2);
+
+ return task;
+}
#include <TChain.h>
#include <AliCFContainer.h>
+#include <AliInputEventHandler.h>
+#include <AliESDInputHandler.h>
+#include <AliAnalysisManager.h>
#include <AliVEvent.h>
#include "AliDielectron.h"
#include "AliDielectronHistos.h"
#include "AliDielectronCF.h"
+#include "AliDielectronMC.h"
#include "AliAnalysisTaskMultiDielectron.h"
ClassImp(AliAnalysisTaskMultiDielectron)
AliAnalysisTaskSE(),
fListDielectron(),
fListHistos(),
- fListCF()
+ fListCF(),
+ fSelectPhysics(kFALSE)
{
//
// Constructor
AliAnalysisTaskSE(name),
fListDielectron(),
fListHistos(),
- fListCF()
+ fListCF(),
+ fSelectPhysics(kFALSE)
{
//
// Constructor
// Add all histogram manager histogram lists to the output TList
//
- if (!fListHistos.IsEmpty()) return; //already initialised
+ if (!fListHistos.IsEmpty()||!fListCF.IsEmpty()) return; //already initialised
TIter nextDie(&fListDielectron);
AliDielectron *die=0;
// Main loop. Called for every event
//
- if (fListHistos.IsEmpty()) return;
+ if (fListHistos.IsEmpty()&&fListCF.IsEmpty()) 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());
+ Bool_t isSelected = kTRUE;
+ if( fSelectPhysics && inputHandler && inputHandler->GetEventSelection() ) {
+ isSelected = inputHandler->IsEventSelected();
+ }
+
+ if (!isSelected) return;
//bz for AliKF
Double_t bz = InputEvent()->GetMagneticField();
PostData(2, &fListCF);
}
+//_________________________________________________________________________________
+void AliAnalysisTaskMultiDielectron::FinishTaskOutput()
+{
+ //
+ // Write debug tree
+ //
+ TIter nextDie(&fListDielectron);
+ AliDielectron *die=0;
+ while ( (die=static_cast<AliDielectron*>(nextDie())) ){
+ die->SaveDebugTree();
+ }
+}
+
public:
AliAnalysisTaskMultiDielectron();
AliAnalysisTaskMultiDielectron(const char *name);
- virtual ~AliAnalysisTaskMultiDielectron(){}
+ virtual ~AliAnalysisTaskMultiDielectron(){ }
- virtual void UserExec(Option_t *option);
- virtual void UserCreateOutputObjects();
+ virtual void UserExec(Option_t *option);
+ virtual void UserCreateOutputObjects();
+ virtual void FinishTaskOutput();
+ void UsePhysicsSelection(Bool_t phy=kTRUE) {fSelectPhysics=phy;}
void AddDielectron(AliDielectron * const die) { fListDielectron.Add(die); }
TList fListDielectron; // List of dielectron framework instances
TList fListHistos; //! List of histogram manager lists in the framework classes
TList fListCF; //! List with CF Managers
+
+ Bool_t fSelectPhysics; // Whether to use physics selection
AliAnalysisTaskMultiDielectron(const AliAnalysisTaskMultiDielectron &c);
AliAnalysisTaskMultiDielectron& operator= (const AliAnalysisTaskMultiDielectron &c);
#include "AliDielectronCF.h"
#include "AliDielectronMC.h"
#include "AliDielectronVarManager.h"
+#include "AliDielectronDebugTree.h"
+
#include "AliDielectron.h"
ClassImp(AliDielectron)
fPdgMother(443),
fHistos(0x0),
fPairCandidates(new TObjArray(10)),
- fCfManagerPair(0x0)
+ fCfManagerPair(0x0),
+ fDebugTree(0x0)
{
//
// Default constructor
fPdgMother(443),
fHistos(0x0),
fPairCandidates(new TObjArray(10)),
- fCfManagerPair(0x0)
+ fCfManagerPair(0x0),
+ fDebugTree(0x0)
{
//
// Named constructor
//
if (fHistos) delete fHistos;
if (fPairCandidates) delete fPairCandidates;
+ if (fDebugTree) delete fDebugTree;
}
//________________________________________________________________
// Initialise objects
//
if (fCfManagerPair) fCfManagerPair->InitialiseContainer(fPairFilter);
- AliDielectronVarManager::InitESDpid(); //currently this will take the values for MC
-}
+
+}
//________________________________________________________________
void AliDielectron::Process(AliVEvent *ev1, AliVEvent *ev2)
// Process the events
//
- AliDielectronVarManager::SetEvent(0x0);
+ AliDielectronVarManager::SetEvent(ev1);
//in case we have MC load the MC event and process the MC particles
- if (AliDielectronMC::Instance()->ConnectMCEvent()) ProcessMC();
+ if (AliDielectronMC::Instance()->HasMC()) {
+ AliDielectronMC::Instance()->ConnectMCEvent();
+ ProcessMC();
+ }
//if candidate array doesn't exist, create it
if (!fPairCandidates->UncheckedAt(0)) {
//in case there is a histogram manager, fill the QA histograms
if (fHistos) FillHistograms(ev1);
-
+
+ //fill debug tree if a manager is attached
+ if (fDebugTree) FillDebugTree();
}
//________________________________________________________________
delete candidate;
}
+//________________________________________________________________
+void AliDielectron::FillDebugTree()
+{
+ //
+ // Fill Histogram information for tracks and pairs
+ //
+
+ //Fill Debug tree
+ for (Int_t i=0; i<10; ++i){
+ Int_t ntracks=PairArray(i)->GetEntriesFast();
+ for (Int_t ipair=0; ipair<ntracks; ++ipair){
+ fDebugTree->Fill(static_cast<AliDielectronPair*>(PairArray(i)->UncheckedAt(ipair)));
+ }
+ }
+}
-
+//________________________________________________________________
+void AliDielectron::SaveDebugTree()
+{
+ //
+ // delete the debug tree, this will also write the tree
+ //
+ if (fDebugTree) fDebugTree->DeleteStreamer();
+}
class AliVEvent;
class THashList;
class AliDielectronCF;
+class AliDielectronDebugTree;
//________________________________________________________________
class AliDielectron : public TNamed {
void SetCFManagerPair(AliDielectronCF * const cf) { fCfManagerPair=cf; }
AliDielectronCF* GetCFManagerPair() const { return fCfManagerPair; }
+
+ void SetDebugTree(AliDielectronDebugTree * const tree) { fDebugTree=tree; }
static const char* TrackClassName(Int_t i) { return (i>=0&&i<4)?fgkTrackClassNames[i]:""; }
static const char* PairClassName(Int_t i) { return (i>=0&&i<10)?fgkPairClassNames[i]:""; }
+
+ void SaveDebugTree();
private:
AliDielectronCF *fCfManagerPair;//Correction Framework Manager for the Pair
+ AliDielectronDebugTree *fDebugTree; // Debug tree output
+
void FillTrackArrays(AliVEvent * const ev, Int_t eventNr=0);
void FillPairArrays(Int_t arr1, Int_t arr2);
void ProcessMC();
void FillHistograms(const AliVEvent *ev);
+ void FillDebugTree();
AliDielectron(const AliDielectron &c);
AliDielectron &operator=(const AliDielectron &c);
- ClassDef(AliDielectron,1);
+ ClassDef(AliDielectron,2);
};
inline void AliDielectron::InitPairCandidateArrays()
fStepForAfterAllCuts(kTRUE),
fStepsForEachCut(kFALSE),
fStepsForCutsIncreasing(kFALSE),
+ fStepsForSignal(kTRUE),
+ fStepsForBackground(kFALSE),
fNStepMasks(0),
fPdgMother(-1),
- fCfContainer(0x0)
+ fCfContainer(0x0),
+ fHasMC(kFALSE),
+ fNAddSteps(0)
{
//
// Default constructor
fStepForAfterAllCuts(kTRUE),
fStepsForEachCut(kFALSE),
fStepsForCutsIncreasing(kFALSE),
+ fStepsForSignal(kTRUE),
+ fStepsForBackground(kFALSE),
fNStepMasks(0),
fPdgMother(-1),
- fCfContainer(0x0)
+ fCfContainer(0x0),
+ fHasMC(kFALSE),
+ fNAddSteps(0)
{
//
// Named constructor
fNCuts=filter.GetCuts()->GetEntries();
+ fHasMC=AliDielectronMC::Instance()->HasMC();
+ fNAddSteps=1;
+ if (fHasMC){
+ if (fStepsForSignal) ++fNAddSteps;
+ if (fStepsForBackground) ++fNAddSteps;
+ } else {
+ //if
+ fStepForMCtruth=kFALSE;
+ fStepForNoCutsMCmotherPid=kFALSE;
+ fStepsForSignal=kFALSE;
+ fStepsForBackground=kFALSE;
+ }
+
fNSteps=0;
if (fStepForMCtruth) ++fNSteps;
if (fStepForNoCutsMCmotherPid) ++fNSteps;
- if (fStepForAfterAllCuts) fNSteps+=2;
-
- if (fStepsForEachCut&&fNCuts>1) fNSteps+=(2*fNCuts); //one step for each cut + Signal (MC)
- if (fStepsForCutsIncreasing&&fNCuts>2) fNSteps+=(2*(fNCuts-2)); //one step for the increasing cuts + Signal (MC)
+ if (fStepForAfterAllCuts) fNSteps+=fNAddSteps;
+
+ if (fStepsForEachCut&&fNCuts>1) fNSteps+=(fNAddSteps*fNCuts); //one step for each cut + Signal (MC)
+ if (fStepsForCutsIncreasing&&fNCuts>2) fNSteps+=(fNAddSteps*(fNCuts-2)); //one step for the increasing cuts + Signal (MC)
// e.g. cut2&cut3, cut2&cut3&cut4
- fNSteps+=(2*fNStepMasks); // cuts for the additional cut masks
+ fNSteps+=(fNAddSteps*fNStepMasks); // cuts for the additional cut masks
// create the container
Int_t *nbins=new Int_t[fNVars+2*fNVarsLeg];
for (Int_t i=0;i<fNVars;++i) nbins[i]=fNBins[i];
fCfContainer->SetStepTitle(step++,"No cuts (Signal)");
}
- //After All cuts
TString cutName;
- if (fStepForAfterAllCuts){
- cutName="All Cuts"; //TODO: User GetTitle???
- fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
- cutName+=" (Signal)";
- fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut with MC truth
- }
-
//Steps for each of the cuts
if (fStepsForEachCut&&fNCuts>1){
for (Int_t iCut=0; iCut<fNCuts;++iCut) {
cutName=filter.GetCuts()->At(iCut)->GetName(); //TODO: User GetTitle???
+
fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
- cutName+=" (Signal)";
- fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut with MC truth
+
+ if (fHasMC){
+ if (fStepsForSignal)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Signal)").Data()); //Step for the cut with MC truth
+ if (fStepsForBackground)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
+ }
}
}
for (Int_t iCut=1; iCut<fNCuts-1;++iCut) {
cutName+="&";
cutName+=filter.GetCuts()->At(iCut)->GetName();
+
fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
- cutName+=" (Signal)";
- fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut with MC truth
+
+ if (fHasMC){
+ if (fStepsForSignal)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Signal)").Data()); //Step for the cut with MC truth
+ if (fStepsForBackground)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
+ }
}
}
}
}
}
+
fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
- cutName+=" (Signal)";
- fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut with MC truth
+
+ if (fHasMC){
+ if (fStepsForSignal)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Signal)").Data()); //Step for the cut with MC truth
+ if (fStepsForBackground)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
+ }
+ }
+
+ //After All cuts
+ if (fStepForAfterAllCuts){
+ cutName="No pair cuts";
+ if (filter.GetCuts()->At(0)){
+ cutName=filter.GetCuts()->At(0)->GetName(); //TODO: User GetTitle???
+ for (Int_t iCut=1; iCut<fNCuts;++iCut) {
+ cutName+="&";
+ cutName+=filter.GetCuts()->At(iCut)->GetName();
+ }
+ fCfContainer->SetStepTitle(step++, cutName.Data()); //Step for the cut
+ }
+ if (fHasMC){
+ if (fStepsForSignal)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Signal)").Data()); //Step for the cut with MC truth
+ if (fStepsForBackground)
+ fCfContainer->SetStepTitle(step++, (cutName+" (Background)").Data()); //Step for the cut with MC truth
+ }
}
if (step!=fNSteps) {
//
Bool_t isMCTruth=kFALSE;
- if (fPdgMother>=0) isMCTruth=AliDielectronMC::Instance()->IsMotherPdg(particle,fPdgMother);
+ if (fHasMC) isMCTruth=AliDielectronMC::Instance()->IsMotherPdg(particle,fPdgMother);
Double_t valuesPair[AliDielectronVarManager::kNMaxValues];
AliDielectronVarManager::Fill(particle,valuesPair);
++step;
}
- //All cuts
- if (fStepForAfterAllCuts){
- if (mask == selectedMask){
- fCfContainer->Fill(fValues,step);
- ++step;
- if (isMCTruth) fCfContainer->Fill(fValues,step);
- ++step;
- } else {
- step+=2;
- }
- }
-
//Steps for each of the cuts
if (fStepsForEachCut&&fNCuts>1){
for (Int_t iCut=0; iCut<fNCuts;++iCut) {
if (mask&(1<<iCut)) {
fCfContainer->Fill(fValues,step);
++step;
- if (isMCTruth) fCfContainer->Fill(fValues,step);
- ++step;
+
+ if (fHasMC){
+ if ( fStepsForSignal){
+ if (isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ if ( fStepsForBackground ){
+ if (!isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ }
} else {
- step+=2;
+ step+=fNAddSteps;
}
}
}
if (mask&(1<<((iCut+1)-1))) {
fCfContainer->Fill(fValues,step);
++step;
- if (isMCTruth) fCfContainer->Fill(fValues,step);
- ++step;
+
+ if (fHasMC){
+ if ( fStepsForSignal){
+ if (isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ if ( fStepsForBackground ){
+ if (!isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ }
} else {
- step+=2;
+ step+=fNAddSteps;
}
}
}
if (mask&userMask) {
fCfContainer->Fill(fValues,step);
++step;
- if (isMCTruth) fCfContainer->Fill(fValues,step);
+
+ if (fHasMC){
+ if ( fStepsForSignal){
+ if (isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ if ( fStepsForBackground ){
+ if (!isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ }
+ } else {
+ step+=fNAddSteps;
+ }
+ }
+
+ //All cuts
+ if (fStepForAfterAllCuts){
+ if (mask == selectedMask){
+ fCfContainer->Fill(fValues,step);
++step;
+
+ if (fHasMC){
+ if ( fStepsForSignal){
+ if (isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ if ( fStepsForBackground ){
+ if (!isMCTruth) fCfContainer->Fill(fValues,step);
+ ++step;
+ }
+ }
} else {
- step+=2;
+ step+=fNAddSteps;
}
}
+ if (step!=fNSteps) {
+ AliError("Something went wrong in the step filling!!!");
+ }
+
}
//________________________________________________________________
void SetStepForAfterAllCuts(Bool_t steps=kTRUE) { fStepForAfterAllCuts=steps; }
void SetStepsForEachCut(Bool_t steps=kTRUE) { fStepsForEachCut=steps; }
void SetStepsForCutsIncreasing(Bool_t steps=kTRUE) { fStepsForCutsIncreasing=steps; }
+ void SetStepsForSignal(Bool_t steps=kTRUE) { fStepsForSignal=steps; }
+ void SetStepsForBackground(Bool_t steps=kTRUE) { fStepsForBackground=steps; }
void SetPdgMother(Int_t pdg) { fPdgMother=pdg; }
Int_t fNCuts; // Number of cuts in the filter concerned
- Double_t *fValues; // Value array for filling the container
+ Double_t *fValues; //! Value array for filling the container
Bool_t fStepForMCtruth; //create a step for the MC truth
Bool_t fStepForNoCutsMCmotherPid; //create a step for before cuts, but with MC truth of the mother
Bool_t fStepsForEachCut; //create steps for each cut?
Bool_t fStepsForCutsIncreasing; //create steps for increasing cut combinatons?
//e.g. cut1&cut2, cut1&cut2&cut3 ...
-
+ Bool_t fStepsForSignal; //steps for pure signal
+ Bool_t fStepsForBackground; //steps for pure background
+
UInt_t fStepMasks[kNmaxAddSteps]; //steps for additional cut combinatons
UInt_t fNStepMasks; //number of configured step masks
Int_t fPdgMother; //Pdg code of MCtruth validation
AliCFContainer* fCfContainer; //the CF container
+
+ Bool_t fHasMC; //if MC info is available
+ Int_t fNAddSteps; //number of additional MC related steps per cut step
AliDielectronCF(const AliDielectronCF &c);
AliDielectronCF &operator=(const AliDielectronCF &c);
#include <TList.h>
#include <TClass.h>
#include <TObject.h>
+#include <TVirtualPS.h>
#include <TFile.h>
#include <TString.h>
#include <TObjString.h>
+#include <TVectorD.h>
#include <TMath.h>
#include <TH1.h>
#include <TH2.h>
AliDielectronCFdraw::AliDielectronCFdraw() :
TNamed(),
fCfContainer(0x0),
- fEffGrid(0x0)
+ fEffGrid(0x0),
+ fVdata(0)
{
//
// Ctor
AliDielectronCFdraw::AliDielectronCFdraw(const char*name, const char* title) :
TNamed(name,title),
fCfContainer(0x0),
- fEffGrid(0x0)
+ fEffGrid(0x0),
+ fVdata(0)
{
//
// Named Ctor
AliDielectronCFdraw::AliDielectronCFdraw(AliCFContainer *cont) :
TNamed(cont->GetName(), cont->GetTitle()),
fCfContainer(cont),
- fEffGrid(new AliCFEffGrid("eff","eff",*cont))
+ fEffGrid(new AliCFEffGrid("eff","eff",*cont)),
+ fVdata(0)
{
//
// directly set the CF container
AliDielectronCFdraw::AliDielectronCFdraw(const char* filename) :
TNamed(),
fCfContainer(0x0),
- fEffGrid(0x0)
+ fEffGrid(0x0),
+ fVdata(0)
{
//
// get CF container(s) from file 'filename'
fCfContainer=new AliCFContainer(GetName(), GetTitle(), nstep, 1, nbins);
//delete unneeded steps
- for (Int_t istep=0; istep<nstep; ++istep) delete fCfContainer->GetGrid(istep);
+// for (Int_t istep=0; istep<nstep; ++istep) delete fCfContainer->GetGrid(istep);
//add step to the new container
Int_t istep=0;
if (arr->GetEntriesFast()==0){
// all slices in case of 0 entries
for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
- fCfContainer->SetRangeUser(ivar,fCfContainer->GetAxis(ivar,istep)->GetXmin(),
- fCfContainer->GetAxis(ivar,istep)->GetXmax(),istep);
+ fCfContainer->GetAxis(ivar,istep)->SetRange(0,0);
}
} else {
TIter next(arr);
TObjString *ostr=0x0;
while ( (ostr=static_cast<TObjString*>(next())) ) {
Int_t istep=ostr->GetString().Atoi();
- fCfContainer->SetRangeUser(ivar,fCfContainer->GetAxis(ivar,istep)->GetXmin(),
- fCfContainer->GetAxis(ivar,istep)->GetXmax(),istep);
+ fCfContainer->GetAxis(ivar,istep)->SetRange(0,0);
}
}
delete arr;
}
//________________________________________________________________
-void AliDielectronCFdraw::Draw(const Option_t* varnames, const char* slices, const char* opt)
+void AliDielectronCFdraw::Draw(const Option_t* varnames, const char* opt, const char* slices)
{
//
// Draw 'variables' of 'slices'
switch (entries){
case 1:
- Draw(ivar[0],slices,opt);
+ Draw(ivar[0],opt,slices);
break;
case 2:
- Draw(ivar[1],ivar[0],slices,opt);
+ Draw(ivar[1],ivar[0],opt,slices);
break;
case 3:
- Draw(ivar[2],ivar[1],ivar[0],slices,opt);
+ Draw(ivar[2],ivar[1],ivar[0],opt,slices);
break;
}
delete arrVars;
arrHists=new TObjArray(fCfContainer->GetNStep());
for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
TH1 *hproj=Project(dim,vars,istep);
+ if (!hproj) continue;
hproj->SetName(Form("proj_%02d",istep));
hproj->SetTitle(fCfContainer->GetStepTitle(istep));
arrHists->Add(hproj);
while ( (ostr=static_cast<TObjString*>(next())) ) {
Int_t istep=ostr->GetString().Atoi();
TH1 *hproj=Project(dim,vars,istep);
+ if (!hproj) continue;
hproj->SetName(Form("proj_%02d",istep));
hproj->SetTitle(fCfContainer->GetStepTitle(istep));
arrHists->Add(hproj);
const Int_t ndim=1;
Int_t vars[ndim]={var};
TObjArray *arr=CollectHistosEff(ndim,vars,nominators,denominator,type);
- Draw(arr,opt);
+ TString drawOpt=opt;
+ drawOpt+="eff";
+ Draw(arr,drawOpt);
delete arr;
}
const Int_t ndim=2;
Int_t vars[ndim]={var0,var1};
TObjArray *arr=CollectHistosEff(ndim,vars,nominators,denominator,type);
- Draw(arr,opt);
+ TString drawOpt=opt;
+ drawOpt+="eff";
+ Draw(arr,drawOpt);
delete arr;
}
const Int_t ndim=3;
Int_t vars[ndim]={var0,var1,var2};
TObjArray *arr=CollectHistosEff(ndim,vars,nominators,denominator,type);
- Draw(arr,opt);
+ TString drawOpt=opt;
+ drawOpt+="eff";
+ Draw(arr,drawOpt);
delete arr;
}
if (arr->GetEntriesFast()==0){
// all slices in case of 0 entries
arrHists=new TObjArray(fCfContainer->GetNStep());
+ fVdata.ResizeTo(arrHists->GetSize());
for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
fEffGrid->CalculateEfficiency(istep,denominator);
TH1 *hproj=ProjectEff(dim,vars);
+ if (!hproj) continue;
Float_t eff=fEffGrid->GetAverage();
+ fVdata(istep)=eff;
hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
- hproj->SetTitle(Form("%s (%f)",fCfContainer->GetStepTitle(istep),eff));
+ hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
arrHists->Add(hproj);
}
} else {
arrHists=new TObjArray(arr->GetEntriesFast());
TIter next(arr);
TObjString *ostr=0x0;
+ fVdata.ResizeTo(arrHists->GetSize());
+ Int_t count=0;
while ( (ostr=static_cast<TObjString*>(next())) ) {
Int_t istep=ostr->GetString().Atoi();
fEffGrid->CalculateEfficiency(istep,denominator);
TH1 *hproj=ProjectEff(dim,vars);
+ if (!hproj) continue;
Float_t eff=fEffGrid->GetAverage();
+ fVdata(count++)=eff;
hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
- hproj->SetTitle(Form("%s (%f)",fCfContainer->GetStepTitle(istep),eff));
+ hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
arrHists->Add(hproj);
}
}
if (arr->GetEntriesFast()==0){
// all slices in case of 0 entries
arrHists=new TObjArray(fCfContainer->GetNStep());
+ fVdata.ResizeTo(arrHists->GetSize());
for (Int_t istep=0; istep<fCfContainer->GetNStep(); ++istep){
TH1 *hproj=Project(dim,vars,istep);
+ if (!hproj) continue;
Float_t eff=0;
if (entriesDen>0) eff=hproj->GetEffectiveEntries()/entriesDen;
+ fVdata(istep)=eff;
hproj->Divide(hDen);
hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
- hproj->SetTitle(Form("%s (%f)",fCfContainer->GetStepTitle(istep),eff));
+ hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
arrHists->Add(hproj);
}
} else {
arrHists=new TObjArray(arr->GetEntriesFast());
+ fVdata.ResizeTo(arrHists->GetSize());
TIter next(arr);
TObjString *ostr=0x0;
+ Int_t count=0;
while ( (ostr=static_cast<TObjString*>(next())) ) {
Int_t istep=ostr->GetString().Atoi();
TH1 *hproj=Project(dim,vars,istep);
+ if (!hproj) continue;
Float_t eff=0;
if (entriesDen>0) eff=hproj->GetEffectiveEntries()/entriesDen;
+ fVdata(count++)=eff;
hproj->Divide(hDen);
hproj->SetName(Form("eff_%02d/%02d",istep,denominator));
- hproj->SetTitle(Form("%s (%f)",fCfContainer->GetStepTitle(istep),eff));
+ hproj->SetTitle(Form("%s (%.3f)",fCfContainer->GetStepTitle(istep),eff));
arrHists->Add(hproj);
}
}
//
TString optStr(opt);
optStr.ToLower();
- Bool_t drawSame=optStr.Contains("same");
- Bool_t optLeg=optStr.Contains("leg");
- optStr.ReplaceAll("same","");
+ Bool_t drawSame = optStr.Contains("same");
+ Bool_t drawSamePlus = optStr.Contains("same+");
+ Bool_t drawEff = optStr.Contains("eff");
+ Bool_t optLeg = optStr.Contains("leg");
+ Bool_t optScaleMax = optStr.Contains("max");
+
+ if (!drawSamePlus) optStr.ReplaceAll("same","");
+
+ optStr.ReplaceAll("+","");
+ optStr.ReplaceAll("eff","");
+ optStr.ReplaceAll("leg","");
+ optStr.ReplaceAll("max","");
if (!gPad) new TCanvas;
Int_t nPads = arr->GetEntriesFast();
if (nPads==0) return;
- if (nPads==1){
- arr->UncheckedAt(0)->Draw(optStr.Data());
- return;
- }
+// if (nPads==1){
+// arr->UncheckedAt(0)->Draw(optStr.Data());
+// return;
+// }
TCanvas *c=gPad->GetCanvas();
- c->Clear();
+ if (!gVirtualPS&&!drawSamePlus) c->Clear();
if (!drawSame){
}
} else {
TLegend *leg=0;
- if (optLeg) leg=new TLegend(.8,.3,.99,.9);
+ if (drawSamePlus){
+ //find already existing legend to attach entries to it
+ TIter nextPrimitive(gPad->GetListOfPrimitives());
+ TObject *o=0x0;
+ while ((o=nextPrimitive())) if (o->IsA()==TLegend::Class()) leg=(TLegend*)o;
+ }
+
+ if (optLeg&&!leg) leg=new TLegend(.2,.3,.99,.9);
+ Int_t addColor=0;
+ if (leg) addColor=leg->GetNRows();
Int_t offset=20;
if (nPads<7) offset=24;
for (Int_t i=0; i<nPads; ++i){
- if (i==1) optStr+="same";
+ if (i==1&&!drawSamePlus) optStr+="same";
TH1 *hist=static_cast<TH1*>(arr->UncheckedAt(i));
- hist->SetLineColor(i+1);
+ hist->SetLineColor(i+1+addColor);
hist->SetLineWidth(2);
- hist->SetMarkerColor(i+1);
- hist->SetMarkerStyle(offset+i);
+ hist->SetMarkerColor(i+1+addColor);
+ hist->SetMarkerStyle(offset+i+addColor);
hist->Draw(optStr.Data());
+
+ if (drawEff&&i==0&&!drawSamePlus) {
+ hist->GetYaxis()->SetRangeUser(0,2);
+ hist->SetYTitle("Rec. Signal / Gen. Signal");
+ }
+
if (leg) leg->AddEntry(hist,hist->GetTitle(),"lp");
}
if (leg){
leg->SetFillColor(10);
- leg->SetY1(.9-nPads*.05);
+ leg->SetY1(.9-leg->GetNRows()*.05);
+ leg->SetY1NDC(.9-leg->GetNRows()*.05);
leg->SetMargin(.1);
- leg->Draw();
+ if (!drawSamePlus) leg->Draw();
+ }
+ if (optScaleMax){
+ TIter nextPrimitive(gPad->GetListOfPrimitives());
+ TObject *o=0x0;
+ TH1 *hfirst=0x0;
+ Float_t max=0;
+ while ((o=nextPrimitive())) if (o->InheritsFrom(TH1::Class())){
+ TH1 *h=(TH1*)o;
+ if (!hfirst) hfirst=h;
+ if (h->GetMaximum()>max) max=h->GetMaximum();
+ }
+ max*=1.1;
+ hfirst->SetMaximum(max);
}
}
#include <TNamed.h>
+#include <TVectorD.h>
#include <AliCFContainer.h>
class TObjArray;
TObjArray* CollectHistosEff(Int_t dim, Int_t *vars, const char* nominators, Int_t denominator, Int_t type=0);
TH1* ProjectEff(Int_t ndim, Int_t *vars);
-
+ const TVectorD& GetData() const {return fVdata;}
void Draw(const TObjArray *arr, const char* opt="");
private:
AliCFContainer *fCfContainer; // CF container
AliCFEffGrid *fEffGrid; // Efficiency calculation
+
+ TVectorD fVdata; // vector with data, like mean efficiencies
AliDielectronCFdraw(const AliDielectronCFdraw &c);
AliDielectronCFdraw &operator=(const AliDielectronCFdraw &c);
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// Dielectron DebugTree //
+// //
+// //
+/*
+register variables from the variable manager. The output will be written
+to a tree
+
+NOTE: Please use with extream care! Only for debugging and test purposes!!!
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include <TFile.h>
+#include <TTreeStream.h>
+#include <TString.h>
+
+#include <AliAnalysisManager.h>
+
+#include "AliDielectronPair.h"
+
+#include "AliDielectronDebugTree.h"
+
+ClassImp(AliDielectronDebugTree)
+
+AliDielectronDebugTree::AliDielectronDebugTree() :
+ TNamed(),
+ fFileName("jpsi_debug.root"),
+ fNVars(0),
+ fNVarsLeg(0),
+ fStreamer(0x0)
+{
+ //
+ // Default Constructor
+ //
+ for (Int_t i=0; i<AliDielectronVarManager::kNMaxValues;++i){
+ fVariables[i]=0;
+ fVariablesLeg[i]=0;
+ }
+}
+
+//______________________________________________
+AliDielectronDebugTree::AliDielectronDebugTree(const char* name, const char* title) :
+ TNamed(name, title),
+ fFileName("jpsi_debug.root"),
+ fNVars(0),
+ fNVarsLeg(0),
+ fStreamer(0x0)
+{
+ //
+ // Named Constructor
+ //
+ for (Int_t i=0; i<AliDielectronVarManager::kNMaxValues;++i){
+ fVariables[i]=0;
+ fVariablesLeg[i]=0;
+ }
+}
+
+//______________________________________________
+AliDielectronDebugTree::~AliDielectronDebugTree()
+{
+ //
+ // Default Destructor
+ //
+ if (fStreamer){
+ fStreamer->GetFile()->Write();
+ delete fStreamer;
+ }
+}
+
+//______________________________________________
+void AliDielectronDebugTree::Fill(AliDielectronPair *pair)
+{
+ //
+ // Fill configured variables to the tree
+ //
+
+ //is there anything to fill
+ if (fNVars==0&&fNVarsLeg==0) return;
+
+ //only in local mode!!!
+ AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
+ if (man && man->GetAnalysisType()!=AliAnalysisManager::kLocalAnalysis) return;
+
+ if (!fStreamer) fStreamer=new TTreeSRedirector(fFileName.Data());
+
+ Int_t var=0;
+ Double_t values[AliDielectronVarManager::kNMaxValues];
+ Double_t valuesLeg1[AliDielectronVarManager::kNMaxValues];
+ Double_t valuesLeg2[AliDielectronVarManager::kNMaxValues];
+// fill pair values
+ if (fNVars>0){
+ AliDielectronVarManager::Fill(pair,values);
+
+ for (Int_t i=0; i<fNVars; ++i){
+ var=fVariables[i];
+ (*fStreamer) << "Pair"
+ << Form("%s=",AliDielectronVarManager::GetValueName(var))
+ << values[var];
+ }
+ }
+
+ if (fNVarsLeg>0){
+ //leg1
+ AliDielectronVarManager::Fill(pair->GetFirstDaughter(),valuesLeg1);
+ //leg2
+ AliDielectronVarManager::Fill(pair->GetSecondDaughter(),valuesLeg2);
+
+ for (Int_t i=0; i<fNVarsLeg; ++i){
+ var=fVariablesLeg[i];
+ (*fStreamer) << "Pair"
+ << Form("Leg1_%s=",AliDielectronVarManager::GetValueName(var))
+ << valuesLeg1[var]
+ << Form("Leg2_%s=",AliDielectronVarManager::GetValueName(var))
+ << valuesLeg2[var];
+ }
+
+ }
+
+ (*fStreamer) << "Pair" << "\n";
+
+
+}
+
+//______________________________________________
+void AliDielectronDebugTree::DeleteStreamer()
+{
+ //
+ // delete the streamer
+ //
+ if (!fStreamer) return;
+ delete fStreamer;
+ fStreamer=0x0;
+
+}
+
+//______________________________________________
+void AliDielectronDebugTree::WriteTree()
+{
+ //
+ // Write file
+ //
+ if (!fStreamer) return;
+ fStreamer->GetFile()->Write();
+}
--- /dev/null
+#ifndef ALIDIELECTRONDEBUGTREE_H
+#define ALIDIELECTRONDEBUGTREE_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronDebugTree #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <TNamed.h>
+#include <TString.h>
+
+#include "AliDielectronVarManager.h"
+
+class TTreeSRedirector;
+class AliDielectronPair;
+
+class AliDielectronDebugTree : public TNamed {
+public:
+ AliDielectronDebugTree();
+ AliDielectronDebugTree(const char*name, const char* title);
+
+ virtual ~AliDielectronDebugTree();
+
+ void SetOutputFileName(const char* file) { fFileName=file; }
+
+ void AddPairVariable(AliDielectronVarManager::ValueTypes type) { fVariables[fNVars++]=(Int_t)type; }
+ void AddLegVariable(AliDielectronVarManager::ValueTypes type) { fVariablesLeg[fNVarsLeg++]=(Int_t)type; }
+
+ void Fill(AliDielectronPair *pair);
+
+ void DeleteStreamer();
+ void WriteTree();
+private:
+ TString fFileName; //output file name
+
+ Int_t fNVars; //number of configured variables
+ Int_t fVariables[AliDielectronVarManager::kNMaxValues]; //configured variables
+ Int_t fNVarsLeg; //number of configured variables
+ Int_t fVariablesLeg[AliDielectronVarManager::kNMaxValues]; //configured variables for the legs
+
+ TTreeSRedirector *fStreamer; //! Tree Redirector
+
+ AliDielectronDebugTree(const AliDielectronDebugTree &c);
+ AliDielectronDebugTree &operator=(const AliDielectronDebugTree &c);
+
+
+ ClassDef(AliDielectronDebugTree,1) // Dielectron DebugTree
+};
+
+
+
+#endif
#include <TROOT.h>
#include <TLegend.h>
#include <TKey.h>
+#include <TAxis.h>
+#include <TVirtualPS.h>
// #include <TVectorD.h>
#include "AliDielectronHistos.h"
//_____________________________________________________________________________
void AliDielectronHistos::UserHistogram(const char* histClass,const char *name, const char* title,
Int_t nbinsX, Double_t xmin, Double_t xmax,
- UInt_t valTypeX)
+ UInt_t valTypeX, Bool_t logBinX)
{
//
// Default histogram creation 1D case
//
if (!IsHistogramOk(histClass,name)) return;
+
+ Double_t *binLimX=0x0;
+
+ if (logBinX) {
+ binLimX=MakeLogBinning(nbinsX, xmin, xmax);
+ } else {
+ binLimX=MakeLinBinning(nbinsX, xmin, xmax);
+ }
+
+ TH1* hist=new TH1F(name,title,nbinsX,binLimX);
+
+ delete [] binLimX;
- TH1* hist=new TH1F(name,title,nbinsX,xmin,xmax);
Bool_t isReserved=fReservedWords->Contains(histClass);
if (isReserved)
UserHistogramReservedWords(histClass, hist, valTypeX);
void AliDielectronHistos::UserHistogram(const char* histClass,const char *name, const char* title,
Int_t nbinsX, Double_t xmin, Double_t xmax,
Int_t nbinsY, Double_t ymin, Double_t ymax,
- UInt_t valTypeX, UInt_t valTypeY)
+ UInt_t valTypeX, UInt_t valTypeY,
+ Bool_t logBinX, Bool_t logBinY)
{
//
// Default histogram creation 2D case
//
if (!IsHistogramOk(histClass,name)) return;
- TH1* hist=new TH2F(name,title,nbinsX,xmin,xmax,nbinsY,ymin,ymax);
+ Double_t *binLimX=0x0;
+ Double_t *binLimY=0x0;
+
+ if (logBinX) {
+ binLimX=MakeLogBinning(nbinsX, xmin, xmax);
+ } else {
+ binLimX=MakeLinBinning(nbinsX, xmin, xmax);
+ }
+ if (logBinY) {
+ binLimY=MakeLogBinning(nbinsY, ymin, ymax);
+ } else {
+ binLimY=MakeLinBinning(nbinsY, ymin, ymax);
+ }
+
+ TH1* hist=new TH2F(name,title,nbinsX,binLimX,nbinsY,binLimY);
+
+ delete [] binLimX;
+ delete [] binLimY;
+
+
Bool_t isReserved=fReservedWords->Contains(histClass);
if (isReserved)
UserHistogramReservedWords(histClass, hist, valTypeX+100*valTypeY);
Int_t nbinsX, Double_t xmin, Double_t xmax,
Int_t nbinsY, Double_t ymin, Double_t ymax,
Int_t nbinsZ, Double_t zmin, Double_t zmax,
- UInt_t valTypeX, UInt_t valTypeY, UInt_t valTypeZ)
+ UInt_t valTypeX, UInt_t valTypeY, UInt_t valTypeZ,
+ Bool_t logBinX, Bool_t logBinY, Bool_t logBinZ)
{
//
// Default histogram creation 3D case
//
if (!IsHistogramOk(histClass,name)) return;
- TH1* hist=new TH3F(name,title,nbinsX,xmin,xmax,nbinsY,ymin,ymax,nbinsZ,zmin,zmax);
+ Double_t *binLimX=0x0;
+ Double_t *binLimY=0x0;
+ Double_t *binLimZ=0x0;
+
+ if (logBinX) {
+ binLimX=MakeLogBinning(nbinsX, xmin, xmax);
+ } else {
+ binLimX=MakeLinBinning(nbinsX, xmin, xmax);
+ }
+
+ if (logBinY) {
+ binLimY=MakeLogBinning(nbinsY, ymin, ymax);
+ } else {
+ binLimY=MakeLinBinning(nbinsY, ymin, ymax);
+ }
+
+ if (logBinZ) {
+ binLimZ=MakeLogBinning(nbinsZ, zmin, zmax);
+ } else {
+ binLimZ=MakeLinBinning(nbinsZ, zmin, zmax);
+ }
+
+ TH1* hist=new TH3F(name,title,nbinsX,binLimX,nbinsY,binLimY,nbinsZ,binLimZ);
+
+ delete [] binLimX;
+ delete [] binLimY;
+ delete [] binLimZ;
+
Bool_t isReserved=fReservedWords->Contains(histClass);
if (isReserved)
UserHistogramReservedWords(histClass, hist, valTypeX+100*valTypeY+10000*valTypeZ);
drawStr.ToLower();
//options
// Bool_t same=drawOpt.Contains("same"); //FIXME not yet implemented
+
+ TCanvas *c=0x0;
+ if (gVirtualPS) {
+ if (!gPad){
+ Error("Draw","When writing to a file you have to create a canvas before opening the file!!!");
+ return;
+ }
+ c=gPad->GetCanvas();
+// c=new TCanvas;
+ }
TIter nextClass(&fHistoList);
THashList *classTable=0;
+ Bool_t first=kTRUE;
while ( (classTable=(THashList*)nextClass()) ){
//optimised division
Int_t nPads = classTable->GetEntries();
Int_t nRows = (Int_t)TMath::Ceil( (Double_t)nPads/(Double_t)nCols );
//create canvas
- TString canvasName;
- canvasName.Form("c%s",classTable->GetName());
- TCanvas *c=(TCanvas*)gROOT->FindObject(canvasName.Data());
- if (!c) c=new TCanvas(canvasName.Data(),classTable->GetName());
- c->Clear();
- c->Divide(nCols,nRows);
-
+ if (!gVirtualPS){
+ TString canvasName;
+ canvasName.Form("c%s_%s",GetName(),classTable->GetName());
+ c=(TCanvas*)gROOT->FindObject(canvasName.Data());
+ if (!c) c=new TCanvas(canvasName.Data(),Form("%s: %s",GetName(),classTable->GetName()));
+ c->Clear();
+ } else {
+ if (first){
+ first=kFALSE;
+ } else {
+ c->Clear();
+ }
+ }
+ if (nCols>1||nRows>1) c->Divide(nCols,nRows);
+
//loop over histograms and draw them
TIter nextHist(classTable);
Int_t iPad=0;
while ( (h=(TH1*)nextHist()) ){
TString drawOpt;
if ( (h->InheritsFrom(TH2::Class())) ) drawOpt="colz";
- c->cd(++iPad);
+ if (nCols>1||nRows>1) c->cd(++iPad);
+ if ( TMath::Abs(h->GetXaxis()->GetBinWidth(1)-h->GetXaxis()->GetBinWidth(2))>1e-10 ) gPad->SetLogx();
+ if ( TMath::Abs(h->GetYaxis()->GetBinWidth(1)-h->GetYaxis()->GetBinWidth(2))>1e-10 ) gPad->SetLogy();
+ if ( TMath::Abs(h->GetZaxis()->GetBinWidth(1)-h->GetZaxis()->GetBinWidth(2))>1e-10 ) gPad->SetLogz();
h->Draw(drawOpt.Data());
}
+ if (gVirtualPS) c->Update();
}
+// if (gVirtualPS) delete c;
}
//_____________________________________________________________________________
//
// set histogram classes and histograms to this instance. It will take onwnership!
//
+ ResetHistogramList();
+ SetName(list.GetName());
TIter next(&list);
TObject *o;
while ( (o=next()) ){
f.Close();
}
+//_____________________________________________________________________________
void AliDielectronHistos::DrawSame(const char* histName, const Option_t *opt)
{
//
(*fReservedWords)=words;
}
+//_____________________________________________________________________________
+Double_t* AliDielectronHistos::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax) const
+{
+ //
+ // Make logarithmic binning
+ // the user has to delete the array afterwards!!!
+ //
+
+ //check limits
+ if (xmin<1e-20 || xmax<1e-20){
+ Error("MakeLogBinning","For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
+ return MakeLinBinning(nbinsX, xmin, xmax);
+ }
+ Double_t *binLim=new Double_t[nbinsX+1];
+ Double_t first=xmin;
+ Double_t last=xmax;
+ Double_t expMax=TMath::Log(last/first);
+ for (Int_t i=0; i<nbinsX+1; ++i){
+ binLim[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
+ }
+ return binLim;
+}
+
+//_____________________________________________________________________________
+Double_t* AliDielectronHistos::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax) const
+{
+ //
+ // Make logarithmic binning
+ // the user has to delete the array afterwards!!!
+ //
+ Double_t *binLim=new Double_t[nbinsX+1];
+ Double_t first=xmin;
+ Double_t last=xmax;
+ Double_t binWidth=(last-first)/nbinsX;
+ for (Int_t i=0; i<nbinsX+1; ++i){
+ binLim[i]=first+binWidth*(Double_t)i;
+ }
+ return binLim;
+}
+
void UserHistogram(const char* histClass,const char *name, const char* title,
Int_t nbinsX, Double_t xmin, Double_t xmax,
- UInt_t valTypeX=kNoAutoFill);
+ UInt_t valTypeX=kNoAutoFill, Bool_t logBinX=kFALSE);
void UserHistogram(const char* histClass,const char *name, const char* title,
Int_t nbinsX, Double_t xmin, Double_t xmax,
Int_t nbinsY, Double_t ymin, Double_t ymax,
- UInt_t valTypeX=kNoAutoFill, UInt_t valTypeY=0);
+ UInt_t valTypeX=kNoAutoFill, UInt_t valTypeY=0,
+ Bool_t logBinX=kFALSE, Bool_t logBinY=kFALSE);
void UserHistogram(const char* histClass,const char *name, const char* title,
Int_t nbinsX, Double_t xmin, Double_t xmax,
Int_t nbinsY, Double_t ymin, Double_t ymax,
Int_t nbinsZ, Double_t zmin, Double_t zmax,
- UInt_t valTypeX=kNoAutoFill, UInt_t valTypeY=0, UInt_t valTypeZ=0);
+ UInt_t valTypeX=kNoAutoFill, UInt_t valTypeY=0, UInt_t valTypeZ=0,
+ Bool_t logBinX=kFALSE, Bool_t logBinY=kFALSE, Bool_t logBinZ=kFALSE);
void UserHistogram(const char* histClass, TH1* hist, UInt_t valTypes=kNoAutoFill);
TH1* GetHistogram(const char* histClass, const char* name) const;
void SetHistogramList(THashList &list);
+ void ResetHistogramList(){fHistoList.Clear();}
const THashList* GetHistogramList() const {return &fHistoList;}
void AddClass(const char* histClass);
// virtual TObject **GetObjectRef(const TObject *obj) const { return 0; }
// virtual TIterator *MakeIterator(Bool_t dir = kIterForward) const ;
// virtual TObject *Remove(TObject *obj) { return 0; }
-
+
private:
THashList fHistoList; //-> list of histograms
Bool_t IsHistogramOk(const char* classTable, const char* name);
+ Double_t* MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax) const;
+ Double_t* MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax) const;
+
enum {kNoAutoFill=1000000000};
AliDielectronHistos(const AliDielectronHistos &hist);
AnalysisType type=kUNSET;
if (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliESDInputHandler::Class()) type=kESD;
else if (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()->IsA()==AliAODHandler::Class()) type=kAOD;
+
+ AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
fgInstance=new AliDielectronMC(type);
+ fgInstance->SetHasMC(mcHandler!=0x0);
return fgInstance;
}
AliDielectronMC::AliDielectronMC(AnalysisType type):
fMCEvent(0x0),
fStack(0x0),
- fAnaType(type)
+ fAnaType(type),
+ fHasMC(kTRUE)
{
//
// default constructor
if (!mcMother1) return -1;
if (lblMother1!=lblMother2) return -1;
if (TMath::Abs(mcPart1->PdgCode())!=11) return -1;
+ if (mcPart1->PdgCode()!=-mcPart2->PdgCode()) return -1;
if (mcMother1->PdgCode()!=pdgMother) return -1;
return lblMother1;
if (!mcMother1) return -1;
if (lblMother1!=lblMother2) return -1;
if (TMath::Abs(mcPart1->GetPdgCode())!=11) return -1;
+ if (mcPart1->GetPdgCode()!=-mcPart2->GetPdgCode()) return -1;
if (mcMother1->GetPdgCode()!=pdgMother) return -1;
return lblMother1;
AliDielectronMC(AnalysisType type=kUNSET);
virtual ~AliDielectronMC();
+ void SetHasMC(Bool_t hasMC) { fHasMC=hasMC; }
+ Bool_t HasMC() const { return fHasMC; }
+
static AliDielectronMC* Instance();
void Initialize(); // initialization
AliStack *fStack; // MC stack
AnalysisType fAnaType; // Analysis type
+ Bool_t fHasMC; // Do we have an MC handler?
static AliDielectronMC* fgInstance; //! singleton pointer
ClassImp(AliDielectronPair)
AliDielectronPair::AliDielectronPair() :
- fOpeningAngle(-1),
fType(-1),
fLabel(-1),
fPair(),
+ fD1(),
+ fD2(),
fRefD1(),
fRefD2()
{
//______________________________________________
AliDielectronPair::AliDielectronPair(AliVTrack * const particle1, Int_t pid1,
AliVTrack * const particle2, Int_t pid2, Char_t type) :
- fOpeningAngle(-1),
fType(type),
fLabel(-1),
fPair(),
+ fD1(),
+ fD2(),
fRefD1(),
fRefD2()
{
AliVTrack * const particle2, Int_t pid2)
{
//
- // Set the tracks to the pair KF particle
+ // Sort particles by pt, first particle larget Pt
+ // set AliKF daughters and pair
//
fPair.Initialize();
+ fD1.Initialize();
+ fD2.Initialize();
AliKFParticle kf1(*particle1,pid1);
AliKFParticle kf2(*particle2,pid2);
-
+
fPair.AddDaughter(kf1);
fPair.AddDaughter(kf2);
if (particle1->Pt()>particle2->Pt()){
fRefD1 = particle1;
fRefD2 = particle2;
+ fD1+=kf1;
+ fD2+=kf2;
} else {
fRefD1 = particle2;
fRefD2 = particle1;
+ fD1+=kf2;
+ fD2+=kf1;
}
-
- fOpeningAngle=kf1.GetAngle(kf2);
}
virtual const Double_t *PID() const { return 0;} //TODO: check
// Dummy
Int_t PdgCode() const {return 0;}
-
- Double_t OpeningAngle() const { return fOpeningAngle; }
+
UChar_t GetType() const { return fType; }
void SetType(Char_t type) { fType=type; }
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); }
+
// internal KF particle
- const AliKFParticle& GetKFParticle() const { return fPair; }
-
+ const AliKFParticle& GetKFParticle() const { return fPair; }
+ const AliKFParticle& GetKFFirstDaughter() const { return fD1; }
+ const AliKFParticle& GetKFSecondDaughter() const { return fD2; }
+
// daughter references
AliVParticle* GetFirstDaughter() const { return dynamic_cast<AliVParticle*>(fRefD1.GetObject()); }
AliVParticle* GetSecondDaughter() const { return dynamic_cast<AliVParticle*>(fRefD2.GetObject()); }
private:
- Double_t fOpeningAngle; // opening angle of the pair
Char_t fType; // type of the pair e.g. like sign SE, unlike sign SE, ... see AliDielectron
Int_t fLabel; // MC label
AliKFParticle fPair; // KF particle internally used for pair calculation
-
+ AliKFParticle fD1; // KF particle first daughter
+ AliKFParticle fD2; // KF particle1 second daughter
+
TRef fRefD1; // Reference to first daughter
TRef fRefD2; // Reference to second daughter
- ClassDef(AliDielectronPair,2)
+ ClassDef(AliDielectronPair,3)
};
#endif
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// Dielectron SignalBase //
+// //
+// //
+/*
+Base class for signal extraction from a histogram or an array of histograms
+The histogram is assumed to be an inv. mass spectrum,
+the array of histograms is assumed to be an array with inv. mass histograms
+resulting from single and mixed events, as defined in AliDielectron.cxx
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+
+#include "TPaveText.h"
+#include "AliDielectronSignalBase.h"
+
+ClassImp(AliDielectronSignalBase)
+
+AliDielectronSignalBase::AliDielectronSignalBase() :
+ TNamed(),
+ fValues(4),
+ fErrors(4),
+ fIntMin(2.99),
+ fIntMax(3.19)
+{
+ //
+ // Default Constructor
+ //
+
+}
+
+//______________________________________________
+AliDielectronSignalBase::AliDielectronSignalBase(const char* name, const char* title) :
+ TNamed(name, title),
+ fValues(4),
+ fErrors(4),
+ fIntMin(2.99),
+ fIntMax(3.19)
+{
+ //
+ // Named Constructor
+ //
+}
+
+//______________________________________________
+AliDielectronSignalBase::~AliDielectronSignalBase()
+{
+ //
+ // Default Destructor
+ //
+
+}
+
+//______________________________________________
+TPaveText* AliDielectronSignalBase::DrawStats(Double_t x1/*=0.*/, Double_t y1/*=0.*/, Double_t x2/*=0.*/, Double_t y2/*=0.*/)
+{
+ //
+ // Draw extracted values in a TPaveText
+ // with the corners x1,y2,x2,y2
+ //
+ if (TMath::Abs(x1)<1e-20&&TMath::Abs(x2)<1e-20){
+ x1=.6;
+ x2=.9;
+ y1=.7;
+ y2=.9;
+ }
+ TPaveText *t=new TPaveText(x1,y1,x2,y2,"brNDC");
+ t->SetFillColor(kWhite);
+ t->SetBorderSize(1);
+ t->SetTextAlign(12);
+ t->AddText(Form("Singal : %.2f #pm %.2f",GetSignal(),GetSignalError()));
+ t->AddText(Form("Backgnd: %.2f #pm %.2f",GetBackground(),GetBackgroundError()));
+ t->AddText(Form("Signif.: %.2f #pm %.2f",GetSignificance(),GetSignificanceError()));
+ t->AddText(Form("SoB : %.2f #pm %.2f",GetSignalOverBackground(),GetSignalOverBackgroundError()));
+ t->Draw();
+
+ return t;
+}
+
--- /dev/null
+#ifndef ALIDIELECTRONSIGNALBASE_H
+#define ALIDIELECTRONSIGNALBASE_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronSignalBase #
+//# Manage Cuts on the legs of the pair #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <TNamed.h>
+#include <TVectorT.h>
+#include <TMath.h>
+
+class TObjArray;
+class TPaveText;
+class TH1;
+
+class AliDielectronSignalBase : public TNamed {
+public:
+ AliDielectronSignalBase();
+ AliDielectronSignalBase(const char*name, const char* title);
+
+ virtual ~AliDielectronSignalBase();
+
+
+ void SetIntegralRange(Double_t min, Double_t max) {fIntMin=min;fIntMax=max;}
+
+ void SetSignal(Double_t val, Double_t valErr) {fValues(0)=val; fErrors(0)=valErr;}
+ void SetBackground(Double_t val, Double_t valErr) {fValues(1)=val; fErrors(1)=valErr;}
+ void SetSignificance(Double_t val, Double_t valErr) {fValues(2)=val; fErrors(2)=valErr;}
+ void SetSignalOverBackground(Double_t val, Double_t valErr) {fValues(3)=val; fErrors(3)=valErr;}
+
+ const TVectorD& GetValues() const {return fValues;}
+ const TVectorD& GetErrors() const {return fErrors;}
+
+ Double_t GetIntegralMin() const { return fIntMin; }
+ Double_t GetIntegralMax() const { return fIntMax; }
+
+ Double_t GetSignal() {return fValues(0);}
+ Double_t GetBackground() {return fValues(1);}
+ Double_t GetSignificance() {return fValues(2);}
+ Double_t GetSignalOverBackground() {return fValues(3);}
+
+ Double_t GetSignalError() {return fErrors(0);}
+ Double_t GetBackgroundError() {return fErrors(1);}
+ Double_t GetSignificanceError() {return fErrors(2);}
+ Double_t GetSignalOverBackgroundError() {return fErrors(3);}
+
+ void GetSignal(Double_t &val, Double_t &valErr) {val=fValues(0); valErr=fErrors(0);}
+ void GetBackground(Double_t &val, Double_t &valErr) {val=fValues(1); valErr=fErrors(1);}
+ void GetSignificance(Double_t &val, Double_t &valErr) {val=fValues(2); valErr=fErrors(2);}
+ void GetSignalOverBackground(Double_t &val, Double_t &valErr) {val=fValues(3); valErr=fErrors(3);}
+
+ /**
+ This function needs to be implemented by the signal extraction classes.
+
+ The signal extraction is done on the mass spectra.
+ The TObjArray should contain the Inv. Mass spectra of the 10 possible combinations
+ for single and mixed events defined in AliDielectron.cxx
+ In the referece TVectorDs the values and value errors of the result should be stored:
+ Parameter - 0: Signal entries
+ 1: Background entries
+ 2: Significance ( S/sqr(S+B) )
+ 3: Signal/Background
+
+ It is enough to calculate the signal and background and then call
+ SetSignificanceAndSOB(TVectorD &values, TVectorD &errors)
+ Please also read the description there!!!
+ */
+ virtual void Process(TObjArray * const /*arrhist*/) = 0;
+
+protected:
+
+ void SetSignificanceAndSOB();
+ TPaveText* DrawStats(Double_t x1=0., Double_t y1=0., Double_t x2=0., Double_t y2=0.);
+ void Reset() {fValues.Zero(); fErrors.Zero();}
+
+private:
+ TVectorD fValues; // values
+ TVectorD fErrors; // value errors
+
+ Double_t fIntMin; // signal extraction range min
+ Double_t fIntMax; // signal extraction range max
+
+ AliDielectronSignalBase(const AliDielectronSignalBase &c);
+ AliDielectronSignalBase &operator=(const AliDielectronSignalBase &c);
+
+ ClassDef(AliDielectronSignalBase,2) // Dielectron SignalBase
+};
+
+//
+// Inline functions
+//
+
+inline void AliDielectronSignalBase::SetSignificanceAndSOB()
+{
+ //
+ // calculate significance and signal over background from values
+ // it is expected that:
+ // - 'values' and 'errors' have the size 4
+ // - Parameters 0 are given: signal and signal error
+ // - Parameters 1 are given: background and background error
+ // - Optionally parameter 2 can be given: signal+background and its error
+ //
+ // The calculated significance and S/B and the corresponding errors
+ // are written in parameters 2 and 3, the first two parameters (signal and background)
+ // stay untouched
+
+ Double_t signal=fValues(0), signal_err=fErrors(0);
+ Double_t background=fValues(1), background_err=fErrors(1);
+ Double_t sigPlusBack=fValues(2), sigPlusBack_err=fErrors(2);
+ Double_t significance=0, significance_err=0;
+ Double_t sob=0, sob_err=0;
+
+ if (sigPlusBack<1e-20){
+ //calculate signal plus background
+ sigPlusBack=signal+background;
+ sigPlusBack_err=TMath::Sqrt( signal_err*signal_err + background_err*background_err );
+ }
+
+ if (sigPlusBack>1e-30){
+
+ significance=signal/TMath::Sqrt(sigPlusBack);
+
+ significance_err=TMath::Sqrt((signal_err/signal)*(signal_err/signal)+
+ (0.5*sigPlusBack_err/sigPlusBack)*(0.5*sigPlusBack_err/sigPlusBack)
+ )*significance;
+
+ }
+
+ if (background>1e-30){
+ sob=signal/background;
+ sob_err=TMath::Sqrt((signal_err/signal)*(signal_err/signal)+
+ (background_err/background)*(background_err/background))*sob;
+ }
+
+ fValues(2)=significance;
+ fErrors(2)=significance_err;
+
+ fValues(3)=sob;
+ fErrors(3)=sob_err;
+}
+
+#endif
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// //
+// Dielectron SignalExt //
+// //
+/*
+Dielectron signal extraction class using functions as input.
+
+A function to describe the signal as well as one to describe the background
+has to be deployed by the user. Alternatively on of the default implementaions
+can be used.
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include <TF1.h>
+#include <TH1.h>
+#include <TCanvas.h>
+#include <TMath.h>
+#include <TString.h>
+
+#include <AliLog.h>
+
+#include "AliDielectronSignalExt.h"
+
+ClassImp(AliDielectronSignalExt)
+
+AliDielectronSignalExt::AliDielectronSignalExt() :
+ AliDielectronSignalBase(),
+ fSignPM(0x0),
+ fSignPP(0x0),
+ fSignMM(0x0),
+ fBackground(0x0),
+ fSignal(0x0),
+ fMethod(1),
+ fRebin(1),
+ fBins(0),
+ fDrawMin(0.),
+ fDrawMax(0.),
+ fFitMin(2.5),
+ fFitMax(4)
+{
+ //
+ // Default Constructor
+ //
+
+}
+
+//______________________________________________
+AliDielectronSignalExt::AliDielectronSignalExt(const char* name, const char* title) :
+ AliDielectronSignalBase(name, title),
+ fSignPM(0x0),
+ fSignPP(0x0),
+ fSignMM(0x0),
+ fBackground(0x0),
+ fSignal(0x0),
+ fMethod(1),
+ fRebin(1),
+ fBins(0),
+ fDrawMin(0.),
+ fDrawMax(0.),
+ fFitMin(2.5),
+ fFitMax(4)
+{
+ //
+ // Named Constructor
+ //
+}
+
+//______________________________________________
+AliDielectronSignalExt::~AliDielectronSignalExt()
+{
+ //
+ // Default Destructor
+ //
+ if (fSignPM) delete fSignPM;
+ if (fSignPP) delete fSignPP;
+ if (fSignMM) delete fSignMM;
+ if (fBackground) delete fBackground;
+ if (fSignal) delete fSignal;
+}
+
+//______________________________________________
+void AliDielectronSignalExt::SetHistograms(TH1F* const unlike, TH1F* const backg, TH1F* const signal)
+{
+ //
+ // set histograms
+ //
+ fSignPM = (TH1F*)unlike->Clone("fSignPM");
+ fBackground = backg;
+ fSignal = signal;
+}
+
+//______________________________________________
+void AliDielectronSignalExt::Process(TObjArray* const arrhist)
+{
+ //
+ // signal subtraction. support like-sign subtraction and event mixing method
+ //
+ switch ( fMethod ){
+ case 1 :
+ ProcessLS(arrhist); // process like-sign subtraction method
+ break;
+
+ case 2 :
+ ProcessEM(arrhist); // process event mixing method
+ break;
+
+ default :
+ AliWarning("Subtraction method not supported. Please check SetMethod() function.");
+ }
+}
+
+//______________________________________________
+void AliDielectronSignalExt::ProcessLS(TObjArray* const arrhist)
+{
+ //
+ // signal subtraction
+ //
+ fSignPP = (TH1F*)arrhist->At(0); // like sign : plus-plus
+ fSignPM = (TH1F*)arrhist->At(1); // unlike-sign : plus-minus
+ fSignMM = (TH1F*)arrhist->At(2); // like sign : minus-minus
+ fSignPP->Sumw2();
+ fSignPM->Sumw2();
+ fSignMM->Sumw2();
+
+ if ( fRebin>1 ){ Rebin(fRebin); } // rebinning of histogram
+
+ fBins = fSignPM->GetNbinsX(); // number of bins
+ Double_t minX = fSignPM->GetBinLowEdge(1); // minimum X value in axis
+ Double_t maxX = fSignPM->GetBinLowEdge(fBins+1); // maximum X value in axis
+
+ AliDebug(10,Form("histogram #bin = %d , min = %f , max = %f\n",fBins,minX,maxX));
+ TH1F* hBackground = new TH1F("hBackground","Like-sign background",fBins,minX,maxX);
+ TH1F* hSubtracted = new TH1F("hSubtracted","Background subtracted",fBins,minX,maxX);
+
+ // fill out background and subtracted histogram
+ for ( Int_t ibin=1; ibin<fBins+1; ibin++ ){
+ Double_t mass = ibin*(maxX-minX)/(Double_t)fBins;
+ Double_t backgr = 2.*sqrt( fSignPP->GetBinContent(ibin) * fSignMM->GetBinContent(ibin) );
+ Double_t signal = fSignPM->GetBinContent(ibin) - backgr;
+
+ hBackground->Fill(mass, backgr);
+ hSubtracted->Fill(mass, signal);
+ }
+ SetHistograms(fSignPM, hBackground, hSubtracted);
+
+
+ Double_t signal=0, signal_err=0;
+ Double_t background=0, background_err=0;
+
+ signal = fSignal->IntegralAndError(fSignal->FindBin(GetIntegralMin()),
+ fSignal->FindBin(GetIntegralMax()), signal_err);
+ background = fBackground->IntegralAndError(fBackground->FindBin(GetIntegralMin()),
+ fBackground->FindBin(GetIntegralMax()), background_err);
+
+ //reset result arrays
+ Reset();
+ //set values
+ SetSignal(signal,signal_err);
+ SetBackground(background,background_err);
+ SetSignificanceAndSOB();
+
+ // cleanup
+ //delete hBackground;
+ //delete hSubtracted;
+}
+
+//______________________________________________
+void AliDielectronSignalExt::ProcessEM(TObjArray* const arrhist)
+{
+ //
+ // event mixing. not implemented yet.
+ //
+ printf("event mixing method is not implemented yet. Like-sign method will be used.\n");
+ ProcessLS(arrhist);
+}
+
+//______________________________________________
+void AliDielectronSignalExt::Rebin(Int_t rebin)
+{
+ //
+ // rebinning of histograms
+ //
+ fSignPM->Rebin(rebin);
+ fSignPP->Rebin(rebin);
+ fSignMM->Rebin(rebin);
+}
+
+//______________________________________________
+void AliDielectronSignalExt::Draw(const Option_t* option)
+{
+ //
+ // Draw the fitted function
+ //
+ TString drawOpt(option);
+ drawOpt.ToLower();
+
+ TCanvas *cSub = new TCanvas("cSub","signal, background subtracted",1400,1000);
+ cSub->Divide(2,2);
+ cSub->Draw();
+
+ Double_t minX = fSignPM->GetBinLowEdge(1); // minimum X value in axis
+ Double_t maxX = fSignPM->GetBinLowEdge(fBins+1); // maximum X value in axis
+ if ( TMath::Abs(fDrawMin)<1.e-30 ) fDrawMin = minX;
+ if ( TMath::Abs(fDrawMax)<1.e-30 ) fDrawMax = maxX;
+
+ cSub->cd(4);
+ fSignPM->GetXaxis()->SetRangeUser(fDrawMin, fDrawMax);
+ fSignPM->SetLineColor(1);
+ fSignPM->SetLineWidth(2);
+ fSignPM->SetMarkerStyle(6);
+ fSignPM->DrawCopy("P");
+
+ fBackground->SetLineColor(4);
+ fBackground->SetMarkerColor(4);
+ fBackground->SetMarkerStyle(6);
+ fBackground->DrawCopy("Psame");
+
+ fSignal->SetMarkerStyle(20);
+ fSignal->SetMarkerColor(2);
+ fSignal->DrawCopy("Psame");
+
+ cSub->cd(1);
+ fSignal->DrawCopy("P");
+
+ cSub->cd(2);
+ fSignPM->DrawCopy("P");
+
+ cSub->cd(3);
+ fSignPP->DrawCopy("P");
+}
+
--- /dev/null
+#ifndef ALIDIELECTRONSIGNALEXT_H
+#define ALIDIELECTRONSIGNALEXT_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronSignalExt #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <TVectorT.h>
+#include <TString.h>
+#include <TH1.h>
+
+#include "AliDielectronSignalBase.h"
+
+class AliDielectronSignalExt : public AliDielectronSignalBase {
+
+public:
+ AliDielectronSignalExt();
+ AliDielectronSignalExt(const char*name, const char* title);
+
+ virtual ~AliDielectronSignalExt();
+
+ virtual void Process(TObjArray* const arrhist);
+ void ProcessLS(TObjArray* const arrhist); // like-sign method
+ void ProcessEM(TObjArray* const arrhist); // event mixing method
+
+ void SetMethod(Int_t method){ fMethod=method; }
+ void SetHistograms(TH1F* const unlike, TH1F* const backg, TH1F* const signal);
+ void SetRebin(Int_t rebin) {fRebin = rebin;}
+ void SetDrawRange(Double_t min, Double_t max){fDrawMin=min; fDrawMax=max;}
+ void SetFitRange(Double_t min, Double_t max) {fFitMin=min;fFitMax=max;}
+
+ Double_t GetFitMin() const { return fFitMin; }
+ Double_t GetFitMax() const { return fFitMax; }
+ void Rebin(Int_t rebin);
+
+ virtual void Draw(const Option_t* option = "");
+
+
+private:
+
+ TH1F *fSignPM; // histogram of unlike sign (plus-minus)
+ TH1F *fSignPP; // histogram of like sign (plus-plus)
+ TH1F *fSignMM; // histogram of like sign (minus-minus)
+ TH1F *fBackground; // histogram of like-sign background
+ TH1F *fSignal; // histogram of subtracted signal
+
+ Int_t fMethod; // subtraction method. 1(like-sign), 2(event mixing)
+ Int_t fRebin; // number of histogram rebin iteration
+ Int_t fBins; // number of bins in X axis
+ Double_t fDrawMin; // minimum X when drawing
+ Double_t fDrawMax; // maximum X when drawing
+ Double_t fFitMin; // fit range min
+ Double_t fFitMax; // fit range max
+
+ AliDielectronSignalExt(const AliDielectronSignalExt &c);
+ AliDielectronSignalExt &operator=(const AliDielectronSignalExt &c);
+
+ ClassDef(AliDielectronSignalExt,1) // Dielectron SignalFunc
+};
+
+
+
+#endif
--- /dev/null
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// Dielectron SignalFunc //
+// //
+// //
+/*
+Dielectron signal extraction class using functions as input.
+
+A function to describe the signal as well as one to describe the background
+has to be deployed by the user. Alternatively on of the default implementaions
+can be used.
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include <TF1.h>
+#include <TH1.h>
+#include <TGraph.h>
+#include <TMath.h>
+#include <TString.h>
+#include <TPaveText.h>
+#include <RooRealVar.h>
+#include <RooCBShape.h>
+#include <RooExponential.h>
+
+#include <AliLog.h>
+
+#include "AliDielectronSignalFunc.h"
+
+ClassImp(AliDielectronSignalFunc)
+
+AliDielectronSignalFunc::AliDielectronSignalFunc() :
+ AliDielectronSignalBase(),
+ fSignal(0x0),
+ fBackground(0x0),
+ fSigBack(0x0),
+ fFitFunc(0x0),
+ fVInitParam(0),
+ fFitOpt("MNQE"),
+ fUseIntegral(kFALSE),
+ fFitMin(2.5),
+ fFitMax(4),
+ fParM(-1),
+ fParMres(-1)
+{
+ //
+ // Default Constructor
+ //
+
+}
+
+//______________________________________________
+AliDielectronSignalFunc::AliDielectronSignalFunc(const char* name, const char* title) :
+ AliDielectronSignalBase(name, title),
+ fSignal(0x0),
+ fBackground(0x0),
+ fSigBack(0x0),
+ fFitFunc(0x0),
+ fVInitParam(0),
+ fFitOpt("MNQ"),
+ fUseIntegral(kFALSE),
+ fFitMin(2.5),
+ fFitMax(4),
+ fParM(-1),
+ fParMres(-1)
+{
+ //
+ // Named Constructor
+ //
+}
+
+//______________________________________________
+AliDielectronSignalFunc::~AliDielectronSignalFunc()
+{
+ //
+ // Default Destructor
+ //
+ if (fSigBack) delete fSigBack;
+}
+
+
+//______________________________________________
+void AliDielectronSignalFunc::Process(TH1 * const hist)
+{
+ //
+ // Fit the mass histogram with the function and retrieve the parameters
+ //
+
+ //reset result arrays
+ Reset();
+
+ //sanity check
+ if (!fSigBack) {
+ AliError("Use 'SetFunctions(TF1*,TF1*)' or 'SetDefaults(Int_t)' to setup the signal functions first'!");
+ return;
+ }
+
+ //initial the fit function to its default parameters
+ if (fVInitParam.GetMatrixArray()) fSigBack->SetParameters(fVInitParam.GetMatrixArray());
+
+ //Perform fit and retrieve values
+ hist->Fit(fSigBack,fFitOpt.Data(),"",fFitMin,fFitMax);
+
+ //retrieve values and errors
+ //TODO: perhpas implement different methods to retrieve the valus
+ fSignal->SetParameters(fSigBack->GetParameters());
+ fBackground->SetParameters(fSigBack->GetParameters()+fSignal->GetNpar());
+
+ //TODO: proper error estimate?!?
+ // perhaps this is not possible in a general way?
+
+ // signal and background histograms
+ TH1 *hSignal=(TH1*)hist->Clone("DieleSignalHist");
+ TH1 *hBackground=(TH1*)hist->Clone("DieleBackgroundHist");
+
+ fSignal->SetRange(hist->GetXaxis()->GetXmin(),hist->GetXaxis()->GetXmax());
+ fBackground->SetRange(hist->GetXaxis()->GetXmin(),hist->GetXaxis()->GetXmax());
+
+ hSignal->Add(fBackground,-1);
+ hBackground->Add(fSignal,-1);
+
+
+ //get values and errors signal and background
+ Double_t signal=0, signal_err=0;
+ Double_t background=0, background_err=0;
+ Double_t sigPlusBack=0, sigPlusBack_err=0;
+
+ if (!fUseIntegral){
+ signal=hSignal->IntegralAndError(hSignal->FindBin(GetIntegralMin()),
+ hSignal->FindBin(GetIntegralMax()), signal_err);
+ background=hBackground->IntegralAndError(hBackground->FindBin(GetIntegralMin()),
+ hBackground->FindBin(GetIntegralMax()), background_err);
+ sigPlusBack=hist->IntegralAndError(hist->FindBin(GetIntegralMin()),
+ hist->FindBin(GetIntegralMax()), sigPlusBack_err);
+ } else {
+ signal=fSignal->Integral(GetIntegralMin(),GetIntegralMax())/hist->GetBinWidth(hist->FindBin((GetIntegralMax()-GetIntegralMin())/2.));
+ background=fBackground->Integral(GetIntegralMin(),GetIntegralMax())/hist->GetBinWidth(hist->FindBin((GetIntegralMax()-GetIntegralMin())/2.));
+ sigPlusBack=fSigBack->Integral(GetIntegralMin(),GetIntegralMax())/hist->GetBinWidth(hist->FindBin((GetIntegralMax()-GetIntegralMin())/2.));
+ //TODO: Error calculation
+ }
+
+ //set values
+ SetSignal(signal,signal_err);
+ SetBackground(background,background_err);
+ SetSignificanceAndSOB();
+
+ //cleanup
+ delete hSignal;
+ delete hBackground;
+}
+
+//______________________________________________
+void AliDielectronSignalFunc::Process(TObjArray * const arrhist)
+{
+ //
+ // Fit the mass histogram with the function and retrieve the parameters
+ //
+
+ TH1 *hist=(TH1*)arrhist->At(1);
+ Process(hist);
+}
+//______________________________________________
+void AliDielectronSignalFunc::SetFunctions(TF1 * const sig, TF1 * const back, TF1 * const combined,
+ Int_t parM, Int_t parMres)
+{
+ //
+ // Set the signal, background functions and combined fit function
+ // Note: The process method assumes that the first n parameters in the
+ // combined fit function correspond to the n parameters of the signal function
+ // and the n+1 to n+m parameters to the m parameters of the background function!!!
+ // if parM and parMres are set this information can be used in the drawing function
+ // to also show in the statistics box the mass and mass resolution
+
+ if (!sig||!back||!combined) {
+ AliError("Both, signal and background function need to be set!");
+ }
+ fSignal=sig;
+ fBackground=back;
+ fSigBack=combined;
+
+ InitParams();
+ fParM=parM;
+ fParMres=parMres;
+}
+
+//______________________________________________
+void AliDielectronSignalFunc::InitParams()
+{
+ //
+ // initilise common fit parameters
+ //
+ fVInitParam.ResizeTo(fSigBack->GetNpar());
+ fVInitParam.SetElements(fSigBack->GetParameters());
+}
+
+//______________________________________________
+void AliDielectronSignalFunc::SetDefaults(Int_t type)
+{
+ //
+ // Setup some default functions:
+ // type = 0: gaus signal + linear background in 2.5 - 4 GeV inv. mass
+ // type = 1: gaus signal + exponential background in 2.5 - 4 GeV inv. mass
+ // type = 2: half gaussian, half exponential signal function
+ // type = 3: Crystal-Ball function
+ // type = 4: Crystal-Ball signal + exponential background
+ //
+
+
+ if (type==0){
+ fSignal=new TF1("DieleSignal","gaus",2.5,4);
+ fBackground=new TF1("DieleBackground","pol1",2.5,4);
+ fSigBack=new TF1("DieleCombined","gaus+pol1(3)",2.5,4);
+
+ fSigBack->SetParameters(1,3.1,.05,2.5,1);
+ fSigBack->SetParLimits(0,0,10000000);
+ fSigBack->SetParLimits(1,3.05,3.15);
+ fSigBack->SetParLimits(2,.02,.1);
+
+ SetFunctions(fSignal,fBackground,fSigBack,1,2);
+
+ } else if (type==1){
+ fSignal=new TF1("DieleSignal","gaus",2.5,4);
+ fBackground=new TF1("DieleBackground","[0]*exp(-(x-[1])/[2])",2.5,4);
+ fSigBack=new TF1("DieleCombined","gaus+[3]*exp(-(x-[4])/[5])",2.5,4);
+
+ fSigBack->SetParameters(1,3.1,.05,1,2.5,1);
+ fSigBack->SetParLimits(0,0,10000000);
+ fSigBack->SetParLimits(1,3.05,3.15);
+ fSigBack->SetParLimits(2,.02,.1);
+
+ SetFunctions(fSignal,fBackground,fSigBack,1,2);
+
+ } else if (type==2){
+ // half gaussian, half exponential signal function
+ // exponential background
+ fSignal = new TF1("DieleSignal","(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))",2.5,4);
+ fBackground=new TF1("DieleBackground","[0]*exp(-(x-[1])/[2])+[3]",2.5,4);
+ fSigBack=new TF1("DieleCombined","(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))+[4]*exp(-(x-[5])/[6])+[7]",2.5,4);
+ fSigBack->SetParameters(1.,3.1,.05,.1,1,2.5,1,0);
+
+ fSigBack->SetParLimits(0,0,10000000);
+ fSigBack->SetParLimits(1,3.05,3.15);
+ fSigBack->SetParLimits(2,.02,.1);
+ fSigBack->FixParameter(6,2.5);
+ fSigBack->FixParameter(7,0);
+ SetFunctions(fSignal,fBackground,fSigBack,1,2);
+
+ } else if (type==3){
+ // Crystal-Ball function
+ RooRealVar m("m","Invariant mass",2.5,3.8);
+ RooRealVar alpha("alpha","alpha",0.5,0,10);
+ RooRealVar n("n","n",1,0,10);
+ RooRealVar m0("m0","m0",3.1,1,5);
+ RooRealVar sigma("sigma","sigma",0.3,0,3);
+ RooCBShape jpsi("jpsi","jpsi",m,m0,sigma,alpha,n);
+
+ //RooAddPdf model("model","jpsi fitting function", RooArgList(jpsi,expo), fsig);
+ RooRealVar fsig("fsig","signal fraction",0.5,0.,1.);
+ RooAddPdf model("model","jpsi fitting function", RooArgList(jpsi),fsig);
+
+ model.Clone("fFitFunc");
+
+ } else if (type==4){
+ // Crystal-Ball function
+ RooRealVar m("m","Invariant mass",2.5,3.8);
+ RooRealVar alpha("alpha","alpha",0.5,0,10);
+ RooRealVar n("n","n",1,0,10);
+ RooRealVar m0("m0","m0",3.1,1,5);
+ RooRealVar sigma("sigma","sigma",0.3,0,3);
+ RooCBShape jpsi("jpsi","jpsi",m,m0,sigma,alpha,n);
+
+ // Expoenetial function
+ RooRealVar al("al","al",0.5);
+ RooExponential expo("expo","exponential", m, al);
+
+ // add two functions
+ RooRealVar fsig("fsig","signal fraction",0.5,0.,1.);
+ RooAddPdf model("model","jpsi fitting function", RooArgList(jpsi,expo), fsig);
+
+ model.Clone("fFitFunc");
+ }
+}
+
+
+//______________________________________________
+void AliDielectronSignalFunc::Draw(const Option_t* option)
+{
+ //
+ // Draw the fitted function
+ //
+
+ TString drawOpt(option);
+ drawOpt.ToLower();
+
+ Bool_t optStat=drawOpt.Contains("stat");
+
+ fSigBack->SetNpx(200);
+ fSigBack->SetRange(GetIntegralMin(),GetIntegralMax());
+ fBackground->SetNpx(200);
+ fBackground->SetRange(GetIntegralMin(),GetIntegralMax());
+
+ TGraph *grSig=new TGraph(fSigBack);
+ grSig->SetFillColor(kGreen);
+ grSig->SetFillStyle(3001);
+
+ TGraph *grBack=new TGraph(fBackground);
+ grBack->SetFillColor(kRed);
+ grBack->SetFillStyle(3001);
+
+ grSig->SetPoint(0,grBack->GetX()[0],grBack->GetY()[0]);
+ grSig->SetPoint(grSig->GetN()-1,grBack->GetX()[grBack->GetN()-1],grBack->GetY()[grBack->GetN()-1]);
+
+ grBack->SetPoint(0,grBack->GetX()[0],0.);
+ grBack->SetPoint(grBack->GetN()-1,grBack->GetX()[grBack->GetN()-1],0.);
+
+ fSigBack->SetRange(fFitMin,fFitMax);
+ fBackground->SetRange(fFitMin,fFitMax);
+
+ if (!drawOpt.Contains("same")){
+ grSig->Draw("af");
+ } else {
+ grSig->Draw("f");
+ }
+ grBack->Draw("f");
+ fSigBack->Draw("same");
+ fBackground->Draw("same");
+
+ if (optStat) {
+ TPaveText* pave=DrawStats();
+ if (fParM>=0) pave->AddText(Form("Mass: %.3f #pm %.3f", fSigBack->GetParameter(fParM), fSigBack->GetParError(fParM)));
+ if (fParMres>=0) pave->AddText(Form("Mass res.: %.3f #pm %.3f", fSigBack->GetParameter(fParMres), fSigBack->GetParError(fParMres)));
+ }
+}
+
--- /dev/null
+#ifndef ALIDIELECTRONSIGNALFUNC_H
+#define ALIDIELECTRONSIGNALFUNC_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronSignalFunc #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <TVectorT.h>
+#include <TString.h>
+#include <RooFit.h>
+#include <RooAddPdf.h>
+
+#include "AliDielectronSignalBase.h"
+
+class AliDielectronSignalFunc : public AliDielectronSignalBase {
+public:
+ AliDielectronSignalFunc();
+ AliDielectronSignalFunc(const char*name, const char* title);
+
+ virtual ~AliDielectronSignalFunc();
+
+ void Process(TH1 * const hist);
+ virtual void Process(TObjArray * const arrhist);
+
+ void SetFunctions(TF1 * const sig, TF1 * const back, TF1 * const combined, Int_t parM=-1, Int_t parMres=-1);
+ void InitParams();
+ void SetFitOption(const char* opt) {fFitOpt=opt;}
+ void SetDefaults(Int_t type);
+ void SetFitRange(Double_t min, Double_t max) {fFitMin=min;fFitMax=max;}
+
+ void SetUseIntegral(Bool_t integral=kTRUE) { fUseIntegral=integral; }
+
+ const char* GetFitOption() const { return fFitOpt.Data(); }
+ TF1* GetSignalFunction() const { return fSignal; }
+ TF1* GetBackgroundFunction() const { return fBackground; }
+ TF1* GetCombinedFunction() const { return fSigBack; }
+
+ Double_t GetFitMin() const { return fFitMin; }
+ Double_t GetFitMax() const { return fFitMax; }
+
+ virtual void Draw(const Option_t* option = "");
+
+private:
+ TF1 *fSignal; // Function for the signal description
+ TF1 *fBackground; // Function for the background description
+ TF1 *fSigBack; // Combined function signal plus background
+ RooAddPdf fFitFunc; // RooFit fit function
+ TVectorD fVInitParam; // initial fit parameters
+ TString fFitOpt; // fit option used
+
+ Bool_t fUseIntegral; // use integral instead of bin counts
+
+ Double_t fFitMin; // fit range min
+ Double_t fFitMax; // fit range max
+
+ Int_t fParM; // Paramter which defines the mass
+ Int_t fParMres; // Paramter which defines the resolution of the mass
+
+ AliDielectronSignalFunc(const AliDielectronSignalFunc &c);
+ AliDielectronSignalFunc &operator=(const AliDielectronSignalFunc &c);
+
+ ClassDef(AliDielectronSignalFunc,1) // Dielectron SignalFunc
+};
+
+
+#endif
--- /dev/null
+
+/*************************************************************************
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+* *
+* Author: The ALICE Off-line Project. *
+* Contributors are mentioned in the code where appropriate. *
+* *
+* Permission to use, copy, modify and distribute this software and its *
+* documentation strictly for non-commercial purposes is hereby granted *
+* without fee, provided that the above copyright notice appears in all *
+* copies and that both the copyright notice and this permission notice *
+* appear in the supporting documentation. The authors make no claims *
+* about the suitability of this software for any purpose. It is *
+* provided "as is" without express or implied warranty. *
+**************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+// Add general description //
+// //
+// //
+/*
+Detailed description
+
+
+*/
+// //
+///////////////////////////////////////////////////////////////////////////
+
+#include <TMath.h>
+#include <TVectorT.h>
+#include <TH1.h>
+
+#include <AliCFContainer.h>
+#include <AliCFGridSparse.h>
+
+#include "AliDielectronSignalBase.h"
+#include "AliDielectronSignalFunc.h"
+
+#include "AliDielectronSpectrum.h"
+
+ClassImp(AliDielectronSpectrum)
+
+AliDielectronSpectrum::AliDielectronSpectrum() :
+ TNamed(),
+ fCFSignal(0x0),
+ fCFCorrection(0x0),
+ fCFSpectrum(0x0),
+ fCFCorrMatrix(0x0),
+ fStepSignal(kTRUE),
+ fStepSignificance(kFALSE),
+ fStepSOB(kFALSE),
+ fSignalStep(-1),
+ fCorrNom(-1),
+ fCorrDenom(-1),
+ fSignalMethods(0),
+ fVariables("Pt"),
+ fOwnerSpectrum(kTRUE),
+ fNvars(0),
+ fVars(0x0),
+ fNbins(0x0),
+ fCurrentBins(0x0),
+ fCurrentPositions(0x0)
+{
+ //
+ // Default Constructor
+ //
+ fSignalMethods.SetOwner(kFALSE);
+}
+
+//______________________________________________
+AliDielectronSpectrum::AliDielectronSpectrum(const char* name, const char* title) :
+ TNamed(name, title),
+ fCFSignal(0x0),
+ fCFCorrection(0x0),
+ fCFSpectrum(0x0),
+ fCFCorrMatrix(0x0),
+ fStepSignal(kTRUE),
+ fStepSignificance(kFALSE),
+ fStepSOB(kFALSE),
+ fSignalStep(-1),
+ fCorrNom(-1),
+ fCorrDenom(-1),
+ fSignalMethods(0),
+ fVariables("Pt"),
+ fOwnerSpectrum(kTRUE),
+ fNvars(0),
+ fVars(0x0),
+ fNbins(0x0),
+ fCurrentBins(0x0),
+ fCurrentPositions(0x0)
+{
+ //
+ // Named Constructor
+ //
+ fSignalMethods.SetOwner(kFALSE);
+}
+
+//______________________________________________
+AliDielectronSpectrum::~AliDielectronSpectrum()
+{
+ //
+ // Default Destructor
+ //
+
+ if (fNbins) delete [] fNbins;
+ if (fVars) delete [] fVars;
+ if (fCFCorrMatrix) delete fCFCorrMatrix;
+ if ( fOwnerSpectrum && fCFSpectrum ) delete fCFSpectrum;
+}
+
+//______________________________________________
+void AliDielectronSpectrum::Process()
+{
+ //
+ // Extract signal and perform correction in the specified bins
+ //
+
+ //sanity checks
+ if (!fCFSignal){
+ AliError("Cannot perform signal extraction, no signal container set");
+ return;
+ }
+
+ if (fSignalMethods.GetEntries()==0){
+ AliWarning("No Signal extraction method specified, using a default one");
+ AliDielectronSignalFunc *func=new AliDielectronSignalFunc("gaus+exp","Gauss for Signal and Exponential background");
+ func->SetDefaults(1);
+ fSignalMethods.Add(func);
+ fSignalMethods.SetOwner();
+ }
+
+ //setup configured variables
+ if (!SetupVariables()) return;
+
+ //create container for the spectrum
+ CreateCFSpectrum();
+
+ if (!fCFSpectrum){
+ AliError("Could not create the Spectrum container");
+ return;
+ }
+
+ //get efficiency map if correction container is available
+ if (fCFCorrection&&!fCFCorrMatrix){
+ CreateCorrectionMatrix();
+ }
+
+ //loop over all configured bins and extract the signal
+ fCurrentBins=new Int_t[fNvars];
+ fCurrentPositions=new Double_t[fNvars];
+ ExtractSignalInBins();
+ delete [] fCurrentBins;
+ delete [] fCurrentPositions;
+ if (fSignalMethods.IsOwner()) {
+ fSignalMethods.Delete();
+ fSignalMethods.SetOwner(kFALSE);
+ }
+
+}
+
+//______________________________________________
+Bool_t AliDielectronSpectrum::SetupVariables()
+{
+ //
+ // Setup the variables arrays
+ //
+
+ TObjArray *arr=fVariables.Tokenize(":");
+ fNvars=arr->GetEntries();
+ fVars=new Int_t[fNvars];
+ fNbins=new Int_t[fNvars];
+
+ for (Int_t iVar=0; iVar<fNvars; ++iVar){
+ fVars[iVar]=fCFSignal->GetVar(arr->UncheckedAt(iVar)->GetName());
+ if (fVars[iVar]==-1){
+ AliError(Form("Variable '%s' not found in Signal container!",arr->UncheckedAt(iVar)->GetName()));
+ delete [] fVars;
+ fVars=0x0;
+ delete [] fNbins;
+ fNbins=0x0;
+ delete arr;
+ return kFALSE;
+ }
+
+ fNbins[iVar]=fCFSignal->GetNBins(fVars[iVar]);
+ }
+ delete arr;
+ return kTRUE;
+}
+
+//______________________________________________
+void AliDielectronSpectrum::CreateCFSpectrum()
+{
+ //
+ // Create CF container for the spectrum
+ //
+
+ Int_t nAddStep=0;
+ if (fStepSignal) nAddStep+=2;
+ if (fStepSignificance) ++nAddStep;
+ if (fStepSOB) ++nAddStep;
+
+ Int_t nStep=nAddStep*(fSignalMethods.GetEntries());
+ if (fSignalMethods.GetEntries()>1) nStep+=nAddStep;
+
+ fCFSpectrum = new AliCFContainer(GetName(), GetTitle(), nStep, fNvars, fNbins);
+
+ // initialize the variables and their bin limits
+ for (Int_t iVar=0; iVar<fNvars; iVar++) {
+ fCFSpectrum->SetBinLimits(iVar, fCFSignal->GetBinLimits(fVars[iVar]));
+ fCFSpectrum->SetVarTitle(iVar, fCFSignal->GetVarTitle(fVars[iVar]));
+ }
+
+ // setup step titles
+ Int_t steps=0;
+ for (Int_t iMethod=0; iMethod<fSignalMethods.GetEntries(); ++iMethod){
+ TString name(fSignalMethods.UncheckedAt(iMethod)->GetName());
+ if (fStepSignal){
+ fCFSpectrum->SetStepTitle(steps++,(name+" (Siganl)").Data());
+ fCFSpectrum->SetStepTitle(steps++,(name+" (Corrected Siganl)").Data());
+ }
+ if (fStepSignificance){
+ fCFSpectrum->SetStepTitle(steps++,(name+" (Significance)").Data());
+ }
+ if (fStepSOB){
+ fCFSpectrum->SetStepTitle(steps++,(name+" (S/B)").Data());
+ }
+ }
+
+ if (fSignalMethods.GetEntries()>1){
+ fCFSpectrum->SetStepTitle(steps++,"Mean of methods");
+ fCFSpectrum->SetStepTitle(steps++,"Mean of methods (Corrected)");
+ }
+
+ if (nStep!=steps){
+ AliError("Something went wrong in the step creation");
+ delete fCFSpectrum;
+ fCFSpectrum=0x0;
+ }
+}
+
+//______________________________________________
+void AliDielectronSpectrum::CreateCorrectionMatrix()
+{
+ //
+ // Get the correction matrix for the corresponding variables
+ //
+
+ if (!fCFCorrection) return;
+
+ TObjArray *arr=fVariables.Tokenize(":");
+ Int_t nvars=arr->GetEntries();
+ Int_t *vars=new Int_t[nvars];
+
+ for (Int_t iVar=0; iVar<fNvars; ++iVar){
+ vars[iVar]=fCFCorrection->GetVar(arr->UncheckedAt(iVar)->GetName());
+ if (vars[iVar]==-1){
+ AliError(Form("Variable '%s' not found in Correction container!",arr->UncheckedAt(iVar)->GetName()));
+ delete [] vars;
+ delete arr;
+ return;
+ }
+ }
+ delete arr;
+
+ fCFCorrMatrix =fCFCorrection->GetGrid(fCorrNom)->Project(nvars,vars,0,0);
+ AliCFGridSparse *hnDeNom=fCFCorrection->GetGrid(fCorrDenom)->Project(nvars,vars,0,0);
+ fCFCorrMatrix->Divide(hnDeNom);
+ delete hnDeNom;
+}
+
+//______________________________________________
+void AliDielectronSpectrum::ExtractSignalInBins(Int_t variable)
+{
+ //
+ // recursively loop over bins and extract signal
+ //
+
+ Int_t varPairType=fCFSignal->GetVar("PairType");
+ Int_t varMass=fCFSignal->GetVar("M");
+
+ for (Int_t ibin=0; ibin<fNbins[variable]; ++ibin){
+ Int_t bin=ibin+1;
+ fCFSignal->GetGrid(fSignalStep)->GetGrid()->GetAxis(fVars[variable])->SetRange(bin,bin);
+ fCurrentBins[variable]=bin;
+ fCurrentPositions[variable]=fCFSignal->GetGrid(fSignalStep)->GetBinCenter(fVars[variable],bin);
+
+ if (variable != fNvars-1) ExtractSignalInBins(variable+1);
+
+ TObjArray arrPairTypes(10);
+ arrPairTypes.SetOwner();
+
+ for (Int_t itype=0; itype<3; ++itype){
+// Int_t itype=1;
+ fCFSignal->SetRangeUser(varPairType,itype,itype,fSignalStep);
+ TH1 *h=fCFSignal->Project(varMass,fSignalStep);
+ h->SetDirectory(0);
+ arrPairTypes.AddAt(h,itype);
+ }
+ AliInfo(Form("Processing bin: %d (%.2f)",ibin, fCurrentPositions[variable]));
+ //loop over all signal extraction methods and retrieve signals
+ for (Int_t iMethod=0; iMethod<fSignalMethods.GetEntries(); ++iMethod){
+ AliDielectronSignalBase *signalMethod=(AliDielectronSignalBase*)fSignalMethods.At(iMethod);
+ signalMethod->Process(&arrPairTypes);
+ const TVectorD &val=signalMethod->GetValues();
+ const TVectorD &err=signalMethod->GetErrors();
+
+ //Fill sparse containers
+ Int_t step=0;
+ if (fStepSignal){
+ //signal
+ fCFSpectrum->GetGrid(step)->SetElement(fCurrentBins,val(0));
+ fCFSpectrum->GetGrid(step)->SetElementError(fCurrentBins,err(0));
+ ++step;
+
+ //corrected signal
+ if (fCFCorrMatrix){
+ Float_t corrFactor = fCFCorrMatrix->GetElement(fCurrentPositions);
+ Float_t corrError = fCFCorrMatrix->GetElementError(fCurrentPositions);
+
+ Float_t value=val(0)*corrFactor;
+ fCFSpectrum->GetGrid(step)->SetElement(fCurrentBins,value);
+ Float_t error=TMath::Sqrt( (err(0)/val(0))*(err(0)/val(0)) +
+ (corrError/corrFactor)*(corrError/corrFactor)
+ )*value;
+ fCFSpectrum->GetGrid(step)->SetElementError(fCurrentBins,error);
+// printf("corrFactor: %f+-%f\n",value,error);
+ }
+ ++step;
+ }
+
+ if (fStepSignificance) {
+ fCFSpectrum->GetGrid(step)->SetElement(fCurrentBins,val(2));
+ fCFSpectrum->GetGrid(step)->SetElementError(fCurrentBins,err(2));
+ ++step;
+ }
+
+ if (fStepSOB) {
+ fCFSpectrum->GetGrid(step)->SetElement(fCurrentBins,val(3));
+ fCFSpectrum->GetGrid(step)->SetElementError(fCurrentBins,err(3));
+ ++step;
+ }
+ }// end method loop
+
+ arrPairTypes.Delete();
+
+ }// end bin loop
+}
+
--- /dev/null
+#ifndef ALIDIELECTRONSPECTRUM_H
+#define ALIDIELECTRONSPECTRUM_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice */
+
+//#############################################################
+//# #
+//# Class AliDielectronSpectrum #
+//# Manage Cuts on the legs of the pair #
+//# #
+//# Authors: #
+//# Anton Andronic, GSI / A.Andronic@gsi.de #
+//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
+//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Frederick Kramer, Uni Ffm, / Frederick.Kramer@cern.ch #
+//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
+//# WooJin J. Park, GSI / W.J.Park@gsi.de #
+//# Jens Wiechula, Uni HD / Jens.Wiechula@cern.ch #
+//# #
+//#############################################################
+
+#include <TObjArray.h>
+#include <TString.h>
+
+#include <TNamed.h>
+
+class AliCFGridSparse;
+class AliDielectronSignalBase;
+class AliCFGridSparse;
+class AliCFContainer;
+
+class AliDielectronSpectrum : public TNamed {
+public:
+ AliDielectronSpectrum();
+ AliDielectronSpectrum(const char*name, const char* title);
+
+ virtual ~AliDielectronSpectrum();
+
+ void AddMethod(AliDielectronSignalBase * const method) {fSignalMethods.Add(method);}
+
+ void SetCorrectionContainer(AliCFContainer * const container, Int_t nominator, Int_t denominator);
+ void SetSignalContainer(AliCFContainer * const container, Int_t step);
+
+ void SetVariables(const char* vars) { fVariables=vars; }
+
+ void SetStepForSignal(Bool_t step=kTRUE) { fStepSignal=step; }
+ void SetStepForSignificance(Bool_t step=kTRUE) { fStepSignificance=step; }
+ void SetStepForSignalOverBackground(Bool_t step=kTRUE) { fStepSOB=step; }
+
+ void SetNoOwnerSpectrum(Bool_t noOwner=kTRUE) { fOwnerSpectrum=!noOwner; }
+
+ AliCFContainer* GetSpectrumContainer() const { return fCFSpectrum; }
+ AliCFGridSparse* GetCorrectionMatrix() const { return fCFCorrMatrix; }
+
+ void Process();
+
+
+private:
+ AliCFContainer *fCFSignal; // CF container with from which to extract the Signal
+ AliCFContainer *fCFCorrection; // CF container from which to extract the correction matrix
+ AliCFContainer *fCFSpectrum; // CF container with extracted signal
+ AliCFGridSparse *fCFCorrMatrix; // correction matrix
+
+ Bool_t fStepSignal; // if to create a step for the signal
+ Bool_t fStepSignificance; // if to create a step for the significance
+ Bool_t fStepSOB; // if to create a step for signal over background
+
+ Int_t fSignalStep; // step to use from the signal container
+
+ Int_t fCorrNom; // Nominator to use from corr matrix container
+ Int_t fCorrDenom; // Deominator to use from corr matrix container
+
+ TObjArray fSignalMethods; // array with signal extraction methods
+ TString fVariables; // variable names as a function of which to extract the signal
+
+ Bool_t fOwnerSpectrum; // if we own the creted spectrum
+
+ Int_t fNvars; //! number of variables
+ Int_t *fVars; //! variable numbers translated from fVariables
+ Int_t *fNbins; //! number of bins for each variable
+ Int_t *fCurrentBins; //! bin currently selected for each variable
+ Double_t *fCurrentPositions; //! variables values currently selected
+
+ void Fill(Int_t *bin, Int_t step, Double_t value, Double_t error);
+ Bool_t SetupVariables();
+ void CreateCorrectionMatrix();
+ void CreateCFSpectrum();
+ void ExtractSignalInBins(Int_t variable=0);
+
+ AliDielectronSpectrum(const AliDielectronSpectrum &c);
+ AliDielectronSpectrum &operator=(const AliDielectronSpectrum &c);
+
+ ClassDef(AliDielectronSpectrum,1) //Cut class providing cuts for both legs of a pair
+
+};
+
+//
+// inline functions
+//
+inline void AliDielectronSpectrum::SetCorrectionContainer(AliCFContainer * const container, Int_t nominator, Int_t denominator)
+{
+ fCFCorrection=container;
+ fCorrNom=nominator;
+ fCorrDenom=denominator;
+}
+
+inline void AliDielectronSpectrum::SetSignalContainer(AliCFContainer * const container, Int_t step)
+{
+ fCFSignal=container;
+ fSignalStep=step;
+}
+
+#endif
fActiveCuts[i]=0;
fCutMin[i]=0;
fCutMax[i]=0;
+ fCutExclude[i]=kFALSE;
}
}
for (Int_t iCut=0; iCut<fNActiveCuts; ++iCut){
Int_t cut=fActiveCuts[iCut];
SETBIT(fSelectedCutsMask,iCut);
- if ( (values[cut]<fCutMin[iCut]) || (values[cut]>fCutMax[iCut]) ) CLRBIT(fSelectedCutsMask,iCut);
+ if ( ((values[cut]<fCutMin[iCut]) || (values[cut]>fCutMax[iCut]))^fCutExclude[iCut] ) CLRBIT(fSelectedCutsMask,iCut);
}
Bool_t isSelected=(fSelectedCutsMask==fActiveCutsMask);
}
//________________________________________________________________________
-void AliDielectronVarCuts::AddCut(AliDielectronVarManager::ValueTypes type, Double_t min, Double_t max)
+void AliDielectronVarCuts::AddCut(AliDielectronVarManager::ValueTypes type, Double_t min, Double_t max, Bool_t excludeRange)
{
//
// Set cut range and activate it
}
fCutMin[fNActiveCuts]=min;
fCutMax[fNActiveCuts]=max;
+ fCutExclude[fNActiveCuts]=excludeRange;
SETBIT(fActiveCutsMask,fNActiveCuts);
fActiveCuts[fNActiveCuts]=(UShort_t)type;
++fNActiveCuts;
}
for (Int_t iCut=0; iCut<fNActiveCuts; ++iCut){
Int_t cut=(Int_t)fActiveCuts[iCut];
- printf("Cut %02d: %f < %s < %f\n", iCut,
- fCutMin[iCut], AliDielectronVarManager::GetValueName((Int_t)cut), fCutMax[iCut]);
+ Bool_t inverse=fCutExclude[iCut];
+
+ if (!inverse){
+ printf("Cut %02d: %f < %s < %f\n", iCut,
+ fCutMin[iCut], AliDielectronVarManager::GetValueName((Int_t)cut), fCutMax[iCut]);
+ } else {
+ printf("Cut %02d: !(%f < %s < %f)\n", iCut,
+ fCutMin[iCut], AliDielectronVarManager::GetValueName((Int_t)cut), fCutMax[iCut]);
+ }
}
}
//#############################################################
#include <Rtypes.h>
+
#include <AliAnalysisCuts.h>
-#include <AliDielectronVarManager.h>
+#include "AliDielectronVarManager.h"
class AliDielectronVarCuts : public AliAnalysisCuts {
public:
AliDielectronVarCuts(const char* name, const char* title);
virtual ~AliDielectronVarCuts();
//TODO: make copy constructor and assignment operator public
- void AddCut(AliDielectronVarManager::ValueTypes type, Double_t min, Double_t max);
- void AddCut(AliDielectronVarManager::ValueTypes type, Double_t value);
+ void AddCut(AliDielectronVarManager::ValueTypes type, Double_t min, Double_t max, Bool_t excludeRange=kFALSE);
+ void AddCut(AliDielectronVarManager::ValueTypes type, Double_t value, Bool_t excludeRange=kFALSE);
// setters
void SetCutOnMCtruth(Bool_t mc=kTRUE) { fCutOnMCtruth=mc; }
Double_t fCutMin[AliDielectronVarManager::kNMaxValues]; // minimum values for the cuts
Double_t fCutMax[AliDielectronVarManager::kNMaxValues]; // maximum values for the cuts
+ Bool_t fCutExclude[AliDielectronVarManager::kNMaxValues]; // inverse cut logic?
AliDielectronVarCuts(const AliDielectronVarCuts &c);
AliDielectronVarCuts &operator=(const AliDielectronVarCuts &c);
//
//Inline functions
//
-inline void AliDielectronVarCuts::AddCut(AliDielectronVarManager::ValueTypes type, Double_t value)
+inline void AliDielectronVarCuts::AddCut(AliDielectronVarManager::ValueTypes type, Double_t value, Bool_t excludeRange)
{
//
// Set cut in a small delta around value
//
const Double_t kDelta=1e-20;
- AddCut(type,value-kDelta,value+kDelta);
+ AddCut(type,value-kDelta,value+kDelta, excludeRange);
}
#endif
"P_InnerParam",
"TPC_signal",
"TPC_nSigma_Electrons",
+ "TPC_nSigma_Pions",
+ "TPC_nSigma_Muons",
+ "TPC_nSigma_Kaons",
+ "TPC_nSigma_Protons",
+ "TOF_nSigma_Pions",
+ "TOF_nSigma_Muons",
+ "TOF_nSigma_Kaons",
+ "TOF_nSigma_Protons",
//
"Chi2NDF",
"DecayLength",
"R",
"OpeningAngle",
+ "LegDistance",
+ "LegDistanceXY",
"Merr",
"DCA",
"PairType",
"Nevents"
};
-AliESDpid* AliDielectronVarManager::fgESDpid = new AliESDpid;
+AliESDpid* AliDielectronVarManager::fgESDpid = 0x0;
AliVEvent* AliDielectronVarManager::fgEvent = 0x0;
+AliKFVertex* AliDielectronVarManager::fgKFVertex = 0x0;
//________________________________________________________________
AliDielectronVarManager::AliDielectronVarManager() :
TNamed("AliDielectronVarManager","AliDielectronVarManager")
//# Anton Andronic, GSI / A.Andronic@gsi.de #
//# Ionut C. Arsene, GSI / I.C.Arsene@gsi.de #
//# Julian Book, Uni Ffm / Julian.Book@cern.ch #
+//# Markus Köhler, GSI / M.Koehler@gsi.de #
//# Frederick Kramer, Uni Ffm / Frederick.Kramer@cern.ch #
//# Magnus Mager, CERN / Magnus.Mager@cern.ch #
//# WooJin J. Park, GSI / W.J.Park@gsi.de #
#include <AliVParticle.h>
#include <AliESDtrack.h>
#include <AliAODTrack.h>
+#include <AliAODPid.h>
#include <AliKFParticle.h>
+#include <AliKFVertex.h>
#include <AliMCParticle.h>
#include <AliAODMCParticle.h>
#include <AliVTrack.h> // ?
kPdgCode, // PDG code
kPIn, // momentum at inner wall of TPC (if available), used for PID
kTPCsignal, // TPC dE/dx signal
+
kTPCnSigmaEle, // number of sigmas to the dE/dx electron line in the TPC
+ kTPCnSigmaPio, // number of sigmas to the dE/dx pion line in the TPC
+ kTPCnSigmaMuo, // number of sigmas to the dE/dx muon line in the TPC
+ kTPCnSigmaKao, // number of sigmas to the dE/dx kaon line in the TPC
+ kTPCnSigmaPro, // number of sigmas to the dE/dx proton line in the TPC
+
+ kTOFnSigmaPio, // number of sigmas to the pion line in the TOF
+ kTOFnSigmaMuo, // number of sigmas to the muon line in the TOF
+ kTOFnSigmaKao, // number of sigmas to the kaon line in the TOF
+ kTOFnSigmaPro, // number of sigmas to the proton line in the TOF
+
kParticleMax, //
// TODO: kRNClusters ??
// AliDielectronPair specific variables
kDecayLength, // decay length
kR, // distance to the origin
kOpeningAngle, // opening angle
+ kLegDist, // distance of the legs
+ kLegDistXY, // distance of the legs in XY
kMerr, // error of mass calculation
kDCA, // distance of closest approach TODO: not implemented yet
kPairType, // type of the pair, like like sign ++ unlikesign ...
static void Fill(const TObject* particle, Double_t * const values);
static void InitESDpid(Int_t type=0);
- static void SetEvent(AliVEvent * const ev) { fgEvent = ev; }
+ static void SetESDpid(AliESDpid * const pid) {fgESDpid=pid;}
+ static AliESDpid* GetESDpid() {return fgESDpid;}
+ static void SetEvent(AliVEvent * const ev);
static const char* GetValueName(Int_t i) { return (i>=0&&i<kNMaxValues)?fgkParticleNames[i]:""; }
private:
static void FillVarVParticle(const AliVParticle *particle, Double_t * const values);
static void FillVarESDtrack(const AliESDtrack *particle, Double_t * const values);
static void FillVarAODTrack(const AliAODTrack *particle, Double_t * const values);
- static void FillVarMCParticle(const AliMCParticle *particle, Double_t * const values);
+ static void FillVarMCParticle(const AliMCParticle *particle, Double_t * const values);
static void FillVarAODMCParticle(const AliAODMCParticle *particle, Double_t * const values);
static void FillVarDielectronPair(const AliDielectronPair *pair, Double_t * const values);
static void FillVarVEvent(const AliVEvent *event, Double_t * const values);
static AliESDpid* fgESDpid; // ESD pid object
static AliVEvent* fgEvent; // current event pointer
-
+ static AliKFVertex *fgKFVertex; // kf vertex
+
AliDielectronVarManager(const AliDielectronVarManager &c);
AliDielectronVarManager &operator=(const AliDielectronVarManager &c);
// TODO: for the moment we set the bethe bloch parameters manually
// this should be changed in future!
values[AliDielectronVarManager::kTPCnSigmaEle]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kElectron);
+ values[AliDielectronVarManager::kTPCnSigmaPio]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kPion);
+ values[AliDielectronVarManager::kTPCnSigmaMuo]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kMuon);
+ values[AliDielectronVarManager::kTPCnSigmaKao]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kKaon);
+ values[AliDielectronVarManager::kTPCnSigmaPro]=fgESDpid->NumberOfSigmasTPC(particle,AliPID::kProton);
+
+ Double_t t0=fgESDpid->GetTOFResponse().GetTimeZero();
+ values[AliDielectronVarManager::kTOFnSigmaPio]=fgESDpid->NumberOfSigmasTOF(particle,AliPID::kPion,t0);
+ values[AliDielectronVarManager::kTOFnSigmaMuo]=fgESDpid->NumberOfSigmasTOF(particle,AliPID::kMuon,t0);
+ values[AliDielectronVarManager::kTOFnSigmaKao]=fgESDpid->NumberOfSigmasTOF(particle,AliPID::kKaon,t0);
+ values[AliDielectronVarManager::kTOFnSigmaPro]=fgESDpid->NumberOfSigmasTOF(particle,AliPID::kProton,t0);
}
inline void AliDielectronVarManager::FillVarAODTrack(const AliAODTrack *particle, Double_t * const values)
// Fill common AliVParticle interface information
FillVarVParticle(particle, values);
- // Fill AliAODTrack interface information
+
+ // Reset AliESDtrack interface specific information
+ values[AliDielectronVarManager::kNclsITS] = 0;
+ values[AliDielectronVarManager::kNclsTPC] = 0;
+ values[AliDielectronVarManager::kNFclsTPC] = 0;
+ values[AliDielectronVarManager::kNclsTRD] = 0;
+ values[AliDielectronVarManager::kTRDntracklets] = 0;
+ values[AliDielectronVarManager::kTRDpidQuality] = 0;
+
+ //TODO: This is only an approximation!!!
+ values[AliDielectronVarManager::kTPCsignalN] = particle->GetTPCClusterMap().CountBits();
+
+// Fill AliAODTrack interface information
// ...
+ values[AliDielectronVarManager::kImpactParXY] = particle->DCA();
+ values[AliDielectronVarManager::kImpactParZ] = particle->ZAtDCA();
+
+ values[AliDielectronVarManager::kPIn]=0;
+ values[AliDielectronVarManager::kTPCsignal]=0;
+
+ values[AliDielectronVarManager::kTPCnSigmaEle]=0;
+ values[AliDielectronVarManager::kTPCnSigmaPio]=0;
+ values[AliDielectronVarManager::kTPCnSigmaMuo]=0;
+ values[AliDielectronVarManager::kTPCnSigmaKao]=0;
+ values[AliDielectronVarManager::kTPCnSigmaPro]=0;
+
+ AliAODPid *pid=particle->GetDetPid();
+ if (pid){
+ Double_t mom =pid->GetTPCmomentum();
+ //TODO: kTPCsignalN is only an approximation (see above)!!
+
+ Double_t tpcNsigmaEle=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
+ TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]) ,AliPID::kElectron);
+
+ Double_t tpcNsigmaPio=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
+ TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kPion);
+
+ Double_t tpcNsigmaMuo=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
+ TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kMuon);
+
+ Double_t tpcNsigmaKao=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
+ TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kKaon);
+
+ Double_t tpcNsigmaPro=fgESDpid->GetTPCResponse().GetNumberOfSigmas(mom,pid->GetTPCsignal(),
+ TMath::Nint(values[AliDielectronVarManager::kTPCsignalN]),AliPID::kProton);
+
+ values[AliDielectronVarManager::kPIn]=pid->GetTPCmomentum();
+ values[AliDielectronVarManager::kTPCsignal]=pid->GetTPCsignal();
+
+ values[AliDielectronVarManager::kTPCnSigmaEle]=tpcNsigmaEle;
+ values[AliDielectronVarManager::kTPCnSigmaPio]=tpcNsigmaPio;
+ values[AliDielectronVarManager::kTPCnSigmaMuo]=tpcNsigmaMuo;
+ values[AliDielectronVarManager::kTPCnSigmaKao]=tpcNsigmaKao;
+ values[AliDielectronVarManager::kTPCnSigmaPro]=tpcNsigmaPro;
+
+ values[AliDielectronVarManager::kTRDntracklets] = 0;
+ values[AliDielectronVarManager::kTRDpidQuality] = 0;
+
+ }
}
inline void AliDielectronVarManager::FillVarMCParticle(const AliMCParticle *particle, Double_t * const values)
values[AliDielectronVarManager::kDecayLength] = kfPair.GetDecayLength();
values[AliDielectronVarManager::kR] = kfPair.GetR();
values[AliDielectronVarManager::kOpeningAngle] = pair->OpeningAngle();
- values[AliDielectronVarManager::kMerr] = kfPair.GetErrMass()>0?kfPair.GetErrMass()/kfPair.GetMass():1000000;
+ values[AliDielectronVarManager::kLegDist] = pair->DistanceDaughters();
+ values[AliDielectronVarManager::kLegDistXY] = pair->DistanceDaughtersXY();
+ values[AliDielectronVarManager::kMerr] = kfPair.GetErrMass()>1e-30&&kfPair.GetErrMass()>1e-30?kfPair.GetErrMass()/kfPair.GetMass():1000000;
values[AliDielectronVarManager::kPairType] = pair->GetType();
}
// type=0 is simulation
// type=1 is data
+ if (!fgESDpid) fgESDpid=new AliESDpid;
Double_t alephParameters[5];
// simulation
alephParameters[0] = 2.15898e+00/50.;
}
+
+inline void AliDielectronVarManager::SetEvent(AliVEvent * const ev)
+{
+
+ fgEvent = ev;
+ if (fgKFVertex) delete fgKFVertex;
+ fgKFVertex=0x0;
+ if (ev) fgKFVertex=new AliKFVertex(*ev->GetPrimaryVertex());
+}
/*
inline void AliDielectronVarManager::FillValues(const TParticle *particle, Double_t *values)
{
--- /dev/null
+
+void InitHistograms(AliDielectron *die, Int_t cutDefinition);
+void InitCF(AliDielectron* die, Int_t cutDefinition);
+
+void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition);
+void SetupPairCuts(AliDielectron *die, Int_t cutDefinition);
+
+AliESDtrackCuts *SetupESDtrackCuts(Int_t cutDefinition);
+
+TString names=("trackQ+Pt>1GeV+loose_PID;trackQ+Pt>1GeV+tight_PID;trackQ+Pt>.5GeV");
+TObjArray *arrNames=names.Tokenize(";");
+
+const Int_t nDie=arrNames->GetEntries();
+
+AliDielectron* ConfigJpsi2ee(Int_t cutDefinition)
+{
+ //
+ // Setup the instance of AliDielectron
+ //
+
+ // create the actual framework object
+ TString name=Form("%02d",cutDefinition);
+ if (cutDefinition<arrNames->GetEntriesFast()){
+ name=arrNames->At(cutDefinition)->GetName();
+ }
+ AliDielectron *die =
+ new AliDielectron(Form("%s",name.Data()),
+ Form("Track cuts: %s",name.Data()));
+
+ // cut setup
+ SetupTrackCuts(die,cutDefinition);
+ SetupPairCuts(die,cutDefinition);
+
+ //
+ // histogram setup
+ // only if an AliDielectronHistos object is attached to the
+ // dielectron framework the QA 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);
+
+ return die;
+}
+
+//______________________________________________________________________________________
+void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition)
+{
+ //
+ // Setup the track cuts
+ //
+
+ //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);
+ die->GetTrackFilter().AddCuts(pt);
+
+ return;
+ }
+
+
+ AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>1","Pt>1");
+ pt->AddCut(AliDielectronVarManager::kPt,1.,1e30);
+ die->GetTrackFilter().AddCuts(pt);
+
+ //loose PID
+ if (cutDefinition==0){
+ //ESD pid cuts (TPC nSigma Electrons)
+ AliDielectronVarCuts *pid = new AliDielectronVarCuts("TPCpid","TPC nSigma cut |e|<3 pi>2 |P|>2 |K|>2");
+ //include
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaEle, -3., 3.);
+ //exclude
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPio, -20., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPro, -2., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaKao, -2., 2., kTRUE);
+ die->GetTrackFilter().AddCuts(pid);
+ }
+
+ //tight PID
+ if (cutDefinition==1){
+ //ESD pid cuts (TPC nSigma Electrons)
+ AliDielectronVarCuts *pid = new AliDielectronVarCuts("TPCpid","TPC nSigma cut |e|<2 pi>2 |P|>2 |K|>2; |dXY|<200um");
+ //include
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaEle, -2., 2.);
+ //exclude
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPio, -2., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPro, -2., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaKao, -2., 2., kTRUE);
+ //tighten Impact parameter to 200um
+ pid->AddCut(AliDielectronVarManager::kImpactParXY, -0.02, 0.02);
+ die->GetTrackFilter().AddCuts(pid);
+ }
+
+}
+
+//______________________________________________________________________________________
+void SetupPairCuts(AliDielectron *die, Int_t cutDefinition)
+{
+ //
+ // 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);
+}
+
+//______________________________________________________________________________________
+AliESDtrackCuts *SetupESDtrackCuts(Int_t cutDefinition)
+{
+ //
+ // Setup default AliESDtrackCuts
+ //
+ AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts;
+
+ esdTrackCuts->SetMaxDCAToVertexZ(3.0);
+ esdTrackCuts->SetMaxDCAToVertexXY(.07);
+ esdTrackCuts->SetRequireTPCRefit(kTRUE);
+ esdTrackCuts->SetRequireITSRefit(kTRUE);
+ esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
+ esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
+
+ esdTrackCuts->SetMinNClustersTPC(75);
+ esdTrackCuts->SetMaxChi2PerClusterTPC(4);
+
+ return esdTrackCuts;
+}
+
+//______________________________________________________________________________________
+void InitHistograms(AliDielectron *die, Int_t cutDefinition)
+{
+ //
+ // Initialise the QA histograms
+ //
+
+ //Setup QA histograms
+ 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
+ for (Int_t i=0; i<2; ++i){
+ histos->AddClass(Form("Track_%s",AliDielectron::TrackClassName(i)));
+ }
+
+ //Pair classes, only first event
+ 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);
+ }
+
+ //add histograms to Track classes
+ histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",200,0,20.,AliDielectronVarManager::kPt);
+ 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);
+
+ //add histograms to Pair classes
+ histos->UserHistogram("Pair","InvMass","Inv.Mass;Inv. Mass [GeV];#pairs",
+ 500,0.,4.,AliDielectronVarManager::kM);
+ histos->UserHistogram("Pair","Rapidity","Rapidity;Rapidity;#pairs",
+ 100,-2.,2.,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);
+}
+
+
+void InitCF(AliDielectron* die, Int_t cutDefinition)
+{
+ //
+ // Setupd the CF Manager if needed
+ //
+ AliDielectronCF *cf=new AliDielectronCF(die->GetName(),die->GetTitle());
+
+ //pair variables
+ cf->AddVariable(AliDielectronVarManager::kPt,20,0,10);
+ cf->AddVariable(AliDielectronVarManager::kY,40,-2,2);
+ cf->AddVariable(AliDielectronVarManager::kM,200,-0.01,3.99);
+ 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();
+ }
+
+ die->SetCFManagerPair(cf);
+}
--- /dev/null
+
+void InitHistograms(AliDielectron *die, Int_t cutDefinition);
+void InitCF(AliDielectron* die, Int_t cutDefinition);
+
+void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition);
+void SetupPairCuts(AliDielectron *die, Int_t cutDefinition);
+
+TString names=("trackQ+Pt>1GeV+loose_PID;trackQ+Pt>1GeV+tight_PID;trackQ+Pt>.5GeV");
+TObjArray *arrNames=names.Tokenize(";");
+
+const Int_t nDie=arrNames->GetEntries();
+
+AliDielectron* ConfigJpsi2ee(Int_t cutDefinition)
+{
+ //
+ // Setup the instance of AliDielectron
+ //
+
+ // create the actual framework object
+ TString name=Form("%02d",cutDefinition);
+ if (cutDefinition<arrNames->GetEntriesFast()){
+ name=arrNames->At(cutDefinition)->GetName();
+ }
+ AliDielectron *die =
+ new AliDielectron(Form("%s",name.Data()),
+ Form("Track cuts: %s",name.Data()));
+
+ // cut setup
+ SetupTrackCuts(die,cutDefinition);
+ SetupPairCuts(die,cutDefinition);
+
+ //
+ // histogram setup
+ // only if an AliDielectronHistos object is attached to the
+ // dielectron framework the QA 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);
+
+ return die;
+}
+
+//______________________________________________________________________________________
+void SetupTrackCuts(AliDielectron *die, Int_t cutDefinition)
+{
+ //
+ // Setup the track cuts
+ //
+
+ //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);
+ die->GetTrackFilter().AddCuts(pt);
+
+ return;
+ }
+
+
+ AliDielectronVarCuts *pt = new AliDielectronVarCuts("Pt>1","Pt>1");
+ pt->AddCut(AliDielectronVarManager::kPt,1.,1e30);
+ pt->AddCut(AliDielectronVarManager::kImpactParZ,-3.,3.);
+ pt->AddCut(AliDielectronVarManager::kImpactParXY,-.07,.07);
+ die->GetTrackFilter().AddCuts(pt);
+
+ //loose PID
+ if (cutDefinition==0){
+ //ESD pid cuts (TPC nSigma Electrons)
+ AliDielectronVarCuts *pid = new AliDielectronVarCuts("TPCpid","TPC nSigma cut |e|<3 pi>2 |P|>2 |K|>2");
+ //include
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaEle, -3., 3.);
+ //exclude
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPio, -20., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPro, -2., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaKao, -2., 2., kTRUE);
+ die->GetTrackFilter().AddCuts(pid);
+ }
+
+ //tight PID
+ if (cutDefinition==1){
+ //ESD pid cuts (TPC nSigma Electrons)
+ AliDielectronVarCuts *pid = new AliDielectronVarCuts("TPCpid","TPC nSigma cut |e|<2 pi>2 |P|>2 |K|>2; |dXY|<200um");
+ //include
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaEle, -2., 2.);
+ //exclude
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPio, -2., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaPro, -2., 2., kTRUE);
+ pid->AddCut(AliDielectronVarManager::kTPCnSigmaKao, -2., 2., kTRUE);
+ //tighten Impact parameter to 200um
+ pid->AddCut(AliDielectronVarManager::kImpactParXY, -0.02, 0.02);
+ die->GetTrackFilter().AddCuts(pid);
+ }
+
+}
+
+//______________________________________________________________________________________
+void SetupPairCuts(AliDielectron *die, Int_t cutDefinition)
+{
+ //
+ // 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);
+}
+
+//______________________________________________________________________________________
+void InitHistograms(AliDielectron *die, Int_t cutDefinition)
+{
+ //
+ // Initialise the QA histograms
+ //
+
+ //Setup QA histograms
+ 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
+ for (Int_t i=0; i<2; ++i){
+ histos->AddClass(Form("Track_%s",AliDielectron::TrackClassName(i)));
+ }
+
+ //Pair classes, only first event
+ 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);
+ }
+
+ //add histograms to Track classes
+ histos->UserHistogram("Track","Pt","Pt;Pt [GeV];#tracks",200,0,20.,AliDielectronVarManager::kPt);
+ 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);
+
+ //add histograms to Pair classes
+ histos->UserHistogram("Pair","InvMass","Inv.Mass;Inv. Mass [GeV];#pairs",
+ 500,0.,4.,AliDielectronVarManager::kM);
+ histos->UserHistogram("Pair","Rapidity","Rapidity;Rapidity;#pairs",
+ 100,-2.,2.,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);
+}
+
+
+void InitCF(AliDielectron* die, Int_t cutDefinition)
+{
+ //
+ // Setupd the CF Manager if needed
+ //
+ AliDielectronCF *cf=new AliDielectronCF(die->GetName(),die->GetTitle());
+
+ //pair variables
+ cf->AddVariable(AliDielectronVarManager::kPt,20,0,10);
+ cf->AddVariable(AliDielectronVarManager::kY,40,-2,2);
+ cf->AddVariable(AliDielectronVarManager::kM,200,-0.01,3.99);
+ 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();
+ }
+
+ die->SetCFManagerPair(cf);
+}