]>
Commit | Line | Data |
---|---|---|
4c904d93 | 1 | /************************************************************************** |
2 | * This file is property of and copyright by the ALICE HLT Project * | |
3 | * All rights reserved. * | |
4 | * * | |
5 | * Primary Authors: * | |
6 | * Artur Szostak <artursz@iafrica.com> * | |
7 | * * | |
8 | * Permission to use, copy, modify and distribute this software and its * | |
9 | * documentation strictly for non-commercial purposes is hereby granted * | |
10 | * without fee, provided that the above copyright notice appears in all * | |
11 | * copies and that both the copyright notice and this permission notice * | |
12 | * appear in the supporting documentation. The authors make no claims * | |
13 | * about the suitability of this software for any purpose. It is * | |
14 | * provided "as is" without express or implied warranty. * | |
15 | **************************************************************************/ | |
16 | ||
17 | /// \ingroup macros | |
18 | /// \file PlotEfficiency.C | |
19 | /// \brief Generates efficiency plots from the histograms generated by AliAnalysisTaskLinkToMC. | |
20 | /// \author Artur Szostak <artursz@iafrica.com> | |
21 | /// | |
22 | /// This macro is used to generate efficiency plots, fake track ratio plots | |
23 | /// and also calculate the total integrated efficiency and fake track ratio. | |
24 | /// The histograms generated by the AliAnalysisTaskLinkToMC analysis task is | |
25 | /// used as input. The macro can be run with root as follows: | |
26 | /// | |
27 | /// $ root PlotEfficiency.C\(\"hists.root\"\) | |
28 | /// | |
29 | /// where hists.root is the name of the output file containing the histograms | |
30 | /// generated by the analysis task. Note that the '\' character before the '(' | |
31 | /// and '"' characters is required. Alternatively run the macros as follows from | |
32 | /// the root command prompt: | |
33 | /// | |
34 | /// $ root | |
35 | /// .L PlotEfficiency.C | |
36 | /// PlotEfficiency("hists.root") | |
37 | ||
38 | ||
39 | #if !defined(__CINT__) || defined(__MAKECINT__) | |
40 | #include "Riostream.h" | |
41 | #include "TH1.h" | |
42 | #include "TH2.h" | |
43 | #include "TFile.h" | |
44 | #include "TList.h" | |
45 | #include "TGraphAsymmErrors.h" | |
46 | #include "TCanvas.h" | |
47 | #endif | |
48 | ||
49 | ||
50 | void EfficiencyIntegral(const char* prefix, TH1* pass, TH1* total, double& eff, double& error, bool includeOverflow = false) | |
51 | { | |
52 | /// This routine calculates the integrated efficiency for a 1D histogram. | |
53 | /// [in] \param prefix The string to print before the calculated efficiency result. | |
54 | /// [in] \param pass The histogram of reconstructed tracks. | |
55 | /// [in] \param total The histogram of findable (reconstructable) tracks. | |
56 | /// [out] \param eff Filled by the total integrated efficiency calculated. | |
57 | /// [out] \param error Filled by the uncertainty in the efficiency value. | |
58 | /// [in] \param includeOverflow If set then the overflow bins are also used | |
59 | /// in the calculation. | |
60 | ||
61 | TH1D passhist("temppasshist", "", 1, 0, 1); | |
62 | TH1D totalhist("temptotalhist", "", 1, 0, 1); | |
63 | if (includeOverflow) | |
64 | { | |
65 | passhist.SetBinContent(1, pass->Integral(0, pass->GetNbinsX()+1)); | |
66 | totalhist.SetBinContent(1, total->Integral(0, total->GetNbinsX()+1)); | |
67 | } | |
68 | else | |
69 | { | |
70 | passhist.SetBinContent(1, pass->Integral(1, pass->GetNbinsX())); | |
71 | totalhist.SetBinContent(1, total->Integral(1, total->GetNbinsX())); | |
72 | } | |
73 | TGraphAsymmErrors* effgraph = new TGraphAsymmErrors(&passhist, &totalhist); | |
74 | eff = effgraph->GetY()[0]; | |
75 | double effhigh = effgraph->GetEYhigh()[0]; | |
76 | double efflow = effgraph->GetEYlow()[0]; | |
77 | error = effhigh > efflow ? effhigh : efflow; | |
78 | cout << prefix << eff << " + " << effhigh << " - " << efflow << endl; | |
79 | } | |
80 | ||
81 | ||
82 | void PlotEfficiency(const char* histfilename = "hists.root", bool includeOverflow = true) | |
83 | { | |
84 | /// Opens the file containing the histograms generated by the AliAnalysisTaskLinkToMC | |
85 | /// analysis task and tries to find the first TList in the file. | |
86 | /// From the TList we try to find TH2 classes with the following names: | |
87 | /// "findableTracksHist", "foundTracksHistMC" and "fakeTracksHist" | |
88 | /// From these we plot the efficiency plots for the pT and rapidity dimension | |
89 | /// by projecting the 2D histograms. Also, the fake track ratio histograms are | |
90 | /// created. Finally the total integrated efficiency and fake track ratio is | |
91 | /// calculated and printed. | |
92 | /// \param histfilename The name of the root file containing the histograms | |
93 | /// generated by the AliAnalysisTaskLinkToMC analysis task. | |
94 | /// \param includeOverflow Indicates if the overflow bins should be used for | |
95 | /// for the total integral efficiency and fake track ratio calculations. | |
96 | /// The default is not to use the overflow bins. | |
97 | ||
98 | TFile file(histfilename); | |
99 | ||
100 | TList* histolist = NULL; | |
101 | TIter next(file.GetListOfKeys()); | |
102 | for (TObject* key = next(); key != NULL; key = next()) | |
103 | { | |
104 | TList* list = dynamic_cast<TList*>( file.Get(key->GetName()) ); | |
105 | if (list != NULL) | |
106 | { | |
107 | histolist = list; | |
108 | break; | |
109 | } | |
110 | } | |
111 | ||
112 | if (histolist == NULL) | |
113 | { | |
114 | cout << "ERROR: Could not find the histogram list in file '" << histfilename << "'" << endl; | |
115 | return; | |
116 | } | |
117 | ||
118 | TH2* findableTracksHist = dynamic_cast<TH2*>( histolist->FindObject("findableTracksHist") ); | |
119 | TH2* foundTracksHistMC = dynamic_cast<TH2*>( histolist->FindObject("foundTracksHistMC") ); | |
120 | TH2* fakeTracksHist = dynamic_cast<TH2*>( histolist->FindObject("fakeTracksHist") ); | |
121 | ||
122 | if (findableTracksHist == NULL) | |
123 | { | |
124 | cout << "ERROR: Could not find the histogram findableTracksHist in TList '" << histolist->GetName() << "'" << endl; | |
125 | return; | |
126 | } | |
127 | if (foundTracksHistMC == NULL) | |
128 | { | |
129 | cout << "ERROR: Could not find the histogram foundTracksHistMC in TList '" << histolist->GetName() << "'" << endl; | |
130 | return; | |
131 | } | |
132 | if (fakeTracksHist == NULL) | |
133 | { | |
134 | cout << "ERROR: Could not find the histogram fakeTracksHist in TList '" << histolist->GetName() << "'" << endl; | |
135 | return; | |
136 | } | |
137 | ||
138 | TH1* projeffx = foundTracksHistMC->ProjectionX(); | |
139 | TH1* projfindx = findableTracksHist->ProjectionX(); | |
140 | TGraphAsymmErrors* effpt = new TGraphAsymmErrors(projeffx, projfindx); | |
141 | effpt->Draw("ap"); | |
142 | effpt->GetHistogram()->SetTitle("Efficiency of track reconstruction."); | |
143 | effpt->GetHistogram()->SetXTitle(projeffx->GetXaxis()->GetTitle()); | |
144 | effpt->GetHistogram()->SetYTitle("#epsilon"); | |
145 | effpt->GetYaxis()->SetRangeUser(0, 1.1); | |
146 | ||
147 | new TCanvas; | |
148 | TH1* projfakex = fakeTracksHist->ProjectionX(); | |
149 | TGraphAsymmErrors* fakept = new TGraphAsymmErrors(projfakex, projfindx); | |
150 | fakept->Draw("ap"); | |
151 | fakept->GetHistogram()->SetTitle("Fake reconstructed tracks ratio."); | |
152 | fakept->GetHistogram()->SetXTitle(projfakex->GetXaxis()->GetTitle()); | |
153 | fakept->GetHistogram()->SetYTitle("fake tracks / findable tracks"); | |
154 | ||
155 | new TCanvas; | |
156 | TH1* projeffy = foundTracksHistMC->ProjectionY(); | |
157 | TH1* projfindy = findableTracksHist->ProjectionY(); | |
158 | TGraphAsymmErrors* effy = new TGraphAsymmErrors(projeffy, projfindy); | |
159 | effy->Draw("ap"); | |
160 | effy->GetHistogram()->SetTitle("Efficiency of track reconstruction."); | |
161 | effy->GetHistogram()->SetXTitle(projeffy->GetXaxis()->GetTitle()); | |
162 | effy->GetHistogram()->SetYTitle("#epsilon"); | |
163 | effy->GetYaxis()->SetRangeUser(0, 1.1); | |
164 | ||
165 | new TCanvas; | |
166 | TH1* projfakey = fakeTracksHist->ProjectionY(); | |
167 | TGraphAsymmErrors* fakey = new TGraphAsymmErrors(projfakey, projfindy); | |
168 | fakey->Draw("ap"); | |
169 | fakey->GetHistogram()->SetTitle("Fake reconstructed tracks ratio."); | |
170 | fakey->GetHistogram()->SetXTitle(projfakey->GetXaxis()->GetTitle()); | |
171 | fakey->GetHistogram()->SetYTitle("fake tracks / findable tracks"); | |
172 | ||
173 | if (includeOverflow) | |
174 | { | |
175 | cout << "findable tracks = " << projfindx->Integral(0, projfindx->GetNbinsX()+1) << endl; | |
176 | cout << "found tracks = " << projeffx->Integral(0, projeffx->GetNbinsX()+1) << endl; | |
177 | cout << "fake tracks = " << projfakex->Integral(0, projfakex->GetNbinsX()+1) << endl; | |
178 | } | |
179 | else | |
180 | { | |
181 | cout << "findable tracks = " << projfindx->Integral(0, projfindx->GetNbinsX()+1) << endl; | |
182 | cout << "found tracks = " << projeffx->Integral(0, projeffx->GetNbinsX()+1) << endl; | |
183 | cout << "fake tracks = " << projfakex->Integral(0, projfakex->GetNbinsX()+1) << endl; | |
184 | } | |
185 | double eff, fake, deff, dfake; | |
186 | EfficiencyIntegral("Total efficiency = ", projeffx, projfindx, eff, deff, includeOverflow); | |
187 | EfficiencyIntegral("Total fake ratio = ", projfakex, projfindx, fake, dfake, includeOverflow); | |
188 | cout << "Sum = " << eff + fake << " +/- " << sqrt(deff*deff + dfake*dfake) << endl; | |
189 | } | |
190 |