From 8074859b452bd6cabe69f5c8306d56ed96420f33 Mon Sep 17 00:00:00 2001 From: jgrosseo Date: Tue, 22 Jan 2008 14:10:30 +0000 Subject: [PATCH] AliCutTask uses global cut initialization --- PWG0/esdTrackCuts/AliCutTask.cxx | 95 +++++++++++++++++++++++++++----- PWG0/esdTrackCuts/AliCutTask.h | 5 ++ PWG0/esdTrackCuts/CreateCuts.C | 22 -------- PWG0/esdTrackCuts/draw.C | 68 +++++++++++++++++++++++ PWG0/esdTrackCuts/run.C | 28 ++++++---- 5 files changed, 172 insertions(+), 46 deletions(-) delete mode 100644 PWG0/esdTrackCuts/CreateCuts.C create mode 100644 PWG0/esdTrackCuts/draw.C diff --git a/PWG0/esdTrackCuts/AliCutTask.cxx b/PWG0/esdTrackCuts/AliCutTask.cxx index f64b40d1d76..b408fa4539b 100644 --- a/PWG0/esdTrackCuts/AliCutTask.cxx +++ b/PWG0/esdTrackCuts/AliCutTask.cxx @@ -28,9 +28,12 @@ ClassImp(AliCutTask) //________________________________________________________________________ -AliCutTask::AliCutTask(const char *name) - : AliAnalysisTask(name, ""), fESD(0), fTrackCuts(0), fTrackCutsPrimaries(0), - fTrackCutsSecondaries(0), fVertex(0), fTriggerStats(0), fOutput(0) +AliCutTask::AliCutTask(const char *name) : + AliAnalysisTask(name, ""), fESD(0), fTrackCuts(0), + fAnalysisMode(AliPWG0Helper::kTPCITS), + fTrackCutsPrimaries(0), + fTrackCutsSecondaries(0), fVertex(0), fTriggerStats(0), fOutput(0), + fPrimStats(0) { // Constructor @@ -50,7 +53,7 @@ void AliCutTask::ConnectInputData(Option_t *) Printf("ERROR: Could not read chain from input slot 0"); } else { // Disable all branches and enable only the needed ones - tree->SetBranchStatus("*", kFALSE); + //tree->SetBranchStatus("*", kFALSE); // old esd tree->SetBranchStatus("fTriggerMask", kTRUE); @@ -60,8 +63,8 @@ void AliCutTask::ConnectInputData(Option_t *) // new esd tree->SetBranchStatus("TriggerMask", kTRUE); tree->SetBranchStatus("AliESDHeader", kTRUE); - AliPWG0Helper::SetBranchStatusRecursive(tree, "Tracks", kTRUE); - AliPWG0Helper::SetBranchStatusRecursive(tree, "SPDVertex", kTRUE); + //AliPWG0Helper::SetBranchStatusRecursive(tree, "Tracks", kTRUE); + //AliPWG0Helper::SetBranchStatusRecursive(tree, "SPDVertex", kTRUE); AliESDInputHandler *esdH = dynamic_cast (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); @@ -98,6 +101,16 @@ void AliCutTask::CreateOutputObjects() fTriggerStats->GetXaxis()->SetBinLabel(5, "VZERO_OR_LEFT | VZERO_OR_RIGHT"); fOutput->Add(fTriggerStats); + fPrimStats = new TH1F("fPrimStats", "fPrimStats", 8, 0.5, 8.5); + fPrimStats->GetXaxis()->SetBinLabel(1, "Primary accepted once"); + fPrimStats->GetXaxis()->SetBinLabel(2, "Primary accepted more than once"); + fPrimStats->GetXaxis()->SetBinLabel(3, "Primary accepted more than once (count)"); + fPrimStats->GetXaxis()->SetBinLabel(4, "Primary track rejected"); + fPrimStats->GetXaxis()->SetBinLabel(5, "Primary track rejected (count)"); + fPrimStats->GetXaxis()->SetBinLabel(6, "Primary track rejected, but other track accepted"); + fPrimStats->GetXaxis()->SetBinLabel(7, "Primary track rejected, but other track accepted (count)"); + fOutput->Add(fPrimStats); + // disable info messages of AliMCEvent (per event) AliLog::SetClassDebugLevel("AliMCEvent", AliLog::kWarning - AliLog::kDebug + 1); } @@ -133,15 +146,22 @@ void AliCutTask::Exec(Option_t *) if (fESD->GetTriggerMask() & 1 || fESD->GetTriggerMask() & 2) fTriggerStats->Fill(4); - //if (!AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask()), AliPWG0Helper::kMB1) - if (fESD->GetTriggerMask() & (0x1 << 14) == 0) + //if (fESD->GetTriggerMask() & (0x1 << 14) == 0) + if (!AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), AliPWG0Helper::kMB2)) return; - if (!AliPWG0Helper::IsVertexReconstructed(fESD->GetVertex())) + // get the ESD vertex + const AliESDVertex* vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode); + + if (!vtxESD) + { + Printf("No vertex. Skipping event"); return; + } - //Printf("There are %d tracks in this event", fESD->GetNumberOfTracks()); - fTrackCuts->CountAcceptedTracks(fESD); + Printf("There are %d tracks in this event", fESD->GetNumberOfTracks()); + Int_t acceptedTracks = fTrackCuts->CountAcceptedTracks(fESD); + Printf("Accepted %d tracks", acceptedTracks); if (fTrackCutsPrimaries && fTrackCutsSecondaries) { @@ -163,7 +183,17 @@ void AliCutTask::Exec(Option_t *) Printf("ERROR: Stack not available"); return; } - + + // count if primaries are counted several times + Int_t max = stack->GetNprimary(); + Int_t* primAcc = new Int_t[max]; + Int_t* primRej = new Int_t[max]; + for (Int_t i=0; iGetNumberOfTracks(); trackId++) { AliESDtrack* track = fESD->GetTrack(trackId); @@ -182,7 +212,19 @@ void AliCutTask::Exec(Option_t *) Int_t label = TMath::Abs(track->GetLabel()); if (stack->IsPhysicalPrimary(label) == kTRUE) { - fTrackCutsPrimaries->AcceptTrack(track); + if (label >= max) + { + Printf("Warning label %d is higher than number of primaries %d", label, max); + continue; + } + + if (fTrackCutsPrimaries->AcceptTrack(track)) + { + primAcc[label]++; + } + else + primRej[label]++; + } else { @@ -208,10 +250,32 @@ void AliCutTask::Exec(Option_t *) fTrackCutsSecondaries->AcceptTrack(track); } } + + for (Int_t i=0; iFill(1); + } else if (primAcc[i] > 1) { + fPrimStats->Fill(2); + fPrimStats->Fill(3, primAcc[i]); + } + + if (primRej[i] > 0) { + if (primAcc[i] == 0) { + fPrimStats->Fill(4); + fPrimStats->Fill(5, primRej[i]); + } else if (primAcc[i] > 0) { + fPrimStats->Fill(6); + fPrimStats->Fill(7, primRej[i]); + } + } + } + + delete[] primAcc; + delete[] primRej; } // get the ESD vertex - fVertex->Fill(fESD->GetVertex()->GetZv()); + fVertex->Fill(vtxESD->GetZv()); } //________________________________________________________________________ @@ -258,6 +322,9 @@ void AliCutTask::Terminate(Option_t *) if (fTrackCutsSecondaries) fTrackCutsSecondaries->SaveHistograms(); + + if (fPrimStats) + fPrimStats->Write(); file->Close(); diff --git a/PWG0/esdTrackCuts/AliCutTask.h b/PWG0/esdTrackCuts/AliCutTask.h index 238225e1470..a7b6c0c381e 100644 --- a/PWG0/esdTrackCuts/AliCutTask.h +++ b/PWG0/esdTrackCuts/AliCutTask.h @@ -9,6 +9,7 @@ class AliESDEvent; class TList; #include "AliAnalysisTask.h" +#include "AliPWG0Helper.h" class AliCutTask : public AliAnalysisTask { public: @@ -21,11 +22,13 @@ class AliCutTask : public AliAnalysisTask { virtual void Terminate(Option_t *); void SetTrackCuts(AliESDtrackCuts* cuts) { fTrackCuts = cuts; } + void SetAnalysisMode(AliPWG0Helper::AnalysisMode mode) { fAnalysisMode = mode; } void EnableSecondaryStudy(); private: AliESDEvent *fESD; //! ESD object AliESDtrackCuts* fTrackCuts; // track cuts + AliPWG0Helper::AnalysisMode fAnalysisMode; // detector that is used for analysis AliESDtrackCuts* fTrackCutsPrimaries; // cuts for tracks from primary particles AliESDtrackCuts* fTrackCutsSecondaries; // cuts for tracks from secondary particles @@ -35,6 +38,8 @@ class AliCutTask : public AliAnalysisTask { TList* fOutput; //! list send on output slot 0 + TH1F* fPrimStats; //! statistics about primaries, see bin names in CreateOutputData + AliCutTask(const AliCutTask&); // not implemented AliCutTask& operator=(const AliCutTask&); // not implemented diff --git a/PWG0/esdTrackCuts/CreateCuts.C b/PWG0/esdTrackCuts/CreateCuts.C deleted file mode 100644 index e0e44359f17..00000000000 --- a/PWG0/esdTrackCuts/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/esdTrackCuts/draw.C b/PWG0/esdTrackCuts/draw.C new file mode 100644 index 00000000000..d458a72fbe1 --- /dev/null +++ b/PWG0/esdTrackCuts/draw.C @@ -0,0 +1,68 @@ +void draw() +{ + const char* files[] = { "trackCuts_normal.root", "trackCuts_increased.root", "trackCuts_decreased.root" }; + const char* titles[] = { "default geometry", "+ 10% material", "- 10% material" }; + Int_t colors[] = { 1, 2, 4 }; + + TCanvas* c = new TCanvas; + + TLegend* legend = new TLegend(0.7, 0.7, 0.9, 0.9); + + for (Int_t i=0; i<3; i++) { + TFile::Open(files[i]); + + TH1* ptPrim = gFile->Get("fTrackCutsPrimaries/before_cuts/pt"); + TH1* ptPrimCut = gFile->Get("fTrackCutsPrimaries/after_cuts/pt_cut"); + + TH1* ptSec = gFile->Get("fTrackCutsSecondaries/before_cuts/pt"); + TH1* ptSecCut = gFile->Get("fTrackCutsSecondaries/after_cuts/pt_cut"); + + ptPrim->Add(ptSec); + ptPrimCut->Add(ptSecCut); + + ptPrim->Sumw2(); + ptPrimCut->Sumw2(); + ptSec->Sumw2(); + ptSecCut->Sumw2(); + + Printf("%s", titles[i]); + Printf("%.2f %% secondaries before cuts", 100.0 * ptSec->GetEntries() / ptPrim->GetEntries()); + Printf("%.2f %% secondaries after cuts", 100.0 * ptSecCut->GetEntries() / ptPrimCut->GetEntries()); + Printf(""); + + ptSec->Divide(ptSec, ptPrim, 1, 1, "B"); + ptSecCut->Divide(ptSecCut, ptPrimCut, 1, 1, "B"); + + ptSec->SetLineColor(colors[i]); + ptSecCut->SetLineColor(colors[i]); + ptSec->SetStats(kFALSE); + + ptSec->GetXaxis()->SetRangeUser(0, 2); + ptSec->GetYaxis()->SetRangeUser(0, 1); + + ptSec->SetTitle(""); + ptSec->GetYaxis()->SetTitle("N_{Secondaries} / N_{All}"); + + ptSec->DrawCopy((i > 0) ? "SAME" : ""); + ptSecCut->DrawCopy("SAME"); + + legend->AddEntry(ptSec, titles[i]); + } + + legend->Draw(); + + c->SaveAs("secondaries_changedmaterial.gif"); +} + +void drawStats(const char* fileName = "trackCuts.root") +{ + TFile::Open(fileName); + + TH1* stat1 = gFile->Get("fTrackCutsPrimaries/cut_statistics"); + TH1* stat2 = gFile->Get("fTrackCutsSecondaries/cut_statistics"); + + new TCanvas; + stat1->Draw(); + stat2->SetLineColor(2); + stat2->Draw("SAME"); +} diff --git a/PWG0/esdTrackCuts/run.C b/PWG0/esdTrackCuts/run.C index f4b201a7e61..70cef463ae9 100644 --- a/PWG0/esdTrackCuts/run.C +++ b/PWG0/esdTrackCuts/run.C @@ -34,15 +34,6 @@ void run(Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, B // Create the analysis manager mgr = new AliAnalysisManager("testAnalysis"); - // selection of esd tracks - gROOT->ProcessLine(".L CreateCuts.C"); - AliESDtrackCuts* esdTrackCuts = CreateTrackCuts(); - if (!esdTrackCuts) - { - printf("ERROR: esdTrackCuts could not be created\n"); - return; - } - TString taskName("AliCutTask.cxx+"); if (aDebug) taskName += "+g"; @@ -54,7 +45,24 @@ void run(Char_t* data, Int_t nRuns=20, Int_t offset=0, Bool_t aDebug = kFALSE, B gROOT->Macro(taskName); task = new AliCutTask; - task->SetTrackCuts(esdTrackCuts); + + AliPWG0Helper::AnalysisMode analysisMode = AliPWG0Helper::kTPC; + task->SetAnalysisMode(analysisMode); + + if (analysisMode != AliPWG0Helper::kSPD) + { + // selection of esd tracks + gROOT->ProcessLine(".L ../CreateStandardCuts.C"); + AliESDtrackCuts* esdTrackCuts = CreateTrackCuts(analysisMode); + if (!esdTrackCuts) + { + printf("ERROR: esdTrackCuts could not be created\n"); + return; + } + + task->SetTrackCuts(esdTrackCuts); + } + mgr->AddTask(task); // Add ESD handler -- 2.43.0