From 745d6088cbdf350ac6aad7a84b8b28db30d5515f Mon Sep 17 00:00:00 2001 From: jgrosseo Date: Wed, 21 May 2008 14:00:33 +0000 Subject: [PATCH] fixing warnings removing selectors removing AliESDtrackCuts (now in ANALYSISalice) adding acceptance range to SPD analysis --- PWG0/AliCorrection.cxx | 11 +- PWG0/AliCorrectionMatrix2D.cxx | 10 + PWG0/AliCorrectionMatrix2D.h | 2 + PWG0/AliPWG0Helper.h | 2 +- PWG0/PWG0Helper.C | 201 --- PWG0/PWG0baseLinkDef.h | 2 - PWG0/dNdEta/AliMultiplicityCorrection.cxx | 2 +- PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.cxx | 229 ---- PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.h | 39 - PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx | 345 ----- PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h | 48 - PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx | 7 +- PWG0/dNdEta/AlidNdEtaSystematicsSelector.cxx | 2 +- PWG0/dNdEta/AlidNdEtaTask.cxx | 7 +- PWG0/dNdEta/dNdEtaAnalysis.cxx | 91 +- PWG0/dNdEta/dNdEtaAnalysis.h | 8 +- PWG0/dNdEta/drawPlots.C | 114 +- PWG0/dNdEta/plotsMultiplicity.C | 10 +- PWG0/dNdEta/run.C | 5 +- PWG0/esdTrackCuts/AliCutTask.cxx | 2 +- PWG0/esdTrackCuts/AliESDtrackCuts.cxx | 1186 ----------------- PWG0/esdTrackCuts/AliESDtrackCuts.h | 174 --- .../AliTestESDtrackCutsSelector.cxx | 2 +- .../AliHighMultiplicitySelector.cxx | 42 +- PWG0/libPWG0base.pkg | 1 - PWG0/multiplicity/AliMultiplicityTask.cxx | 26 +- PWG0/multiplicity/AliMultiplicityTask.h | 8 +- PWG0/multiplicity/CreateCuts.C | 22 - PWG0/multiplicity/run.C | 52 +- 29 files changed, 298 insertions(+), 2352 deletions(-) delete mode 100644 PWG0/PWG0Helper.C delete mode 100644 PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.cxx delete mode 100644 PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.h delete mode 100644 PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx delete mode 100644 PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h delete mode 100644 PWG0/esdTrackCuts/AliESDtrackCuts.cxx delete mode 100644 PWG0/esdTrackCuts/AliESDtrackCuts.h delete mode 100644 PWG0/multiplicity/CreateCuts.C diff --git a/PWG0/AliCorrection.cxx b/PWG0/AliCorrection.cxx index c4739e17365..9b821edc699 100644 --- a/PWG0/AliCorrection.cxx +++ b/PWG0/AliCorrection.cxx @@ -45,7 +45,7 @@ AliCorrection::AliCorrection(const Char_t* name, const Char_t* title, AliPWG0Hel static Float_t binLimitsPtTmp[] = {0.0, 0.05, 0.1, 0.125, 0.15, 0.175, 0.2, 0.225, 0.25, 0.275, 0.3, 0.325, 0.35, 0.375, 0.4, 0.425, 0.45, 0.475, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.5, 2.0, 5.0, 10.0, 100.0}; binLimitsPt = (Float_t*) binLimitsPtTmp; nBinsPt = 28; - } + } else if (analysisMode == AliPWG0Helper::kSPD) { static Float_t binLimitsPtTmp[] = {-0.5, 100.0}; @@ -63,10 +63,11 @@ AliCorrection::AliCorrection(const Char_t* name, const Char_t* title, AliPWG0Hel //Float_t binLimitsVtx[] = {-20,-15,-10,-6,-3,0,3,6,10,15,20}; //Float_t binLimitsVtx[] = {-20,-19,-18,-17,-16,-15,-14,-13,-12,-11,-10,-9,-8,-7,-6,-5,-4,-3,-2,-1,0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20}; Float_t binLimitsVtx[] = {-20,-15,-10,-8,-6,-4,-2,0,2,4,6,8,10,15,20}; - Float_t binLimitsEta[] = {-2,-1.9,-1.8,-1.7,-1.6,-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0}; + Float_t binLimitsEta[] = {-2,-1.8,-1.6,-1.4,-1.2,-1.0,-0.8,-0.6,-0.4,-0.2,0.0,0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,1.8,2.0}; +// Float_t binLimitsEta[] = {-2,-1.9,-1.8,-1.7,-1.6,-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0}; // Float_t binLimitsEta[] = { -7.0, -6.9, -6.8, -6.7, -6.6, -6.5, -6.4, -6.3, -6.2, -6.1, -6.0, -5.9, -5.8, -5.7, -5.6, -5.5, -5.4, -5.3, -5.2, -5.1, -5.0, -4.9, -4.8, -4.7, -4.6, -4.5, -4.4, -4.3, -4.2, -4.1, -4.0, -3.9, -3.8, -3.7, -3.6, -3.5, -3.4, -3.3, -3.2, -3.1, -3.0, -2.9, -2.8, -2.7, -2.6, -2.5, -2.4, -2.3, -2.2, -2.1, -2.0, -1.9, -1.8, -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1, -1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, -0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9, 3.0, 3.1, 3.2, 3.3, 3.4, 3.5, 3.6, 3.7, 3.8, 3.9, 4.0, 4.1, 4.2, 4.3, 4.4, 4.5, 4.6, 4.7, 4.8, 4.9, 5.0, 5.1, 5.2, 5.3, 5.4, 5.5, 5.6, 5.7, 5.8, 5.9, 6.0, 6.1, 6.2, 6.3, 6.4, 6.5, 6.6, 6.7, 6.8, 6.9, 7.0 }; - TH3F* dummyBinning = new TH3F("dummyBinning","dummyBinning",14, binLimitsVtx, 40, binLimitsEta , nBinsPt, binLimitsPt); + TH3F* dummyBinning = new TH3F("dummyBinning","dummyBinning",14, binLimitsVtx, 20, binLimitsEta , nBinsPt, binLimitsPt); fEventCorr = new AliCorrectionMatrix2D("EventCorrection", Form("%s EventCorrection", fTitle.Data()), 14, binLimitsVtx, 22, binLimitsN); fTrackCorr = new AliCorrectionMatrix3D("TrackCorrection", Form("%s TrackCorrection", fTitle.Data()), dummyBinning); @@ -370,7 +371,7 @@ void AliCorrection::PrintInfo(Float_t ptCut) Printf("AliCorrection::PrintInfo: Whole phasespace:"); - Printf("tracks measured: %.1f tracks generated: %.1f, events measured: %1.f, events generated %1.f", measured->Integral(), generated->Integral(), measuredEvents->Integral(), generatedEvents->Integral()); + Printf("tracks measured: %.1f tracks generated: %.1f, events measured: %.1f, events generated %.1f", measured->Integral(), generated->Integral(), measuredEvents->Integral(), generatedEvents->Integral()); Printf("AliCorrection::PrintInfo: Values in |eta| < 0.8, |vtx-z| < 10 and pt > %.2f:", ptCut); @@ -380,7 +381,7 @@ void AliCorrection::PrintInfo(Float_t ptCut) Float_t eventsM = measuredEvents->Integral(measuredEvents->GetXaxis()->FindBin(-9.9), measuredEvents->GetXaxis()->FindBin(9.9), 1, measuredEvents->GetNbinsY()); Float_t eventsG = generatedEvents->Integral(generatedEvents->GetXaxis()->FindBin(-9.9), generatedEvents->GetXaxis()->FindBin(9.9), 1, generatedEvents->GetNbinsY()); - Printf("tracks measured: %.1f tracks generated: %.1f, events measured: %1.f, events generated %1.f", tracksM, tracksG, eventsM, eventsG); + Printf("tracks measured: %.1f tracks generated: %.1f, events measured: %.1f, events generated %.1f", tracksM, tracksG, eventsM, eventsG); if (tracksM > 0) Printf("Effective average correction factor for TRACKS: %.3f", tracksG / tracksM); diff --git a/PWG0/AliCorrectionMatrix2D.cxx b/PWG0/AliCorrectionMatrix2D.cxx index e208f4357b6..4d9a0dcc24b 100644 --- a/PWG0/AliCorrectionMatrix2D.cxx +++ b/PWG0/AliCorrectionMatrix2D.cxx @@ -245,3 +245,13 @@ void AliCorrectionMatrix2D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBi } +//____________________________________________________________________ +void AliCorrectionMatrix2D::Rebin(Int_t x, Int_t y) +{ + // rebins the histograms, recalculates the correction + + GetGeneratedHistogram()->Rebin2D(x, y); + GetMeasuredHistogram()->Rebin2D(x, y); + GetCorrectionHistogram()->Rebin2D(x, y); + Divide(); +} diff --git a/PWG0/AliCorrectionMatrix2D.h b/PWG0/AliCorrectionMatrix2D.h index 9ee513f85b4..4726a897491 100644 --- a/PWG0/AliCorrectionMatrix2D.h +++ b/PWG0/AliCorrectionMatrix2D.h @@ -39,6 +39,8 @@ public: TH1F* Get1DCorrection(Char_t* opt="x", Float_t min=0, Float_t max=0) {return Get1DCorrectionHistogram(opt,min,max);} TH1F* Get1DCorrectionHistogram(Char_t* opt="x", Float_t min=0, Float_t max=0); + void Rebin(Int_t x = 1, Int_t y = 1); + void FillMeas(Float_t ax, Float_t ay); void FillGene(Float_t ax, Float_t ay); Float_t GetCorrection(Float_t ax, Float_t ay) const; diff --git a/PWG0/AliPWG0Helper.h b/PWG0/AliPWG0Helper.h index 1a50fe46c3d..583e346ce5f 100644 --- a/PWG0/AliPWG0Helper.h +++ b/PWG0/AliPWG0Helper.h @@ -22,7 +22,7 @@ class AliPWG0Helper : public TObject { public: enum Trigger { kMB1 = 0, kMB2 }; // definition from ALICE-INT-2005-025 - enum AnalysisMode { kSPD = 0, kTPC, kTPCITS }; + enum AnalysisMode { kInvalid = -1, kSPD = 0, kTPC, kTPCITS }; static Bool_t IsEventTriggered(const AliESD* aEsd, Trigger trigger = kMB2); static Bool_t IsEventTriggered(ULong64_t triggerMask, Trigger trigger = kMB2); diff --git a/PWG0/PWG0Helper.C b/PWG0/PWG0Helper.C deleted file mode 100644 index 02b0bdb1d40..00000000000 --- a/PWG0/PWG0Helper.C +++ /dev/null @@ -1,201 +0,0 @@ -/* $Id$ */ - -// Helper macros can be found in this file -// A set of them can be used to connect to proof and execute selectors. - -TProof* connectProof(const char* proofServer) -{ - TProof* proof = TProof::Open(proofServer); - - if (!proof) - { - printf("ERROR: PROOF connection not established.\n"); - return 0; - } - - // enable the new packetizer - //proof->AddInput(new TNamed("PROOF_Packetizer", "TPacketizerProgressive")); - - proof->ClearInput(); - - return proof; -} - -Bool_t prepareQuery(TString libraries, TString packages, Int_t useAliRoot) -{ - // if not proof load libraries - if (!gProof) - { - TObjArray* librariesList = libraries.Tokenize(";"); - for (Int_t i=0; iGetEntries(); ++i) - { - TObjString* str = dynamic_cast (librariesList->At(i)); - if (!str) - continue; - - printf("Loading %s...", str->String().Data()); - Int_t result = CheckLoadLibrary(str->String()); - if (result < 0) - { - printf("failed\n"); - //return kFALSE; - } - else - printf("succeeded\n"); - } - } - else - { - if (useAliRoot > 0) - ProofEnableAliRoot(useAliRoot); - - TObjArray* packagesList = packages.Tokenize(";"); - for (Int_t i=0; iGetEntries(); ++i) - { - TObjString* str = dynamic_cast (packagesList->At(i)); - if (!str) - continue; - - /*if (!EnablePackageLocal(str->String())) - { - printf("Loading of package %s locally failed\n", str->String().Data()); - return kFALSE; - }*/ - - if (gProof->UploadPackage(Form("%s.par", str->String().Data()))) - { - printf("Uploading of package %s failed\n", str->String().Data()); - return kFALSE; - } - - if (gProof->EnablePackage(str->String())) - { - printf("Loading of package %s failed\n", str->String().Data()); - return kFALSE; - } - } - } - - return kTRUE; -} - -Int_t executeQuery(TChain* chain, TList* inputList, TString selectorName, const char* option = "", Long64_t entries = TChain::kBigNumber) -{ - if (!gProof) - chain->GetUserInfo()->AddAll(inputList); - else - { - for (Int_t i=0; iGetEntries(); ++i) - gProof->AddInput(inputList->At(i)); - } - - TStopwatch timer; - timer.Start(); - - Long64_t result = -1; - - if (gProof) - { - chain->SetProof(); - result = chain->Process(selectorName, option, entries); - } - else - result = chain->Process(selectorName, option, entries); - - if (result < 0) - printf("ERROR: Executing process failed with %d.\n", result); - - timer.Stop(); - timer.Print(); - - return result; -} - -const char* GetAliRootLocation(Int_t aliroot) -{ - switch (aliroot) - { - case 1: return "/afs/cern.ch/alice/caf/sw/ALICE/$ROOTVERSIONTAG/v4-04-Release/slc4_ia32_gcc34/aliroot"; break; - default: return 0; - } -} - -void ProofEnableAliRoot(Int_t aliroot) -{ - // enables a locally deployed AliRoot in a PROOF cluster - - /* executes the following commands on each node: - gSystem->Setenv("ALICE_ROOT", "/home/alicecaf/ALICE/aliroot-head") - gSystem->AddIncludePath("/home/alicecaf/ALICE/aliroot-head/include"); - gSystem->SetDynamicPath(Form("%s:%s", gSystem->GetDynamicPath(), "/home/alicecaf/ALICE/aliroot-head/lib/tgt_linux")) - gSystem->Load("libMinuit"); - gROOT->Macro("$ALICE_ROOT/macros/loadlibs.C"); - */ - - const char* location = GetAliRootLocation(aliroot); - const char* target = "tgt_linux"; - - gProof->Exec(Form("TString str(gSystem->ExpandPathName(\"%s\")); gSystem->Setenv(\"ALICE_ROOT\", str);", location), kTRUE); - - gProof->AddIncludePath(Form("%s/include", location)); - gProof->AddDynamicPath(Form("%s/lib/%s", location, target)); - - // load all libraries - gProof->Exec("gROOT->Macro(\"$ALICE_ROOT/macros/loadlibs.C\")"); -} - -void ProofAddAliRootIncludePath(Int_t aliroot, const char* dir) -{ - // adds an include path inside the aliroot structure - - const char* location = GetAliRootLocation(aliroot); - gProof->AddIncludePath(Form("%s/%s", location, dir)); -} - -Bool_t EnablePackageLocal(const char* package) -{ - printf("Enabling package %s locally...\n", package); - - TString currentDir(gSystem->pwd()); - if (!gSystem->cd(package)) - return kFALSE; - - gROOT->ProcessLine(".x PROOF-INF/SETUP.C"); - gSystem->cd(currentDir); - - return kTRUE; -} - -Int_t CheckLoadLibrary(const char* library) -{ - // checks if a library is already loaded, if not loads the library - - if (strlen(gSystem->GetLibraries(Form("%s.so", library), "", kFALSE)) > 0) - return 1; - - return gSystem->Load(library); -} - -void redeployPackages(const char* proofServer, Bool_t localAliRoot = kTRUE) -{ - // deploys PWG0base and PWG0dep (the latter only when localAliRoot is true) that are expected in $ALICE_ROOT - // when localAliRoot is false ESD.par is also deployed - - TProof::Reset(proofServer); - TVirtualProof* proof = TProof::Open(proofServer); - proof->ClearPackages(); - - if (localAliRoot) - ProofEnableAliRoot(); - else - { - proof->UploadPackage("$ALICE_ROOT/ESD.par"); - proof->EnablePackage("ESD"); - } - - proof->UploadPackage("$ALICE_ROOT/PWG0base.par"); - proof->EnablePackage("PWG0base"); - - proof->UploadPackage("$ALICE_ROOT/PWG0dep.par"); - proof->EnablePackage("PWG0dep"); -} diff --git a/PWG0/PWG0baseLinkDef.h b/PWG0/PWG0baseLinkDef.h index 3f4137c1b54..ef2767c54eb 100644 --- a/PWG0/PWG0baseLinkDef.h +++ b/PWG0/PWG0baseLinkDef.h @@ -11,8 +11,6 @@ #pragma link C++ class dNdEtaAnalysis+; #pragma link C++ class AlidNdEtaCorrection+; -#pragma link C++ class AliESDtrackCuts+; - #pragma link C++ class AliCorrectionMatrix+; #pragma link C++ class AliCorrectionMatrix2D+; #pragma link C++ class AliCorrectionMatrix3D+; diff --git a/PWG0/dNdEta/AliMultiplicityCorrection.cxx b/PWG0/dNdEta/AliMultiplicityCorrection.cxx index 31fd28d5539..0c62c7ff536 100644 --- a/PWG0/dNdEta/AliMultiplicityCorrection.cxx +++ b/PWG0/dNdEta/AliMultiplicityCorrection.cxx @@ -2622,7 +2622,7 @@ void AliMultiplicityCorrection::ApplyLaszloMethod(Int_t inputRange, Bool_t fullP Double_t chi2 = 0; // use smoothed true (normalized) as new prior - Float_t norm = 1; //hTrueSmoothed->Integral(); + norm = 1; //hTrueSmoothed->Integral(); for (Int_t t=1; tGetNbinsX(); t++) { diff --git a/PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.cxx b/PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.cxx deleted file mode 100644 index 852f581f7f0..00000000000 --- a/PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.cxx +++ /dev/null @@ -1,229 +0,0 @@ -/* $Id$ */ - -#include "AlidNdEtaAnalysisESDSelector.h" - -#include -#include -#include -#include -#include -#include -#include - -#include -#include -#include - -#include "esdTrackCuts/AliESDtrackCuts.h" -#include "dNdEta/dNdEtaAnalysis.h" -#include "AliPWG0Helper.h" - -#include "AliGenEventHeader.h" -#include "AliHeader.h" -#include "AliStack.h" -#include "TParticle.h" - -ClassImp(AlidNdEtaAnalysisESDSelector) - -AlidNdEtaAnalysisESDSelector::AlidNdEtaAnalysisESDSelector() : - AliSelector(), - fdNdEtaAnalysis(0), - fMult(0), - fEsdTrackCuts(0) -{ - // - // Constructor. Initialization of pointers - // - - AliLog::SetClassDebugLevel("AlidNdEtaAnalysisESDSelector", AliLog::kDebug); -} - -AlidNdEtaAnalysisESDSelector::~AlidNdEtaAnalysisESDSelector() -{ - // - // Destructor - // - - // histograms are in the output list and deleted when the output - // list is deleted by the TSelector dtor -} - -void AlidNdEtaAnalysisESDSelector::Begin(TTree* tree) -{ - // Begin function - - ReadUserObjects(tree); -} - -void AlidNdEtaAnalysisESDSelector::ReadUserObjects(TTree* tree) -{ - // read the user objects, called from slavebegin and begin - - if (!fEsdTrackCuts && fInput) - fEsdTrackCuts = dynamic_cast (fInput->FindObject("AliESDtrackCuts")); - - if (!fEsdTrackCuts && tree) - fEsdTrackCuts = dynamic_cast (tree->GetUserInfo()->FindObject("AliESDtrackCuts")); - - if (!fEsdTrackCuts) - AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list."); -} - -void AlidNdEtaAnalysisESDSelector::SlaveBegin(TTree* tree) -{ - // The SlaveBegin() function is called after the Begin() function. - // When running with PROOF SlaveBegin() is called on each slave server. - // The tree argument is deprecated (on PROOF 0 is passed). - - AliSelector::SlaveBegin(tree); - - ReadUserObjects(tree); - - fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta"); - fMult = new TH1F("fMult", "fMult;Ntracks;Count", 201, -0.5, 200.5); -} - -void AlidNdEtaAnalysisESDSelector::Init(TTree* tree) -{ - // read the user objects - - AliSelector::Init(tree); - - // Enable only the needed branches - if (tree) - { - tree->SetBranchStatus("*", 0); - tree->SetBranchStatus("fTriggerMask", 1); - tree->SetBranchStatus("fSPDVertex*", 1); - tree->SetBranchStatus("fTracks.fLabel", 1); - - AliESDtrackCuts::EnableNeededBranches(tree); - } -} - -Bool_t AlidNdEtaAnalysisESDSelector::Process(Long64_t entry) -{ - // loop over all events - - if (AliSelector::Process(entry) == kFALSE) - return kFALSE; - - // Check prerequisites - if (!fESD) - { - AliDebug(AliLog::kError, "ESD branch not available"); - return kFALSE; - } - - if (!fEsdTrackCuts) - { - AliDebug(AliLog::kError, "fESDTrackCuts not available"); - return kFALSE; - } - - if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE) - { - AliDebug(AliLog::kDebug+1, Form("Skipping event %d because it was not triggered", (Int_t) entry)); - return kTRUE; - } - - if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE) - { - AliDebug(AliLog::kDebug+1, Form("Skipping event %d because its vertex was not reconstructed", (Int_t) entry)); - return kTRUE; - } - - // get the ESD vertex - const AliESDVertex* vtxESD = fESD->GetVertex(); - Double_t vtx[3]; - vtxESD->GetXYZ(vtx); - - // get number of "good" tracks - TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD); - Int_t nGoodTracks = list->GetEntries(); - - // loop over esd tracks - for (Int_t t=0; t (list->At(t)); - if (!esdTrack) - { - AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t)); - continue; - } - - Double_t p[3]; - esdTrack->GetConstrainedPxPyPz(p); // ### TODO should be okay because we have a vertex, however GetInnerPxPyPy / GetOuterPxPyPy also exist - TVector3 vector(p); - - Float_t theta = vector.Theta(); - Float_t eta = -TMath::Log(TMath::Tan(theta/2.)); - Float_t pt = vector.Pt(); - - fdNdEtaAnalysis->FillTrack(vtx[2], eta, pt); - } // end of track loop - - delete list; - list = 0; - - // for event count per vertex - fdNdEtaAnalysis->FillEvent(vtx[2], nGoodTracks); - fMult->Fill(nGoodTracks); - - return kTRUE; -} - -void AlidNdEtaAnalysisESDSelector::SlaveTerminate() -{ - // The SlaveTerminate() function is called after all entries or objects - // have been processed. When running with PROOF SlaveTerminate() is called - // on each slave server. - - AliSelector::SlaveTerminate(); - - // Add the histograms to the output on each slave server - if (!fOutput) - { - AliDebug(AliLog::kError, Form("ERROR: Output list not initialized.")); - return; - } - - // Add the objects to the output list and set them to 0, so that the destructor does not delete them. - - fOutput->Add(fdNdEtaAnalysis); - fOutput->Add(fMult); - - fdNdEtaAnalysis = 0; -} - -void AlidNdEtaAnalysisESDSelector::Terminate() -{ - // The Terminate() function is the last function to be called during - // a query. It always runs on the client, it can be used to present - // the results graphically or save the results to file. - - AliSelector::Terminate(); - - fdNdEtaAnalysis = dynamic_cast (fOutput->FindObject("dndeta")); - fMult = dynamic_cast (fOutput->FindObject("fMult")); - - if (!fdNdEtaAnalysis) - { - AliDebug(AliLog::kError, "ERROR: Histograms not available"); - return; - } - - TFile* fout = new TFile("analysis_esd_raw.root", "RECREATE"); - - if (fdNdEtaAnalysis) - fdNdEtaAnalysis->SaveHistograms(); - - if (fEsdTrackCuts) - fEsdTrackCuts->SaveHistograms("esd_tracks_cuts"); - - if (fMult) - fMult->Write(); - - fout->Write(); - fout->Close(); -} diff --git a/PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.h b/PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.h deleted file mode 100644 index 517f984f5ca..00000000000 --- a/PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.h +++ /dev/null @@ -1,39 +0,0 @@ -/* $Id$ */ - -#ifndef ALIDNDETAANALYSISESDSELECTOR_H -#define ALIDNDETAANALYSISESDSELECTOR_H - -#include "AliSelector.h" - -class AliESDtrackCuts; -class dNdEtaAnalysis; -class TH1F; - -class AlidNdEtaAnalysisESDSelector : public AliSelector { - public: - AlidNdEtaAnalysisESDSelector(); - virtual ~AlidNdEtaAnalysisESDSelector(); - - virtual void Begin(TTree* tree); - virtual void SlaveBegin(TTree *tree); - virtual void Init(TTree *tree); - virtual Bool_t Process(Long64_t entry); - virtual void SlaveTerminate(); - virtual void Terminate(); - - protected: - void ReadUserObjects(TTree* tree); - - dNdEtaAnalysis* fdNdEtaAnalysis; // contains the uncorrected histograms - TH1F* fMult; // raw multiplicity histogram (control histogram) - - AliESDtrackCuts* fEsdTrackCuts; // Object containing the parameters of the esd track cuts - - private: - AlidNdEtaAnalysisESDSelector(const AlidNdEtaAnalysisESDSelector&); - AlidNdEtaAnalysisESDSelector& operator=(const AlidNdEtaAnalysisESDSelector&); - - ClassDef(AlidNdEtaAnalysisESDSelector, 0); -}; - -#endif diff --git a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx b/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx deleted file mode 100644 index 4b57713ffbf..00000000000 --- a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx +++ /dev/null @@ -1,345 +0,0 @@ -/* $Id$ */ - -#include "AlidNdEtaAnalysisMCSelector.h" - -#include -#include -#include -#include -#include -#include -#include -#include -#include -#include - -#include -#include -#include -#include -#include - -#include "dNdEta/dNdEtaAnalysis.h" -#include -#include -#include -#include "esdTrackCuts/AliESDtrackCuts.h" - - -ClassImp(AlidNdEtaAnalysisMCSelector) - -AlidNdEtaAnalysisMCSelector::AlidNdEtaAnalysisMCSelector() : - AliSelectorRL(), - fEsdTrackCuts(0), - fdNdEtaAnalysis(0), - fdNdEtaAnalysisTr(0), - fdNdEtaAnalysisTrVtx(0), - fdNdEtaAnalysisTracks(0), - fVertex(0), - fPartPt(0), - fEvents(0) -{ - // - // Constructor. Initialization of pointers - // -} - -AlidNdEtaAnalysisMCSelector::~AlidNdEtaAnalysisMCSelector() -{ - // - // Destructor - // -} - -void AlidNdEtaAnalysisMCSelector::Begin(TTree* tree) -{ - // Begin function - - ReadUserObjects(tree); -} - -void AlidNdEtaAnalysisMCSelector::ReadUserObjects(TTree* tree) -{ - // read the user objects, called from slavebegin and begin - - if (!fEsdTrackCuts && fInput) - fEsdTrackCuts = dynamic_cast (fInput->FindObject("AliESDtrackCuts")); - - if (!fEsdTrackCuts && tree) - fEsdTrackCuts = dynamic_cast (tree->GetUserInfo()->FindObject("AliESDtrackCuts")); - - if (!fEsdTrackCuts) - AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list."); -} - -void AlidNdEtaAnalysisMCSelector::SlaveBegin(TTree * tree) -{ - // The SlaveBegin() function is called after the Begin() function. - // When running with PROOF SlaveBegin() is called on each slave server. - // The tree argument is deprecated (on PROOF 0 is passed). - - AliSelectorRL::SlaveBegin(tree); - - ReadUserObjects(tree); - - fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta"); - fdNdEtaAnalysisTr = new dNdEtaAnalysis("dndetaTr", "dndetaTr"); - fdNdEtaAnalysisTrVtx = new dNdEtaAnalysis("dndetaTrVtx", "dndetaTrVtx"); - fdNdEtaAnalysisTracks = new dNdEtaAnalysis("dndetaTracks", "dndetaTracks"); - fVertex = new TH3F("vertex_check", "vertex_check", 50, -50, 50, 50, -50, 50, 50, -50, 50); - for (Int_t i=0; i<3; ++i) - { - fPartEta[i] = new TH1F(Form("dndeta_check_%d", i), Form("dndeta_check_%d", i), 60, -6, 6); - fPartEta[i]->Sumw2(); - } - fPartPt = new TH1F("dndeta_check_pt", "dndeta_check_pt", 1000, 0, 10); - fPartPt->Sumw2(); - fEvents = new TH1F("dndeta_check_vertex", "dndeta_check_vertex", 40, -20, 20); -} - -void AlidNdEtaAnalysisMCSelector::Init(TTree *tree) -{ - AliSelectorRL::Init(tree); - - if (!tree) - { - AliDebug(AliLog::kError, "tree not available"); - return; - } - - tree->SetBranchStatus("*", 0); - tree->SetBranchStatus("fTriggerMask", 1); - tree->SetBranchStatus("fSPDVertex*", 1); - tree->SetBranchStatus("fTracks.fLabel", 1); - - AliESDtrackCuts::EnableNeededBranches(tree); -} - -Bool_t AlidNdEtaAnalysisMCSelector::Process(Long64_t entry) -{ - // fill the dNdEtaAnalysis class from the monte carlo - - if (AliSelectorRL::Process(entry) == kFALSE) - return kFALSE; - - // Check prerequisites - if (!fESD) - { - AliDebug(AliLog::kError, "ESD branch not available"); - return kFALSE; - } - - AliStack* stack = GetStack(); - if (!stack) - { - AliDebug(AliLog::kError, "Stack not available"); - return kFALSE; - } - - AliHeader* header = GetHeader(); - if (!header) - { - AliDebug(AliLog::kError, "Header not available"); - return kFALSE; - } - - Bool_t vertexReconstructed = AliPWG0Helper::IsVertexReconstructed(fESD); - Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD); - - // get the MC vertex - AliGenEventHeader* genHeader = header->GenEventHeader(); - - TArrayF vtxMC(3); - genHeader->PrimaryVertex(vtxMC); - - // loop over mc particles - Int_t nPrim = stack->GetNprimary(); - - Int_t nAcceptedParticles = 0; - - for (Int_t iMc = 0; iMc < nPrim; ++iMc) - { - TParticle* particle = stack->Particle(iMc); - - if (!particle) - continue; - - if (AliPWG0Helper::IsPrimaryCharged(particle, nPrim) == kFALSE) - continue; - - AliDebug(AliLog::kDebug+1, Form("Accepted primary %d, unique ID: %d", iMc, particle->GetUniqueID())); - ++nAcceptedParticles; - - fdNdEtaAnalysis->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); - fVertex->Fill(particle->Vx(), particle->Vy(), particle->Vz()); - - if (eventTriggered) - { - fdNdEtaAnalysisTr->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); - if (vertexReconstructed) - fdNdEtaAnalysisTrVtx->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); - } - - if (TMath::Abs(vtxMC[2]) < 20) - { - fPartEta[0]->Fill(particle->Eta()); - - if (vtxMC[2] < 0) - fPartEta[1]->Fill(particle->Eta()); - else - fPartEta[2]->Fill(particle->Eta()); - } - - if (TMath::Abs(particle->Eta()) < 0.8) - fPartPt->Fill(particle->Pt()); - } - - fEvents->Fill(vtxMC[2]); - - fdNdEtaAnalysis->FillEvent(vtxMC[2], nAcceptedParticles); - if (eventTriggered) - { - fdNdEtaAnalysisTr->FillEvent(vtxMC[2], nAcceptedParticles); - if (vertexReconstructed) - fdNdEtaAnalysisTrVtx->FillEvent(vtxMC[2], nAcceptedParticles); - } - - if (!eventTriggered || !vertexReconstructed) - return kTRUE; - - // from tracks is only done for triggered and vertex reconstructed events - - // get number of "good" tracks - TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD); - Int_t nGoodTracks = list->GetEntries(); - - // loop over esd tracks - for (Int_t t=0; t (list->At(t)); - if (!esdTrack) - { - AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t)); - continue; - } - - Int_t label = TMath::Abs(esdTrack->GetLabel()); - - TParticle* particle = stack->Particle(label); - if (!particle) - { - AliDebug(AliLog::kError, Form("ERROR: Could not retrieve particle %d.", esdTrack->GetLabel())); - continue; - } - - fdNdEtaAnalysisTracks->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); - } // end of track loop - - delete list; - list = 0; - - // for event count per vertex - fdNdEtaAnalysisTracks->FillEvent(vtxMC[2], nGoodTracks); - - return kTRUE; -} - -void AlidNdEtaAnalysisMCSelector::SlaveTerminate() -{ - // The SlaveTerminate() function is called after all entries or objects - // have been processed. When running with PROOF SlaveTerminate() is called - // on each slave server. - - AliSelectorRL::SlaveTerminate(); - - // Add the histograms to the output on each slave server - if (!fOutput) - { - AliDebug(AliLog::kError, Form("ERROR: Output list not initialized.")); - return; - } - - fOutput->Add(fdNdEtaAnalysis); - fOutput->Add(fdNdEtaAnalysisTr); - fOutput->Add(fdNdEtaAnalysisTrVtx); - fOutput->Add(fdNdEtaAnalysisTracks); - - fOutput->Add(fPartPt); - fOutput->Add(fEvents); - for (Int_t i=0; i<3; ++i) - fOutput->Add(fPartEta[i]); -} - -void AlidNdEtaAnalysisMCSelector::Terminate() -{ - // - - AliSelectorRL::Terminate(); - - fdNdEtaAnalysis = dynamic_cast (fOutput->FindObject("dndeta")); - fdNdEtaAnalysisTr = dynamic_cast (fOutput->FindObject("dndetaTr")); - fdNdEtaAnalysisTrVtx = dynamic_cast (fOutput->FindObject("dndetaTrVtx")); - fdNdEtaAnalysisTracks = dynamic_cast (fOutput->FindObject("dndetaTracks")); - fPartPt = dynamic_cast (fOutput->FindObject("dndeta_check_pt")); - fEvents = dynamic_cast (fOutput->FindObject("dndeta_check_vertex")); - for (Int_t i=0; i<3; ++i) - fPartEta[i] = dynamic_cast (fOutput->FindObject(Form("dndeta_check_%d", i))); - - if (!fdNdEtaAnalysis || !fdNdEtaAnalysisTr || !fdNdEtaAnalysisTrVtx || !fPartPt || !fEvents) - { - AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p", (void*) fdNdEtaAnalysis, (void*) fPartPt)); - return; - } - - fdNdEtaAnalysis->Finish(0, -1, AlidNdEtaCorrection::kNone); - fdNdEtaAnalysisTr->Finish(0, -1, AlidNdEtaCorrection::kNone); - fdNdEtaAnalysisTrVtx->Finish(0, -1, AlidNdEtaCorrection::kNone); - fdNdEtaAnalysisTracks->Finish(0, -1, AlidNdEtaCorrection::kNone); - - Int_t events = (Int_t) fdNdEtaAnalysis->GetData()->GetEventCorrection()->GetMeasuredHistogram()->Integral(); - fPartPt->Scale(1.0/events); - fPartPt->Scale(1.0/fPartPt->GetBinWidth(1)); - - if (fPartEta[0]) - { - Int_t events1 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(-19.9), fEvents->GetXaxis()->FindBin(-0.001)); - Int_t events2 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(0.001), fEvents->GetXaxis()->FindBin(19.9)); - - fPartEta[0]->Scale(1.0 / (events1 + events2)); - fPartEta[1]->Scale(1.0 / events1); - fPartEta[2]->Scale(1.0 / events2); - - for (Int_t i=0; i<3; ++i) - fPartEta[i]->Scale(1.0 / fPartEta[i]->GetBinWidth(1)); - - new TCanvas("control", "control", 500, 500); - for (Int_t i=0; i<3; ++i) - { - fPartEta[i]->SetLineColor(i+1); - fPartEta[i]->Draw((i==0) ? "" : "SAME"); - } - } - - TFile* fout = new TFile("analysis_mc.root","RECREATE"); - - fdNdEtaAnalysis->SaveHistograms(); - fdNdEtaAnalysisTr->SaveHistograms(); - fdNdEtaAnalysisTrVtx->SaveHistograms(); - fdNdEtaAnalysisTracks->SaveHistograms(); - fPartPt->Write(); - - fout->Write(); - fout->Close(); - - if (fPartPt) - { - new TCanvas("control2", "control2", 500, 500); - fPartPt->DrawCopy(); - } - - if (fEvents) - { - new TCanvas("control3", "control3", 500, 500); - fEvents->DrawCopy(); - } -} diff --git a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h b/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h deleted file mode 100644 index aade0a83fca..00000000000 --- a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h +++ /dev/null @@ -1,48 +0,0 @@ -/* $Id$ */ - -#ifndef ALIDNDETAANALYSISSELECTORMC_H -#define ALIDNDETAANALYSISSELECTORMC_H - -#include "AliSelectorRL.h" - -class TH3F; -class TH1F; -class dNdEtaAnalysis; -class AliESDtrackCuts; - -class AlidNdEtaAnalysisMCSelector : public AliSelectorRL { - public: - AlidNdEtaAnalysisMCSelector(); - virtual ~AlidNdEtaAnalysisMCSelector(); - - virtual void Begin(TTree *tree); - virtual void SlaveBegin(TTree *tree); - virtual void SlaveTerminate(); - virtual void Init(TTree *tree); - virtual Bool_t Process(Long64_t entry); - virtual void Terminate(); - - protected: - void ReadUserObjects(TTree* tree); - - private: - AliESDtrackCuts* fEsdTrackCuts; // Object containing the parameters of the esd track cuts, needed for fdNdEtaAnalysisTracks - - dNdEtaAnalysis* fdNdEtaAnalysis; // contains the dndeta from the full sample - dNdEtaAnalysis* fdNdEtaAnalysisTr; // contains the dndeta from the triggered events - dNdEtaAnalysis* fdNdEtaAnalysisTrVtx; // contains the dndeta from the triggered events with vertex - dNdEtaAnalysis* fdNdEtaAnalysisTracks; // contains the dndeta from the triggered events with vertex counted from the mc particles associated to the tracks (comparing this to the raw values from the esd shows the effect of the detector resolution) - - // the following are control histograms to check the dNdEtaAnalysis class - TH3F* fVertex; // vertex of counted particles - TH1F* fPartEta[3]; // counted particles as function of eta (full vertex range, below 0 range, above 0 range) - TH1F* fPartPt; // counted particles as function of pt - TH1F* fEvents; // events counted as function of vtx - - AlidNdEtaAnalysisMCSelector(const AlidNdEtaAnalysisMCSelector&); - AlidNdEtaAnalysisMCSelector& operator=(const AlidNdEtaAnalysisMCSelector&); - - ClassDef(AlidNdEtaAnalysisMCSelector, 0); -}; - -#endif diff --git a/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx b/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx index 2048da79c9a..11b771769a8 100644 --- a/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx +++ b/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx @@ -23,7 +23,7 @@ #include #include -#include "esdTrackCuts/AliESDtrackCuts.h" +#include "AliESDtrackCuts.h" #include "AliPWG0Helper.h" //#include "AliCorrection.h" //#include "AliCorrectionMatrix3D.h" @@ -517,9 +517,9 @@ void AlidNdEtaCorrectionTask::Exec(Option_t*) else { fdNdEtaCorrection->FillTrackedParticle(vtxMC[2], etaArr[i], ptArr[i]); - fEtaProfile->Fill(particle->Eta(), particle->Eta() - etaArr[i]); } - + + fEtaProfile->Fill(particle->Eta(), particle->Eta() - etaArr[i]); fdNdEtaAnalysisESD->FillTrack(vtxMC[2], particle->Eta(), particle->Pt()); fEtaCorrelation->Fill(etaArr[i], particle->Eta()); @@ -544,6 +544,7 @@ void AlidNdEtaCorrectionTask::Exec(Option_t*) fdNdEtaCorrectionProcessType[2]->FillTrackedParticle(vtxMC[2], particle->Eta(), particle->Pt()); } + // control histograms Int_t hist = -1; if (label == label2) { diff --git a/PWG0/dNdEta/AlidNdEtaSystematicsSelector.cxx b/PWG0/dNdEta/AlidNdEtaSystematicsSelector.cxx index 276191927f0..0a52ae7cb85 100644 --- a/PWG0/dNdEta/AlidNdEtaSystematicsSelector.cxx +++ b/PWG0/dNdEta/AlidNdEtaSystematicsSelector.cxx @@ -21,7 +21,7 @@ #include <../STEER/AliGenPythiaEventHeader.h> #include <../STEER/AliGenCocktailEventHeader.h> -#include "esdTrackCuts/AliESDtrackCuts.h" +#include "AliESDtrackCuts.h" #include "AliPWG0Helper.h" #include "dNdEta/AlidNdEtaCorrection.h" diff --git a/PWG0/dNdEta/AlidNdEtaTask.cxx b/PWG0/dNdEta/AlidNdEtaTask.cxx index 71fbb0e1d6e..e6ae1748a17 100644 --- a/PWG0/dNdEta/AlidNdEtaTask.cxx +++ b/PWG0/dNdEta/AlidNdEtaTask.cxx @@ -29,7 +29,7 @@ #include #include -#include "esdTrackCuts/AliESDtrackCuts.h" +#include "AliESDtrackCuts.h" #include "AliPWG0Helper.h" #include "AliCorrection.h" #include "AliCorrectionMatrix3D.h" @@ -259,7 +259,6 @@ void AlidNdEtaTask::Exec(Option_t*) if (fUseMCKine) { - stack = mcEvent->Stack(); if (!stack) { @@ -411,7 +410,7 @@ void AlidNdEtaTask::Exec(Option_t*) return; } - AliStack* stack = mcEvent->Stack(); + stack = mcEvent->Stack(); if (!stack) { AliDebug(AliLog::kError, "Stack not available"); @@ -644,7 +643,7 @@ void AlidNdEtaTask::Terminate(Option_t *) fPartPt->Scale(1.0/events); fPartPt->Scale(1.0/fPartPt->GetBinWidth(1)); - TFile* fout = new TFile("analysis_mc.root","RECREATE"); + fout = new TFile("analysis_mc.root","RECREATE"); fdNdEtaAnalysis->SaveHistograms(); fdNdEtaAnalysisNSD->SaveHistograms(); diff --git a/PWG0/dNdEta/dNdEtaAnalysis.cxx b/PWG0/dNdEta/dNdEtaAnalysis.cxx index 53f50eb366d..2c629b0c7d1 100644 --- a/PWG0/dNdEta/dNdEtaAnalysis.cxx +++ b/PWG0/dNdEta/dNdEtaAnalysis.cxx @@ -13,6 +13,7 @@ #include #include #include +#include #include "AlidNdEtaCorrection.h" #include @@ -28,7 +29,9 @@ dNdEtaAnalysis::dNdEtaAnalysis() : TNamed(), fData(0), fMult(0), - fPtDist(0) + fPtDist(0), + fAnalysisMode(AliPWG0Helper::kInvalid), + fTag() { // default constructor @@ -44,7 +47,9 @@ dNdEtaAnalysis::dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::Analy TNamed(name, title), fData(0), fMult(0), - fPtDist(0) + fPtDist(0), + fAnalysisMode(analysisMode), + fTag() { // constructor @@ -59,7 +64,8 @@ dNdEtaAnalysis::dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::Analy fMult = new TH1F("TriggeredMultiplicity", "Triggered Events;raw multiplicity;entries", 1000, -0.5, 999.5); - fdNdEta[0] = new TH1F("dNdEta", "dN_{ch}/d#eta;#eta;dN_{ch}/d#eta", 40, -2, 2); + // TODO get binning from AliCorrection + fdNdEta[0] = new TH1F("dNdEta", "dN_{ch}/d#eta;#eta;dN_{ch}/d#eta", 20, -2, 2); fdNdEtaPtCutOffCorrected[0] = dynamic_cast (fdNdEta[0]->Clone(Form("%s_corrected", fdNdEta[0]->GetName()))); @@ -117,7 +123,9 @@ dNdEtaAnalysis::dNdEtaAnalysis(const dNdEtaAnalysis &c) : TNamed(c), fData(0), fMult(0), - fPtDist(0) + fPtDist(0), + fAnalysisMode(AliPWG0Helper::kInvalid), + fTag() { // // dNdEtaAnalysis copy constructor @@ -157,6 +165,9 @@ void dNdEtaAnalysis::Copy(TObject &c) const target.fPtDist = dynamic_cast (fPtDist->Clone()); + target.fAnalysisMode = fAnalysisMode; + target.fTag = fTag; + TNamed::Copy((TNamed &) c); } @@ -193,21 +204,14 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid // the measured result is not used up to a multiplicity of multCut (the bin at multCut is the first that is used!) // - // TODO put tag somewhere which corrections have been applied - - Printf("\n\nCorrecting dN/deta spectrum >>> %s <<<. Correction type: %d, pt cut: %.2f.", tag, (Int_t) correctionType, ptCut); + fTag.Form("Correcting dN/deta spectrum >>> %s <<<. Correction type: %d, pt cut: %.2f.", tag, (Int_t) correctionType, ptCut); + Printf("\n\n%s", fTag.Data()); // set corrections to 1 fData->SetCorrectionToUnity(); if (correction && correctionType != AlidNdEtaCorrection::kNone) { - /* Printf("FAKE: Rebinning. For test only"); - correction->GetVertexRecoCorrection()->GetEventCorrection()->Rebin(2, 1); - correction->GetTriggerBiasCorrectionINEL()->GetEventCorrection()->Rebin(2, 1); - correction->GetTriggerBiasCorrectionNSD()->GetEventCorrection()->Rebin(2, 1); - fData->GetEventCorrection()->Rebin(2, 1);*/ - TH3F* trackCorr = fData->GetTrackCorrection()->GetCorrectionHistogram(); TH2F* eventCorr = fData->GetEventCorrection()->GetCorrectionHistogram(); @@ -266,7 +270,7 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid //new TCanvas; correctedEvents->DrawCopy("TEXT"); // start above 0 mult. bin with integration - TH1* vertexDist = correctedEvents->ProjectionX("vertexdist_measured", 2); + TH1* vertexDist = measuredEvents->ProjectionX("vertexdist_measured", 2, measuredEvents->GetNbinsY()+1); //new TCanvas; vertexDist->DrawCopy(); Int_t allEventsWithVertex = (Int_t) vertexDist->Integral(0, vertexDist->GetNbinsX()+1); // include under/overflow! @@ -276,13 +280,13 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid for (Int_t i = 1; i <= measuredEvents->GetNbinsX(); i++) { - Double_t alpha = vertexDist->GetBinContent(i) / allEventsWithVertex; + Double_t alpha = (Double_t) vertexDist->GetBinContent(i) / allEventsWithVertex; Double_t events = alpha * triggeredEventsWith0Mult; // multiply with trigger correction if set above events *= fData->GetEventCorrection()->GetCorrectionHistogram()->GetBinContent(i, 1); - Printf("Bin %d, alpha is %.2f, number of events with 0 mult.: %.2f", i, alpha, events); + Printf("Bin %d, alpha is %.2f, number of events with 0 mult.: %.2f", i, alpha * 100., events); correctedEvents->SetBinContent(i, 1, events); } @@ -357,32 +361,35 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid return; } - //new TCanvas; vtxVsEta->DrawCopy("COLZ"); - //vtxVsEta->Rebin2D(1, 4); + const Float_t vertexRangeBegin[kVertexBinning] = { -9.99, -9.99, 0 }; + const Float_t vertexRangeEnd[kVertexBinning] = { 9.99, 0, 9.99 }; - const Float_t vertexRange = 9.99; + // TODO acceptance range binning-independent; only for SPD + const Int_t binBegin[20] = { 11, 10, 9, 7, 5, 4, 3, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1 }; for (Int_t iEta=1; iEta<=vtxVsEta->GetNbinsY(); iEta++) { - // do we have several histograms for different vertex positions? - Int_t vertexBinGlobalBegin = vertexHist->GetXaxis()->FindBin(-vertexRange); - Int_t vertexBinWidth = (vertexHist->GetXaxis()->FindBin(vertexRange) - vertexBinGlobalBegin + 1) / (kVertexBinning-1); - //printf("vertexBinGlobalBegin = %d, vertexBinWidth = %d\n", vertexBinGlobalBegin, vertexBinWidth); + // loop over vertex ranges for (Int_t vertexPos=0; vertexPosGetXaxis()->FindBin(vertexRangeBegin[vertexPos]); + Int_t vertexBinEnd = vertexHist->GetXaxis()->FindBin(vertexRangeEnd[vertexPos]); - // the first histogram is always for the whole vertex range - if (vertexPos > 0) + // adjust acceptance range + + if (fAnalysisMode == AliPWG0Helper::kSPD) { - vertexBinBegin = vertexBinGlobalBegin + vertexBinWidth * (vertexPos-1); - vertexBinEnd = vertexBinBegin + vertexBinWidth; + //printf("Eta bin: %d Vertex range: %d; before: %d %d ", iEta, vertexPos, vertexBinBegin, vertexBinEnd); + vertexBinBegin = TMath::Max(vertexBinBegin, binBegin[iEta - 1]); + vertexBinEnd = TMath::Min(vertexBinEnd, 15 - binBegin[20 - iEta]); + //Printf("after: %d %d", vertexBinBegin, vertexBinEnd); } - //printf("vertexBinBegin = %d, vertexBinEnd = %d\n", vertexBinBegin, vertexBinEnd); + // no data for this bin + if (vertexBinBegin > vertexBinEnd) + continue; - Float_t totalEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd - 1); + Float_t totalEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd); if (totalEvents == 0) { printf("WARNING: No events for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd); @@ -391,14 +398,14 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid Float_t sum = 0; Float_t sumError2 = 0; - for (Int_t iVtx = vertexBinBegin; iVtx < vertexBinEnd; iVtx++) + for (Int_t iVtx = vertexBinBegin; iVtx <= vertexBinEnd; iVtx++) { if (vtxVsEta->GetBinContent(iVtx, iEta) != 0) { sum = sum + vtxVsEta->GetBinContent(iVtx, iEta); - if (sumError2 > 10e30) - printf("WARNING: sum of error2 is dangerously large - be prepared for crash... "); + if (sumError2 > 10e30) + Printf("WARNING: sum of error2 is dangerously large - be prepared for crash... "); sumError2 = sumError2 + TMath::Power(vtxVsEta->GetBinError(iVtx, iEta),2); } @@ -417,7 +424,7 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid //printf("Eta: %d Vertex Range: %d %d, Event Count %f, Track Sum: %f, Track Sum corrected: %f\n", iEta, vertexBinBegin, vertexBinEnd, totalEvents, sum, sum / ptCutOffCorrection); Int_t bin = fdNdEta[vertexPos]->FindBin(vtxVsEta->GetYaxis()->GetBinCenter(iEta)); - if (bin > 0 && bin < fdNdEta[vertexPos]->GetNbinsX()) + if (bin > 0 && bin <= fdNdEta[vertexPos]->GetNbinsX()) { Float_t dndeta = sum / totalEvents; Float_t error = TMath::Sqrt(sumError2) / totalEvents; @@ -438,8 +445,6 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid } } } - - //new TCanvas; fdNdEta[0]->DrawCopy(); } //____________________________________________________________________ @@ -482,6 +487,12 @@ void dNdEtaAnalysis::SaveHistograms() printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fdNdEtaPtCutOffCorrected[%d] is 0\n", i); } + TNamed named("fTag", fTag.Data()); + named.Write(); + + TParameter param("fAnalysisMode", fAnalysisMode); + param.Write(); + gDirectory->cd("../"); } @@ -505,6 +516,12 @@ void dNdEtaAnalysis::LoadHistograms(const Char_t* dir) fPtDist = dynamic_cast (gDirectory->Get(fPtDist->GetName())); + if (dynamic_cast (gFile->Get("fTag"))) + fTag = (dynamic_cast (gFile->Get("fTag")))->GetTitle(); + + if (dynamic_cast*> (gFile->Get("fAnalysisMode"))) + fAnalysisMode = (AliPWG0Helper::AnalysisMode) (dynamic_cast*> (gFile->Get("fAnalysisMode")))->GetVal(); + gDirectory->cd("../"); } diff --git a/PWG0/dNdEta/dNdEtaAnalysis.h b/PWG0/dNdEta/dNdEtaAnalysis.h index 94aee8527c3..027f13881da 100644 --- a/PWG0/dNdEta/dNdEtaAnalysis.h +++ b/PWG0/dNdEta/dNdEtaAnalysis.h @@ -20,6 +20,7 @@ #include #include "AlidNdEtaCorrection.h" #include "AliPWG0Helper.h" +#include class TH1F; class TCollection; @@ -32,7 +33,7 @@ public: enum { kVertexBinning = 1+2 }; // the first is for the whole vertex range, the others divide the vertex range dNdEtaAnalysis(); - dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::AnalysisMode analysisMode = AliPWG0Helper::kTPC); + dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::AnalysisMode analysisMode = AliPWG0Helper::kInvalid); virtual ~dNdEtaAnalysis(); dNdEtaAnalysis(const dNdEtaAnalysis &c); @@ -67,7 +68,10 @@ protected: TH1F* fdNdEta[kVertexBinning]; // dndeta results for different vertex bins (0 = full range) TH1F* fdNdEtaPtCutOffCorrected[kVertexBinning]; // dndeta results for different vertex bins (0 = full range), pt cut off corrected - ClassDef(dNdEtaAnalysis, 1) + AliPWG0Helper::AnalysisMode fAnalysisMode; // detector (combination) used for the analysis + TString fTag; // tag saved that describes the applied correction + + ClassDef(dNdEtaAnalysis, 2) }; #endif diff --git a/PWG0/dNdEta/drawPlots.C b/PWG0/dNdEta/drawPlots.C index e188a69be9e..b08607b9aa6 100644 --- a/PWG0/dNdEta/drawPlots.C +++ b/PWG0/dNdEta/drawPlots.C @@ -303,7 +303,9 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) histESDMBVtxNoPt->SetMarkerStyle(22); histESDMBTracksNoPt->SetMarkerStyle(23); - Float_t etaLimit = 1.2999; + //Float_t etaLimit = 1.2999; + Float_t etaLimit = 2.01; + Float_t etaPlotLimit = 2.2; histESDMBVtx->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); histESDMB->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); @@ -318,7 +320,7 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) Float_t max = TMath::Max(histESDMBVtx->GetMaximum(), histESDMB->GetMaximum()); max = TMath::Max(max, histESD->GetMaximum()); - TH2F* dummy = new TH2F("dummy", "", 100, -1.5, 1.5, 1000, 0, max * 1.1); + TH2F* dummy = new TH2F("dummy", "", 100, -etaPlotLimit, etaPlotLimit, 1000, 0, max * 1.1); Prepare1DPlot(dummy); dummy->SetStats(kFALSE); dummy->SetXTitle("#eta"); @@ -346,9 +348,13 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) TFile* file2 = TFile::Open("analysis_mc.root"); TH1* histMC = (TH1*) file2->Get("dndeta/dNdEta_corrected")->Clone("cloned"); TH1* histMCnsd = (TH1*) file2->Get("dndetaNSD/dNdEta_corrected"); - // FAKE - if (!histMCnsd) - histMCnsd = histMC; + if (histMCnsd) + { + histMCnsd = (TH1*) histMCnsd->Clone("clonedNSD"); + } + else + histMCnsd = new TH1F; + TH1* histMCTr = (TH1*) file2->Get("dndetaTr/dNdEta_corrected")->Clone("cloned2"); TH1* histMCTrVtx = (TH1*) file2->Get("dndetaTrVtx/dNdEta_corrected")->Clone("cloned3"); @@ -374,8 +380,6 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) fdNdEtaAnalysis->Finish(0, 0.3, AlidNdEtaCorrection::kNone, "MC: Tracks w/o resolution effect"); TH1* histMCTracksPtCut = fdNdEtaAnalysis->GetdNdEtaHistogram(0); - TCanvas* canvas2 = new TCanvas("dNdEta2", "dNdEta2", 500, 500); - Prepare1DPlot(histMC); Prepare1DPlot(histMCnsd); Prepare1DPlot(histMCTr); @@ -384,8 +388,17 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) Prepare1DPlot(histMCPtCut); Prepare1DPlot(histMCTrPtCut); Prepare1DPlot(histMCTrVtxPtCut); - if (histMCTracksPtCut) - Prepare1DPlot(histMCTracksPtCut); + Prepare1DPlot(histMCTracksPtCut); + + histMC->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + histMCnsd->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + histMCTr->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + histMCTrVtx->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + + histMCPtCut->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + histMCTrPtCut->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + histMCTrVtxPtCut->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); + histMCTracksPtCut->GetXaxis()->SetRangeUser(-etaLimit, etaLimit); histMC->SetLineColor(1); histMCnsd->SetLineColor(6); @@ -398,8 +411,9 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) if (histMCTracksPtCut) histMCTracksPtCut->SetLineColor(4); + TCanvas* canvas2 = new TCanvas("dNdEta2", "dNdEta2", 500, 500); + TH2* dummy2 = (TH2F*) dummy->Clone("dummy2"); - Prepare1DPlot(dummy2); dummy2->GetYaxis()->SetRangeUser(0, max * 1.1); dummy2->DrawCopy(); @@ -427,6 +441,10 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) canvas2->SaveAs("dNdEta2.eps"); } + DrawdNdEtaRatio(histESD, histMC, "full_inelastic", etaPlotLimit); + DrawdNdEtaRatio(histESDMB, histMCTr, "triggered", etaPlotLimit); + DrawdNdEtaRatio(histESDMBVtx, histMCTrVtx, "triggered_vertex", etaPlotLimit); + new TCanvas; dummy2->DrawCopy(); histMCnsd->Draw("SAME"); @@ -503,9 +521,10 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) Float_t minR = TMath::Min(0.961, ratio->GetMinimum() * 0.95); Float_t maxR = TMath::Max(1.049, ratio->GetMaximum() * 1.05); - TH1F dummy3("dummy3", ";#eta;Ratio: MC / ESD", 1, -1.5, 1.5); + TH1F dummy3("dummy3", ";#eta;Ratio: MC / ESD", 100, -etaPlotLimit, etaPlotLimit); dummy3.SetStats(kFALSE); - dummy3.SetBinContent(1, 1); + for (Int_t i=1; i<=100; ++i) + dummy3.SetBinContent(i, 1); dummy3.GetYaxis()->SetRangeUser(minR, maxR); dummy3.SetLineWidth(2); dummy3.GetXaxis()->SetLabelSize(0.06); @@ -539,6 +558,77 @@ void dNdEta(Bool_t onlyESD = kFALSE, Bool_t save = kTRUE) legend->Draw(); } +void DrawdNdEtaRatio(TH1* corr, TH1* mc, const char* name, Float_t etaPlotLimit) +{ + TCanvas* canvas3 = new TCanvas(name, name, 700, 600); + canvas3->Range(0, 0, 1, 1); + + TPad* pad1 = new TPad(Form("%s_1", name), "", 0, 0.5, 0.98, 0.98); + pad1->Draw(); + + TPad* pad2 = new TPad(Form("%s_2", name), "", 0, 0.02, 0.98, 0.5); + pad2->Draw(); + + pad1->SetRightMargin(0.05); + pad2->SetRightMargin(0.05); + + // no border between them + pad1->SetBottomMargin(0); + pad2->SetTopMargin(0); + + pad1->cd(); + + TLegend* legend = new TLegend(0.4, 0.05, 0.65, 0.3); + legend->SetFillColor(0); + legend->AddEntry(corr, "corrected"); + legend->AddEntry(mc, "MC prediction"); + + TH2F* dummy = new TH2F("dummy", "", 100, -etaPlotLimit, etaPlotLimit, 1000, 0, corr->GetMaximum() * 1.1); + Prepare1DPlot(dummy); + dummy->SetStats(kFALSE); + dummy->SetXTitle("#eta"); + dummy->SetYTitle("dN_{ch}/d#eta"); + dummy->GetYaxis()->SetTitleOffset(1); + + dummy->GetXaxis()->SetLabelSize(0.06); + dummy->GetYaxis()->SetLabelSize(0.06); + dummy->GetXaxis()->SetTitleSize(0.06); + dummy->GetYaxis()->SetTitleSize(0.06); + dummy->GetYaxis()->SetTitleOffset(0.7); + dummy->DrawCopy(); + + corr->Draw("SAME"); + mc->Draw("SAME"); + + legend->Draw(); + + pad2->cd(); + pad2->SetBottomMargin(0.15); + + TH1* ratio = (TH1*) mc->Clone("ratio"); + ratio->Divide(corr); + + Float_t minR = TMath::Min(0.961, ratio->GetMinimum() * 0.95); + Float_t maxR = TMath::Max(1.049, ratio->GetMaximum() * 1.05); + + TH1F dummy3("dummy3", ";#eta;Ratio: MC / corr", 100, -etaPlotLimit, etaPlotLimit); + dummy3.SetStats(kFALSE); + for (Int_t i=1; i<=100; ++i) + dummy3.SetBinContent(i, 1); + dummy3.GetYaxis()->SetRangeUser(minR, maxR); + dummy3.SetLineWidth(2); + dummy3.GetXaxis()->SetLabelSize(0.06); + dummy3.GetYaxis()->SetLabelSize(0.06); + dummy3.GetXaxis()->SetTitleSize(0.06); + dummy3.GetYaxis()->SetTitleSize(0.06); + dummy3.GetYaxis()->SetTitleOffset(0.7); + dummy3.DrawCopy(); + + ratio->Draw("SAME"); + + canvas3->Modified(); +} + void ptSpectrum() { TFile* file = TFile::Open("analysis_esd.root"); diff --git a/PWG0/dNdEta/plotsMultiplicity.C b/PWG0/dNdEta/plotsMultiplicity.C index b2c8acb23d9..2f6a61bc3ad 100644 --- a/PWG0/dNdEta/plotsMultiplicity.C +++ b/PWG0/dNdEta/plotsMultiplicity.C @@ -1225,7 +1225,7 @@ void EfficiencySpecies() TCanvas* canvas = new TCanvas("EfficiencySpecies", "EfficiencySpecies", 1000, 500); canvas->Divide(2, 1); - for (Int_t loop=0; loop<2; ++loop) + for (Int_t loop=0; loop<1; ++loop) { Printf("%s", fileName[loop]); @@ -1252,6 +1252,9 @@ void EfficiencySpecies() return; } + Float_t sumGen = 0; + Float_t sumMeas = 0; + for (Int_t i=0; i<3; ++i) { Printf("correction %d", i); @@ -1287,6 +1290,9 @@ void EfficiencySpecies() genePt->Sumw2(); measPt->Sumw2(); + sumGen += genePt->Integral(); + sumMeas += measPt->Integral(); + TH1* effPt = (TH1*) genePt->Clone(Form("effPt_%d", i)); effPt->Reset(); effPt->Divide(measPt, genePt, 1, 1, "B"); @@ -1326,6 +1332,8 @@ void EfficiencySpecies() Printf("In total %.4f of the particles are below their effective pt cut off", (Float_t) below / total); + Printf("%f measured, %f generated, effiency: %f", sumGen, sumMeas, sumMeas / sumGen); + legend->Draw(); } diff --git a/PWG0/dNdEta/run.C b/PWG0/dNdEta/run.C index d967bc918b6..22758e10ae0 100644 --- a/PWG0/dNdEta/run.C +++ b/PWG0/dNdEta/run.C @@ -1,4 +1,4 @@ -void run(Int_t runWhat, Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, Int_t aProof = kFALSE, Bool_t mc = kTRUE, const char* option = "") +void run(Int_t runWhat, const Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, Int_t aProof = kFALSE, Bool_t mc = kTRUE, const char* option = "") { // runWhat options: 0 = AlidNdEtaTask // 1 = AlidNdEtaCorrectionTask @@ -60,6 +60,7 @@ void run(Int_t runWhat, Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDe gSystem->Load("libSTEERBase"); gSystem->Load("libESD"); gSystem->Load("libANALYSIS"); + gSystem->Load("libANALYSISalice"); gSystem->Load("libPWG0base"); } @@ -218,7 +219,7 @@ void FinishAnalysisAll(const char* dataInput = "analysis_esd_raw.root", const ch file->cd(); fdNdEtaAnalysis = new dNdEtaAnalysis("dndetaTrVtx", "dndetaTrVtx"); fdNdEtaAnalysis->LoadHistograms("fdNdEtaAnalysisESD"); - fdNdEtaAnalysis->Finish(dNdEtaCorrection, 0.3, AlidNdEtaCorrection::kTrack2Particle, "ESD -> MB with trigger"); + fdNdEtaAnalysis->Finish(dNdEtaCorrection, 0.3, AlidNdEtaCorrection::kTrack2Particle, "ESD -> MB with vertex"); //fdNdEtaAnalysis->DrawHistograms(kTRUE); file2->cd(); fdNdEtaAnalysis->SaveHistograms(); diff --git a/PWG0/esdTrackCuts/AliCutTask.cxx b/PWG0/esdTrackCuts/AliCutTask.cxx index 26b07942661..79977c92766 100644 --- a/PWG0/esdTrackCuts/AliCutTask.cxx +++ b/PWG0/esdTrackCuts/AliCutTask.cxx @@ -7,7 +7,7 @@ #include "TParticle.h" #include "TParticlePDG.h" -#include "esdTrackCuts/AliESDtrackCuts.h" +#include "AliESDtrackCuts.h" #include "AliPWG0Helper.h" #include "AliAnalysisTask.h" diff --git a/PWG0/esdTrackCuts/AliESDtrackCuts.cxx b/PWG0/esdTrackCuts/AliESDtrackCuts.cxx deleted file mode 100644 index 17f856d33ad..00000000000 --- a/PWG0/esdTrackCuts/AliESDtrackCuts.cxx +++ /dev/null @@ -1,1186 +0,0 @@ -/************************************************************************** - * Copyright(c) 1998-1999, 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. * - **************************************************************************/ - -/* $Id$ */ - -#include "AliESDtrackCuts.h" - -#include -#include -#include -#include - -#include -#include -#include - -//____________________________________________________________________ -ClassImp(AliESDtrackCuts) - -// Cut names -const Char_t* AliESDtrackCuts::fgkCutNames[kNCuts] = { - "require TPC refit", - "require ITS refit", - "n clusters TPC", - "n clusters ITS", - "#Chi^{2}/clusters TPC", - "#Chi^{2}/clusters ITS", - "cov 11", - "cov 22", - "cov 33", - "cov 44", - "cov 55", - "trk-to-vtx", - "trk-to-vtx failed", - "kink daughters", - "p", - "p_{T}", - "p_{x}", - "p_{y}", - "p_{z}", - "y", - "eta" -}; - -//____________________________________________________________________ -AliESDtrackCuts::AliESDtrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title), - fCutMinNClusterTPC(0), - fCutMinNClusterITS(0), - fCutMaxChi2PerClusterTPC(0), - fCutMaxChi2PerClusterITS(0), - fCutMaxC11(0), - fCutMaxC22(0), - fCutMaxC33(0), - fCutMaxC44(0), - fCutMaxC55(0), - fCutAcceptKinkDaughters(0), - fCutRequireTPCRefit(0), - fCutRequireITSRefit(0), - fCutNsigmaToVertex(0), - fCutSigmaToVertexRequired(0), - fPMin(0), - fPMax(0), - fPtMin(0), - fPtMax(0), - fPxMin(0), - fPxMax(0), - fPyMin(0), - fPyMax(0), - fPzMin(0), - fPzMax(0), - fEtaMin(0), - fEtaMax(0), - fRapMin(0), - fRapMax(0), - fHistogramsOn(0), - ffDTheoretical(0), - fhCutStatistics(0), - fhCutCorrelation(0) -{ - // - // constructor - // - - Init(); - - //############################################################################## - // setting default cuts - SetMinNClustersTPC(); - SetMinNClustersITS(); - SetMaxChi2PerClusterTPC(); - SetMaxChi2PerClusterITS(); - SetMaxCovDiagonalElements(); - SetRequireTPCRefit(); - SetRequireITSRefit(); - SetAcceptKingDaughters(); - SetMinNsigmaToVertex(); - SetRequireSigmaToVertex(); - SetPRange(); - SetPtRange(); - SetPxRange(); - SetPyRange(); - SetPzRange(); - SetEtaRange(); - SetRapRange(); - - SetHistogramsOn(); -} - -//_____________________________________________________________________________ -AliESDtrackCuts::AliESDtrackCuts(const AliESDtrackCuts &c) : AliAnalysisCuts(c), - fCutMinNClusterTPC(0), - fCutMinNClusterITS(0), - fCutMaxChi2PerClusterTPC(0), - fCutMaxChi2PerClusterITS(0), - fCutMaxC11(0), - fCutMaxC22(0), - fCutMaxC33(0), - fCutMaxC44(0), - fCutMaxC55(0), - fCutAcceptKinkDaughters(0), - fCutRequireTPCRefit(0), - fCutRequireITSRefit(0), - fCutNsigmaToVertex(0), - fCutSigmaToVertexRequired(0), - fPMin(0), - fPMax(0), - fPtMin(0), - fPtMax(0), - fPxMin(0), - fPxMax(0), - fPyMin(0), - fPyMax(0), - fPzMin(0), - fPzMax(0), - fEtaMin(0), - fEtaMax(0), - fRapMin(0), - fRapMax(0), - fHistogramsOn(0), - ffDTheoretical(0), - fhCutStatistics(0), - fhCutCorrelation(0) -{ - // - // copy constructor - // - - ((AliESDtrackCuts &) c).Copy(*this); -} - -AliESDtrackCuts::~AliESDtrackCuts() -{ - // - // destructor - // - - for (Int_t i=0; i<2; i++) { - - if (fhNClustersITS[i]) - delete fhNClustersITS[i]; - if (fhNClustersTPC[i]) - delete fhNClustersTPC[i]; - if (fhChi2PerClusterITS[i]) - delete fhChi2PerClusterITS[i]; - if (fhChi2PerClusterTPC[i]) - delete fhChi2PerClusterTPC[i]; - if (fhC11[i]) - delete fhC11[i]; - if (fhC22[i]) - delete fhC22[i]; - if (fhC33[i]) - delete fhC33[i]; - if (fhC44[i]) - delete fhC44[i]; - if (fhC55[i]) - delete fhC55[i]; - - if (fhDXY[i]) - delete fhDXY[i]; - if (fhDZ[i]) - delete fhDZ[i]; - if (fhDXYvsDZ[i]) - delete fhDXYvsDZ[i]; - - if (fhDXYNormalized[i]) - delete fhDXYNormalized[i]; - if (fhDZNormalized[i]) - delete fhDZNormalized[i]; - if (fhDXYvsDZNormalized[i]) - delete fhDXYvsDZNormalized[i]; - if (fhNSigmaToVertex[i]) - delete fhNSigmaToVertex[i]; - if (fhPt[i]) - delete fhPt[i]; - if (fhEta[i]) - delete fhEta[i]; - } - - if (ffDTheoretical) - delete ffDTheoretical; - - if (fhCutStatistics) - delete fhCutStatistics; - if (fhCutCorrelation) - delete fhCutCorrelation; -} - -void AliESDtrackCuts::Init() -{ - // - // sets everything to zero - // - - fCutMinNClusterTPC = 0; - fCutMinNClusterITS = 0; - - fCutMaxChi2PerClusterTPC = 0; - fCutMaxChi2PerClusterITS = 0; - - fCutMaxC11 = 0; - fCutMaxC22 = 0; - fCutMaxC33 = 0; - fCutMaxC44 = 0; - fCutMaxC55 = 0; - - fCutAcceptKinkDaughters = 0; - fCutRequireTPCRefit = 0; - fCutRequireITSRefit = 0; - - fCutNsigmaToVertex = 0; - fCutSigmaToVertexRequired = 0; - - fPMin = 0; - fPMax = 0; - fPtMin = 0; - fPtMax = 0; - fPxMin = 0; - fPxMax = 0; - fPyMin = 0; - fPyMax = 0; - fPzMin = 0; - fPzMax = 0; - fEtaMin = 0; - fEtaMax = 0; - fRapMin = 0; - fRapMax = 0; - - fHistogramsOn = kFALSE; - - for (Int_t i=0; i<2; ++i) - { - fhNClustersITS[i] = 0; - fhNClustersTPC[i] = 0; - - fhChi2PerClusterITS[i] = 0; - fhChi2PerClusterTPC[i] = 0; - - fhC11[i] = 0; - fhC22[i] = 0; - fhC33[i] = 0; - fhC44[i] = 0; - fhC55[i] = 0; - - fhDXY[i] = 0; - fhDZ[i] = 0; - fhDXYvsDZ[i] = 0; - - fhDXYNormalized[i] = 0; - fhDZNormalized[i] = 0; - fhDXYvsDZNormalized[i] = 0; - fhNSigmaToVertex[i] = 0; - - fhPt[i] = 0; - fhEta[i] = 0; - } - ffDTheoretical = 0; - - fhCutStatistics = 0; - fhCutCorrelation = 0; -} - -//_____________________________________________________________________________ -AliESDtrackCuts &AliESDtrackCuts::operator=(const AliESDtrackCuts &c) -{ - // - // Assignment operator - // - - if (this != &c) ((AliESDtrackCuts &) c).Copy(*this); - return *this; -} - -//_____________________________________________________________________________ -void AliESDtrackCuts::Copy(TObject &c) const -{ - // - // Copy function - // - - AliESDtrackCuts& target = (AliESDtrackCuts &) c; - - target.Init(); - - target.fCutMinNClusterTPC = fCutMinNClusterTPC; - target.fCutMinNClusterITS = fCutMinNClusterITS; - - target.fCutMaxChi2PerClusterTPC = fCutMaxChi2PerClusterTPC; - target.fCutMaxChi2PerClusterITS = fCutMaxChi2PerClusterITS; - - target.fCutMaxC11 = fCutMaxC11; - target.fCutMaxC22 = fCutMaxC22; - target.fCutMaxC33 = fCutMaxC33; - target.fCutMaxC44 = fCutMaxC44; - target.fCutMaxC55 = fCutMaxC55; - - target.fCutAcceptKinkDaughters = fCutAcceptKinkDaughters; - target.fCutRequireTPCRefit = fCutRequireTPCRefit; - target.fCutRequireITSRefit = fCutRequireITSRefit; - - target.fCutNsigmaToVertex = fCutNsigmaToVertex; - target.fCutSigmaToVertexRequired = fCutSigmaToVertexRequired; - - target.fPMin = fPMin; - target.fPMax = fPMax; - target.fPtMin = fPtMin; - target.fPtMax = fPtMax; - target.fPxMin = fPxMin; - target.fPxMax = fPxMax; - target.fPyMin = fPyMin; - target.fPyMax = fPyMax; - target.fPzMin = fPzMin; - target.fPzMax = fPzMax; - target.fEtaMin = fEtaMin; - target.fEtaMax = fEtaMax; - target.fRapMin = fRapMin; - target.fRapMax = fRapMax; - - target.fHistogramsOn = fHistogramsOn; - - for (Int_t i=0; i<2; ++i) - { - if (fhNClustersITS[i]) target.fhNClustersITS[i] = (TH1F*) fhNClustersITS[i]->Clone(); - if (fhNClustersTPC[i]) target.fhNClustersTPC[i] = (TH1F*) fhNClustersTPC[i]->Clone(); - - if (fhChi2PerClusterITS[i]) target.fhChi2PerClusterITS[i] = (TH1F*) fhChi2PerClusterITS[i]->Clone(); - if (fhChi2PerClusterTPC[i]) target.fhChi2PerClusterTPC[i] = (TH1F*) fhChi2PerClusterTPC[i]->Clone(); - - if (fhC11[i]) target.fhC11[i] = (TH1F*) fhC11[i]->Clone(); - if (fhC22[i]) target.fhC22[i] = (TH1F*) fhC22[i]->Clone(); - if (fhC33[i]) target.fhC33[i] = (TH1F*) fhC33[i]->Clone(); - if (fhC44[i]) target.fhC44[i] = (TH1F*) fhC44[i]->Clone(); - if (fhC55[i]) target.fhC55[i] = (TH1F*) fhC55[i]->Clone(); - - if (fhDXY[i]) target.fhDXY[i] = (TH1F*) fhDXY[i]->Clone(); - if (fhDZ[i]) target.fhDZ[i] = (TH1F*) fhDZ[i]->Clone(); - if (fhDXYvsDZ[i]) target.fhDXYvsDZ[i] = (TH2F*) fhDXYvsDZ[i]->Clone(); - - if (fhDXYNormalized[i]) target.fhDXYNormalized[i] = (TH1F*) fhDXYNormalized[i]->Clone(); - if (fhDZNormalized[i]) target.fhDZNormalized[i] = (TH1F*) fhDZNormalized[i]->Clone(); - if (fhDXYvsDZNormalized[i]) target.fhDXYvsDZNormalized[i] = (TH2F*) fhDXYvsDZNormalized[i]->Clone(); - if (fhNSigmaToVertex[i]) target.fhNSigmaToVertex[i] = (TH1F*) fhNSigmaToVertex[i]->Clone(); - - if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone(); - if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone(); - } - if (ffDTheoretical) target.ffDTheoretical = (TF1*) ffDTheoretical->Clone(); - - if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone(); - if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone(); - - TNamed::Copy(c); -} - -//_____________________________________________________________________________ -Long64_t AliESDtrackCuts::Merge(TCollection* list) { - // Merge a list of AliESDtrackCuts objects with this (needed for PROOF) - // Returns the number of merged objects (including this) - - if (!list) - return 0; - - if (list->IsEmpty()) - return 1; - - if (!fHistogramsOn) - return 0; - - TIterator* iter = list->MakeIterator(); - TObject* obj; - - - // collection of measured and generated histograms - Int_t count = 0; - while ((obj = iter->Next())) { - - AliESDtrackCuts* entry = dynamic_cast(obj); - if (entry == 0) - continue; - - if (!entry->fHistogramsOn) - continue; - - for (Int_t i=0; i<2; i++) { - - fhNClustersITS[i] ->Add(entry->fhNClustersITS[i] ); - fhNClustersTPC[i] ->Add(entry->fhNClustersTPC[i] ); - - fhChi2PerClusterITS[i] ->Add(entry->fhChi2PerClusterITS[i]); - fhChi2PerClusterTPC[i] ->Add(entry->fhChi2PerClusterTPC[i]); - - fhC11[i] ->Add(entry->fhC11[i] ); - fhC22[i] ->Add(entry->fhC22[i] ); - fhC33[i] ->Add(entry->fhC33[i] ); - fhC44[i] ->Add(entry->fhC44[i] ); - fhC55[i] ->Add(entry->fhC55[i] ); - - fhDXY[i] ->Add(entry->fhDXY[i] ); - fhDZ[i] ->Add(entry->fhDZ[i] ); - fhDXYvsDZ[i] ->Add(entry->fhDXYvsDZ[i] ); - - fhDXYNormalized[i] ->Add(entry->fhDXYNormalized[i] ); - fhDZNormalized[i] ->Add(entry->fhDZNormalized[i] ); - fhDXYvsDZNormalized[i] ->Add(entry->fhDXYvsDZNormalized[i]); - fhNSigmaToVertex[i] ->Add(entry->fhNSigmaToVertex[i]); - - fhPt[i] ->Add(entry->fhPt[i]); - fhEta[i] ->Add(entry->fhEta[i]); - } - - fhCutStatistics ->Add(entry->fhCutStatistics); - fhCutCorrelation ->Add(entry->fhCutCorrelation); - - count++; - } - - return count+1; -} - - -//____________________________________________________________________ -Float_t AliESDtrackCuts::GetSigmaToVertex(AliESDtrack* esdTrack) -{ - // Calculates the number of sigma to the vertex. - - Float_t b[2]; - Float_t bRes[2]; - Float_t bCov[3]; - esdTrack->GetImpactParameters(b,bCov); - if (bCov[0]<=0 || bCov[2]<=0) { - AliDebug(1, "Estimated b resolution lower or equal zero!"); - bCov[0]=0; bCov[2]=0; - } - bRes[0] = TMath::Sqrt(bCov[0]); - bRes[1] = TMath::Sqrt(bCov[2]); - - // ----------------------------------- - // How to get to a n-sigma cut? - // - // The accumulated statistics from 0 to d is - // - // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma) - // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma) - // - // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2) - // Can this be expressed in a different way? - - if (bRes[0] == 0 || bRes[1] ==0) - return -1; - - Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2)); - - // stupid rounding problem screws up everything: - // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :( - if (TMath::Exp(-d * d / 2) < 1e-10) - return 1000; - - d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2); - return d; -} - -void AliESDtrackCuts::EnableNeededBranches(TTree* tree) -{ - // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack - - tree->SetBranchStatus("fTracks.fFlags", 1); - tree->SetBranchStatus("fTracks.fITSncls", 1); - tree->SetBranchStatus("fTracks.fTPCncls", 1); - tree->SetBranchStatus("fTracks.fITSchi2", 1); - tree->SetBranchStatus("fTracks.fTPCchi2", 1); - tree->SetBranchStatus("fTracks.fC*", 1); - tree->SetBranchStatus("fTracks.fD", 1); - tree->SetBranchStatus("fTracks.fZ", 1); - tree->SetBranchStatus("fTracks.fCdd", 1); - tree->SetBranchStatus("fTracks.fCdz", 1); - tree->SetBranchStatus("fTracks.fCzz", 1); - tree->SetBranchStatus("fTracks.fP*", 1); - tree->SetBranchStatus("fTracks.fR*", 1); - tree->SetBranchStatus("fTracks.fKinkIndexes*", 1); -} - -//____________________________________________________________________ -Bool_t -AliESDtrackCuts::AcceptTrack(AliESDtrack* esdTrack) { - // - // figure out if the tracks survives all the track cuts defined - // - // the different quality parameter and kinematic values are first - // retrieved from the track. then it is found out what cuts the - // track did not survive and finally the cuts are imposed. - - // this function needs the following branches: - // fTracks.fFlags - // fTracks.fITSncls - // fTracks.fTPCncls - // fTracks.fITSchi2 - // fTracks.fTPCchi2 - // fTracks.fC //GetExternalCovariance - // fTracks.fD //GetImpactParameters - // fTracks.fZ //GetImpactParameters - // fTracks.fCdd //GetImpactParameters - // fTracks.fCdz //GetImpactParameters - // fTracks.fCzz //GetImpactParameters - // fTracks.fP //GetPxPyPz - // fTracks.fR //GetMass - // fTracks.fP //GetMass - // fTracks.fKinkIndexes - - UInt_t status = esdTrack->GetStatus(); - - // dummy array - Int_t fIdxInt[200]; - - // getting quality parameters from the ESD track - Int_t nClustersITS = esdTrack->GetITSclusters(fIdxInt); - Int_t nClustersTPC = esdTrack->GetTPCclusters(fIdxInt); - - - - Float_t chi2PerClusterITS = -1; - Float_t chi2PerClusterTPC = -1; - if (nClustersITS!=0) - chi2PerClusterITS = esdTrack->GetITSchi2()/Float_t(nClustersITS); - if (nClustersTPC!=0) - chi2PerClusterTPC = esdTrack->GetTPCchi2()/Float_t(nClustersTPC); - - Double_t extCov[15]; - esdTrack->GetExternalCovariance(extCov); - - // getting the track to vertex parameters - Float_t nSigmaToVertex = GetSigmaToVertex(esdTrack); - - // getting the kinematic variables of the track - // (assuming the mass is known) - Double_t p[3]; - esdTrack->GetPxPyPz(p); - Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2)); - Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2)); - Float_t energy = TMath::Sqrt(TMath::Power(esdTrack->GetMass(),2) + TMath::Power(momentum,2)); - - - //y-eta related calculations - Float_t eta = -100.; - Float_t y = -100.; - if((momentum != TMath::Abs(p[2]))&&(momentum != 0)) - eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2])); - if((energy != TMath::Abs(p[2]))&&(momentum != 0)) - y = 0.5*TMath::Log((energy + p[2])/(energy - p[2])); - - - //######################################################################## - // cut the track? - - Bool_t cuts[kNCuts]; - for (Int_t i=0; ifCutMaxChi2PerClusterTPC) - cuts[4]=kTRUE; - if (chi2PerClusterITS>fCutMaxChi2PerClusterITS) - cuts[5]=kTRUE; - if (extCov[0] > fCutMaxC11) - cuts[6]=kTRUE; - if (extCov[2] > fCutMaxC22) - cuts[7]=kTRUE; - if (extCov[5] > fCutMaxC33) - cuts[8]=kTRUE; - if (extCov[9] > fCutMaxC44) - cuts[9]=kTRUE; - if (extCov[14] > fCutMaxC55) - cuts[10]=kTRUE; - if (nSigmaToVertex > fCutNsigmaToVertex && fCutSigmaToVertexRequired) - cuts[11] = kTRUE; - // if n sigma could not be calculated - if (nSigmaToVertex<0 && fCutSigmaToVertexRequired) - cuts[12]=kTRUE; - if (!fCutAcceptKinkDaughters && esdTrack->GetKinkIndex(0)>0) - cuts[13]=kTRUE; - // track kinematics cut - if((momentum < fPMin) || (momentum > fPMax)) - cuts[14]=kTRUE; - if((pt < fPtMin) || (pt > fPtMax)) - cuts[15] = kTRUE; - if((p[0] < fPxMin) || (p[0] > fPxMax)) - cuts[16] = kTRUE; - if((p[1] < fPyMin) || (p[1] > fPyMax)) - cuts[17] = kTRUE; - if((p[2] < fPzMin) || (p[2] > fPzMax)) - cuts[18] = kTRUE; - if((eta < fEtaMin) || (eta > fEtaMax)) - cuts[19] = kTRUE; - if((y < fRapMin) || (y > fRapMax)) - cuts[20] = kTRUE; - - Bool_t cut=kFALSE; - for (Int_t i=0; iFill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks"))); - - if (cut) - fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks"))); - - for (Int_t i=0; iFill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i]))); - - for (Int_t j=i; jGetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i])); - Float_t y = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j])); - fhCutCorrelation->Fill(x,y); - } - } - } - - fhNClustersITS[0]->Fill(nClustersITS); - fhNClustersTPC[0]->Fill(nClustersTPC); - fhChi2PerClusterITS[0]->Fill(chi2PerClusterITS); - fhChi2PerClusterTPC[0]->Fill(chi2PerClusterTPC); - - fhC11[0]->Fill(extCov[0]); - fhC22[0]->Fill(extCov[2]); - fhC33[0]->Fill(extCov[5]); - fhC44[0]->Fill(extCov[9]); - fhC55[0]->Fill(extCov[14]); - - fhPt[0]->Fill(pt); - fhEta[0]->Fill(eta); - - Float_t b[2]; - Float_t bRes[2]; - Float_t bCov[3]; - esdTrack->GetImpactParameters(b,bCov); - if (bCov[0]<=0 || bCov[2]<=0) { - AliDebug(1, "Estimated b resolution lower or equal zero!"); - bCov[0]=0; bCov[2]=0; - } - bRes[0] = TMath::Sqrt(bCov[0]); - bRes[1] = TMath::Sqrt(bCov[2]); - - fhDZ[0]->Fill(b[1]); - fhDXY[0]->Fill(b[0]); - fhDXYvsDZ[0]->Fill(b[1],b[0]); - - if (bRes[0]!=0 && bRes[1]!=0) { - fhDZNormalized[0]->Fill(b[1]/bRes[1]); - fhDXYNormalized[0]->Fill(b[0]/bRes[0]); - fhDXYvsDZNormalized[0]->Fill(b[1]/bRes[1], b[0]/bRes[0]); - fhNSigmaToVertex[0]->Fill(nSigmaToVertex); - } - } - - //######################################################################## - // cut the track! - if (cut) return kFALSE; - - //######################################################################## - // filling histograms after cut - if (fHistogramsOn) { - fhNClustersITS[1]->Fill(nClustersITS); - fhNClustersTPC[1]->Fill(nClustersTPC); - fhChi2PerClusterITS[1]->Fill(chi2PerClusterITS); - fhChi2PerClusterTPC[1]->Fill(chi2PerClusterTPC); - - fhC11[1]->Fill(extCov[0]); - fhC22[1]->Fill(extCov[2]); - fhC33[1]->Fill(extCov[5]); - fhC44[1]->Fill(extCov[9]); - fhC55[1]->Fill(extCov[14]); - - fhPt[1]->Fill(pt); - fhEta[1]->Fill(eta); - - Float_t b[2]; - Float_t bRes[2]; - Float_t bCov[3]; - esdTrack->GetImpactParameters(b,bCov); - if (bCov[0]<=0 || bCov[2]<=0) { - AliDebug(1, "Estimated b resolution lower or equal zero!"); - bCov[0]=0; bCov[2]=0; - } - bRes[0] = TMath::Sqrt(bCov[0]); - bRes[1] = TMath::Sqrt(bCov[2]); - - fhDZ[1]->Fill(b[1]); - fhDXY[1]->Fill(b[0]); - fhDXYvsDZ[1]->Fill(b[1],b[0]); - - if (bRes[0]!=0 && bRes[1]!=0) - { - fhDZNormalized[1]->Fill(b[1]/bRes[1]); - fhDXYNormalized[1]->Fill(b[0]/bRes[0]); - fhDXYvsDZNormalized[1]->Fill(b[1]/bRes[1], b[0]/bRes[0]); - fhNSigmaToVertex[1]->Fill(nSigmaToVertex); - } - } - - return kTRUE; -} - -//____________________________________________________________________ -TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESD* esd) -{ - // - // returns an array of all tracks that pass the cuts - // - - TObjArray* acceptedTracks = new TObjArray(); - - // loop over esd tracks - for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) { - AliESDtrack* track = esd->GetTrack(iTrack); - - if (AcceptTrack(track)) - acceptedTracks->Add(track); - } - - return acceptedTracks; -} - -//____________________________________________________________________ -Int_t AliESDtrackCuts::CountAcceptedTracks(AliESD* esd) -{ - // - // returns an the number of tracks that pass the cuts - // - - Int_t count = 0; - - // loop over esd tracks - for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) { - AliESDtrack* track = esd->GetTrack(iTrack); - - if (AcceptTrack(track)) - count++; - } - - return count; -} - -//____________________________________________________________________ -TObjArray* AliESDtrackCuts::GetAcceptedTracks(AliESDEvent* esd) -{ - // - // returns an array of all tracks that pass the cuts - // - - TObjArray* acceptedTracks = new TObjArray(); - - // loop over esd tracks - for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) { - AliESDtrack* track = esd->GetTrack(iTrack); - - if (AcceptTrack(track)) - acceptedTracks->Add(track); - } - - return acceptedTracks; -} - -//____________________________________________________________________ -Int_t AliESDtrackCuts::CountAcceptedTracks(AliESDEvent* esd) -{ - // - // returns an the number of tracks that pass the cuts - // - - Int_t count = 0; - - // loop over esd tracks - for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) { - AliESDtrack* track = esd->GetTrack(iTrack); - - if (AcceptTrack(track)) - count++; - } - - return count; -} - -//____________________________________________________________________ - void AliESDtrackCuts::DefineHistograms(Int_t color) { - // - // diagnostics histograms are defined - // - - fHistogramsOn=kTRUE; - - Bool_t oldStatus = TH1::AddDirectoryStatus(); - TH1::AddDirectory(kFALSE); - - //################################################################################### - // defining histograms - - fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5); - - fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks"); - fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks"); - - fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);; - - for (Int_t i=0; iGetXaxis()->SetBinLabel(i+4,fgkCutNames[i]); - fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]); - fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]); - } - - fhCutStatistics ->SetLineColor(color); - fhCutCorrelation ->SetLineColor(color); - fhCutStatistics ->SetLineWidth(2); - fhCutCorrelation ->SetLineWidth(2); - - Char_t str[256]; - for (Int_t i=0; i<2; i++) { - if (i==0) sprintf(str," "); - else sprintf(str,"_cut"); - - fhNClustersITS[i] = new TH1F(Form("nClustersITS%s",str) ,"",8,-0.5,7.5); - fhNClustersTPC[i] = new TH1F(Form("nClustersTPC%s",str) ,"",165,-0.5,164.5); - fhChi2PerClusterITS[i] = new TH1F(Form("chi2PerClusterITS%s",str),"",500,0,10); - fhChi2PerClusterTPC[i] = new TH1F(Form("chi2PerClusterTPC%s",str),"",500,0,10); - - fhC11[i] = new TH1F(Form("covMatrixDiagonal11%s",str),"",2000,0,20); - fhC22[i] = new TH1F(Form("covMatrixDiagonal22%s",str),"",2000,0,20); - fhC33[i] = new TH1F(Form("covMatrixDiagonal33%s",str),"",1000,0,1); - fhC44[i] = new TH1F(Form("covMatrixDiagonal44%s",str),"",1000,0,5); - fhC55[i] = new TH1F(Form("covMatrixDiagonal55%s",str),"",1000,0,5); - - fhDXY[i] = new TH1F(Form("dXY%s",str) ,"",500,-10,10); - fhDZ[i] = new TH1F(Form("dZ%s",str) ,"",500,-10,10); - fhDXYvsDZ[i] = new TH2F(Form("dXYvsDZ%s",str),"",200,-10,10,200,-10,10); - - fhDXYNormalized[i] = new TH1F(Form("dXYNormalized%s",str) ,"",500,-10,10); - fhDZNormalized[i] = new TH1F(Form("dZNormalized%s",str) ,"",500,-10,10); - fhDXYvsDZNormalized[i] = new TH2F(Form("dXYvsDZNormalized%s",str),"",200,-10,10,200,-10,10); - - fhNSigmaToVertex[i] = new TH1F(Form("nSigmaToVertex%s",str),"",500,0,50); - - fhPt[i] = new TH1F(Form("pt%s",str) ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0); - fhEta[i] = new TH1F(Form("eta%s",str) ,"#eta distribution;#eta",40,-2.0,2.0); - - fhNClustersITS[i]->SetTitle("n ITS clusters"); - fhNClustersTPC[i]->SetTitle("n TPC clusters"); - fhChi2PerClusterITS[i]->SetTitle("#Chi^{2} per ITS cluster"); - fhChi2PerClusterTPC[i]->SetTitle("#Chi^{2} per TPC cluster"); - - fhC11[i]->SetTitle("cov 11 : #sigma_{y}^{2} [cm^{2}]"); - fhC22[i]->SetTitle("cov 22 : #sigma_{z}^{2} [cm^{2}]"); - fhC33[i]->SetTitle("cov 33 : #sigma_{sin(#phi)}^{2}"); - fhC44[i]->SetTitle("cov 44 : #sigma_{tan(#theta_{dip})}^{2}"); - fhC55[i]->SetTitle("cov 55 : #sigma_{1/p_{T}}^{2} [(c/GeV)^2]"); - - fhDXY[i]->SetTitle("transverse impact parameter"); - fhDZ[i]->SetTitle("longitudinal impact parameter"); - fhDXYvsDZ[i]->SetTitle("longitudinal impact parameter"); - fhDXYvsDZ[i]->SetYTitle("transverse impact parameter"); - - fhDXYNormalized[i]->SetTitle("normalized trans impact par"); - fhDZNormalized[i]->SetTitle("normalized long impact par"); - fhDXYvsDZNormalized[i]->SetTitle("normalized long impact par"); - fhDXYvsDZNormalized[i]->SetYTitle("normalized trans impact par"); - fhNSigmaToVertex[i]->SetTitle("n #sigma to vertex"); - - fhNClustersITS[i]->SetLineColor(color); fhNClustersITS[i]->SetLineWidth(2); - fhNClustersTPC[i]->SetLineColor(color); fhNClustersTPC[i]->SetLineWidth(2); - fhChi2PerClusterITS[i]->SetLineColor(color); fhChi2PerClusterITS[i]->SetLineWidth(2); - fhChi2PerClusterTPC[i]->SetLineColor(color); fhChi2PerClusterTPC[i]->SetLineWidth(2); - - fhC11[i]->SetLineColor(color); fhC11[i]->SetLineWidth(2); - fhC22[i]->SetLineColor(color); fhC22[i]->SetLineWidth(2); - fhC33[i]->SetLineColor(color); fhC33[i]->SetLineWidth(2); - fhC44[i]->SetLineColor(color); fhC44[i]->SetLineWidth(2); - fhC55[i]->SetLineColor(color); fhC55[i]->SetLineWidth(2); - - fhDXY[i]->SetLineColor(color); fhDXY[i]->SetLineWidth(2); - fhDZ[i]->SetLineColor(color); fhDZ[i]->SetLineWidth(2); - - fhDXYNormalized[i]->SetLineColor(color); fhDXYNormalized[i]->SetLineWidth(2); - fhDZNormalized[i]->SetLineColor(color); fhDZNormalized[i]->SetLineWidth(2); - fhNSigmaToVertex[i]->SetLineColor(color); fhNSigmaToVertex[i]->SetLineWidth(2); - } - - // The number of sigmas to the vertex is per definition gaussian - ffDTheoretical = new TF1("nSigmaToVertexTheoretical","([0]/2.506628274)*exp(-(x**2)/2)",0,50); - ffDTheoretical->SetParameter(0,1); - - TH1::AddDirectory(oldStatus); -} - -//____________________________________________________________________ -Bool_t AliESDtrackCuts::LoadHistograms(const Char_t* dir) -{ - // - // loads the histograms from a file - // if dir is empty a directory with the name of this object is taken (like in SaveHistogram) - // - - if (!dir) - dir = GetName(); - - if (!gDirectory->cd(dir)) - return kFALSE; - - ffDTheoretical = dynamic_cast (gDirectory->Get("nSigmaToVertexTheory")); - - fhCutStatistics = dynamic_cast (gDirectory->Get("cut_statistics")); - fhCutCorrelation = dynamic_cast (gDirectory->Get("cut_correlation")); - - Char_t str[5]; - for (Int_t i=0; i<2; i++) { - if (i==0) - { - gDirectory->cd("before_cuts"); - str[0] = 0; - } - else - { - gDirectory->cd("after_cuts"); - sprintf(str,"_cut"); - } - - fhNClustersITS[i] = dynamic_cast (gDirectory->Get(Form("nClustersITS%s",str) )); - fhNClustersTPC[i] = dynamic_cast (gDirectory->Get(Form("nClustersTPC%s",str) )); - fhChi2PerClusterITS[i] = dynamic_cast (gDirectory->Get(Form("chi2PerClusterITS%s",str))); - fhChi2PerClusterTPC[i] = dynamic_cast (gDirectory->Get(Form("chi2PerClusterTPC%s",str))); - - fhC11[i] = dynamic_cast (gDirectory->Get(Form("covMatrixDiagonal11%s",str))); - fhC22[i] = dynamic_cast (gDirectory->Get(Form("covMatrixDiagonal22%s",str))); - fhC33[i] = dynamic_cast (gDirectory->Get(Form("covMatrixDiagonal33%s",str))); - fhC44[i] = dynamic_cast (gDirectory->Get(Form("covMatrixDiagonal44%s",str))); - fhC55[i] = dynamic_cast (gDirectory->Get(Form("covMatrixDiagonal55%s",str))); - - fhDXY[i] = dynamic_cast (gDirectory->Get(Form("dXY%s",str) )); - fhDZ[i] = dynamic_cast (gDirectory->Get(Form("dZ%s",str) )); - fhDXYvsDZ[i] = dynamic_cast (gDirectory->Get(Form("dXYvsDZ%s",str))); - - fhDXYNormalized[i] = dynamic_cast (gDirectory->Get(Form("dXYNormalized%s",str) )); - fhDZNormalized[i] = dynamic_cast (gDirectory->Get(Form("dZNormalized%s",str) )); - fhDXYvsDZNormalized[i] = dynamic_cast (gDirectory->Get(Form("dXYvsDZNormalized%s",str))); - fhNSigmaToVertex[i] = dynamic_cast (gDirectory->Get(Form("nSigmaToVertex%s",str))); - - fhPt[i] = dynamic_cast (gDirectory->Get(Form("pt%s",str))); - fhEta[i] = dynamic_cast (gDirectory->Get(Form("eta%s",str))); - - gDirectory->cd("../"); - } - - gDirectory->cd(".."); - - return kTRUE; -} - -//____________________________________________________________________ -void AliESDtrackCuts::SaveHistograms(const Char_t* dir) { - // - // saves the histograms in a directory (dir) - // - - if (!fHistogramsOn) { - AliDebug(0, "Histograms not on - cannot save histograms!!!"); - return; - } - - if (!dir) - dir = GetName(); - - gDirectory->mkdir(dir); - gDirectory->cd(dir); - - gDirectory->mkdir("before_cuts"); - gDirectory->mkdir("after_cuts"); - - // a factor of 2 is needed since n sigma is positive - ffDTheoretical->SetParameter(0,2*fhNSigmaToVertex[0]->Integral("width")); - ffDTheoretical->Write("nSigmaToVertexTheory"); - - fhCutStatistics->Write(); - fhCutCorrelation->Write(); - - for (Int_t i=0; i<2; i++) { - if (i==0) - gDirectory->cd("before_cuts"); - else - gDirectory->cd("after_cuts"); - - fhNClustersITS[i] ->Write(); - fhNClustersTPC[i] ->Write(); - fhChi2PerClusterITS[i] ->Write(); - fhChi2PerClusterTPC[i] ->Write(); - - fhC11[i] ->Write(); - fhC22[i] ->Write(); - fhC33[i] ->Write(); - fhC44[i] ->Write(); - fhC55[i] ->Write(); - - fhDXY[i] ->Write(); - fhDZ[i] ->Write(); - fhDXYvsDZ[i] ->Write(); - - fhDXYNormalized[i] ->Write(); - fhDZNormalized[i] ->Write(); - fhDXYvsDZNormalized[i] ->Write(); - fhNSigmaToVertex[i] ->Write(); - - fhPt[i] ->Write(); - fhEta[i] ->Write(); - - gDirectory->cd("../"); - } - - gDirectory->cd("../"); -} - -//____________________________________________________________________ -void AliESDtrackCuts::DrawHistograms() -{ - // draws some histograms - - TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Quality Results1", 800, 800); - canvas1->Divide(2, 2); - - canvas1->cd(1); - fhNClustersTPC[0]->SetStats(kFALSE); - fhNClustersTPC[0]->Draw(); - - canvas1->cd(2); - fhChi2PerClusterTPC[0]->SetStats(kFALSE); - fhChi2PerClusterTPC[0]->Draw(); - - canvas1->cd(3); - fhNSigmaToVertex[0]->SetStats(kFALSE); - fhNSigmaToVertex[0]->GetXaxis()->SetRangeUser(0, 10); - fhNSigmaToVertex[0]->Draw(); - - canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName())); - - TCanvas* canvas2 = new TCanvas(Form("%s_2", GetName()), "Track Quality Results2", 1200, 800); - canvas2->Divide(3, 2); - - canvas2->cd(1); - fhC11[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhC11[0]->Draw(); - - canvas2->cd(2); - fhC22[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhC22[0]->Draw(); - - canvas2->cd(3); - fhC33[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhC33[0]->Draw(); - - canvas2->cd(4); - fhC44[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhC44[0]->Draw(); - - canvas2->cd(5); - fhC55[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhC55[0]->Draw(); - - canvas2->SaveAs(Form("%s_%s.gif", GetName(), canvas2->GetName())); - - TCanvas* canvas3 = new TCanvas(Form("%s_3", GetName()), "Track Quality Results3", 1200, 800); - canvas3->Divide(3, 2); - - canvas3->cd(1); - fhDXY[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhDXY[0]->Draw(); - - canvas3->cd(2); - fhDZ[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhDZ[0]->Draw(); - - canvas3->cd(3); - fhDXYvsDZ[0]->SetStats(kFALSE); - gPad->SetLogz(); - gPad->SetRightMargin(0.15); - fhDXYvsDZ[0]->Draw("COLZ"); - - canvas3->cd(4); - fhDXYNormalized[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhDXYNormalized[0]->Draw(); - - canvas3->cd(5); - fhDZNormalized[0]->SetStats(kFALSE); - gPad->SetLogy(); - fhDZNormalized[0]->Draw(); - - canvas3->cd(6); - fhDXYvsDZNormalized[0]->SetStats(kFALSE); - gPad->SetLogz(); - gPad->SetRightMargin(0.15); - fhDXYvsDZNormalized[0]->Draw("COLZ"); - - canvas3->SaveAs(Form("%s_%s.gif", GetName(), canvas3->GetName())); - - TCanvas* canvas4 = new TCanvas(Form("%s_4", GetName()), "Track Quality Results4", 800, 500); - canvas4->Divide(2, 1); - - canvas4->cd(1); - fhCutStatistics->SetStats(kFALSE); - fhCutStatistics->LabelsOption("v"); - gPad->SetBottomMargin(0.3); - fhCutStatistics->Draw(); - - canvas4->cd(2); - fhCutCorrelation->SetStats(kFALSE); - fhCutCorrelation->LabelsOption("v"); - gPad->SetBottomMargin(0.3); - gPad->SetLeftMargin(0.3); - fhCutCorrelation->Draw("COLZ"); - - canvas4->SaveAs(Form("%s_%s.gif", GetName(), canvas4->GetName())); - - /*canvas->cd(1); - fhDXYvsDZNormalized[0]->SetStats(kFALSE); - fhDXYvsDZNormalized[0]->DrawCopy("COLZ"); - - canvas->cd(2); - fhNClustersTPC[0]->SetStats(kFALSE); - fhNClustersTPC[0]->DrawCopy(); - - canvas->cd(3); - fhChi2PerClusterITS[0]->SetStats(kFALSE); - fhChi2PerClusterITS[0]->DrawCopy(); - fhChi2PerClusterITS[1]->SetLineColor(2); - fhChi2PerClusterITS[1]->DrawCopy("SAME"); - - canvas->cd(4); - fhChi2PerClusterTPC[0]->SetStats(kFALSE); - fhChi2PerClusterTPC[0]->DrawCopy(); - fhChi2PerClusterTPC[1]->SetLineColor(2); - fhChi2PerClusterTPC[1]->DrawCopy("SAME");*/ -} - diff --git a/PWG0/esdTrackCuts/AliESDtrackCuts.h b/PWG0/esdTrackCuts/AliESDtrackCuts.h deleted file mode 100644 index 9d2e3085cdc..00000000000 --- a/PWG0/esdTrackCuts/AliESDtrackCuts.h +++ /dev/null @@ -1,174 +0,0 @@ -// -// Class for handling of ESD track cuts. -// -// The class manages a number of track quality cuts, a -// track-to-vertex cut and a number of kinematic cuts. Two methods -// can be used to figure out if an ESD track survives the cuts: -// AcceptTrack which takes a single AliESDtrack as argument and -// returns kTRUE/kFALSE or GetAcceptedTracks which takes an AliESD -// object and returns an TObjArray (of AliESDtracks) with the tracks -// in the ESD that survived the cuts. -// -// -// TODO: -// - add functionality to save and load cuts -// - add different ways to make track to vertex cut -// - add histograms for kinematic cut variables? -// - upper and lower cuts for all (non-boolean) cuts -// - update print method -// - is there a smarter way to manage the cuts? -// - put comments to each variable -// - -#ifndef ALIESDTRACKCUTS_H -#define ALIESDTRACKCUTS_H - -#include -#include -#include "AliAnalysisCuts.h" - -class AliESD; -class AliESDEvent; -class AliESDtrack; -class AliLog; -class TTree; - -class AliESDtrackCuts : public AliAnalysisCuts -{ -public: - AliESDtrackCuts(const Char_t* name = "AliESDtrackCuts", const Char_t* title = ""); - virtual ~AliESDtrackCuts(); - Bool_t IsSelected(TObject* obj) - {return AcceptTrack((AliESDtrack*)obj);} - Bool_t AcceptTrack(AliESDtrack* esdTrack); - TObjArray* GetAcceptedTracks(AliESD* esd); - Int_t CountAcceptedTracks(AliESD* esd); - TObjArray* GetAcceptedTracks(AliESDEvent* esd); - Int_t CountAcceptedTracks(AliESDEvent* esd); - - virtual Long64_t Merge(TCollection* list); - virtual void Copy(TObject &c) const; - AliESDtrackCuts(const AliESDtrackCuts& pd); // Copy Constructor - AliESDtrackCuts &operator=(const AliESDtrackCuts &c); - - //###################################################### - // track quality cut setters - void SetMinNClustersTPC(Int_t min=-1) {fCutMinNClusterTPC=min;} - void SetMinNClustersITS(Int_t min=-1) {fCutMinNClusterITS=min;} - void SetMaxChi2PerClusterTPC(Float_t max=1e10) {fCutMaxChi2PerClusterTPC=max;} - void SetMaxChi2PerClusterITS(Float_t max=1e10) {fCutMaxChi2PerClusterITS=max;} - void SetRequireTPCRefit(Bool_t b=kFALSE) {fCutRequireTPCRefit=b;} - void SetRequireITSRefit(Bool_t b=kFALSE) {fCutRequireITSRefit=b;} - void SetAcceptKingDaughters(Bool_t b=kFALSE) {fCutAcceptKinkDaughters=b;} - void SetMaxCovDiagonalElements(Float_t c1=1e10, Float_t c2=1e10, Float_t c3=1e10, Float_t c4=1e10, Float_t c5=1e10) - {fCutMaxC11=c1; fCutMaxC22=c2; fCutMaxC33=c3; fCutMaxC44=c4; fCutMaxC55=c5;} - - // track to vertex cut setters - void SetMinNsigmaToVertex(Float_t sigma=1e10) {fCutNsigmaToVertex = sigma;} - void SetRequireSigmaToVertex(Bool_t b=kTRUE ) {fCutSigmaToVertexRequired = b;} - - // getters - Float_t GetMinNsigmaToVertex() { return fCutNsigmaToVertex;} - Bool_t GetRequireSigmaToVertex( ) { return fCutSigmaToVertexRequired;} - - // track kinmatic cut setters - void SetPRange(Float_t r1=0, Float_t r2=1e10) {fPMin=r1; fPMax=r2;} - void SetPtRange(Float_t r1=0, Float_t r2=1e10) {fPtMin=r1; fPtMax=r2;} - void SetPxRange(Float_t r1=-1e10, Float_t r2=1e10) {fPxMin=r1; fPxMax=r2;} - void SetPyRange(Float_t r1=-1e10, Float_t r2=1e10) {fPyMin=r1; fPyMax=r2;} - void SetPzRange(Float_t r1=-1e10, Float_t r2=1e10) {fPzMin=r1; fPzMax=r2;} - void SetEtaRange(Float_t r1=-1e10, Float_t r2=1e10) {fEtaMin=r1; fEtaMax=r2;} - void SetRapRange(Float_t r1=-1e10, Float_t r2=1e10) {fRapMin=r1; fRapMax=r2;} - - //###################################################### - void SetHistogramsOn(Bool_t b=kFALSE) {fHistogramsOn = b;} - void DefineHistograms(Int_t color=1); - virtual Bool_t LoadHistograms(const Char_t* dir = 0); - void SaveHistograms(const Char_t* dir = 0); - void DrawHistograms(); - - Float_t GetSigmaToVertex(AliESDtrack* esdTrack); - - static void EnableNeededBranches(TTree* tree); - - // void SaveQualityCuts(Char_t* file) - // void LoadQualityCuts(Char_t* file) - - TH1* GetDZNormalized(Int_t i) const { return fhDZNormalized[i]; } - -protected: - void Init(); // sets everything to 0 - - enum { kNCuts = 21 }; - - //###################################################### - // esd track quality cuts - static const Char_t* fgkCutNames[kNCuts]; //! names of cuts (for internal use) - - Int_t fCutMinNClusterTPC; // min number of tpc clusters - Int_t fCutMinNClusterITS; // min number of its clusters - - Float_t fCutMaxChi2PerClusterTPC; // max tpc fit chi2 per tpc cluster - Float_t fCutMaxChi2PerClusterITS; // max its fit chi2 per its cluster - - Float_t fCutMaxC11; // max cov. matrix diag. elements (res. y^2) - Float_t fCutMaxC22; // max cov. matrix diag. elements (res. z^2) - Float_t fCutMaxC33; // max cov. matrix diag. elements (res. sin(phi)^2) - Float_t fCutMaxC44; // max cov. matrix diag. elements (res. tan(theta_dip)^2) - Float_t fCutMaxC55; // max cov. matrix diag. elements (res. 1/pt^2) - - Bool_t fCutAcceptKinkDaughters; // accepting kink daughters? - Bool_t fCutRequireTPCRefit; // require TPC refit - Bool_t fCutRequireITSRefit; // require ITS refit - - // track to vertex cut - Float_t fCutNsigmaToVertex; // max number of estimated sigma from track-to-vertex - Bool_t fCutSigmaToVertexRequired; // cut track if sigma from track-to-vertex could not be calculated - - // esd kinematics cuts - Float_t fPMin, fPMax; // definition of the range of the P - Float_t fPtMin, fPtMax; // definition of the range of the Pt - Float_t fPxMin, fPxMax; // definition of the range of the Px - Float_t fPyMin, fPyMax; // definition of the range of the Py - Float_t fPzMin, fPzMax; // definition of the range of the Pz - Float_t fEtaMin, fEtaMax; // definition of the range of the eta - Float_t fRapMin, fRapMax; // definition of the range of the y - - //###################################################### - // diagnostics histograms - Bool_t fHistogramsOn; // histograms on/off - - TH1F* fhNClustersITS[2]; //-> - TH1F* fhNClustersTPC[2]; //-> - - TH1F* fhChi2PerClusterITS[2]; //-> - TH1F* fhChi2PerClusterTPC[2]; //-> - - TH1F* fhC11[2]; //-> - TH1F* fhC22[2]; //-> - TH1F* fhC33[2]; //-> - TH1F* fhC44[2]; //-> - TH1F* fhC55[2]; //-> - - TH1F* fhDXY[2]; //-> - TH1F* fhDZ[2]; //-> - TH2F* fhDXYvsDZ[2]; //-> - - TH1F* fhDXYNormalized[2]; //-> - TH1F* fhDZNormalized[2]; //-> - TH2F* fhDXYvsDZNormalized[2]; //-> - TH1F* fhNSigmaToVertex[2]; //-> - - TH1F* fhPt[2]; //-> pt of esd tracks - TH1F* fhEta[2]; //-> eta of esd tracks - - TF1* ffDTheoretical; //-> theoretical distance to vertex normalized (2d gauss) - - TH1F* fhCutStatistics; //-> statistics of what cuts the tracks did not survive - TH2F* fhCutCorrelation; //-> 2d statistics plot - - ClassDef(AliESDtrackCuts, 2) -}; - - -#endif diff --git a/PWG0/esdTrackCuts/AliTestESDtrackCutsSelector.cxx b/PWG0/esdTrackCuts/AliTestESDtrackCutsSelector.cxx index f23f8210ba1..3a49e5db778 100644 --- a/PWG0/esdTrackCuts/AliTestESDtrackCutsSelector.cxx +++ b/PWG0/esdTrackCuts/AliTestESDtrackCutsSelector.cxx @@ -21,7 +21,7 @@ #include #include -#include "esdTrackCuts/AliESDtrackCuts.h" +#include "AliESDtrackCuts.h" #include "AliPWG0Helper.h" diff --git a/PWG0/highMultiplicity/AliHighMultiplicitySelector.cxx b/PWG0/highMultiplicity/AliHighMultiplicitySelector.cxx index 559b225b4fd..c778f410ffc 100644 --- a/PWG0/highMultiplicity/AliHighMultiplicitySelector.cxx +++ b/PWG0/highMultiplicity/AliHighMultiplicitySelector.cxx @@ -1003,11 +1003,34 @@ void AliHighMultiplicitySelector::Ntrigger(Bool_t relative) } else proj2->DrawCopy(" P"); - legend2->AddEntry(proj2, Form("%lld evts, FO > %d chips", nTrigger, cut)); + TString eventStr; + if (nTrigger > 1e6) + { + eventStr.Form("%lld M", nTrigger / 1000 / 1000); + } + else if (nTrigger > 1e3) + { + eventStr.Form("%lld K", nTrigger / 1000); + } + else + eventStr.Form("%lld", nTrigger); + + TString triggerStr; + if (cut == 0) + { + triggerStr = "minimum bias"; + } + else + triggerStr.Form("FO > %d chips", cut); - TLine* line = new TLine(triggerLimit, proj2->GetMinimum(), triggerLimit, proj2->GetMaximum()); - line->SetLineColor(colors[currentCut]); - line->Draw(); + legend2->AddEntry(proj2, Form("%s evts, %s", eventStr.Data(), triggerStr.Data())); + + if (triggerLimit > 1) + { + TLine* line = new TLine(triggerLimit, proj2->GetMinimum(), triggerLimit, proj2->GetMaximum()); + line->SetLineColor(colors[currentCut]); + line->Draw(); + } } legend2->Draw(); @@ -1305,6 +1328,7 @@ void AliHighMultiplicitySelector::DrawHistograms() x->ReadHistograms("highmult_central.root"); x->DrawHistograms(); + gSystem->Load("libANALYSIS"); gSystem->Load("libPWG0base"); .L AliHighMultiplicitySelector.cxx+ x = new AliHighMultiplicitySelector(); @@ -1339,6 +1363,16 @@ void AliHighMultiplicitySelector::DrawHistograms() canvas->SaveAs("L1NoCurve.gif"); canvas->SaveAs("L1NoCurve.eps"); + TLine* line = new TLine(fMvsL1->GetXaxis()->GetXmin(), 150, fMvsL1->GetXaxis()->GetXmax(), 150); + line->SetLineWidth(2); + line->SetLineColor(kRed); + line->Draw(); + + canvas->SaveAs("L1NoCurveCut.gif"); + canvas->SaveAs("L1NoCurveCut.eps"); + + return; + // draw corresponding theoretical curve TF1* func = new TF1("func", "[0]*(1-(1-1/[0])**x)", 1, 1000); func->SetParameter(0, 400-5*2); diff --git a/PWG0/libPWG0base.pkg b/PWG0/libPWG0base.pkg index 0677322a403..44a87d5784a 100644 --- a/PWG0/libPWG0base.pkg +++ b/PWG0/libPWG0base.pkg @@ -4,7 +4,6 @@ HDRS = dNdEta/dNdEtaAnalysis.h \ dNdEta/AlidNdEtaCorrection.h \ - esdTrackCuts/AliESDtrackCuts.h \ AliCorrectionMatrix.h \ AliCorrectionMatrix2D.h \ AliCorrectionMatrix3D.h \ diff --git a/PWG0/multiplicity/AliMultiplicityTask.cxx b/PWG0/multiplicity/AliMultiplicityTask.cxx index 0be9afb38a2..3b4bceb0c40 100644 --- a/PWG0/multiplicity/AliMultiplicityTask.cxx +++ b/PWG0/multiplicity/AliMultiplicityTask.cxx @@ -28,7 +28,7 @@ #include #include -#include "esdTrackCuts/AliESDtrackCuts.h" +#include "AliESDtrackCuts.h" #include "AliPWG0Helper.h" #include "dNdEta/AliMultiplicityCorrection.h" #include "AliCorrection.h" @@ -40,7 +40,7 @@ AliMultiplicityTask::AliMultiplicityTask(const char* opt) : AliAnalysisTask("AliMultiplicityTask", ""), fESD(0), fOption(opt), - fAnalysisMode(kSPD), + fAnalysisMode(AliPWG0Helper::kSPD), fReadMC(kFALSE), fMultiplicity(0), fEsdTrackCuts(0), @@ -95,10 +95,10 @@ void AliMultiplicityTask::ConnectInputData(Option_t *) tree->SetBranchStatus("fTriggerMask", 1); tree->SetBranchStatus("fSPDVertex*", 1); - if (fAnalysisMode == kSPD) + if (fAnalysisMode == AliPWG0Helper::kSPD) tree->SetBranchStatus("fSPDMult*", 1); - if (fAnalysisMode == kTPC) { + if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) { AliESDtrackCuts::EnableNeededBranches(tree); tree->SetBranchStatus("fTracks.fLabel", 1); } @@ -197,14 +197,14 @@ void AliMultiplicityTask::Exec(Option_t*) // Check prerequisites if (!fESD) { - AliDebug(AliLog::kError, "ESD branch not available"); + AliDebug(AliLog::kError, "ESD not available"); return; } - // MB1 definition - //Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), AliPWG0Helper::kMB1); + //MB1 definition + Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), AliPWG0Helper::kMB1); // only FASTOR - Bool_t eventTriggered = fESD->GetTriggerMask() & 32; + //Bool_t eventTriggered = fESD->GetTriggerMask() & 32; Bool_t eventVertex = AliPWG0Helper::IsVertexReconstructed(fESD->GetVertex()); @@ -222,7 +222,7 @@ void AliMultiplicityTask::Exec(Option_t*) Int_t inputCount = 0; Int_t* labelArr = 0; Float_t* etaArr = 0; - if (fAnalysisMode == kSPD) + if (fAnalysisMode == AliPWG0Helper::kSPD) { // get tracklets const AliMultiplicity* mult = fESD->GetMultiplicity(); @@ -250,7 +250,7 @@ void AliMultiplicityTask::Exec(Option_t*) ++inputCount; } } - else if (fAnalysisMode == kTPC) + else if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS) { if (!fEsdTrackCuts) { @@ -286,8 +286,10 @@ void AliMultiplicityTask::Exec(Option_t*) // eta range for nMCTracksSpecies, nESDTracksSpecies Float_t etaRange = -1; switch (fAnalysisMode) { - case kTPC: etaRange = 0.9; break; - case kSPD: etaRange = 2.0; break; + case AliPWG0Helper::kTPC: + case AliPWG0Helper::kTPCITS: + etaRange = 0.9; break; + case AliPWG0Helper::kSPD: etaRange = 2.0; break; } if (!fReadMC) // Processing of ESD information diff --git a/PWG0/multiplicity/AliMultiplicityTask.h b/PWG0/multiplicity/AliMultiplicityTask.h index fc63aaa38d7..3fa8ee4dc5e 100644 --- a/PWG0/multiplicity/AliMultiplicityTask.h +++ b/PWG0/multiplicity/AliMultiplicityTask.h @@ -4,7 +4,9 @@ #define AliMultiplicityTask_H #include "AliAnalysisTask.h" + #include +#include "AliPWG0Helper.h" class AliESDtrackCuts; class AliMultiplicityCorrection; @@ -15,8 +17,6 @@ class AliESDEvent; class AliMultiplicityTask : public AliAnalysisTask { public: - enum AnalysisMethod { kSPD = 0, kTPC }; - AliMultiplicityTask(const char* opt = ""); virtual ~AliMultiplicityTask(); @@ -28,14 +28,14 @@ class AliMultiplicityTask : public AliAnalysisTask { void SetTrackCuts(AliESDtrackCuts* cuts) { fEsdTrackCuts = cuts; } void SetPtSpectrum(TH1* hist) { fPtSpectrum = hist; } - void SetAnalysisMode(AnalysisMethod mode) { fAnalysisMode = mode; } + void SetAnalysisMode(AliPWG0Helper::AnalysisMode mode) { fAnalysisMode = mode; } void SetReadMC(Bool_t flag = kTRUE) { fReadMC = flag; } protected: AliESDEvent *fESD; //! ESD object TString fOption; // option string - AnalysisMethod fAnalysisMode; // detector that is used for analysis + AliPWG0Helper::AnalysisMode fAnalysisMode; // detector that is used for analysis Bool_t fReadMC; // if true reads MC data (to build correlation maps) AliMultiplicityCorrection* fMultiplicity; //! object containing the extracted data diff --git a/PWG0/multiplicity/CreateCuts.C b/PWG0/multiplicity/CreateCuts.C deleted file mode 100644 index e0e44359f17..00000000000 --- a/PWG0/multiplicity/CreateCuts.C +++ /dev/null @@ -1,22 +0,0 @@ -/* $Id$ */ - -// this macro creates the track and event cuts used in this analysis - -AliESDtrackCuts* CreateTrackCuts(Bool_t hists = kTRUE) -{ - AliESDtrackCuts* esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts"); - - if (hists) - esdTrackCuts->DefineHistograms(1); - - esdTrackCuts->SetMinNClustersTPC(50); - esdTrackCuts->SetMaxChi2PerClusterTPC(3.5); - esdTrackCuts->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); - esdTrackCuts->SetRequireTPCRefit(kTRUE); - - esdTrackCuts->SetMinNsigmaToVertex(3); - esdTrackCuts->SetRequireSigmaToVertex(kTRUE); - esdTrackCuts->SetAcceptKingDaughters(kFALSE); - - return esdTrackCuts; -} diff --git a/PWG0/multiplicity/run.C b/PWG0/multiplicity/run.C index fc7b08fbb89..faa4074ad67 100644 --- a/PWG0/multiplicity/run.C +++ b/PWG0/multiplicity/run.C @@ -1,20 +1,31 @@ -void run(Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, Bool_t aProof = kFALSE, Bool_t mc = kTRUE, const char* option = "") +void run(Char_t* data, Long64_t nRuns = -1, Long64_t offset = 0, Bool_t aDebug = kFALSE, Int_t aProof = 0, Bool_t mc = kTRUE, const char* option = "") { + // aProof option: 0 no proof + // 1 proof with chain + // 2 proof with dataset + + if (nRuns < 0) + nRuns = 1234567890; + if (aProof) { TProof::Open("lxb6046"); // Enable the needed package - gProof->UploadPackage("STEERBase"); + /*gProof->UploadPackage("STEERBase"); gProof->EnablePackage("STEERBase"); gProof->UploadPackage("ESD"); gProof->EnablePackage("ESD"); gProof->UploadPackage("AOD"); gProof->EnablePackage("AOD"); gProof->UploadPackage("ANALYSIS"); - gProof->EnablePackage("ANALYSIS"); - gProof->UploadPackage("PWG0base"); - gProof->EnablePackage("PWG0base"); + gProof->EnablePackage("ANALYSIS");*/ + + gProof->UploadPackage("$ALICE_ROOT/AF-v4-12"); + gProof->EnablePackage("$ALICE_ROOT/AF-v4-12"); + + gProof->UploadPackage("$ALICE_ROOT/PWG0base"); + gProof->EnablePackage("$ALICE_ROOT/PWG0base"); } else { @@ -26,16 +37,14 @@ void run(Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, B gSystem->Load("libPWG0base"); } - // Create chain of input files - gROOT->LoadMacro("../CreateESDChain.C"); - chain = CreateESDChain(data, nRuns, offset); - // Create the analysis manager mgr = new AliAnalysisManager("testAnalysis"); + AliPWG0Helper::AnalysisMode analysisMode = AliPWG0Helper::kSPD; + // selection of esd tracks - gROOT->ProcessLine(".L CreateCuts.C"); - AliESDtrackCuts* esdTrackCuts = CreateTrackCuts(); + gROOT->ProcessLine(".L ../CreateStandardCuts.C"); + AliESDtrackCuts* esdTrackCuts = CreateTrackCuts(analysisMode); if (!esdTrackCuts) { printf("ERROR: esdTrackCuts could not be created\n"); @@ -47,14 +56,14 @@ void run(Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, B taskName += "+g"; // Create, add task - if (aProof) { + if (aProof > 0) { gProof->Load(taskName); } else gROOT->Macro(taskName); task = new AliMultiplicityTask(option); task->SetTrackCuts(esdTrackCuts); - task->SetAnalysisMode(AliMultiplicityTask::kSPD); + task->SetAnalysisMode(analysisMode); if (mc) task->SetReadMC(); @@ -88,5 +97,20 @@ void run(Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, B // Run analysis mgr->InitAnalysis(); mgr->PrintStatus(); - mgr->StartAnalysis((aProof) ? "proof" : "local", chain); + + if (aProof == 2) + { + // process dataset + + mgr->StartAnalysis("proof", data, nRuns, offset); + } + else + { + // Create chain of input files + gROOT->LoadMacro("../CreateESDChain.C"); + chain = CreateESDChain(data, nRuns, offset); + + mgr->StartAnalysis((aProof > 0) ? "proof" : "local", chain); + } + } -- 2.43.0