Logo ROOT   6.10/09
Reference Guide
TEfficiency.cxx
Go to the documentation of this file.
1 #ifndef ROOT_TEfficiency_cxx
2 #define ROOT_TEfficiency_cxx
3 
4 //standard header
5 #include <vector>
6 #include <string>
7 #include <cmath>
8 #include <stdlib.h>
9 #include <cassert>
10 
11 //ROOT headers
12 #include "Math/DistFuncMathCore.h"
14 #include "TDirectory.h"
15 #include "TF1.h"
16 #include "TGraphAsymmErrors.h"
17 #include "TH1.h"
18 #include "TH2.h"
19 #include "TH3.h"
20 #include "TList.h"
21 #include "TMath.h"
22 #include "TROOT.h"
23 #include "TStyle.h"
24 #include "TVirtualPad.h"
25 #include "TError.h"
26 #include "Math/BrentMinimizer1D.h"
27 #include "Math/WrappedFunction.h"
28 
29 //custom headers
30 #include "TEfficiency.h"
31 
32 // file with extra class for FC method
33 #include "TEfficiencyHelper.h"
34 
35 //default values
38 const Double_t kDefConfLevel = 0.682689492137; // 1 sigma
40 const Double_t kDefWeight = 1;
41 
43 
44 ////////////////////////////////////////////////////////////////////////////////
45 /** \class TEfficiency
46  \ingroup Hist
47  \brief Class to handle efficiency histograms
48 
49 ## I. Overview
50 This class handles the calculation of efficiencies and their uncertainties. It
51 provides several statistical methods for calculating frequentist and bayesian
52 confidence intervals as well as a function for combining several efficiencies.
53 
54 Efficiencies have a lot of applications and meanings but in principle they can
55 be described by the fraction of good/passed events k out of sample containing
56 N events. One is usually interested in the dependency of the efficiency on other
57 (binned) variables. The number of passed and total events is therefore stored
58 internally in two histograms (TEfficiency::fTotalHistogram and TEfficiency::fPassedHistogram).
59 Then the efficiency as well as its upper and lower error an be calculated for each bin
60 individually.
61 
62 As the efficiency can be regarded as a parameter of a binomial distribution, the
63 number of pass ed and total events must always be integer numbers. Therefore a
64 filling with weights is not possible however you can assign a global weight to each
65 TEfficiency object (TEfficiency::SetWeight).
66 It is necessary to create one TEfficiency object
67 for each weight if you investigate a process involving different weights. This
68 procedure needs more effort but enables you to re-use the filled object in cases
69 where you want to change one or more weights. This would not be possible if all
70 events with different weights were filled in the same histogram.
71 
72 ## II. Creating a TEfficiency object
73 If you start a new analysis, it is highly recommended to use the TEfficiency class
74 from the beginning. You can then use one of the constructors for fixed or
75 variable bin size and your desired dimension. These constructors append the
76 created TEfficiency object to the current directory. So it will be written
77 automatically to a file during the next TFile:Write command.
78 
79 Example: create a twodimensional TEfficiency object with
80 - name = "eff"
81 - title = "my efficiency"
82 - axistitles: x, y and LaTeX formated epsilon as label for Z axis
83 - 10 bins with constant bin width (= 1) along X axis starting at 0 (lower edge
84  from first bin) upto 10 (upper edge of last bin)
85 - 20 bins with constant bin width (= 0.5) along Y axis starting at -5 (lower
86  edge from first bin) upto 5 (upper edge of last bin)
87 
88  TEfficiency* pEff = new TEfficiency("eff","my efficiency;x;y;#epsilon",10,0,10,20,-5,5);
89 
90 If you already have two histograms filled with the number of passed and total
91 events, you will use the constructor TEfficiency(const TH1& passed,const TH1& total)
92 to construct the TEfficiency object. The histograms "passed" and "total" have
93 to fulfill the conditions mentioned in TEfficiency::CheckConsistency, otherwise the construction will fail.
94 As the histograms already exist, the new TEfficiency is by default **not** attached
95 to the current directory to avoid duplication of data. If you want to store the
96 new object anyway, you can either write it directly by calling TObject::Write or attach it to a directory using TEfficiency::SetDirectory.
97 This also applies for TEfficiency objects created by the copy constructor TEfficiency::TEfficiency(const TEfficiency& rEff).
98 
99 
100 ### Example 1
101 
102 ~~~~~~~~~~~~~~~{.cpp}
103 TEfficiency* pEff = 0;
104 TFile* pFile = new TFile("myfile.root","recreate");
105 
106 //h_pass and h_total are valid and consistent histograms
107 if(TEfficiency::CheckConsistency(h_pass,h_total))
108 {
109  pEff = new TEfficiency(h_pass,h_total);
110  // this will write the TEfficiency object to "myfile.root"
111  // AND pEff will be attached to the current directory
112  pEff->Write();
113 }
114 ~~~~~~~~~~~~~~~
115 
116 ### Example 2
117 
118 ~~~~~~~~~~~~~~~{.cpp}
119 TEfficiency* pEff = 0;
120 TFile* pFile = new TFile("myfile.root","recreate");
121 
122 //h_pass and h_total are valid and consistent histograms
123 if(TEfficiency::CheckConsistency(h_pass,h_total))
124 {
125  pEff = new TEfficiency(h_pass,h_total);
126  //this will attach the TEfficiency object to the current directory
127  pEff->SetDirectory(gDirectory);
128  //now all objects in gDirectory will be written to "myfile.root"
129  pFile->Write();
130 }
131 ~~~~~~~~~~~~~~~
132 
133 In the case that you already have two filled histograms and you only want to
134 plot them as a graph, you should rather use TGraphAsymmErrors::TGraphAsymmErrors(const TH1* pass,const TH1* total,Option_t* opt)
135 to create a graph object.
136 
137 ## III. Filling with events
138 You can fill the TEfficiency object by calling the TEfficiency::Fill(Bool_t bPassed,Double_t x,Double_t y,Double_t z) method.
139 The boolean flag "bPassed" indicates whether the current event is a good
140  (both histograms are filled) or not (only TEfficiency::fTotalHistogram is filled).
141 The variables x,y and z determine the bin which is filled. For lower dimensions the z- or even the y-value may be omitted.
142 
143 Begin_Macro(source)
144 {
145  //canvas only needed for this documentation
146  TCanvas* c1 = new TCanvas("example","",600,400);
147  c1->SetFillStyle(1001);
148  c1->SetFillColor(kWhite);
149 
150  //create one-dimensional TEfficiency object with fixed bin size
151  TEfficiency* pEff = new TEfficiency("eff","my efficiency;x;#epsilon",20,0,10);
152  TRandom3 rand3;
153 
154  bool bPassed;
155  double x;
156  for(int i=0; i<10000; ++i)
157  {
158  //simulate events with variable under investigation
159  x = rand3.Uniform(10);
160  //check selection: bPassed = DoesEventPassSelection(x)
161  bPassed = rand3.Rndm() < TMath::Gaus(x,5,4);
162  pEff->Fill(bPassed,x);
163  }
164 
165  pEff->Draw("AP");
166 
167  //only for this documentation
168  return c1;
169 }
170 End_Macro
171 
172 You can also set the number of passed or total events for a bin directly by
173 using the TEfficiency::SetPassedEvents or TEfficiency::SetTotalEvents method.
174 
175 ## IV. Statistic options
176 The calculation of the estimated efficiency depends on the chosen statistic
177 option. Let k denotes the number of passed events and N the number of total
178 events.
179 
180 ###Frequentist methods
181 The expectation value of the number of passed events is given by the true
182 efficiency times the total number of events. One can estimate the efficiency
183 by replacing the expected number of passed events by the observed number of
184 passed events.
185 
186 \f[
187  k = \epsilon \times N \Rightarrow \hat{\varepsilon} = \frac{k}{N}
188 \f]
189 
190 ### Bayesian methods
191 In bayesian statistics a likelihood-function (how probable is it to get the
192 observed data assuming a true efficiency) and a prior probability (what is the
193 probability that a certain true efficiency is actually realised) are used to
194 determine a posterior probability by using Bayes theorem. At the moment, only
195 beta distributions (have 2 free parameters) are supported as prior
196 probabilities.
197 
198 \f{eqnarray*}{
199  P(\epsilon | k ; N) &=& \frac{1}{norm} \times P(k | \epsilon ; N) \times Prior(\epsilon) \\
200  P(k | \epsilon ; N) &=& Binomial(N,k) \times \epsilon^{k} \times (1 - \epsilon)^{N - k} ...\ binomial\ distribution \\
201  Prior(\epsilon) &=& \frac{1}{B(\alpha,\beta)} \times \epsilon ^{\alpha - 1} \times (1 - \epsilon)^{\beta - 1} \equiv Beta(\epsilon; \alpha,\beta) \\
202  \Rightarrow P(\epsilon | k ; N) &=& \frac{1}{norm'} \times \epsilon^{k + \alpha - 1} \times (1 - \epsilon)^{N - k + \beta - 1} \equiv Beta(\epsilon; k + \alpha, N - k + \beta)
203 \f}
204 
205 By default the expectation value of this posterior distribution is used as estimator for the efficiency:
206 
207 \f[
208  \hat{\varepsilon} = \frac{k + \alpha}{N + \alpha + \beta}
209 \f]
210 
211 Optionally the mode can also be used as value for the estimated efficiency. This can be done by calling
212 SetBit(kPosteriorMode) or TEfficiency::SetPosteriorMode. In this case the estimated efficiency is:
213 
214 \f[
215  \hat{\varepsilon} = \frac{k + \alpha -1}{N + \alpha + \beta - 2}
216 \f]
217 
218 In the case of a uniform prior distribution, B(x,1,1), the posterior mode is k/n, equivalent to the frequentist
219 estimate (the maximum likelihood value).
220 
221 The statistic options also specify which confidence interval is used for calculating
222 the uncertainties of the efficiency. The following properties define the error
223 calculation:
224 - **fConfLevel:** desired confidence level: 0 < fConfLevel < 1 (TEfficiency::GetConfidenceLevel / TEfficiency::SetConfidenceLevel)
225 - **fStatisticOption** defines which method is used to calculate the boundaries of the confidence interval (TEfficiency::SetStatisticOption)
226 - **fBeta_alpha, fBeta_beta:** parameters for the prior distribution which is only used in the bayesian case (TEfficiency::GetBetaAlpha / TEfficiency::GetBetaBeta / TEfficiency::SetBetaAlpha / TEfficiency::SetBetaBeta)
227 - **kIsBayesian:** flag whether bayesian statistics are used or not (TEfficiency::UsesBayesianStat)
228 - **kShortestInterval:** flag whether shortest interval (instead of central one) are used in case of Bayesian statistics (TEfficiency::UsesShortestInterval). Normally shortest interval should be used in combination with the mode (see TEfficiency::UsesPosteriorMode)
229 - **fWeight:** global weight for this TEfficiency object which is used during combining or merging with other TEfficiency objects(TEfficiency::GetWeight / TEfficiency::SetWeight)
230 
231 In the following table the implemented confidence intervals are listed
232 with their corresponding statistic option. For more details on the calculation,
233 please have a look at the the mentioned functions.
234 
235 
236 | name | statistic option | function | kIsBayesian | parameters |
237 |------------------|------------------|---------------------|-------------|------------|
238 | Clopper-Pearson | kFCP | TEfficiency::ClopperPearson |false |total events, passed events, confidence level |
239 | normal approximation | kFNormal | TEfficiency::Normal | false | total events, passed events, confidence level |
240 | Wilson | kFWilson | TEfficiency::Wilson | false | total events, passed events, confidence level |
241 | Agresti-Coull | kFAC | TEfficiency::AgrestiCoull | false | total events, passed events. confidence level |
242 | Feldman-Cousins | kFFC | TEfficiency::FeldmanCousins | false | total events, passed events, confidence level |
243 | Mid-P Lancaster | kMidP | TEfficiency::MidPInterval | false | total events, passed events, confidence level |
244 | Jeffrey | kBJeffrey | TEfficiency::Bayesian | true | total events, passed events, confidence level, fBeta_alpha = 0.5, fBeta_beta = 0.5 |
245 | Uniform prior | kBUniform |TEfficiency::Bayesian | true |total events, passed events, confidence level, fBeta_alpha = 1, fBeta_beta = 1 |
246 | custom prior | kBBayesian |TEfficiency::Bayesian | true |total events, passed events, confidence level, fBeta_alpha, fBeta_beta |
247 
248 The following example demonstrates the effect of different statistic options and
249 confidence levels.
250 
251 Begin_Macro(source)
252 {
253  //canvas only needed for the documentation
254  TCanvas* c1 = new TCanvas("c1","",600,400);
255  c1->Divide(2);
256  c1->SetFillStyle(1001);
257  c1->SetFillColor(kWhite);
258 
259  //create one-dimensional TEfficiency object with fixed bin size
260  TEfficiency* pEff = new TEfficiency("eff","different confidence levels;x;#epsilon",20,0,10);
261  TRandom3 rand3;
262 
263  bool bPassed;
264  double x;
265  for(int i=0; i<1000; ++i)
266  {
267  //simulate events with variable under investigation
268  x = rand3.Uniform(10);
269  //check selection: bPassed = DoesEventPassSelection(x)
270  bPassed = rand3.Rndm() < TMath::Gaus(x,5,4);
271  pEff->Fill(bPassed,x);
272  }
273 
274  //set style attributes
275  pEff->SetFillStyle(3004);
276  pEff->SetFillColor(kRed);
277 
278  //copy current TEfficiency object and set new confidence level
279  TEfficiency* pCopy = new TEfficiency(*pEff);
280  pCopy->SetConfidenceLevel(0.90);
281 
282  //set style attributes
283  pCopy->SetFillStyle(3005);
284  pCopy->SetFillColor(kBlue);
285 
286  c1->cd(1);
287 
288  //add legend
289  TLegend* leg1 = new TLegend(0.3,0.1,0.7,0.5);
290  leg1->AddEntry(pEff,"95%","F");
291  leg1->AddEntry(pCopy,"68.3%","F");
292 
293  pEff->Draw("A4");
294  pCopy->Draw("same4");
295  leg1->Draw("same");
296 
297  //use same confidence level but different statistic methods
298  TEfficiency* pEff2 = new TEfficiency(*pEff);
299  TEfficiency* pCopy2 = new TEfficiency(*pEff);
300 
301  pEff2->SetStatisticOption(TEfficiency::kFNormal);
302  pCopy2->SetStatisticOption(TEfficiency::kFAC);
303 
304  pEff2->SetTitle("different statistic options;x;#epsilon");
305 
306  //set style attributes
307  pCopy2->SetFillStyle(3005);
308  pCopy2->SetFillColor(kBlue);
309 
310  c1->cd(2);
311 
312  //add legend
313  TLegend* leg2 = new TLegend(0.3,0.1,0.7,0.5);
314  leg2->AddEntry(pEff2,"kFNormal","F");
315  leg2->AddEntry(pCopy2,"kFAC","F");
316 
317  pEff2->Draw("a4");
318  pCopy2->Draw("same4");
319  leg2->Draw("same");
320 
321  //only for this documentation
322  c1->cd(0);
323  return c1;
324 }
325 End_Macro
326 
327 The prior probability of the efficiency in bayesian statistics can be given
328 in terms of a beta distribution. The beta distribution has to positive shape
329 parameters. The resulting priors for different combinations of these shape
330 parameters are shown in the plot below.
331 
332 Begin_Macro(source)
333 {
334  //canvas only needed for the documentation
335  TCanvas* c1 = new TCanvas("c1","",600,400);
336  c1->SetFillStyle(1001);
337  c1->SetFillColor(kWhite);
338 
339  //create different beta distributions
340  TF1* f1 = new TF1("f1","TMath::BetaDist(x,1,1)",0,1);
341  f1->SetLineColor(kBlue);
342  TF1* f2 = new TF1("f2","TMath::BetaDist(x,0.5,0.5)",0,1);
343  f2->SetLineColor(kRed);
344  TF1* f3 = new TF1("f3","TMath::BetaDist(x,1,5)",0,1);
345  f3->SetLineColor(kGreen+3);
346  f3->SetTitle("Beta distributions as priors;#epsilon;P(#epsilon)");
347  TF1* f4 = new TF1("f4","TMath::BetaDist(x,4,3)",0,1);
348  f4->SetLineColor(kViolet);
349 
350  //add legend
351  TLegend* leg = new TLegend(0.25,0.5,0.85,0.89);
352  leg->SetFillColor(kWhite);
353  leg->SetFillStyle(1001);
354  leg->AddEntry(f1,"a=1, b=1","L");
355  leg->AddEntry(f2,"a=0.5, b=0.5","L");
356  leg->AddEntry(f3,"a=1, b=5","L");
357  leg->AddEntry(f4,"a=4, b=3","L");
358 
359  f3->Draw();
360  f1->Draw("same");
361  f2->Draw("Same");
362  f4->Draw("same");
363  leg->Draw("same");
364 
365  //only for this documentation
366  return c1;
367 }
368 End_Macro
369 
370 
371 ## IV.1 Coverage probabilities for different methods
372 The following pictures illustrate the actual coverage probability for the
373 different values of the true efficiency and the total number of events when a
374 confidence level of 95% is desired.
375 
376 \image html normal95.gif "Normal Approximation"
377 
378 
379 \image html wilson95.gif "Wilson"
380 
381 
382 \image html ac95.gif "Agresti Coull"
383 
384 
385 \image html cp95.gif "Clopper Pearson"
386 
387 
388 \image html uni95.gif "Bayesian with Uniform Prior"
389 
390 
391 \image html jeffrey95.gif "Bayesian with Jeffrey Prior"
392 
393 The average (over all possible true efficiencies) coverage probability for
394 different number of total events is shown in the next picture.
395 \image html av_cov.png "Average Coverage"
396 
397 ## V. Merging and combining TEfficiency objects
398 In many applications the efficiency should be calculated for an inhomogeneous
399 sample in the sense that it contains events with different weights. In order
400 to be able to determine the correct overall efficiency, it is necessary to
401 use for each subsample (= all events with the same weight) a different
402 TEfficiency object. After finishing your analysis you can then construct the
403 overall efficiency with its uncertainty.
404 
405 This procedure has the advantage that you can change the weight of one
406 subsample easily without rerunning the whole analysis. On the other hand more
407 effort is needed to handle several TEfficiency objects instead of one
408 histogram. In the case of many different or even continuously distributed
409 weights this approach becomes cumbersome. One possibility to overcome this
410 problem is the usage of binned weights.
411 
412 ### Example
413 In particle physics weights arises from the fact that you want to
414 normalise your results to a certain reference value. A very common formula for
415 calculating weights is
416 
417 \f{eqnarray*}{
418  w &=& \frac{\sigma L}{N_{gen} \epsilon_{trig}} \\
419  &-& \sigma ...\ cross\ section \\
420  &-& L ...\ luminosity \\
421  &-& N_{gen}\ ... number\ of\ generated\ events \\
422  &-& \epsilon_{trig}\ ...\ (known)\ trigger\ efficiency \\
423 \f}
424 
425 The reason for different weights can therefore be:
426 - different processes
427 - other integrated luminosity
428 - varying trigger efficiency
429 - different sample sizes
430 - ...
431 - or even combination of them
432 
433 Depending on the actual meaning of different weights in your case, you
434 should either merge or combine them to get the overall efficiency.
435 
436 ### V.1 When should I use merging?
437 If the weights are artificial and do not represent real alternative hypotheses,
438 you should merge the different TEfficiency objects. That means especially for
439 the bayesian case that the prior probability should be the same for all merged
440 TEfficiency objects. The merging can be done by invoking one of the following
441 operations:
442 - eff1.Add(eff2)
443 - eff1 += eff2
444 - eff1 = eff1 + eff2
445 
446 The result of the merging is stored in the TEfficiency object which is marked
447 bold above. The contents of the internal histograms of both TEfficiency
448 objects are added and a new weight is assigned. The statistic options are not
449 changed.
450 
451 \f[
452  \frac{1}{w_{new}} = \frac{1}{w_{1}} + \frac{1}{w_{2}}
453 \f]
454 
455 ### Example:
456 If you use two samples with different numbers of generated events for the same
457 process and you want to normalise both to the same integrated luminosity and
458 trigger efficiency, the different weights then arise just from the fact that
459 you have different numbers of events. The TEfficiency objects should be merged
460 because the samples do not represent true alternatives. You expect the same
461 result as if you would have a big sample with all events in it.
462 
463 \f[
464  w_{1} = \frac{\sigma L}{\epsilon N_{1}}, w_{2} = \frac{\sigma L}{\epsilon N_{2}} \Rightarrow w_{new} = \frac{\sigma L}{\epsilon (N_{1} + N_{2})} = \frac{1}{\frac{1}{w_{1}} + \frac{1}{w_{2}}}
465 \f]
466 
467 ### V.2 When should I use combining?
468 You should combine TEfficiency objects whenever the weights represent
469 alternatives processes for the efficiency. As the combination of two TEfficiency
470 objects is not always consistent with the representation by two internal
471 histograms, the result is not stored in a TEfficiency object but a TGraphAsymmErrors
472 is returned which shows the estimated combined efficiency and its uncertainty
473 for each bin.
474 At the moment the combination method TEfficiency::Combine only supports combination of 1-dimensional
475 efficiencies in a bayesian approach.
476 
477 
478 For calculating the combined efficiency and its uncertainty for each bin only Bayesian statistics
479 is used. No frequentists methods are presently supported for computing the combined efficiency and
480 its confidence interval.
481 In the case of the Bayesian statistics a combined posterior is constructed taking into account the
482 weight of each TEfficiency object. The same prior is used for all the TEfficiency objects.
483 
484 \f{eqnarray*}{
485  P_{comb}(\epsilon | {w_{i}}, {k_{i}} , {N_{i}}) = \frac{1}{norm} \prod_{i}{L(k_{i} | N_{i}, \epsilon)}^{w_{i}} \Pi( \epsilon )\\
486 L(k_{i} | N_{i}, \epsilon)\ is\ the\ likelihood\ function\ for\ the\ sample\ i\ (a\ Binomial\ distribution)\\
487 \Pi( \epsilon)\ is\ the\ prior,\ a\ beta\ distribution\ B(\epsilon, \alpha, \beta).\\
488 The\ resulting\ combined\ posterior\ is \\
489 P_{comb}(\epsilon |{w_{i}}; {k_{i}}; {N_{i}}) = B(\epsilon, \sum_{i}{ w_{i} k_{i}} + \alpha, \sum_{i}{ w_{i}(n_{i}-k_{i})}+\beta) \\
490 \hat{\varepsilon} = \int_{0}^{1} \epsilon \times P_{comb}(\epsilon | {k_{i}} , {N_{i}}) d\epsilon \\
491 confidence\ level = 1 - \alpha \\
492 \frac{\alpha}{2} = \int_{0}^{\epsilon_{low}} P_{comb}(\epsilon | {k_{i}} , {N_{i}}) d\epsilon ...\ defines\ lower\ boundary \\
493 1- \frac{\alpha}{2} = \int_{0}^{\epsilon_{up}} P_{comb}(\epsilon | {k_{i}} , {N_{i}}) d\epsilon ...\ defines\ upper\ boundary
494 \f}
495 
496 
497 ###Example:
498 If you use cuts to select electrons which can originate from two different
499 processes, you can determine the selection efficiency for each process. The
500 overall selection efficiency is then the combined efficiency. The weights to be used in the
501 combination should be the probability that an
502 electron comes from the corresponding process.
503 
504 \f[
505 p_{1} = \frac{\sigma_{1}}{\sigma_{1} + \sigma_{2}} = \frac{N_{1}w_{1}}{N_{1}w_{1} + N_{2}w_{2}}\\
506 p_{2} = \frac{\sigma_{2}}{\sigma_{1} + \sigma_{2}} = \frac{N_{2}w_{2}}{N_{1}w_{1} + N_{2}w_{2}}
507 \f]
508 
509 ## VI. Further operations
510 
511 ### VI.Information about the internal histograms
512 The methods TEfficiency::GetPassedHistogram and TEfficiency::GetTotalHistogram
513 return a constant pointer to the internal histograms. They can be used to
514 obtain information about the internal histograms (e.g. the binning, number of passed / total events in a bin, mean values...).
515 One can obtain a clone of the internal histograms by calling TEfficiency::GetCopyPassedHisto or TEfficiency::GetCopyTotalHisto.
516 The returned histograms are completely independent from the current
517 TEfficiency object. By default, they are not attached to a directory to
518 avoid the duplication of data and the user is responsible for deleting them.
519 
520 
521 ~~~~~~~~~~~~~~~{.cpp}
522 //open a root file which contains a TEfficiency object
523 TFile* pFile = new TFile("myfile.root","update");
524 
525 //get TEfficiency object with name "my_eff"
526 TEfficiency* pEff = (TEfficiency*)pFile->Get("my_eff");
527 
528 //get clone of total histogram
529 TH1* clone = pEff->GetCopyTotalHisto();
530 
531 //change clone...
532 //save changes of clone directly
533 clone->Write();
534 //or append it to the current directory and write the file
535 //clone->SetDirectory(gDirectory);
536 //pFile->Write();
537 
538 //delete histogram object
539 delete clone;
540 clone = 0;
541 ~~~~~~~~~~~~~~~
542 
543 It is also possible to set the internal total or passed histogram by using the
544 methods TEfficiency::SetPassedHistogram or TEfficiency::SetTotalHistogram.
545 
546 In order to ensure the validity of the TEfficiency object, the consistency of the
547 new histogram and the stored histogram is checked. It sometimes might be
548 impossible to change the histograms in a consistent way. Therefore one can force
549 the replacement by passing the option "f". Then the user has to ensure that the
550 other internal histogram is replaced as well and that the TEfficiency object is
551 in a valid state.
552 
553 ### VI.2 Fitting
554 The efficiency can be fitted using the TEfficiency::Fit function which uses
555 internally the TBinomialEfficiencyFitter::Fit method.
556 As this method is using a maximum-likelihood-fit, it is necessary to initialise
557 the given fit function with reasonable start values.
558 The resulting fit function is attached to the list of associated functions and
559 will be drawn automatically during the next TEfficiency::Draw command.
560 The list of associated function can be modified by using the pointer returned
561 by TEfficiency::GetListOfFunctions.
562 
563 Begin_Macro(source)
564 {
565  //canvas only needed for this documentation
566  TCanvas* c1 = new TCanvas("example","",600,400);
567  c1->SetFillStyle(1001);
568  c1->SetFillColor(kWhite);
569 
570  //create one-dimensional TEfficiency object with fixed bin size
571  TEfficiency* pEff = new TEfficiency("eff","my efficiency;x;#epsilon",20,0,10);
572  TRandom3 rand3;
573 
574  bool bPassed;
575  double x;
576  for(int i=0; i<10000; ++i)
577  {
578  //simulate events with variable under investigation
579  x = rand3.Uniform(10);
580  //check selection: bPassed = DoesEventPassSelection(x)
581  bPassed = rand3.Rndm() < TMath::Gaus(x,5,4);
582  pEff->Fill(bPassed,x);
583  }
584 
585  //create a function for fitting and do the fit
586  TF1* f1 = new TF1("f1","gaus",0,10);
587  f1->SetParameters(1,5,2);
588  pEff->Fit(f1);
589 
590  //create a threshold function
591  TF1* f2 = new TF1("thres","0.8",0,10);
592  f2->SetLineColor(kRed);
593  //add it to the list of functions
594  //use add first because the parameters of the last function will be displayed
595  pEff->GetListOfFunctions()->AddFirst(f2);
596 
597  pEff->Draw("AP");
598 
599  //only for this documentation
600  return c1;
601 }
602 End_Macro
603 
604 ### VI.3 Draw a TEfficiency object
605 A TEfficiency object can be drawn by calling the usual TEfficiency::Draw method.
606 At the moment drawing is only supported for 1- and 2-dimensional TEfficiency objects.
607 In the 1-dimensional case you can use the same options as for the TGraphAsymmErrors::Draw
608 method. For 2-dimensional TEfficiency objects you can pass the same options as
609 for a TH2::Draw object.
610 
611 ********************************************************************************/
612 //______________________________________________________________________________
613 
614 ////////////////////////////////////////////////////////////////////////////////
615 ///default constructor
616 ///
617 ///should not be used explicitly
618 
620 fBeta_alpha(kDefBetaAlpha),
621 fBeta_beta(kDefBetaBeta),
622 fBoundary(0),
623 fConfLevel(kDefConfLevel),
624 fDirectory(0),
625 fFunctions(0),
626 fPaintGraph(0),
627 fPaintHisto(0),
628 fPassedHistogram(0),
629 fTotalHistogram(0),
630 fWeight(kDefWeight)
631 {
632  SetStatisticOption(kDefStatOpt);
633 
634  // create 2 dummy histograms
635  fPassedHistogram = new TH1F("h_passed","passed",10,0,10);
636  fTotalHistogram = new TH1F("h_total","total",10,0,10);
637 }
638 
639 ////////////////////////////////////////////////////////////////////////////////
640 ///constructor using two existing histograms as input
641 ///
642 ///Input: passed - contains the events fulfilling some criteria
643 /// total - contains all investigated events
644 ///
645 ///Notes: - both histograms have to fulfill the conditions of CheckConsistency
646 /// - dimension of the resulting efficiency object depends
647 /// on the dimension of the given histograms
648 /// - Clones of both histograms are stored internally
649 /// - The function SetName(total.GetName() + "_clone") is called to set
650 /// the names of the new object and the internal histograms..
651 /// - The created TEfficiency object is NOT appended to a directory. It
652 /// will not be written to disk during the next TFile::Write() command
653 /// in order to prevent duplication of data. If you want to save this
654 /// TEfficiency object anyway, you can either append it to a
655 /// directory by calling SetDirectory(TDirectory*) or write it
656 /// explicitly to disk by calling Write().
657 
658 TEfficiency::TEfficiency(const TH1& passed,const TH1& total):
659 fBeta_alpha(kDefBetaAlpha),
660 fBeta_beta(kDefBetaBeta),
661 fConfLevel(kDefConfLevel),
662 fDirectory(0),
663 fFunctions(0),
664 fPaintGraph(0),
665 fPaintHisto(0),
666 fWeight(kDefWeight)
667 {
668  //check consistency of histograms
669  if(CheckConsistency(passed,total)) {
670  Bool_t bStatus = TH1::AddDirectoryStatus();
672  fTotalHistogram = (TH1*)total.Clone();
673  fPassedHistogram = (TH1*)passed.Clone();
674  TH1::AddDirectory(bStatus);
675 
676  TString newName = total.GetName();
677  newName += TString("_clone");
678  SetName(newName);
679 
680  // are the histograms filled with weights?
681  if(CheckWeights(passed,total))
682  {
683  Info("TEfficiency","given histograms are filled with weights");
685  }
686  }
687  else {
688  Error("TEfficiency(const TH1&,const TH1&)","histograms are not consistent -> results are useless");
689  Warning("TEfficiency(const TH1&,const TH1&)","using two empty TH1D('h1','h1',10,0,10)");
690 
691  Bool_t bStatus = TH1::AddDirectoryStatus();
693  fTotalHistogram = new TH1D("h1_total","h1 (total)",10,0,10);
694  fPassedHistogram = new TH1D("h1_passed","h1 (passed)",10,0,10);
695  TH1::AddDirectory(bStatus);
696  }
697 
698  SetBit(kPosteriorMode,false);
699  SetBit(kShortestInterval,false);
700 
702  SetDirectory(0);
703 }
704 
705 ////////////////////////////////////////////////////////////////////////////////
706 /// Create 1-dimensional TEfficiency object with variable bin size.
707 ///
708 /// Constructor creates two new and empty histograms with a given binning
709 ///
710 /// Input:
711 ///
712 /// - `name`: the common part of the name for both histograms (no blanks)
713 /// fTotalHistogram has name: name + "_total"
714 /// fPassedHistogram has name: name + "_passed"
715 /// - `title`: the common part of the title for both histogram
716 /// fTotalHistogram has title: title + " (total)"
717 /// fPassedHistogram has title: title + " (passed)"
718 /// It is possible to label the axis by passing a title with
719 /// the following format: "title;xlabel;ylabel".
720 /// - `nbins`: number of bins on the x-axis
721 /// - `xbins`: array of length (nbins + 1) with low-edges for each bin
722 /// xbins[nbinsx] ... lower edge for overflow bin
723 
724 TEfficiency::TEfficiency(const char* name,const char* title,Int_t nbins,
725  const Double_t* xbins):
729 fDirectory(0),
730 fFunctions(0),
731 fPaintGraph(0),
732 fPaintHisto(0),
734 {
735  Bool_t bStatus = TH1::AddDirectoryStatus();
737  fTotalHistogram = new TH1D("total","total",nbins,xbins);
738  fPassedHistogram = new TH1D("passed","passed",nbins,xbins);
739  TH1::AddDirectory(bStatus);
740 
741  Build(name,title);
742 }
743 
744 ////////////////////////////////////////////////////////////////////////////////
745 /// Create 1-dimensional TEfficiency object with fixed bins size.
746 ///
747 /// Constructor creates two new and empty histograms with a fixed binning.
748 ///
749 /// Input:
750 ///
751 /// - `name`: the common part of the name for both histograms(no blanks)
752 /// fTotalHistogram has name: name + "_total"
753 /// fPassedHistogram has name: name + "_passed"
754 /// - `title`: the common part of the title for both histogram
755 /// fTotalHistogram has title: title + " (total)"
756 /// fPassedHistogram has title: title + " (passed)"
757 /// It is possible to label the axis by passing a title with
758 /// the following format: "title;xlabel;ylabel".
759 /// - `nbinsx`: number of bins on the x-axis
760 /// - `xlow`: lower edge of first bin
761 /// - `xup`: upper edge of last bin
762 
763 TEfficiency::TEfficiency(const char* name,const char* title,Int_t nbinsx,
764  Double_t xlow,Double_t xup):
768 fDirectory(0),
769 fFunctions(0),
770 fPaintGraph(0),
771 fPaintHisto(0),
773 {
774  Bool_t bStatus = TH1::AddDirectoryStatus();
776  fTotalHistogram = new TH1D("total","total",nbinsx,xlow,xup);
777  fPassedHistogram = new TH1D("passed","passed",nbinsx,xlow,xup);
778  TH1::AddDirectory(bStatus);
779 
780  Build(name,title);
781 }
782 
783 ////////////////////////////////////////////////////////////////////////////////
784 /// Create 2-dimensional TEfficiency object with fixed bin size.
785 ///
786 /// Constructor creates two new and empty histograms with a fixed binning.
787 ///
788 /// Input:
789 ///
790 /// - `name`: the common part of the name for both histograms(no blanks)
791 /// fTotalHistogram has name: name + "_total"
792 /// fPassedHistogram has name: name + "_passed"
793 /// - `title`: the common part of the title for both histogram
794 /// fTotalHistogram has title: title + " (total)"
795 /// fPassedHistogram has title: title + " (passed)"
796 /// It is possible to label the axis by passing a title with
797 /// the following format: "title;xlabel;ylabel;zlabel".
798 /// - `nbinsx`: number of bins on the x-axis
799 /// - `xlow`: lower edge of first x-bin
800 /// - `xup`: upper edge of last x-bin
801 /// - `nbinsy`: number of bins on the y-axis
802 /// - `ylow`: lower edge of first y-bin
803 /// - `yup`: upper edge of last y-bin
804 
805 TEfficiency::TEfficiency(const char* name,const char* title,Int_t nbinsx,
806  Double_t xlow,Double_t xup,Int_t nbinsy,
807  Double_t ylow,Double_t yup):
811 fDirectory(0),
812 fFunctions(0),
813 fPaintGraph(0),
814 fPaintHisto(0),
816 {
817  Bool_t bStatus = TH1::AddDirectoryStatus();
819  fTotalHistogram = new TH2D("total","total",nbinsx,xlow,xup,nbinsy,ylow,yup);
820  fPassedHistogram = new TH2D("passed","passed",nbinsx,xlow,xup,nbinsy,ylow,yup);
821  TH1::AddDirectory(bStatus);
822 
823  Build(name,title);
824 }
825 
826 ////////////////////////////////////////////////////////////////////////////////
827 /// Create 2-dimensional TEfficiency object with variable bin size.
828 ///
829 /// Constructor creates two new and empty histograms with a given binning.
830 ///
831 /// Input:
832 ///
833 /// - `name`: the common part of the name for both histograms(no blanks)
834 /// fTotalHistogram has name: name + "_total"
835 /// fPassedHistogram has name: name + "_passed"
836 /// - `title`: the common part of the title for both histogram
837 /// fTotalHistogram has title: title + " (total)"
838 /// fPassedHistogram has title: title + " (passed)"
839 /// It is possible to label the axis by passing a title with
840 /// the following format: "title;xlabel;ylabel;zlabel".
841 /// - `nbinsx`: number of bins on the x-axis
842 /// - `xbins`: array of length (nbins + 1) with low-edges for each bin
843 /// xbins[nbinsx] ... lower edge for overflow x-bin
844 /// - `nbinsy`: number of bins on the y-axis
845 /// - `ybins`: array of length (nbins + 1) with low-edges for each bin
846 /// ybins[nbinsy] ... lower edge for overflow y-bin
847 
848 TEfficiency::TEfficiency(const char* name,const char* title,Int_t nbinsx,
849  const Double_t* xbins,Int_t nbinsy,
850  const Double_t* ybins):
854 fDirectory(0),
855 fFunctions(0),
856 fPaintGraph(0),
857 fPaintHisto(0),
859 {
860  Bool_t bStatus = TH1::AddDirectoryStatus();
862  fTotalHistogram = new TH2D("total","total",nbinsx,xbins,nbinsy,ybins);
863  fPassedHistogram = new TH2D("passed","passed",nbinsx,xbins,nbinsy,ybins);
864  TH1::AddDirectory(bStatus);
865 
866  Build(name,title);
867 }
868 
869 ////////////////////////////////////////////////////////////////////////////////
870 /// Create 3-dimensional TEfficiency object with fixed bin size.
871 ///
872 /// Constructor creates two new and empty histograms with a fixed binning.
873 ///
874 /// Input:
875 ///
876 /// - `name`: the common part of the name for both histograms(no blanks)
877 /// fTotalHistogram has name: name + "_total"
878 /// fPassedHistogram has name: name + "_passed"
879 /// - `title`: the common part of the title for both histogram
880 /// fTotalHistogram has title: title + " (total)"
881 /// fPassedHistogram has title: title + " (passed)"
882 /// It is possible to label the axis by passing a title with
883 /// the following format: "title;xlabel;ylabel;zlabel".
884 /// - `nbinsx`: number of bins on the x-axis
885 /// - `xlow`: lower edge of first x-bin
886 /// - `xup`: upper edge of last x-bin
887 /// - `nbinsy`: number of bins on the y-axis
888 /// - `ylow`: lower edge of first y-bin
889 /// - `yup`: upper edge of last y-bin
890 /// - `nbinsz`: number of bins on the z-axis
891 /// - `zlow`: lower edge of first z-bin
892 /// - `zup`: upper edge of last z-bin
893 
894 TEfficiency::TEfficiency(const char* name,const char* title,Int_t nbinsx,
895  Double_t xlow,Double_t xup,Int_t nbinsy,
896  Double_t ylow,Double_t yup,Int_t nbinsz,
897  Double_t zlow,Double_t zup):
901 fDirectory(0),
902 fFunctions(0),
903 fPaintGraph(0),
904 fPaintHisto(0),
906 {
907  Bool_t bStatus = TH1::AddDirectoryStatus();
909  fTotalHistogram = new TH3D("total","total",nbinsx,xlow,xup,nbinsy,ylow,yup,nbinsz,zlow,zup);
910  fPassedHistogram = new TH3D("passed","passed",nbinsx,xlow,xup,nbinsy,ylow,yup,nbinsz,zlow,zup);
911  TH1::AddDirectory(bStatus);
912 
913  Build(name,title);
914 }
915 
916 ////////////////////////////////////////////////////////////////////////////////
917 /// Create 3-dimensional TEfficiency object with variable bin size.
918 ///
919 /// Constructor creates two new and empty histograms with a given binning.
920 ///
921 /// Input:
922 ///
923 /// - `name`: the common part of the name for both histograms(no blanks)
924 /// fTotalHistogram has name: name + "_total"
925 /// fPassedHistogram has name: name + "_passed"
926 /// - `title`: the common part of the title for both histogram
927 /// fTotalHistogram has title: title + " (total)"
928 /// fPassedHistogram has title: title + " (passed)"
929 /// It is possible to label the axis by passing a title with
930 /// the following format: "title;xlabel;ylabel;zlabel".
931 /// - `nbinsx`: number of bins on the x-axis
932 /// - `xbins`: array of length (nbins + 1) with low-edges for each bin
933 /// xbins[nbinsx] ... lower edge for overflow x-bin
934 /// - `nbinsy`: number of bins on the y-axis
935 /// - `ybins`: array of length (nbins + 1) with low-edges for each bin
936 /// xbins[nbinsx] ... lower edge for overflow y-bin
937 /// - `nbinsz`: number of bins on the z-axis
938 /// - `zbins`: array of length (nbins + 1) with low-edges for each bin
939 /// xbins[nbinsx] ... lower edge for overflow z-bin
940 
941 TEfficiency::TEfficiency(const char* name,const char* title,Int_t nbinsx,
942  const Double_t* xbins,Int_t nbinsy,
943  const Double_t* ybins,Int_t nbinsz,
944  const Double_t* zbins):
948 fDirectory(0),
949 fFunctions(0),
950 fPaintGraph(0),
951 fPaintHisto(0),
953 {
954  Bool_t bStatus = TH1::AddDirectoryStatus();
956  fTotalHistogram = new TH3D("total","total",nbinsx,xbins,nbinsy,ybins,nbinsz,zbins);
957  fPassedHistogram = new TH3D("passed","passed",nbinsx,xbins,nbinsy,ybins,nbinsz,zbins);
958  TH1::AddDirectory(bStatus);
959 
960  Build(name,title);
961 }
962 
963 ////////////////////////////////////////////////////////////////////////////////
964 /// Copy constructor.
965 ///
966 ///The list of associated objects (e.g. fitted functions) is not copied.
967 ///
968 ///Note:
969 ///
970 /// - SetName(rEff.GetName() + "_copy") is called to set the names of the
971 /// object and the histograms.
972 /// - The titles are set by calling SetTitle("[copy] " + rEff.GetTitle()).
973 /// - The copied TEfficiency object is NOT appended to a directory. It
974 /// will not be written to disk during the next TFile::Write() command
975 /// in order to prevent duplication of data. If you want to save this
976 /// TEfficiency object anyway, you can either append it to a directory
977 /// by calling SetDirectory(TDirectory*) or write it explicitly to disk
978 /// by calling Write().
979 
981 TNamed(),
982 TAttLine(),
983 TAttFill(),
984 TAttMarker(),
986 fBeta_beta(rEff.fBeta_beta),
988 fConfLevel(rEff.fConfLevel),
989 fDirectory(0),
990 fFunctions(0),
991 fPaintGraph(0),
992 fPaintHisto(0),
993 fWeight(rEff.fWeight)
994 {
995  // copy TObject bits
996  ((TObject&)rEff).Copy(*this);
997 
998  Bool_t bStatus = TH1::AddDirectoryStatus();
1000  fTotalHistogram = (TH1*)((rEff.fTotalHistogram)->Clone());
1001  fPassedHistogram = (TH1*)((rEff.fPassedHistogram)->Clone());
1002  TH1::AddDirectory(bStatus);
1003 
1004  TString name = rEff.GetName();
1005  name += "_copy";
1006  SetName(name);
1007  TString title = "[copy] ";
1008  title += rEff.GetTitle();
1009  SetTitle(title);
1010 
1012 
1013  SetDirectory(0);
1014 
1015  //copy style
1016  rEff.TAttLine::Copy(*this);
1017  rEff.TAttFill::Copy(*this);
1018  rEff.TAttMarker::Copy(*this);
1019 }
1020 
1021 ////////////////////////////////////////////////////////////////////////////////
1022 ///default destructor
1023 
1025 {
1026  //delete all function in fFunctions
1027  // use same logic as in TH1 destructor
1028  // (see TH1::~TH1 code in TH1.cxx)
1029  if(fFunctions) {
1031  TObject* obj = 0;
1032  while ((obj = fFunctions->First())) {
1033  while(fFunctions->Remove(obj)) { }
1034  if (!obj->TestBit(kNotDeleted)) {
1035  break;
1036  }
1037  delete obj;
1038  obj = 0;
1039  }
1040  delete fFunctions;
1041  fFunctions = 0;
1042  }
1043 
1044  if(fDirectory)
1045  fDirectory->Remove(this);
1046 
1047  delete fTotalHistogram;
1048  delete fPassedHistogram;
1049  delete fPaintGraph;
1050  delete fPaintHisto;
1051 }
1052 
1053 ////////////////////////////////////////////////////////////////////////////////
1054 /**
1055  Calculates the boundaries for the frequentist Agresti-Coull interval
1056 
1057  \param total number of total events
1058  \param passed 0 <= number of passed events <= total
1059  \param level confidence level
1060  \param bUpper true - upper boundary is returned
1061  false - lower boundary is returned
1062 
1063 
1064  \f{eqnarray*}{
1065  \alpha &=& 1 - \frac{level}{2} \\
1066  \kappa &=& \Phi^{-1}(1 - \alpha,1)\ ... normal\ quantile\ function\\
1067  mode &=& \frac{passed + \frac{\kappa^{2}}{2}}{total + \kappa^{2}}\\
1068  \Delta &=& \kappa * \sqrt{\frac{mode * (1 - mode)}{total + \kappa^{2}}}\\
1069  return &=& max(0,mode - \Delta)\ or\ min(1,mode + \Delta)
1070  \f}
1071 
1072 */
1073 
1075 {
1076  Double_t alpha = (1.0 - level)/2;
1077  Double_t kappa = ROOT::Math::normal_quantile(1 - alpha,1);
1078 
1079  Double_t mode = (passed + 0.5 * kappa * kappa) / (total + kappa * kappa);
1080  Double_t delta = kappa * std::sqrt(mode * (1 - mode) / (total + kappa * kappa));
1081 
1082  if(bUpper)
1083  return ((mode + delta) > 1) ? 1.0 : (mode + delta);
1084  else
1085  return ((mode - delta) < 0) ? 0.0 : (mode - delta);
1086 }
1087 
1088 ////////////////////////////////////////////////////////////////////////////////
1089 /// Calculates the boundaries for the frequentist Feldman-Cousins interval
1090 ///
1091 /// \param total number of total events
1092 /// \param passed 0 <= number of passed events <= total
1093 /// \param level confidence level
1094 /// \param bUpper: true - upper boundary is returned
1095 /// false - lower boundary is returned
1096 
1098 {
1099  Double_t lower = 0;
1100  Double_t upper = 1;
1101  if (!FeldmanCousinsInterval(total,passed,level, lower, upper)) {
1102  ::Error("FeldmanCousins","Error running FC method - return 0 or 1");
1103  }
1104  return (bUpper) ? upper : lower;
1105 }
1106 
1107 ////////////////////////////////////////////////////////////////////////////////
1108 /// Calculates the interval boundaries using the frequentist methods of Feldman-Cousins
1109 ///
1110 /// \param[in] total number of total events
1111 /// \param[in] passed 0 <= number of passed events <= total
1112 /// \param[in] level confidence level
1113 /// \param[out] lower lower boundary returned on exit
1114 /// \param[out] upper lower boundary returned on exit
1115 /// \return a flag with the status of the calculation
1116 ///
1117 /// Calculation:
1118 ///
1119 /// The Feldman-Cousins is a frequentist method where the interval is estimated using a Neyman construction where the ordering
1120 /// is based on the likelihood ratio:
1121 /// \f[
1122 /// LR = \frac{Binomial(k | N, \epsilon)}{Binomial(k | N, \hat{\epsilon} ) }
1123 /// \f]
1124 /// See G. J. Feldman and R. D. Cousins, Phys. Rev. D57 (1998) 3873
1125 /// and R. D. Cousins, K. E. Hymes, J. Tucker, Nuclear Instruments and Methods in Physics Research A 612 (2010) 388
1126 ///
1127 /// Implemented using classes developed by Jordan Tucker and Luca Lista
1128 /// See File hist/hist/src/TEfficiencyHelper.h
1129 
1131 {
1133  double alpha = 1.-level;
1134  fc.Init(alpha);
1135  fc.Calculate(passed, total);
1136  lower = fc.Lower();
1137  upper = fc.Upper();
1138  return true;
1139 }
1140 
1141 ////////////////////////////////////////////////////////////////////////////////
1142 /// Calculates the boundaries using the mid-P binomial
1143 /// interval (Lancaster method) from B. Cousing and J. Tucker.
1144 /// See http://arxiv.org/abs/0905.3831 for a description and references for the method
1145 ///
1146 /// Modify equal_tailed to get the kind of interval you want.
1147 /// Can also be converted to interval on ratio of poisson means X/Y by the substitutions
1148 /// ~~~ {.cpp}
1149 /// X = passed
1150 /// total = X + Y
1151 /// lower_poisson = lower/(1 - lower)
1152 /// upper_poisson = upper/(1 - upper)
1153 /// ~~~
1154 
1156 {
1157  const double alpha = 1. - level;
1158  const bool equal_tailed = true; // change if you don;t want equal tailed interval
1159  const double alpha_min = equal_tailed ? alpha/2 : alpha;
1160  const double tol = 1e-9; // tolerance
1161  double pmin = 0;
1162  double pmax = 0;
1163  double p = 0;
1164 
1165  pmin = 0; pmax = 1;
1166 
1167 
1168  // treat special case for 0<passed<1
1169  // do a linear interpolation of the upper limit values
1170  if ( passed > 0 && passed < 1) {
1171  double p0 = MidPInterval(total,0.0,level,bUpper);
1172  double p1 = MidPInterval(total,1.0,level,bUpper);
1173  p = (p1 - p0) * passed + p0;
1174  return p;
1175  }
1176 
1177  while (std::abs(pmax - pmin) > tol) {
1178  p = (pmin + pmax)/2;
1179  //double v = 0.5 * ROOT::Math::binomial_pdf(int(passed), p, int(total));
1180  // make it work for non integer using the binomial - beta relationship
1181  double v = 0.5 * ROOT::Math::beta_pdf(p, passed+1., total-passed+1)/(total+1);
1182  //if (passed > 0) v += ROOT::Math::binomial_cdf(int(passed - 1), p, int(total));
1183  // compute the binomial cdf at passed -1
1184  if ( (passed-1) >= 0) v += ROOT::Math::beta_cdf_c(p, passed, total-passed+1);
1185 
1186  double vmin = (bUpper) ? alpha_min : 1.- alpha_min;
1187  if (v > vmin)
1188  pmin = p;
1189  else
1190  pmax = p;
1191  }
1192 
1193  return p;
1194 }
1195 
1196 
1197 ////////////////////////////////////////////////////////////////////////////////
1198 /**
1199 Calculates the boundaries for a Bayesian confidence interval (shortest or central interval depending on the option)
1200 
1201 \param[in] total number of total events
1202 \param[in] passed 0 <= number of passed events <= total
1203 \param[in] level confidence level
1204 \param[in] alpha shape parameter > 0 for the prior distribution (fBeta_alpha)
1205 \param[in] beta shape parameter > 0 for the prior distribution (fBeta_beta)
1206 \param[in] bUpper
1207  - true - upper boundary is returned
1208  - false - lower boundary is returned
1209 \param[in] bShortest ??
1210 
1211 Note: In the case central confidence interval is calculated.
1212  when passed = 0 (or passed = total) the lower (or upper)
1213  interval values will be larger than 0 (or smaller than 1).
1214 
1215 Calculation:
1216 
1217 The posterior probability in bayesian statistics is given by:
1218 \f[
1219  P(\varepsilon |k,N) \propto L(\varepsilon|k,N) \times Prior(\varepsilon)
1220 \f]
1221 As an efficiency can be interpreted as probability of a positive outcome of
1222 a Bernoullli trial the likelihood function is given by the binomial
1223 distribution:
1224 \f[
1225  L(\varepsilon|k,N) = Binomial(N,k) \varepsilon ^{k} (1 - \varepsilon)^{N-k}
1226 \f]
1227 At the moment only beta distributions are supported as prior probabilities
1228 of the efficiency (\f$ B(\alpha,\beta)\f$ is the beta function):
1229 \f[
1230  Prior(\varepsilon) = \frac{1}{B(\alpha,\beta)} \varepsilon ^{\alpha - 1} (1 - \varepsilon)^{\beta - 1}
1231 \f]
1232 The posterior probability is therefore again given by a beta distribution:
1233 \f[
1234  P(\varepsilon |k,N) \propto \varepsilon ^{k + \alpha - 1} (1 - \varepsilon)^{N - k + \beta - 1}
1235 \f]
1236 In case of central intervals
1237 the lower boundary for the equal-tailed confidence interval is given by the
1238 inverse cumulative (= quantile) function for the quantile \f$ \frac{1 - level}{2} \f$.
1239 The upper boundary for the equal-tailed confidence interval is given by the
1240 inverse cumulative (= quantile) function for the quantile \f$ \frac{1 + level}{2} \f$.
1241 Hence it is the solution \f$ \varepsilon \f$ of the following equation:
1242 \f[
1243  I_{\varepsilon}(k + \alpha,N - k + \beta) = \frac{1}{norm} \int_{0}^{\varepsilon} dt t^{k + \alpha - 1} (1 - t)^{N - k + \beta - 1} = \frac{1 \pm level}{2}
1244 \f]
1245 In the case of shortest interval the minimum interval around the mode is found by minimizing the length of all intervals width the
1246 given probability content. See TEfficiency::BetaShortestInterval
1247 */
1248 
1250 {
1251  Double_t a = double(passed)+alpha;
1252  Double_t b = double(total-passed)+beta;
1253 
1254  if (bShortest) {
1255  double lower = 0;
1256  double upper = 1;
1257  BetaShortestInterval(level,a,b,lower,upper);
1258  return (bUpper) ? upper : lower;
1259  }
1260  else
1261  return BetaCentralInterval(level, a, b, bUpper);
1262 }
1263 
1264 ////////////////////////////////////////////////////////////////////////////////
1265 /// Calculates the boundaries for a central confidence interval for a Beta distribution
1266 ///
1267 /// \param[in] level confidence level
1268 /// \param[in] a parameter > 0 for the beta distribution (for a posterior is passed + prior_alpha
1269 /// \param[in] b parameter > 0 for the beta distribution (for a posterior is (total-passed) + prior_beta
1270 /// \param[in] bUpper true - upper boundary is returned
1271 /// false - lower boundary is returned
1272 
1274 {
1275  if(bUpper) {
1276  if((a > 0) && (b > 0))
1277  return ROOT::Math::beta_quantile((1+level)/2,a,b);
1278  else {
1279  gROOT->Error("TEfficiency::BayesianCentral","Invalid input parameters - return 1");
1280  return 1;
1281  }
1282  }
1283  else {
1284  if((a > 0) && (b > 0))
1285  return ROOT::Math::beta_quantile((1-level)/2,a,b);
1286  else {
1287  gROOT->Error("TEfficiency::BayesianCentral","Invalid input parameters - return 0");
1288  return 0;
1289  }
1290  }
1291 }
1292 
1293 struct Beta_interval_length {
1294  Beta_interval_length(Double_t level,Double_t alpha,Double_t beta ) :
1295  fCL(level), fAlpha(alpha), fBeta(beta)
1296  {}
1297 
1298  Double_t LowerMax() {
1299  // max allowed value of lower given the interval size
1300  return ROOT::Math::beta_quantile_c(fCL, fAlpha,fBeta);
1301  }
1302 
1303  Double_t operator() (double lower) const {
1304  // return length of interval
1305  Double_t plow = ROOT::Math::beta_cdf(lower, fAlpha, fBeta);
1306  Double_t pup = plow + fCL;
1307  double upper = ROOT::Math::beta_quantile(pup, fAlpha,fBeta);
1308  return upper-lower;
1309  }
1310  Double_t fCL; // interval size (confidence level)
1311  Double_t fAlpha; // beta distribution alpha parameter
1312  Double_t fBeta; // beta distribution beta parameter
1313 
1314 };
1315 
1316 ////////////////////////////////////////////////////////////////////////////////
1317 /// Calculates the boundaries for a shortest confidence interval for a Beta distribution
1318 ///
1319 /// \param[in] level confidence level
1320 /// \param[in] a parameter > 0 for the beta distribution (for a posterior is passed + prior_alpha
1321 /// \param[in] b parameter > 0 for the beta distribution (for a posterior is (total-passed) + prior_beta
1322 /// \param[out] upper upper boundary is returned
1323 /// \param[out] lower lower boundary is returned
1324 ///
1325 /// The lower/upper boundary are then obtained by finding the shortest interval of the beta distribution
1326 /// contained the desired probability level.
1327 /// The length of all possible intervals is minimized in order to find the shortest one
1328 
1330 {
1331  if (a <= 0 || b <= 0) {
1332  lower = 0; upper = 1;
1333  gROOT->Error("TEfficiency::BayesianShortest","Invalid input parameters - return [0,1]");
1334  return kFALSE;
1335  }
1336 
1337  // treat here special cases when mode == 0 or 1
1338  double mode = BetaMode(a,b);
1339  if (mode == 0.0) {
1340  lower = 0;
1341  upper = ROOT::Math::beta_quantile(level, a, b);
1342  return kTRUE;
1343  }
1344  if (mode == 1.0) {
1345  lower = ROOT::Math::beta_quantile_c(level, a, b);
1346  upper = 1.0;
1347  return kTRUE;
1348  }
1349  // special case when the shortest interval is undefined return the central interval
1350  // can happen for a posterior when passed=total=0
1351  //
1352  if ( a==b && a<=1.0) {
1353  lower = BetaCentralInterval(level,a,b,kFALSE);
1354  upper = BetaCentralInterval(level,a,b,kTRUE);
1355  return kTRUE;
1356  }
1357 
1358  // for the other case perform a minimization
1359  // make a function of the length of the posterior interval as a function of lower bound
1360  Beta_interval_length intervalLength(level,a,b);
1361  // minimize the interval length
1364  minim.SetFunction(func, 0, intervalLength.LowerMax() );
1365  minim.SetNpx(2); // no need to bracket with many iterations. Just do few times to estimate some better points
1366  bool ret = minim.Minimize(100, 1.E-10,1.E-10);
1367  if (!ret) {
1368  gROOT->Error("TEfficiency::BayesianShortes","Error finding the shortest interval");
1369  return kFALSE;
1370  }
1371  lower = minim.XMinimum();
1372  upper = lower + minim.FValMinimum();
1373  return kTRUE;
1374 }
1375 
1376 ////////////////////////////////////////////////////////////////////////////////
1377 /// Compute the mean (average) of the beta distribution
1378 ///
1379 /// \param[in] a parameter > 0 for the beta distribution (for a posterior is passed + prior_alpha
1380 /// \param[in] b parameter > 0 for the beta distribution (for a posterior is (total-passed) + prior_beta
1381 ///
1382 
1384 {
1385  if (a <= 0 || b <= 0 ) {
1386  gROOT->Error("TEfficiency::BayesianMean","Invalid input parameters - return 0");
1387  return 0;
1388  }
1389 
1390  Double_t mean = a / (a + b);
1391  return mean;
1392 }
1393 
1394 ////////////////////////////////////////////////////////////////////////////////
1395 /// Compute the mode of the beta distribution
1396 ///
1397 /// \param[in] a parameter > 0 for the beta distribution (for a posterior is passed + prior_alpha
1398 /// \param[in] b parameter > 0 for the beta distribution (for a posterior is (total-passed) + prior_beta
1399 ///
1400 /// note the mode is defined for a Beta(a,b) only if (a,b)>1 (a = passed+alpha; b = total-passed+beta)
1401 /// return then the following in case (a,b) < 1:
1402 /// - if (a==b) return 0.5 (it is really undefined)
1403 /// - if (a < b) return 0;
1404 /// - if (a > b) return 1;
1405 
1407 {
1408  if (a <= 0 || b <= 0 ) {
1409  gROOT->Error("TEfficiency::BayesianMode","Invalid input parameters - return 0");
1410  return 0;
1411  }
1412  if ( a <= 1 || b <= 1) {
1413  if ( a < b) return 0;
1414  if ( a > b) return 1;
1415  if (a == b) return 0.5; // cannot do otherwise
1416  }
1417 
1418  // since a and b are > 1 here denominator cannot be 0 or < 0
1419  Double_t mode = (a - 1.0) / (a + b -2.0);
1420  return mode;
1421 }
1422 ////////////////////////////////////////////////////////////////////////////////
1423 /// Building standard data structure of a TEfficiency object
1424 ///
1425 /// Notes:
1426 /// - calls: SetName(name), SetTitle(title)
1427 /// - set the statistic option to the default (kFCP)
1428 /// - appends this object to the current directory SetDirectory(gDirectory)
1429 
1430 void TEfficiency::Build(const char* name,const char* title)
1431 {
1432  SetName(name);
1433  SetTitle(title);
1434 
1437 
1438  SetBit(kPosteriorMode,false);
1439  SetBit(kShortestInterval,false);
1440  SetBit(kUseWeights,false);
1441 
1442  //set normalisation factors to 0, otherwise the += may not work properly
1445 }
1446 
1447 ////////////////////////////////////////////////////////////////////////////////
1448 /// Checks binning for each axis
1449 ///
1450 /// It is assumed that the passed histograms have the same dimension.
1451 
1453 {
1454 
1455  const TAxis* ax1 = 0;
1456  const TAxis* ax2 = 0;
1457 
1458  //check binning along axis
1459  for(Int_t j = 0; j < pass.GetDimension(); ++j) {
1460  switch(j) {
1461  case 0:
1462  ax1 = pass.GetXaxis();
1463  ax2 = total.GetXaxis();
1464  break;
1465  case 1:
1466  ax1 = pass.GetYaxis();
1467  ax2 = total.GetYaxis();
1468  break;
1469  case 2:
1470  ax1 = pass.GetZaxis();
1471  ax2 = total.GetZaxis();
1472  break;
1473  }
1474 
1475  if(ax1->GetNbins() != ax2->GetNbins()) {
1476  gROOT->Info("TEfficiency::CheckBinning","Histograms are not consistent: they have different number of bins");
1477  return false;
1478  }
1479 
1480  for(Int_t i = 1; i <= ax1->GetNbins() + 1; ++i)
1481  if(!TMath::AreEqualRel(ax1->GetBinLowEdge(i), ax2->GetBinLowEdge(i), 1.E-15)) {
1482  gROOT->Info("TEfficiency::CheckBinning","Histograms are not consistent: they have different bin edges");
1483  return false;
1484  }
1485 
1486 
1487  }
1488 
1489  return true;
1490 }
1491 
1492 ////////////////////////////////////////////////////////////////////////////////
1493 /// Checks the consistence of the given histograms
1494 ///
1495 /// The histograms are considered as consistent if:
1496 /// - both have the same dimension
1497 /// - both have the same binning
1498 /// - pass.GetBinContent(i) <= total.GetBinContent(i) for each bin i
1499 ///
1500 
1502 {
1503  if(pass.GetDimension() != total.GetDimension()) {
1504  gROOT->Error("TEfficiency::CheckConsistency","passed TEfficiency objects have different dimensions");
1505  return false;
1506  }
1507 
1508  if(!CheckBinning(pass,total)) {
1509  gROOT->Error("TEfficiency::CheckConsistency","passed TEfficiency objects have different binning");
1510  return false;
1511  }
1512 
1513  if(!CheckEntries(pass,total)) {
1514  gROOT->Error("TEfficiency::CheckConsistency","passed TEfficiency objects do not have consistent bin contents");
1515  return false;
1516  }
1517 
1518  return true;
1519 }
1520 
1521 ////////////////////////////////////////////////////////////////////////////////
1522 /// Checks whether bin contents are compatible with binomial statistics
1523 ///
1524 /// The following inequality has to be valid for each bin i:
1525 /// total.GetBinContent(i) >= pass.GetBinContent(i)
1526 ///
1527 ///
1528 ///
1529 /// Note:
1530 ///
1531 /// - It is assumed that both histograms have the same dimension and binning.
1532 
1534 {
1535 
1536  //check: pass <= total
1537  Int_t nbinsx, nbinsy, nbinsz, nbins;
1538 
1539  nbinsx = pass.GetNbinsX();
1540  nbinsy = pass.GetNbinsY();
1541  nbinsz = pass.GetNbinsZ();
1542 
1543  switch(pass.GetDimension()) {
1544  case 1: nbins = nbinsx + 2; break;
1545  case 2: nbins = (nbinsx + 2) * (nbinsy + 2); break;
1546  case 3: nbins = (nbinsx + 2) * (nbinsy + 2) * (nbinsz + 2); break;
1547  default: nbins = 0;
1548  }
1549 
1550  for(Int_t i = 0; i < nbins; ++i) {
1551  if(pass.GetBinContent(i) > total.GetBinContent(i)) {
1552  gROOT->Info("TEfficiency::CheckEntries","Histograms are not consistent: passed bin content > total bin content");
1553  return false;
1554  }
1555  }
1556 
1557  return true;
1558 }
1559 
1560 ////////////////////////////////////////////////////////////////////////////////
1561 /// Check if both histogram are weighted. If they are weighted a true is returned
1562 ///
1564 {
1565  if (pass.GetSumw2N() == 0 && total.GetSumw2N() == 0) return false;
1566 
1567  // check also that the total sum of weight and weight squares are consistent
1568  Double_t statpass[TH1::kNstat];
1569  Double_t stattotal[TH1::kNstat];
1570 
1571  pass.GetStats(statpass);
1572  total.GetStats(stattotal);
1573 
1574  double tolerance = (total.IsA() == TH1F::Class() ) ? 1.E-5 : 1.E-12;
1575 
1576  //require: sum of weights == sum of weights^2
1577  if(!TMath::AreEqualRel(statpass[0],statpass[1],tolerance) ||
1578  !TMath::AreEqualRel(stattotal[0],stattotal[1],tolerance) ) {
1579  return true;
1580  }
1581 
1582  // histograms are not weighted
1583  return false;
1584 
1585 }
1586 
1587 
1588 ////////////////////////////////////////////////////////////////////////////////
1589 /// Create the graph used be painted (for dim=1 TEfficiency)
1590 /// The return object is managed by the caller
1591 
1593 {
1594  if (GetDimension() != 1) {
1595  Error("CreatePaintingGraph","Call this function only for dimension == 1");
1596  return 0;
1597  }
1598 
1599 
1600  Int_t npoints = fTotalHistogram->GetNbinsX();
1601  TGraphAsymmErrors * graph = new TGraphAsymmErrors(npoints);
1602  graph->SetName("eff_graph");
1603  FillGraph(graph,opt);
1604 
1605  return graph;
1606 }
1607 
1608 
1609 ////////////////////////////////////////////////////////////////////////////////
1610 /// Fill the graph to be painted with information from TEfficiency
1611 /// Internal method called by TEfficiency::Paint or TEfficiency::CreateGraph
1612 
1614 {
1615  TString option = opt;
1616  option.ToLower();
1617 
1618  Bool_t plot0Bins = false;
1619  if (option.Contains("e0") ) plot0Bins = true;
1620 
1621  Double_t x,y,xlow,xup,ylow,yup;
1622  //point i corresponds to bin i+1 in histogram
1623  // point j is point graph index
1624  // LM: cannot use TGraph::SetPoint because it deletes the underlying
1625  // histogram each time (see TGraph::SetPoint)
1626  // so use it only when extra points are added to the graph
1627  Int_t j = 0;
1628  double * px = graph->GetX();
1629  double * py = graph->GetY();
1630  double * exl = graph->GetEXlow();
1631  double * exh = graph->GetEXhigh();
1632  double * eyl = graph->GetEYlow();
1633  double * eyh = graph->GetEYhigh();
1634  Int_t npoints = fTotalHistogram->GetNbinsX();
1635  for (Int_t i = 0; i < npoints; ++i) {
1636  if (!plot0Bins && fTotalHistogram->GetBinContent(i+1) == 0 ) continue;
1637  x = fTotalHistogram->GetBinCenter(i+1);
1638  y = GetEfficiency(i+1);
1640  xup = fTotalHistogram->GetBinWidth(i+1) - xlow;
1641  ylow = GetEfficiencyErrorLow(i+1);
1642  yup = GetEfficiencyErrorUp(i+1);
1643  // in the case the graph already existed and extra points have been added
1644  if (j >= graph->GetN() ) {
1645  graph->SetPoint(j,x,y);
1646  graph->SetPointError(j,xlow,xup,ylow,yup);
1647  }
1648  else {
1649  px[j] = x;
1650  py[j] = y;
1651  exl[j] = xlow;
1652  exh[j] = xup;
1653  eyl[j] = ylow;
1654  eyh[j] = yup;
1655  }
1656  j++;
1657  }
1658 
1659  // tell the graph the effective number of points
1660  graph->Set(j);
1661  //refresh title before painting if changed
1662  TString oldTitle = graph->GetTitle();
1663  TString newTitle = GetTitle();
1664  if (oldTitle != newTitle ) {
1665  graph->SetTitle(newTitle);
1666  }
1667 
1668  // set the axis labels
1669  TString xlabel = fTotalHistogram->GetXaxis()->GetTitle();
1670  TString ylabel = fTotalHistogram->GetYaxis()->GetTitle();
1671  if (xlabel) graph->GetXaxis()->SetTitle(xlabel);
1672  if (ylabel) graph->GetYaxis()->SetTitle(ylabel);
1673 
1674  //copying style information
1675  TAttLine::Copy(*graph);
1676  TAttFill::Copy(*graph);
1677  TAttMarker::Copy(*graph);
1678 
1679  // this method forces the graph to compute correctly the axis
1680  // according to the given points
1681  graph->GetHistogram();
1682 
1683 }
1684 
1685 ////////////////////////////////////////////////////////////////////////////////
1686 /// Create the histogram used to be painted (for dim=2 TEfficiency)
1687 /// The return object is managed by the caller
1688 
1690 {
1691  if (GetDimension() != 2) {
1692  Error("CreatePaintingistogram","Call this function only for dimension == 2");
1693  return 0;
1694  }
1695 
1696  Int_t nbinsx = fTotalHistogram->GetNbinsX();
1697  Int_t nbinsy = fTotalHistogram->GetNbinsY();
1698  TAxis * xaxis = fTotalHistogram->GetXaxis();
1699  TAxis * yaxis = fTotalHistogram->GetYaxis();
1700  TH2 * hist = 0;
1701 
1702  if (xaxis->IsVariableBinSize() && yaxis->IsVariableBinSize() )
1703  hist = new TH2F("eff_histo",GetTitle(),nbinsx,xaxis->GetXbins()->GetArray(),
1704  nbinsy,yaxis->GetXbins()->GetArray());
1705  else if (xaxis->IsVariableBinSize() && ! yaxis->IsVariableBinSize() )
1706  hist = new TH2F("eff_histo",GetTitle(),nbinsx,xaxis->GetXbins()->GetArray(),
1707  nbinsy,yaxis->GetXmin(), yaxis->GetXmax());
1708  else if (!xaxis->IsVariableBinSize() && yaxis->IsVariableBinSize() )
1709  hist = new TH2F("eff_histo",GetTitle(),nbinsx,xaxis->GetXmin(), xaxis->GetXmax(),
1710  nbinsy,yaxis->GetXbins()->GetArray());
1711  else
1712  hist = new TH2F("eff_histo",GetTitle(),nbinsx,xaxis->GetXmin(), xaxis->GetXmax(),
1713  nbinsy,yaxis->GetXmin(), yaxis->GetXmax());
1714 
1715 
1716  hist->SetDirectory(0);
1717 
1718  FillHistogram(hist);
1719 
1720  return hist;
1721 }
1722 
1723 ////////////////////////////////////////////////////////////////////////////////
1724 /// Fill the 2d histogram to be painted with information from TEfficiency 2D
1725 /// Internal method called by TEfficiency::Paint or TEfficiency::CreatePaintingGraph
1726 
1727 void TEfficiency::FillHistogram(TH2 * hist ) const
1728 {
1729  //refresh title before each painting
1730  hist->SetTitle(GetTitle());
1731 
1732  // set the axis labels
1733  TString xlabel = fTotalHistogram->GetXaxis()->GetTitle();
1734  TString ylabel = fTotalHistogram->GetYaxis()->GetTitle();
1735  if (xlabel) hist->GetXaxis()->SetTitle(xlabel);
1736  if (ylabel) hist->GetYaxis()->SetTitle(ylabel);
1737 
1738  Int_t bin;
1739  Int_t nbinsx = hist->GetNbinsX();
1740  Int_t nbinsy = hist->GetNbinsY();
1741  for(Int_t i = 0; i < nbinsx + 2; ++i) {
1742  for(Int_t j = 0; j < nbinsy + 2; ++j) {
1743  bin = GetGlobalBin(i,j);
1744  hist->SetBinContent(bin,GetEfficiency(bin));
1745  }
1746  }
1747 
1748  //copying style information
1749  TAttLine::Copy(*hist);
1750  TAttFill::Copy(*hist);
1751  TAttMarker::Copy(*hist);
1752  hist->SetStats(0);
1753 
1754  return;
1755 
1756 }
1757 ////////////////////////////////////////////////////////////////////////////////
1758 /**
1759 Calculates the boundaries for the frequentist Clopper-Pearson interval
1760 
1761 This interval is recommended by the PDG.
1762 
1763 \param[in] total number of total events
1764 \param[in] passed 0 <= number of passed events <= total
1765 \param[in] level confidence level
1766 \param[in] bUpper true - upper boundary is returned
1767  ;false - lower boundary is returned
1768 
1769 Calculation:
1770 
1771 The lower boundary of the Clopper-Pearson interval is the "exact" inversion
1772 of the test:
1773  \f{eqnarray*}{
1774  P(x \geq passed; total) &=& \frac{1 - level}{2}\\
1775  P(x \geq passed; total) &=& 1 - P(x \leq passed - 1; total)\\
1776  &=& 1 - \frac{1}{norm} * \int_{0}^{1 - \varepsilon} t^{total - passed} (1 - t)^{passed - 1} dt\\
1777  &=& 1 - \frac{1}{norm} * \int_{\varepsilon}^{1} t^{passed - 1} (1 - t)^{total - passed} dt\\
1778  &=& \frac{1}{norm} * \int_{0}^{\varepsilon} t^{passed - 1} (1 - t)^{total - passed} dt\\
1779  &=& I_{\varepsilon}(passed,total - passed + 1)
1780  \f}
1781 The lower boundary is therefore given by the \f$ \frac{1 - level}{2}\f$ quantile
1782 of the beta distribution.
1783 
1784 The upper boundary of the Clopper-Pearson interval is the "exact" inversion
1785 of the test:
1786  \f{eqnarray*}{
1787  P(x \leq passed; total) &=& \frac{1 - level}{2}\\
1788  P(x \leq passed; total) &=& \frac{1}{norm} * \int_{0}^{1 - \varepsilon} t^{total - passed - 1} (1 - t)^{passed} dt\\
1789  &=& \frac{1}{norm} * \int_{\varepsilon}^{1} t^{passed} (1 - t)^{total - passed - 1} dt\\
1790  &=& 1 - \frac{1}{norm} * \int_{0}^{\varepsilon} t^{passed} (1 - t)^{total - passed - 1} dt\\
1791  \Rightarrow 1 - \frac{1 - level}{2} &=& \frac{1}{norm} * \int_{0}^{\varepsilon} t^{passed} (1 - t)^{total - passed -1} dt\\
1792  \frac{1 + level}{2} &=& I_{\varepsilon}(passed + 1,total - passed)
1793  \f}
1794 The upper boundary is therefore given by the \f$\frac{1 + level}{2}\f$ quantile
1795 of the beta distribution.
1796 
1797 Note: The connection between the binomial distribution and the regularized
1798  incomplete beta function \f$ I_{\varepsilon}(\alpha,\beta)\f$ has been used.
1799 */
1800 
1802 {
1803  Double_t alpha = (1.0 - level) / 2;
1804  if(bUpper)
1805  return ((passed == total) ? 1.0 : ROOT::Math::beta_quantile(1 - alpha,passed + 1,total-passed));
1806  else
1807  return ((passed == 0) ? 0.0 : ROOT::Math::beta_quantile(alpha,passed,total-passed+1.0));
1808 }
1809 ////////////////////////////////////////////////////////////////////////////////
1810 /**
1811  Calculates the combined efficiency and its uncertainties
1812 
1813  This method does a bayesian combination of the given samples.
1814 
1815  \param[in] up contains the upper limit of the confidence interval afterwards
1816  \param[in] low contains the lower limit of the confidence interval afterwards
1817  \param[in] n number of samples which are combined
1818  \param[in] pass array of length n containing the number of passed events
1819  \param[in] total array of length n containing the corresponding numbers of total events
1820  \param[in] alpha shape parameters for the beta distribution as prior
1821  \param[in] beta shape parameters for the beta distribution as prior
1822  \param[in] level desired confidence level
1823  \param[in] w weights for each sample; if not given, all samples get the weight 1
1824  The weights do not need to be normalized, since they are internally renormalized
1825  to the number of effective entries.
1826  \param[in] opt
1827  - mode : The mode is returned instead of the mean of the posterior as best value
1828  When using the mode the shortest interval is also computed instead of the central one
1829  - shortest: compute shortest interval (done by default if mode option is set)
1830  - central: compute central interval (done by default if mode option is NOT set)
1831 
1832  Calculation:
1833 
1834  The combined posterior distributions is calculated from the Bayes theorem assuming a common prior Beta distribution.
1835  It is easy to proof that the combined posterior is then:
1836  \f{eqnarray*}{
1837  P_{comb}(\epsilon |{w_{i}}; {k_{i}}; {N_{i}}) &=& B(\epsilon, \sum_{i}{ w_{i} k_{i}} + \alpha, \sum_{i}{ w_{i}(n_{i}-k_{i})}+\beta)\\
1838  w_{i} &=& weight\ for\ each\ sample\ renormalized\ to\ the\ effective\ entries\\
1839  w^{'}_{i} &=& w_{i} \frac{ \sum_{i} {w_{i} } } { \sum_{i} {w_{i}^{2} } }
1840  \f}
1841 
1842  The estimated efficiency is the mode (or the mean) of the obtained posterior distribution
1843 
1844  The boundaries of the confidence interval for a confidence level (1 - a)
1845  are given by the a/2 and 1-a/2 quantiles of the resulting cumulative
1846  distribution.
1847 
1848  Example (uniform prior distribution):
1849 
1850 Begin_Macro(source)
1851 {
1852  TCanvas* c1 = new TCanvas("c1","",600,800);
1853  c1->Divide(1,2);
1854  c1->SetFillStyle(1001);
1855  c1->SetFillColor(kWhite);
1856 
1857  TF1* p1 = new TF1("p1","TMath::BetaDist(x,19,9)",0,1);
1858  TF1* p2 = new TF1("p2","TMath::BetaDist(x,4,8)",0,1);
1859  TF1* comb = new TF1("comb2","TMath::BetaDist(x,[0],[1])",0,1);
1860  double nrm = 1./(0.6*0.6+0.4*0.4); // weight normalization
1861  double a = 0.6*18.0 + 0.4*3.0 + 1.0; // new alpha parameter of combined beta dist.
1862  double b = 0.6*10+0.4*7+1.0; // new beta parameter of combined beta dist.
1863  comb->SetParameters(nrm*a ,nrm *b );
1864  TF1* const1 = new TF1("const1","0.05",0,1);
1865  TF1* const2 = new TF1("const2","0.95",0,1);
1866 
1867  p1->SetLineColor(kRed);
1868  p1->SetTitle("combined posteriors;#epsilon;P(#epsilon|k,N)");
1869  p2->SetLineColor(kBlue);
1870  comb->SetLineColor(kGreen+2);
1871 
1872  TLegend* leg1 = new TLegend(0.12,0.65,0.5,0.85);
1873  leg1->AddEntry(p1,"k1 = 18, N1 = 26","l");
1874  leg1->AddEntry(p2,"k2 = 3, N2 = 10","l");
1875  leg1->AddEntry(comb,"combined: p1 = 0.6, p2=0.4","l");
1876 
1877  c1->cd(1);
1878  comb->Draw();
1879  p1->Draw("same");
1880  p2->Draw("same");
1881  leg1->Draw("same");
1882  c1->cd(2);
1883  const1->SetLineWidth(1);
1884  const2->SetLineWidth(1);
1885  TGraph* gr = (TGraph*)comb->DrawIntegral();
1886  gr->SetTitle("cumulative function of combined posterior with boundaries for cl = 95%;#epsilon;CDF");
1887  const1->Draw("same");
1888  const2->Draw("same");
1889 
1890  c1->cd(0);
1891  return c1;
1892 }
1893 End_Macro
1894 
1895 **/
1896 ////////////////////////////////////////////////////////////////////
1898  const Int_t* pass,const Int_t* total,
1899  Double_t alpha, Double_t beta,
1900  Double_t level,const Double_t* w,Option_t* opt)
1901 {
1902  TString option(opt);
1903  option.ToLower();
1904 
1905  //LM: new formula for combination
1906  // works only if alpha beta are the same always
1907  // the weights are normalized to w(i) -> N_eff w(i)/ Sum w(i)
1908  // i.e. w(i) -> Sum (w(i) / Sum (w(i)^2) * w(i)
1909  // norm = Sum (w(i) / Sum (w(i)^2)
1910  double ntot = 0;
1911  double ktot = 0;
1912  double sumw = 0;
1913  double sumw2 = 0;
1914  for (int i = 0; i < n ; ++i) {
1915  if(pass[i] > total[i]) {
1916  ::Error("TEfficiency::Combine","total events = %i < passed events %i",total[i],pass[i]);
1917  ::Info("TEfficiency::Combine","stop combining");
1918  return -1;
1919  }
1920 
1921  ntot += w[i] * total[i];
1922  ktot += w[i] * pass[i];
1923  sumw += w[i];
1924  sumw2 += w[i]*w[i];
1925  //mean += w[i] * (pass[i] + alpha[i])/(total[i] + alpha[i] + beta[i]);
1926  }
1927  double norm = sumw/sumw2;
1928  ntot *= norm;
1929  ktot *= norm;
1930  if(ktot > ntot) {
1931  ::Error("TEfficiency::Combine","total = %f < passed %f",ntot,ktot);
1932  ::Info("TEfficiency::Combine","stop combining");
1933  return -1;
1934  }
1935 
1936  double a = ktot + alpha;
1937  double b = ntot - ktot + beta;
1938 
1939  double mean = a/(a+b);
1940  double mode = BetaMode(a,b);
1941 
1942 
1943  Bool_t shortestInterval = option.Contains("sh") || ( option.Contains("mode") && !option.Contains("cent") );
1944 
1945  if (shortestInterval)
1946  BetaShortestInterval(level, a, b, low, up);
1947  else {
1948  low = BetaCentralInterval(level, a, b, false);
1949  up = BetaCentralInterval(level, a, b, true);
1950  }
1951 
1952  if (option.Contains("mode")) return mode;
1953  return mean;
1954 
1955 }
1956 ////////////////////////////////////////////////////////////////////////////////
1957 /// Combines a list of 1-dimensional TEfficiency objects
1958 ///
1959 /// A TGraphAsymmErrors object is returned which contains the estimated
1960 /// efficiency and its uncertainty for each bin.
1961 /// If the combination fails, a zero pointer is returned.
1962 ///
1963 /// At the moment the combining is only implemented for bayesian statistics.
1964 ///
1965 /// \param[in] pList list containing TEfficiency objects which should be combined
1966 /// only one-dimensional efficiencies are taken into account
1967 /// \param[in] option
1968 /// - s : strict combining; only TEfficiency objects with the same beta
1969 /// prior and the flag kIsBayesian == true are combined
1970 /// If not specified the prior parameter of the first TEfficiency object is used
1971 /// - v : verbose mode; print information about combining
1972 /// - cl=x : set confidence level (0 < cl < 1). If not specified, the
1973 /// confidence level of the first TEfficiency object is used.
1974 /// - mode Use mode of combined posterior as estimated value for the efficiency
1975 /// - shortest: compute shortest interval (done by default if mode option is set)
1976 /// - central: compute central interval (done by default if mode option is NOT set)
1977 /// \param[in] n number of weights (has to be the number of one-dimensional
1978 /// TEfficiency objects in pList)
1979 /// If no weights are passed, the internal weights GetWeight() of
1980 /// the given TEfficiency objects are used.
1981 /// \param[in] w array of length n with weights for each TEfficiency object in
1982 /// pList (w[0] correspond to pList->First ... w[n-1] -> pList->Last)
1983 /// The weights do not have to be normalised.
1984 ///
1985 /// For each bin the calculation is done by the Combine(double&, double& ...) method.
1986 
1988  Int_t n,const Double_t* w)
1989 {
1990  TString opt = option;
1991  opt.ToLower();
1992 
1993  //parameter of prior distribution, confidence level and normalisation factor
1994  Double_t alpha = -1;
1995  Double_t beta = -1;
1996  Double_t level = 0;
1997 
1998  //flags for combining
1999  Bool_t bStrict = false;
2000  Bool_t bOutput = false;
2001  Bool_t bWeights = false;
2002  //list of all information needed to weight and combine efficiencies
2003  std::vector<TH1*> vTotal; vTotal.reserve(n);
2004  std::vector<TH1*> vPassed; vPassed.reserve(n);
2005  std::vector<Double_t> vWeights; vWeights.reserve(n);
2006  // std::vector<Double_t> vAlpha;
2007  // std::vector<Double_t> vBeta;
2008 
2009  if(opt.Contains("s")) {
2010  opt.ReplaceAll("s","");
2011  bStrict = true;
2012  }
2013 
2014  if(opt.Contains("v")) {
2015  opt.ReplaceAll("v","");
2016  bOutput = true;
2017  }
2018 
2019  if(opt.Contains("cl=")) {
2020  Ssiz_t pos = opt.Index("cl=") + 3;
2021  level = atof( opt(pos,opt.Length() ).Data() );
2022  if((level <= 0) || (level >= 1))
2023  level = 0;
2024  opt.ReplaceAll("cl=","");
2025  }
2026 
2027  //are weights explicitly given
2028  if(n && w) {
2029  bWeights = true;
2030  for(Int_t k = 0; k < n; ++k) {
2031  if(w[k] > 0)
2032  vWeights.push_back(w[k]);
2033  else {
2034  gROOT->Error("TEfficiency::Combine","invalid custom weight found w = %.2lf",w[k]);
2035  gROOT->Info("TEfficiency::Combine","stop combining");
2036  return 0;
2037  }
2038  }
2039  }
2040 
2041  TIter next(pList);
2042  TObject* obj = 0;
2043  TEfficiency* pEff = 0;
2044  while((obj = next())) {
2045  pEff = dynamic_cast<TEfficiency*>(obj);
2046  //is object a TEfficiency object?
2047  if(pEff) {
2048  if(pEff->GetDimension() > 1)
2049  continue;
2050  if(!level) level = pEff->GetConfidenceLevel();
2051 
2052  if(alpha<1) alpha = pEff->GetBetaAlpha();
2053  if(beta<1) beta = pEff->GetBetaBeta();
2054 
2055  //if strict combining, check priors, confidence level and statistic
2056  if(bStrict) {
2057  if(alpha != pEff->GetBetaAlpha())
2058  continue;
2059  if(beta != pEff->GetBetaBeta())
2060  continue;
2061  if(!pEff->UsesBayesianStat())
2062  continue;
2063  }
2064 
2065  vTotal.push_back(pEff->fTotalHistogram);
2066  vPassed.push_back(pEff->fPassedHistogram);
2067 
2068  //no weights given -> use weights of TEfficiency objects
2069  if(!bWeights)
2070  vWeights.push_back(pEff->fWeight);
2071 
2072  //strict combining -> using global prior
2073  // if(bStrict) {
2074  // vAlpha.push_back(alpha);
2075  // vBeta.push_back(beta);
2076  // }
2077  // else {
2078  // vAlpha.push_back(pEff->GetBetaAlpha());
2079  // vBeta.push_back(pEff->GetBetaBeta());
2080  // }
2081  }
2082  }
2083 
2084  //no TEfficiency objects found
2085  if(vTotal.empty()) {
2086  gROOT->Error("TEfficiency::Combine","no TEfficiency objects in given list");
2087  gROOT->Info("TEfficiency::Combine","stop combining");
2088  return 0;
2089  }
2090 
2091  //invalid number of custom weights
2092  if(bWeights && (n != (Int_t)vTotal.size())) {
2093  gROOT->Error("TEfficiency::Combine","number of weights n=%i differs from number of TEfficiency objects k=%i which should be combined",n,(Int_t)vTotal.size());
2094  gROOT->Info("TEfficiency::Combine","stop combining");
2095  return 0;
2096  }
2097 
2098  Int_t nbins_max = vTotal.at(0)->GetNbinsX();
2099  //check binning of all histograms
2100  for(UInt_t i=0; i<vTotal.size(); ++i) {
2101  if (!TEfficiency::CheckBinning(*vTotal.at(0),*vTotal.at(i)) )
2102  gROOT->Warning("TEfficiency::Combine","histograms have not the same binning -> results may be useless");
2103  if(vTotal.at(i)->GetNbinsX() < nbins_max) nbins_max = vTotal.at(i)->GetNbinsX();
2104  }
2105 
2106  //display information about combining
2107  if(bOutput) {
2108  gROOT->Info("TEfficiency::Combine","combining %i TEfficiency objects",(Int_t)vTotal.size());
2109  if(bWeights)
2110  gROOT->Info("TEfficiency::Combine","using custom weights");
2111  if(bStrict) {
2112  gROOT->Info("TEfficiency::Combine","using the following prior probability for the efficiency: P(e) ~ Beta(e,%.3lf,%.3lf)",alpha,beta);
2113  }
2114  else
2115  gROOT->Info("TEfficiency::Combine","using individual priors of each TEfficiency object");
2116  gROOT->Info("TEfficiency::Combine","confidence level = %.2lf",level);
2117  }
2118 
2119  //create TGraphAsymmErrors with efficiency
2120  std::vector<Double_t> x(nbins_max);
2121  std::vector<Double_t> xlow(nbins_max);
2122  std::vector<Double_t> xhigh(nbins_max);
2123  std::vector<Double_t> eff(nbins_max);
2124  std::vector<Double_t> efflow(nbins_max);
2125  std::vector<Double_t> effhigh(nbins_max);
2126 
2127  //parameters for combining:
2128  //number of objects
2129  Int_t num = vTotal.size();
2130  std::vector<Int_t> pass(num);
2131  std::vector<Int_t> total(num);
2132 
2133  //loop over all bins
2134  Double_t low = 0;
2135  Double_t up = 0;
2136  for(Int_t i=1; i <= nbins_max; ++i) {
2137  //the binning of the x-axis is taken from the first total histogram
2138  x[i-1] = vTotal.at(0)->GetBinCenter(i);
2139  xlow[i-1] = x[i-1] - vTotal.at(0)->GetBinLowEdge(i);
2140  xhigh[i-1] = vTotal.at(0)->GetBinWidth(i) - xlow[i-1];
2141 
2142  for(Int_t j = 0; j < num; ++j) {
2143  pass[j] = (Int_t)(vPassed.at(j)->GetBinContent(i) + 0.5);
2144  total[j] = (Int_t)(vTotal.at(j)->GetBinContent(i) + 0.5);
2145  }
2146 
2147  //fill efficiency and errors
2148  eff[i-1] = Combine(up,low,num,&pass[0],&total[0],alpha,beta,level,&vWeights[0],opt.Data());
2149  //did an error occurred ?
2150  if(eff[i-1] == -1) {
2151  gROOT->Error("TEfficiency::Combine","error occurred during combining");
2152  gROOT->Info("TEfficiency::Combine","stop combining");
2153  return 0;
2154  }
2155  efflow[i-1]= eff[i-1] - low;
2156  effhigh[i-1]= up - eff[i-1];
2157  }//loop over all bins
2158 
2159  TGraphAsymmErrors* gr = new TGraphAsymmErrors(nbins_max,&x[0],&eff[0],&xlow[0],&xhigh[0],&efflow[0],&effhigh[0]);
2160 
2161  return gr;
2162 }
2163 
2164 ////////////////////////////////////////////////////////////////////////////////
2165 /// Compute distance from point px,py to a graph.
2166 ///
2167 /// Compute the closest distance of approach from point px,py to this line.
2168 /// The distance is computed in pixels units.
2169 ///
2170 /// Forward the call to the painted graph
2171 
2173 {
2174  if (fPaintGraph) return fPaintGraph->DistancetoPrimitive(px,py);
2175  if (fPaintHisto) return fPaintHisto->DistancetoPrimitive(px,py);
2176  return 0;
2177 }
2178 
2179 
2180 ////////////////////////////////////////////////////////////////////////////////
2181 /// Draws the current TEfficiency object
2182 ///
2183 /// \param[in] opt
2184 /// - 1-dimensional case: same options as TGraphAsymmErrors::Draw()
2185 /// but as default "AP" is used
2186 /// - 2-dimensional case: same options as TH2::Draw()
2187 /// - 3-dimensional case: not yet supported
2188 ///
2189 /// Specific TEfficiency drawing options:
2190 /// - E0 - plot bins where the total number of passed events is zero
2191 /// (the error interval will be [0,1] )
2192 
2194 {
2195  //check options
2196  TString option = opt;
2197  option.ToLower();
2198  // use by default "AP"
2199  if (option.IsNull() ) option = "ap";
2200 
2201  if(gPad && !option.Contains("same"))
2202  gPad->Clear();
2203  else {
2204  // add always "a" if not present
2205  if (!option.Contains("a") ) option += "a";
2206  }
2207 
2208  // add always p to the option
2209  if (!option.Contains("p") ) option += "p";
2210 
2211 
2212  AppendPad(option.Data());
2213 }
2214 
2215 ////////////////////////////////////////////////////////////////////////////////
2216 /// Execute action corresponding to one event.
2217 ///
2218 /// This member function is called when the drawn class is clicked with the locator
2219 /// If Left button clicked on one of the line end points, this point
2220 /// follows the cursor until button is released.
2221 ///
2222 /// if Middle button clicked, the line is moved parallel to itself
2223 /// until the button is released.
2224 /// Forward the call to the underlying graph
2225 
2227 {
2228  if (fPaintGraph) fPaintGraph->ExecuteEvent(event,px,py);
2229  else if (fPaintHisto) fPaintHisto->ExecuteEvent(event,px,py);
2230 }
2231 
2232 ////////////////////////////////////////////////////////////////////////////////
2233 /// This function is used for filling the two histograms.
2234 ///
2235 /// \param[in] bPassed flag whether the current event passed the selection
2236 /// - true: both histograms are filled
2237 /// - false: only the total histogram is filled
2238 /// \param[in] x x-value
2239 /// \param[in] y y-value (use default=0 for 1-D efficiencies)
2240 /// \param[in] z z-value (use default=0 for 2-D or 1-D efficiencies)
2241 
2243 {
2244  switch(GetDimension()) {
2245  case 1:
2246  fTotalHistogram->Fill(x);
2247  if(bPassed)
2248  fPassedHistogram->Fill(x);
2249  break;
2250  case 2:
2251  ((TH2*)(fTotalHistogram))->Fill(x,y);
2252  if(bPassed)
2253  ((TH2*)(fPassedHistogram))->Fill(x,y);
2254  break;
2255  case 3:
2256  ((TH3*)(fTotalHistogram))->Fill(x,y,z);
2257  if(bPassed)
2258  ((TH3*)(fPassedHistogram))->Fill(x,y,z);
2259  break;
2260  }
2261 }
2262 
2263 ////////////////////////////////////////////////////////////////////////////////
2264 ///This function is used for filling the two histograms with a weight.
2265 ///
2266 /// \param[in] bPassed flag whether the current event passed the selection
2267 /// - true: both histograms are filled
2268 /// - false: only the total histogram is filled
2269 /// \param[in] weight weight for the event
2270 /// \param[in] x x-value
2271 /// \param[in] y y-value (use default=0 for 1-D efficiencies)
2272 /// \param[in] z z-value (use default=0 for 2-D or 1-D efficiencies)
2273 ///
2274 /// Note: - this function will call SetUseWeightedEvents if it was not called by the user before
2275 
2277 {
2278  if(!TestBit(kUseWeights))
2279  {
2280  Info("FillWeighted","call SetUseWeightedEvents() manually to ensure correct storage of sum of weights squared");
2282  }
2283 
2284  switch(GetDimension()) {
2285  case 1:
2286  fTotalHistogram->Fill(x,weight);
2287  if(bPassed)
2288  fPassedHistogram->Fill(x,weight);
2289  break;
2290  case 2:
2291  ((TH2*)(fTotalHistogram))->Fill(x,y,weight);
2292  if(bPassed)
2293  ((TH2*)(fPassedHistogram))->Fill(x,y,weight);
2294  break;
2295  case 3:
2296  ((TH3*)(fTotalHistogram))->Fill(x,y,z,weight);
2297  if(bPassed)
2298  ((TH3*)(fPassedHistogram))->Fill(x,y,z,weight);
2299  break;
2300  }
2301 }
2302 
2303 ////////////////////////////////////////////////////////////////////////////////
2304 /// Returns the global bin number containing the given values
2305 ///
2306 /// Note:
2307 ///
2308 /// - values which belong to dimensions higher than the current dimension
2309 /// of the TEfficiency object are ignored (i.e. for 1-dimensional
2310 /// efficiencies only the x-value is considered)
2311 
2313 {
2315  Int_t ny = 0;
2316  Int_t nz = 0;
2317 
2318  switch(GetDimension()) {
2319  case 3: nz = fTotalHistogram->GetZaxis()->FindFixBin(z);
2320  case 2: ny = fTotalHistogram->GetYaxis()->FindFixBin(y);break;
2321  }
2322 
2323  return GetGlobalBin(nx,ny,nz);
2324 }
2325 
2326 ////////////////////////////////////////////////////////////////////////////////
2327 /// Fits the efficiency using the TBinomialEfficiencyFitter class
2328 ///
2329 /// The resulting fit function is added to the list of associated functions.
2330 ///
2331 /// Options:
2332 /// - "+": previous fitted functions in the list are kept, by default
2333 /// all functions in the list are deleted
2334 /// - for more fitting options see TBinomialEfficiencyFitter::Fit
2335 
2337 {
2338  TString option = opt;
2339  option.ToLower();
2340 
2341  //replace existing functions in list with same name
2342  Bool_t bDeleteOld = true;
2343  if(option.Contains("+")) {
2344  option.ReplaceAll("+","");
2345  bDeleteOld = false;
2346  }
2347 
2349 
2350  TFitResultPtr result = Fitter.Fit(f1,option.Data());
2351 
2352  //create copy which is appended to the list
2353  TF1* pFunc = new TF1(*f1);
2354 
2355  if(bDeleteOld) {
2356  TIter next(fFunctions);
2357  TObject* obj = 0;
2358  while((obj = next())) {
2359  if(obj->InheritsFrom(TF1::Class())) {
2360  fFunctions->Remove(obj);
2361  delete obj;
2362  }
2363  }
2364  }
2365 
2366  // create list if necessary
2367  if(!fFunctions)
2368  fFunctions = new TList();
2369 
2370  fFunctions->Add(pFunc);
2371 
2372  return result;
2373 }
2374 
2375 ////////////////////////////////////////////////////////////////////////////////
2376 /// Returns a cloned version of fPassedHistogram
2377 ///
2378 /// Notes:
2379 /// - The histogram is filled with unit weights. You might want to scale
2380 /// it with the global weight GetWeight().
2381 /// - The returned object is owned by the user who has to care about the
2382 /// deletion of the new TH1 object.
2383 /// - This histogram is by default NOT attached to the current directory
2384 /// to avoid duplication of data. If you want to store it automatically
2385 /// during the next TFile::Write() command, you have to attach it to
2386 /// the corresponding directory.
2387 ///
2388 /// ~~~~~~~{.cpp}
2389 /// TFile* pFile = new TFile("passed.root","update");
2390 /// TEfficiency* pEff = (TEfficiency*)gDirectory->Get("my_eff");
2391 /// TH1* copy = pEff->GetCopyPassedHisto();
2392 /// copy->SetDirectory(gDirectory);
2393 /// pFile->Write();
2394 /// ~~~~~~~
2395 
2397 {
2398  Bool_t bStatus = TH1::AddDirectoryStatus();
2400  TH1* tmp = (TH1*)(fPassedHistogram->Clone());
2401  TH1::AddDirectory(bStatus);
2402 
2403  return tmp;
2404 }
2405 
2406 ////////////////////////////////////////////////////////////////////////////////
2407 /// Returns a cloned version of fTotalHistogram
2408 ///
2409 /// Notes:
2410 /// - The histogram is filled with unit weights. You might want to scale
2411 /// it with the global weight GetWeight().
2412 /// - The returned object is owned by the user who has to care about the
2413 /// deletion of the new TH1 object.
2414 /// - This histogram is by default NOT attached to the current directory
2415 /// to avoid duplication of data. If you want to store it automatically
2416 /// during the next TFile::Write() command, you have to attach it to
2417 /// the corresponding directory.
2418 ///
2419 /// ~~~~~~~{.cpp}
2420 /// TFile* pFile = new TFile("total.root","update");
2421 /// TEfficiency* pEff = (TEfficiency*)gDirectory->Get("my_eff");
2422 /// TH1* copy = pEff->GetCopyTotalHisto();
2423 /// copy->SetDirectory(gDirectory);
2424 /// pFile->Write();
2425 /// ~~~~~~~
2426 
2428 {
2429  Bool_t bStatus = TH1::AddDirectoryStatus();
2431  TH1* tmp = (TH1*)(fTotalHistogram->Clone());
2432  TH1::AddDirectory(bStatus);
2433 
2434  return tmp;
2435 }
2436 
2437 ////////////////////////////////////////////////////////////////////////////////
2438 ///returns the dimension of the current TEfficiency object
2439 
2441 {
2442  return fTotalHistogram->GetDimension();
2443 }
2444 
2445 ////////////////////////////////////////////////////////////////////////////////
2446 /// Returns the efficiency in the given global bin
2447 ///
2448 /// Note:
2449 /// - The estimated efficiency depends on the chosen statistic option:
2450 /// for frequentist ones:
2451 /// \f$ \hat{\varepsilon} = \frac{passed}{total} \f$
2452 /// for bayesian ones the expectation value of the resulting posterior
2453 /// distribution is returned:
2454 /// \f$ \hat{\varepsilon} = \frac{passed + \alpha}{total + \alpha + \beta} \f$
2455 /// If the bit kPosteriorMode is set (or the method TEfficiency::UsePosteriorMode() has been called ) the
2456 /// mode (most probable value) of the posterior is returned:
2457 /// \f$ \hat{\varepsilon} = \frac{passed + \alpha -1}{total + \alpha + \beta -2} \f$
2458 /// - If the denominator is equal to 0, an efficiency of 0 is returned.
2459 /// - When \f$ passed + \alpha < 1 \f$ or \f$ total - passed + \beta < 1 \f$ the above
2460 /// formula for the mode is not valid. In these cases values the estimated efficiency is 0 or 1.
2461 
2463 {
2465  Double_t passed = fPassedHistogram->GetBinContent(bin);
2466 
2467  if(TestBit(kIsBayesian)) {
2468 
2469  // parameters for the beta prior distribution
2472 
2473  Double_t aa,bb;
2474  if(TestBit(kUseWeights))
2475  {
2477  Double_t tw2 = fTotalHistogram->GetSumw2()->At(bin);
2479 
2480  if (tw2 <= 0 ) return pw/tw;
2481 
2482  // tw/tw2 renormalize the weights
2483  double norm = tw/tw2;
2484  aa = pw * norm + alpha;
2485  bb = (tw - pw) * norm + beta;
2486  }
2487  else
2488  {
2489  aa = passed + alpha;
2490  bb = total - passed + beta;
2491  }
2492 
2493  if (!TestBit(kPosteriorMode) )
2494  return BetaMean(aa,bb);
2495  else
2496  return BetaMode(aa,bb);
2497 
2498  }
2499  else
2500  return (total)? ((Double_t)passed)/total : 0;
2501 }
2502 
2503 ////////////////////////////////////////////////////////////////////////////////
2504 /// Returns the lower error on the efficiency in the given global bin
2505 ///
2506 /// The result depends on the current confidence level fConfLevel and the
2507 /// chosen statistic option fStatisticOption. See SetStatisticOption(Int_t) for
2508 /// more details.
2509 ///
2510 /// Note: If the histograms are filled with weights, only bayesian methods and the
2511 /// normal approximation are supported.
2512 
2514 {
2516  Int_t passed = (Int_t)fPassedHistogram->GetBinContent(bin);
2517 
2518  Double_t eff = GetEfficiency(bin);
2519 
2520  // check whether weights have been used
2521  if(TestBit(kUseWeights))
2522  {
2524  Double_t tw2 = fTotalHistogram->GetSumw2()->At(bin);
2526  Double_t pw2 = fPassedHistogram->GetSumw2()->At(bin);
2527 
2528  if(TestBit(kIsBayesian))
2529  {
2532 
2533  if (tw2 <= 0) return 0;
2534 
2535  // tw/tw2 renormalize the weights
2536  Double_t norm = tw/tw2;
2537  Double_t aa = pw * norm + alpha;
2538  Double_t bb = (tw - pw) * norm + beta;
2539  Double_t low = 0;
2540  Double_t upper = 1;
2541  if(TestBit(kShortestInterval)) {
2543  }
2544  else {
2545  low = TEfficiency::BetaCentralInterval(fConfLevel,aa,bb,false);
2546  }
2547 
2548  return eff - low;
2549  }
2550  else
2551  {
2552  if(fStatisticOption != kFNormal)
2553  {
2554  Warning("GetEfficiencyErrorLow","frequentist confidence intervals for weights are only supported by the normal approximation");
2555  Info("GetEfficiencyErrorLow","setting statistic option to kFNormal");
2556  const_cast<TEfficiency*>(this)->SetStatisticOption(kFNormal);
2557  }
2558 
2559  Double_t variance = ( pw2 * (1. - 2 * eff) + tw2 * eff *eff ) / ( tw * tw) ;
2560  Double_t sigma = sqrt(variance);
2561 
2562  Double_t prob = 0.5 * (1.- fConfLevel);
2564 
2565  // avoid to return errors which makes eff-err < 0
2566  return (eff - delta < 0) ? eff : delta;
2567  }
2568  }
2569  else
2570  {
2571  if(TestBit(kIsBayesian))
2572  {
2573  // parameters for the beta prior distribution
2576  return (eff - Bayesian(total,passed,fConfLevel,alpha,beta,false,TestBit(kShortestInterval)));
2577  }
2578  else
2579  return (eff - fBoundary(total,passed,fConfLevel,false));
2580  }
2581 }
2582 
2583 ////////////////////////////////////////////////////////////////////////////////
2584 /// Returns the upper error on the efficiency in the given global bin
2585 ///
2586 /// The result depends on the current confidence level fConfLevel and the
2587 /// chosen statistic option fStatisticOption. See SetStatisticOption(Int_t) for
2588 /// more details.
2589 ///
2590 /// Note: If the histograms are filled with weights, only bayesian methods and the
2591 /// normal approximation are supported.
2592 
2594 {
2596  Int_t passed = (Int_t)fPassedHistogram->GetBinContent(bin);
2597 
2598  Double_t eff = GetEfficiency(bin);
2599 
2600  // check whether weights have been used
2601  if(TestBit(kUseWeights))
2602  {
2604  Double_t tw2 = fTotalHistogram->GetSumw2()->At(bin);
2606  Double_t pw2 = fPassedHistogram->GetSumw2()->At(bin);
2607 
2608  if(TestBit(kIsBayesian))
2609  {
2612 
2613  if (tw2 <= 0) return 0;
2614 
2615  // tw/tw2 renormalize the weights
2616  Double_t norm = tw/tw2;
2617  Double_t aa = pw * norm + alpha;
2618  Double_t bb = (tw - pw) * norm + beta;
2619  Double_t low = 0;
2620  Double_t upper = 1;
2621  if(TestBit(kShortestInterval)) {
2623  }
2624  else {
2625  upper = TEfficiency::BetaCentralInterval(fConfLevel,aa,bb,true);
2626  }
2627 
2628  return upper - eff;
2629  }
2630  else
2631  {
2632  if(fStatisticOption != kFNormal)
2633  {
2634  Warning("GetEfficiencyErrorUp","frequentist confidence intervals for weights are only supported by the normal approximation");
2635  Info("GetEfficiencyErrorUp","setting statistic option to kFNormal");
2636  const_cast<TEfficiency*>(this)->SetStatisticOption(kFNormal);
2637  }
2638 
2639  Double_t variance = ( pw2 * (1. - 2 * eff) + tw2 * eff *eff ) / ( tw * tw) ;
2640  Double_t sigma = sqrt(variance);
2641 
2642  Double_t prob = 0.5 * (1.- fConfLevel);
2644 
2645  return (eff + delta > 1) ? 1.-eff : delta;
2646  }
2647  }
2648  else
2649  {
2650  if(TestBit(kIsBayesian))
2651  {
2652  // parameters for the beta prior distribution
2655  return (Bayesian(total,passed,fConfLevel,alpha,beta,true,TestBit(kShortestInterval)) - eff);
2656  }
2657  else
2658  return fBoundary(total,passed,fConfLevel,true) - eff;
2659  }
2660 }
2661 
2662 ////////////////////////////////////////////////////////////////////////////////
2663 /// Returns the global bin number which can be used as argument for the
2664 /// following functions:
2665 ///
2666 /// - GetEfficiency(bin), GetEfficiencyErrorLow(bin), GetEfficiencyErrorUp(bin)
2667 /// - SetPassedEvents(bin), SetTotalEvents(bin)
2668 ///
2669 /// see TH1::GetBin() for conventions on numbering bins
2670 
2672 {
2673  return fTotalHistogram->GetBin(binx,biny,binz);
2674 }
2675 
2676 ////////////////////////////////////////////////////////////////////////////////
2677 
2679 {
2680  return (fFunctions) ? fFunctions : fFunctions = new TList();
2681 }
2682 
2683 ////////////////////////////////////////////////////////////////////////////////
2684 /// Merges the TEfficiency objects in the given list to the given
2685 /// TEfficiency object using the operator+=(TEfficiency&)
2686 ///
2687 /// The merged result is stored in the current object. The statistic options and
2688 /// the confidence level are taken from the current object.
2689 ///
2690 /// This function should be used when all TEfficiency objects correspond to
2691 /// the same process.
2692 ///
2693 /// The new weight is set according to:
2694 /// \f$ \frac{1}{w_{new}} = \sum_{i} \frac{1}{w_{i}} \f$
2695 
2697 {
2698  if(!pList->IsEmpty()) {
2699  TIter next(pList);
2700  TObject* obj = 0;
2701  TEfficiency* pEff = 0;
2702  while((obj = next())) {
2703  pEff = dynamic_cast<TEfficiency*>(obj);
2704  if(pEff) {
2705  *this += *pEff;
2706  }
2707  }
2708  }
2709  return (Long64_t)fTotalHistogram->GetEntries();
2710 }
2711 
2712 ////////////////////////////////////////////////////////////////////////////////
2713 /**
2714 Returns the confidence limits for the efficiency supposing that the
2715 efficiency follows a normal distribution with the rms below
2716 
2717 \param[in] total number of total events
2718 \param[in] passed 0 <= number of passed events <= total
2719 \param[in] level confidence level
2720 \param[in] bUpper
2721  - true - upper boundary is returned
2722  - false - lower boundary is returned
2723 
2724 Calculation:
2725 
2726 \f{eqnarray*}{
2727  \hat{\varepsilon} &=& \frac{passed}{total}\\
2728  \sigma_{\varepsilon} &=& \sqrt{\frac{\hat{\varepsilon} (1 - \hat{\varepsilon})}{total}}\\
2729  \varepsilon_{low} &=& \hat{\varepsilon} \pm \Phi^{-1}(\frac{level}{2},\sigma_{\varepsilon})
2730 \f}
2731 */
2732 
2734 {
2735  Double_t alpha = (1.0 - level)/2;
2736  if (total == 0) return (bUpper) ? 1 : 0;
2737  Double_t average = passed / total;
2738  Double_t sigma = std::sqrt(average * (1 - average) / total);
2739  Double_t delta = ROOT::Math::normal_quantile(1 - alpha,sigma);
2740 
2741  if(bUpper)
2742  return ((average + delta) > 1) ? 1.0 : (average + delta);
2743  else
2744  return ((average - delta) < 0) ? 0.0 : (average - delta);
2745 }
2746 
2747 ////////////////////////////////////////////////////////////////////////////////
2748 /// Adds the histograms of another TEfficiency object to current histograms
2749 ///
2750 /// The statistic options and the confidence level remain unchanged.
2751 ///
2752 /// fTotalHistogram += rhs.fTotalHistogram;
2753 /// fPassedHistogram += rhs.fPassedHistogram;
2754 ///
2755 /// calculates a new weight:
2756 /// current weight of this TEfficiency object = \f$ w_{1} \f$
2757 /// weight of rhs = \f$ w_{2} \f$
2758 /// \f$ w_{new} = \frac{w_{1} \times w_{2}}{w_{1} + w_{2}} \f$
2759 
2761 {
2762 
2763  if (fTotalHistogram == 0 && fPassedHistogram == 0) {
2764  // efficiency is empty just copy it over
2765  *this = rhs;
2766  return *this;
2767  }
2768  else if (fTotalHistogram == 0 || fPassedHistogram == 0) {
2769  Fatal("operator+=","Adding to a non consistent TEfficiency object which has not a total or a passed histogram ");
2770  return *this;
2771  }
2772 
2773  if (rhs.fTotalHistogram == 0 && rhs.fPassedHistogram == 0 ) {
2774  Warning("operator+=","no operation: adding an empty object");
2775  return *this;
2776  }
2777  else if (rhs.fTotalHistogram == 0 || rhs.fPassedHistogram == 0 ) {
2778  Fatal("operator+=","Adding a non consistent TEfficiency object which has not a total or a passed histogram ");
2779  return *this;
2780  }
2781 
2784 
2787 
2788  SetWeight((fWeight * rhs.GetWeight())/(fWeight + rhs.GetWeight()));
2789 
2790  return *this;
2791 }
2792 
2793 ////////////////////////////////////////////////////////////////////////////////
2794 /// Assignment operator
2795 ///
2796 /// The histograms, statistic option, confidence level, weight and paint styles
2797 /// of rhs are copied to the this TEfficiency object.
2798 ///
2799 /// Note: - The list of associated functions is not copied. After this
2800 /// operation the list of associated functions is empty.
2801 
2803 {
2804  if(this != &rhs)
2805  {
2806  //statistic options
2809  SetBetaAlpha(rhs.GetBetaAlpha());
2810  SetBetaBeta(rhs.GetBetaBeta());
2811  SetWeight(rhs.GetWeight());
2812 
2813  //associated list of functions
2814  if(fFunctions)
2815  fFunctions->Delete();
2816 
2817  //copy histograms
2818  delete fTotalHistogram;
2819  delete fPassedHistogram;
2820 
2821  Bool_t bStatus = TH1::AddDirectoryStatus();
2825  TH1::AddDirectory(bStatus);
2826 
2827  //delete temporary paint objects
2828  delete fPaintHisto;
2829  delete fPaintGraph;
2830  fPaintHisto = 0;
2831  fPaintGraph = 0;
2832 
2833  //copy style
2834  rhs.TAttLine::Copy(*this);
2835  rhs.TAttFill::Copy(*this);
2836  rhs.TAttMarker::Copy(*this);
2837  }
2838 
2839  return *this;
2840 }
2841 
2842 ////////////////////////////////////////////////////////////////////////////////
2843 /// Paints this TEfficiency object
2844 ///
2845 /// For details on the possible option see Draw(Option_t*)
2846 ///
2847 /// Note for 1D classes
2848 /// In 1D the TEfficiency uses a TGraphAsymmErrors for drawing
2849 /// The TGraph is created only the first time Paint is used. The user can manipulate the
2850 /// TGraph via the method TEfficiency::GetPaintedGraph()
2851 /// The TGraph creates behing an histogram for the axis. The histogram is created also only the first time.
2852 /// If the axis needs to be updated because in the meantime the class changed use this trick
2853 /// which will trigger a re-calculation of the axis of the graph
2854 /// TEfficiency::GetPaintedGraph()->Set(0)
2855 ///
2856 /// Note that in order to access the painted graph via GetPaintedGraph() you need either to call Paint or better
2857 /// gPad->Update();
2858 ///
2859 
2861 {
2862 
2863 
2864  if(!gPad)
2865  return;
2866 
2867 
2868  //use TGraphAsymmErrors for painting
2869  if(GetDimension() == 1) {
2870  if(!fPaintGraph) {
2871  fPaintGraph = CreateGraph(opt);
2872  }
2873  else
2874  // update existing graph already created
2875  FillGraph(fPaintGraph, opt);
2876 
2877  //paint graph
2878 
2879  fPaintGraph->Paint(opt);
2880 
2881  //paint all associated functions
2882  if(fFunctions) {
2883  //paint box with fit parameters
2884  //the fit statistics will be painted if gStyle->SetOptFit(1) has been
2885  // called by the user
2886  TIter next(fFunctions);
2887  TObject* obj = 0;
2888  while((obj = next())) {
2889  if(obj->InheritsFrom(TF1::Class())) {
2890  fPaintGraph->PaintStats((TF1*)obj);
2891  ((TF1*)obj)->Paint("sameC");
2892  }
2893  }
2894  }
2895 
2896  return;
2897  }
2898 
2899  //use TH2 for painting
2900  if(GetDimension() == 2) {
2901  if(!fPaintHisto) {
2903  }
2904  else
2906 
2907  //paint histogram
2908  fPaintHisto->Paint(opt);
2909  return;
2910  }
2911  Warning("Paint","Painting 3D efficiency is not implemented");
2912 }
2913 
2914 ////////////////////////////////////////////////////////////////////////////////
2915 /// Have histograms fixed bins along each axis?
2916 
2917 void TEfficiency::SavePrimitive(std::ostream& out,Option_t* opt)
2918 {
2919  Bool_t equi_bins = true;
2920 
2921  //indentation
2922  TString indent = " ";
2923  //names for arrays containing the bin edges
2924  //static counter needed if more objects are saved
2925  static Int_t naxis = 0;
2926  TString sxaxis="xAxis",syaxis="yAxis",szaxis="zAxis";
2927 
2928  //note the missing break statements!
2929  switch(GetDimension()) {
2930  case 3:
2931  equi_bins = equi_bins && !fTotalHistogram->GetZaxis()->GetXbins()->fArray
2932  && !fTotalHistogram->GetZaxis()->GetXbins()->fN;
2933  case 2:
2934  equi_bins = equi_bins && !fTotalHistogram->GetYaxis()->GetXbins()->fArray
2935  && !fTotalHistogram->GetYaxis()->GetXbins()->fN;
2936  case 1:
2937  equi_bins = equi_bins && !fTotalHistogram->GetXaxis()->GetXbins()->fArray
2938  && !fTotalHistogram->GetXaxis()->GetXbins()->fN;
2939  }
2940 
2941  //create arrays containing the variable binning
2942  if(!equi_bins) {
2943  Int_t i;
2944  ++naxis;
2945  sxaxis += naxis;
2946  syaxis += naxis;
2947  szaxis += naxis;
2948  //x axis
2949  out << indent << "Double_t " << sxaxis << "["
2950  << fTotalHistogram->GetXaxis()->GetXbins()->fN << "] = {";
2951  for (i = 0; i < fTotalHistogram->GetXaxis()->GetXbins()->fN; ++i) {
2952  if (i != 0) out << ", ";
2953  out << fTotalHistogram->GetXaxis()->GetXbins()->fArray[i];
2954  }
2955  out << "}; " << std::endl;
2956  //y axis
2957  if(GetDimension() > 1) {
2958  out << indent << "Double_t " << syaxis << "["
2959  << fTotalHistogram->GetYaxis()->GetXbins()->fN << "] = {";
2960  for (i = 0; i < fTotalHistogram->GetYaxis()->GetXbins()->fN; ++i) {
2961  if (i != 0) out << ", ";
2962  out << fTotalHistogram->GetYaxis()->GetXbins()->fArray[i];
2963  }
2964  out << "}; " << std::endl;
2965  }
2966  //z axis
2967  if(GetDimension() > 2) {
2968  out << indent << "Double_t " << szaxis << "["
2969  << fTotalHistogram->GetZaxis()->GetXbins()->fN << "] = {";
2970  for (i = 0; i < fTotalHistogram->GetZaxis()->GetXbins()->fN; ++i) {
2971  if (i != 0) out << ", ";
2972  out << fTotalHistogram->GetZaxis()->GetXbins()->fArray[i];
2973  }
2974  out << "}; " << std::endl;
2975  }
2976  }//creating variable binning
2977 
2978  //TEfficiency pointer has efficiency name + counter
2979  static Int_t eff_count = 0;
2980  ++eff_count;
2981  TString eff_name = GetName();
2982  eff_name += eff_count;
2983 
2984  const char* name = eff_name.Data();
2985 
2986  //construct TEfficiency object
2987  const char quote = '"';
2988  out << indent << std::endl;
2989  out << indent << ClassName() << " * " << name << " = new " << ClassName()
2990  << "(" << quote << GetName() << quote << "," << quote
2991  << GetTitle() << quote <<",";
2992  //fixed bin size -> use n,min,max constructor
2993  if(equi_bins) {
2994  out << fTotalHistogram->GetXaxis()->GetNbins() << ","
2995  << fTotalHistogram->GetXaxis()->GetXmin() << ","
2996  << fTotalHistogram->GetXaxis()->GetXmax();
2997  if(GetDimension() > 1) {
2998  out << "," << fTotalHistogram->GetYaxis()->GetNbins() << ","
2999  << fTotalHistogram->GetYaxis()->GetXmin() << ","
3000  << fTotalHistogram->GetYaxis()->GetXmax();
3001  }
3002  if(GetDimension() > 2) {
3003  out << "," << fTotalHistogram->GetZaxis()->GetNbins() << ","
3004  << fTotalHistogram->GetZaxis()->GetXmin() << ","
3005  << fTotalHistogram->GetZaxis()->GetXmax();
3006  }
3007  }
3008  //variable bin size -> use n,*bins constructor
3009  else {
3010  out << fTotalHistogram->GetXaxis()->GetNbins() << "," << sxaxis;
3011  if(GetDimension() > 1)
3012  out << "," << fTotalHistogram->GetYaxis()->GetNbins() << ","
3013  << syaxis;
3014  if(GetDimension() > 2)
3015  out << "," << fTotalHistogram->GetZaxis()->GetNbins() << ","
3016  << szaxis;
3017  }
3018  out << ");" << std::endl;
3019  out << indent << std::endl;
3020 
3021  //set statistic options
3022  out << indent << name << "->SetConfidenceLevel(" << fConfLevel << ");"
3023  << std::endl;
3024  out << indent << name << "->SetBetaAlpha(" << fBeta_alpha << ");"
3025  << std::endl;
3026  out << indent << name << "->SetBetaBeta(" << fBeta_beta << ");" << std::endl;
3027  out << indent << name << "->SetWeight(" << fWeight << ");" << std::endl;
3028  out << indent << name << "->SetStatisticOption(" << fStatisticOption << ");"
3029  << std::endl;
3030  out << indent << name << "->SetPosteriorMode(" << TestBit(kPosteriorMode) << ");" << std::endl;
3031  out << indent << name << "->SetShortestInterval(" << TestBit(kShortestInterval) << ");" << std::endl;
3032  if(TestBit(kUseWeights))
3033  out << indent << name << "->SetUseWeightedEvents();" << std::endl;
3034 
3035  // save bin-by-bin prior parameters
3036  for(unsigned int i = 0; i < fBeta_bin_params.size(); ++i)
3037  {
3038  out << indent << name << "->SetBetaBinParameters(" << i << "," << fBeta_bin_params.at(i).first
3039  << "," << fBeta_bin_params.at(i).second << ");" << std::endl;
3040  }
3041 
3042  //set bin contents
3044  if(GetDimension() > 1)
3045  nbins *= fTotalHistogram->GetNbinsY() + 2;
3046  if(GetDimension() > 2)
3047  nbins *= fTotalHistogram->GetNbinsZ() + 2;
3048 
3049  //important: set first total number than passed number
3050  for(Int_t i = 0; i < nbins; ++i) {
3051  out << indent << name <<"->SetTotalEvents(" << i << "," <<
3052  fTotalHistogram->GetBinContent(i) << ");" << std::endl;
3053  out << indent << name <<"->SetPassedEvents(" << i << "," <<
3054  fPassedHistogram->GetBinContent(i) << ");" << std::endl;
3055  }
3056 
3057  //save list of functions
3058  TIter next(fFunctions);
3059  TObject* obj = 0;
3060  while((obj = next())) {
3061  obj->SavePrimitive(out,"nodraw");
3062  if(obj->InheritsFrom(TF1::Class())) {
3063  out << indent << name << "->GetListOfFunctions()->Add("
3064  << obj->GetName() << ");" << std::endl;
3065  }
3066  }
3067 
3068  //set style
3069  SaveFillAttributes(out,name);
3070  SaveLineAttributes(out,name);
3071  SaveMarkerAttributes(out,name);
3072 
3073  //draw TEfficiency object
3074  TString option = opt;
3075  option.ToLower();
3076  if (!option.Contains("nodraw"))
3077  out<< indent << name<< "->Draw(" << quote << opt << quote << ");"
3078  << std::endl;
3079 }
3080 
3081 ////////////////////////////////////////////////////////////////////////////////
3082 /// Sets the shape parameter &alpha;
3083 ///
3084 /// The prior probability of the efficiency is given by the beta distribution:
3085 /// \f[
3086 /// f(\varepsilon;\alpha;\beta) = \frac{1}{B(\alpha,\beta)} \varepsilon^{\alpha-1} (1 - \varepsilon)^{\beta-1}
3087 /// \f]
3088 ///
3089 /// Note: - both shape parameters have to be positive (i.e. > 0)
3090 
3092 {
3093  if(alpha > 0)
3094  fBeta_alpha = alpha;
3095  else
3096  Warning("SetBetaAlpha(Double_t)","invalid shape parameter %.2lf",alpha);
3097 }
3098 
3099 ////////////////////////////////////////////////////////////////////////////////
3100 /// Sets the shape parameter &beta;
3101 ///
3102 /// The prior probability of the efficiency is given by the beta distribution:
3103 /// \f[
3104 /// f(\varepsilon;\alpha,\beta) = \frac{1}{B(\alpha,\beta)} \varepsilon^{\alpha-1} (1 - \varepsilon)^{\beta-1}
3105 /// \f]
3106 ///
3107 /// Note: - both shape parameters have to be positive (i.e. > 0)
3108 
3110 {
3111  if(beta > 0)
3112  fBeta_beta = beta;
3113  else
3114  Warning("SetBetaBeta(Double_t)","invalid shape parameter %.2lf",beta);
3115 }
3116 
3117 ////////////////////////////////////////////////////////////////////////////////
3118 /// Sets different shape parameter &alpha; and &beta;
3119 /// for the prior distribution for each bin. By default the global parameter are used if they are not set
3120 /// for the specific bin
3121 /// The prior probability of the efficiency is given by the beta distribution:
3122 /// \f[
3123 /// f(\varepsilon;\alpha;\beta) = \frac{1}{B(\alpha,\beta)} \varepsilon^{\alpha-1} (1 - \varepsilon)^{\beta-1}
3124 /// \f]
3125 ///
3126 /// Note:
3127 /// - both shape parameters have to be positive (i.e. > 0)
3128 /// - bin gives the global bin number (cf. GetGlobalBin)
3129 
3131 {
3132  if (!fPassedHistogram || !fTotalHistogram) return;
3133  TH1 * h1 = fTotalHistogram;
3134  // doing this I get h1->fN which is available only for a TH1D
3135  UInt_t n = h1->GetBin(h1->GetNbinsX()+1, h1->GetNbinsY()+1, h1->GetNbinsZ()+1 ) + 1;
3136 
3137  // in case vector is not created do with default alpha, beta params
3138  if (fBeta_bin_params.size() != n )
3139  fBeta_bin_params = std::vector<std::pair<Double_t, Double_t> >(n, std::make_pair(fBeta_alpha, fBeta_beta) );
3140 
3141  // vector contains also values for under/overflows
3142  fBeta_bin_params[bin] = std::make_pair(alpha,beta);
3143  SetBit(kUseBinPrior,true);
3144 
3145 }
3146 
3147 ////////////////////////////////////////////////////////////////////////////////
3148 /// Set the bins for the underlined passed and total histograms
3149 /// If the class have been already filled the previous contents will be lost
3150 
3152 {
3153  if (GetDimension() != 1) {
3154  Error("SetBins","Using wrong SetBins function for a %d-d histogram",GetDimension());
3155  return kFALSE;
3156  }
3157  if (fTotalHistogram->GetEntries() != 0 ) {
3158  Warning("SetBins","Histogram entries will be lost after SetBins");
3161  }
3162  fPassedHistogram->SetBins(nx,xmin,xmax);
3163  fTotalHistogram->SetBins(nx,xmin,xmax);
3164  return kTRUE;
3165 }
3166 
3167 ////////////////////////////////////////////////////////////////////////////////
3168 /// Set the bins for the underlined passed and total histograms
3169 /// If the class have been already filled the previous contents will be lost
3170 
3172 {
3173  if (GetDimension() != 1) {
3174  Error("SetBins","Using wrong SetBins function for a %d-d histogram",GetDimension());
3175  return kFALSE;
3176  }
3177  if (fTotalHistogram->GetEntries() != 0 ) {
3178  Warning("SetBins","Histogram entries will be lost after SetBins");
3181  }
3182  fPassedHistogram->SetBins(nx,xBins);
3183  fTotalHistogram->SetBins(nx,xBins);
3184  return kTRUE;
3185 }
3186 
3187 ////////////////////////////////////////////////////////////////////////////////
3188 /// Set the bins for the underlined passed and total histograms
3189 /// If the class have been already filled the previous contents will be lost
3190 
3192 {
3193  if (GetDimension() != 2) {
3194  Error("SetBins","Using wrong SetBins function for a %d-d histogram",GetDimension());
3195  return kFALSE;
3196  }
3197  if (fTotalHistogram->GetEntries() != 0 ) {
3198  Warning("SetBins","Histogram entries will be lost after SetBins");
3201  }
3202  fPassedHistogram->SetBins(nx,xmin,xmax,ny,ymin,ymax);
3203  fTotalHistogram->SetBins(nx,xmin,xmax,ny,ymin,ymax);
3204  return kTRUE;
3205 }
3206 
3207 ////////////////////////////////////////////////////////////////////////////////
3208 /// Set the bins for the underlined passed and total histograms
3209 /// If the class have been already filled the previous contents will be lost
3210 
3212 {
3213  if (GetDimension() != 2) {
3214  Error("SetBins","Using wrong SetBins function for a %d-d histogram",GetDimension());
3215  return kFALSE;
3216  }
3217  if (fTotalHistogram->GetEntries() != 0 ) {
3218  Warning("SetBins","Histogram entries will be lost after SetBins");
3221  }
3222  fPassedHistogram->SetBins(nx,xBins,ny,yBins);
3223  fTotalHistogram->SetBins(nx,xBins,ny,yBins);
3224  return kTRUE;
3225 }
3226 
3227 ////////////////////////////////////////////////////////////////////////////////
3228 /// Set the bins for the underlined passed and total histograms
3229 /// If the class have been already filled the previous contents will be lost
3230 
3232  Int_t nz, Double_t zmin, Double_t zmax)
3233 {
3234  if (GetDimension() != 3) {
3235  Error("SetBins","Using wrong SetBins function for a %d-d histogram",GetDimension());
3236  return kFALSE;
3237  }
3238  if (fTotalHistogram->GetEntries() != 0 ) {
3239  Warning("SetBins","Histogram entries will be lost after SetBins");
3242  }
3243  fPassedHistogram->SetBins(nx,xmin,xmax,ny,ymin,ymax,nz,zmin,zmax);
3244  fTotalHistogram->SetBins (nx,xmin,xmax,ny,ymin,ymax,nz,zmin,zmax);
3245  return kTRUE;
3246 }
3247 
3248 ////////////////////////////////////////////////////////////////////////////////
3249 /// Set the bins for the underlined passed and total histograms
3250 /// If the class have been already filled the previous contents will be lost
3251 
3252 Bool_t TEfficiency::SetBins(Int_t nx, const Double_t *xBins, Int_t ny, const Double_t *yBins, Int_t nz,
3253  const Double_t *zBins )
3254 {
3255  if (GetDimension() != 3) {
3256  Error("SetBins","Using wrong SetBins function for a %d-d histogram",GetDimension());
3257  return kFALSE;
3258  }
3259  if (fTotalHistogram->GetEntries() != 0 ) {
3260  Warning("SetBins","Histogram entries will be lost after SetBins");
3263  }
3264  fPassedHistogram->SetBins(nx,xBins,ny,yBins,nz,zBins);
3265  fTotalHistogram->SetBins(nx,xBins,ny,yBins,nz,zBins);
3266  return kTRUE;
3267 }
3268 
3269 ////////////////////////////////////////////////////////////////////////////////
3270 /// Sets the confidence level (0 < level < 1)
3271 /// The default value is 1-sigma :~ 0.683
3272 
3274 {
3275  if((level > 0) && (level < 1))
3276  fConfLevel = level;
3277  else
3278  Warning("SetConfidenceLevel(Double_t)","invalid confidence level %.2lf",level);
3279 }
3280 
3281 ////////////////////////////////////////////////////////////////////////////////
3282 /// Sets the directory holding this TEfficiency object
3283 ///
3284 /// A reference to this TEfficiency object is removed from the current
3285 /// directory (if it exists) and a new reference to this TEfficiency object is
3286 /// added to the given directory.
3287 ///
3288 /// Notes: - If the given directory is 0, the TEfficiency object does not
3289 /// belong to any directory and will not be written to file during the
3290 /// next TFile::Write() command.
3291 
3293 {
3294  if(fDirectory == dir)
3295  return;
3296  if(fDirectory)
3297  fDirectory->Remove(this);
3298  fDirectory = dir;
3299  if(fDirectory)
3300  fDirectory->Append(this);
3301 }
3302 
3303 ////////////////////////////////////////////////////////////////////////////////
3304 /// Sets the name
3305 ///
3306 /// Note: The names of the internal histograms are set to "name + _total" and
3307 /// "name + _passed" respectively.
3308 
3309 void TEfficiency::SetName(const char* name)
3310 {
3311  TNamed::SetName(name);
3312 
3313  //setting the names (appending the correct ending)
3314  TString name_total = name + TString("_total");
3315  TString name_passed = name + TString("_passed");
3316  fTotalHistogram->SetName(name_total);
3317  fPassedHistogram->SetName(name_passed);
3318 }
3319 
3320 ////////////////////////////////////////////////////////////////////////////////
3321 /// Sets the number of passed events in the given global bin
3322 ///
3323 /// returns "true" if the number of passed events has been updated
3324 /// otherwise "false" ist returned
3325 ///
3326 /// Note: - requires: 0 <= events <= fTotalHistogram->GetBinContent(bin)
3327 
3329 {
3330  if(events <= fTotalHistogram->GetBinContent(bin)) {
3331  fPassedHistogram->SetBinContent(bin,events);
3332  return true;
3333  }
3334  else {
3335  Error("SetPassedEvents(Int_t,Int_t)","total number of events (%.1lf) in bin %i is less than given number of passed events %i",fTotalHistogram->GetBinContent(bin),bin,events);
3336  return false;
3337  }
3338 }
3339 
3340 ////////////////////////////////////////////////////////////////////////////////
3341 /// Sets the histogram containing the passed events
3342 ///
3343 /// The given histogram is cloned and stored internally as histogram containing
3344 /// the passed events. The given histogram has to be consistent with the current
3345 /// fTotalHistogram (see CheckConsistency(const TH1&,const TH1&)).
3346 /// The method returns whether the fPassedHistogram has been replaced (true) or
3347 /// not (false).
3348 ///
3349 /// Note: The list of associated functions fFunctions is cleared.
3350 ///
3351 /// Option:
3352 /// - "f": force the replacement without checking the consistency
3353 /// This can lead to inconsistent histograms and useless results
3354 /// or unexpected behaviour. But sometimes it might be the only
3355 /// way to change the histograms. If you use this option, you
3356 /// should ensure that the fTotalHistogram is replaced by a
3357 /// consistent one (with respect to rPassed) as well.
3358 
3360 {
3361  TString option = opt;
3362  option.ToLower();
3363 
3364  Bool_t bReplace = option.Contains("f");
3365 
3366  if(!bReplace)
3367  bReplace = CheckConsistency(rPassed,*fTotalHistogram);
3368 
3369  if(bReplace) {
3370  delete fPassedHistogram;
3371  Bool_t bStatus = TH1::AddDirectoryStatus();
3373  fPassedHistogram = (TH1*)(rPassed.Clone());
3375  TH1::AddDirectory(bStatus);
3376 
3377  if(fFunctions)
3378  fFunctions->Delete();
3379 
3380  //check whether both histograms are filled with weights
3381  bool useWeights = CheckWeights(rPassed,*fTotalHistogram);
3382 
3383  SetUseWeightedEvents(useWeights);
3384 
3385  return true;
3386  }
3387  else
3388  return false;
3389 }
3390 
3391 ////////////////////////////////////////////////////////////////////////////////
3392 /// Sets the statistic option which affects the calculation of the confidence interval
3393 ///
3394 /// Options:
3395 /// - kFCP (=0)(default): using the Clopper-Pearson interval (recommended by PDG)
3396 /// sets kIsBayesian = false
3397 /// see also ClopperPearson
3398 /// - kFNormal (=1) : using the normal approximation
3399 /// sets kIsBayesian = false
3400 /// see also Normal
3401 /// - kFWilson (=2) : using the Wilson interval
3402 /// sets kIsBayesian = false
3403 /// see also Wilson
3404 /// - kFAC (=3) : using the Agresti-Coull interval
3405 /// sets kIsBayesian = false
3406 /// see also AgrestiCoull
3407 /// - kFFC (=4) : using the Feldman-Cousins frequentist method
3408 /// sets kIsBayesian = false
3409 /// see also FeldmanCousins
3410 /// - kBJeffrey (=5) : using the Jeffrey interval
3411 /// sets kIsBayesian = true, fBeta_alpha = 0.5 and fBeta_beta = 0.5
3412 /// see also Bayesian
3413 /// - kBUniform (=6) : using a uniform prior
3414 /// sets kIsBayesian = true, fBeta_alpha = 1 and fBeta_beta = 1
3415 /// see also Bayesian
3416 /// - kBBayesian (=7) : using a custom prior defined by fBeta_alpha and fBeta_beta
3417 /// sets kIsBayesian = true
3418 /// see also Bayesian
3419 /// - kMidP (=8) : using the Lancaster Mid-P method
3420 /// sets kIsBayesian = false
3421 
3422 
3424 {
3425  fStatisticOption = option;
3426 
3427  switch(option)
3428  {
3429  case kFCP:
3431  SetBit(kIsBayesian,false);
3432  break;
3433  case kFNormal:
3434  fBoundary = &Normal;
3435  SetBit(kIsBayesian,false);
3436  break;
3437  case kFWilson:
3438  fBoundary = &Wilson;
3439  SetBit(kIsBayesian,false);
3440  break;
3441  case kFAC:
3443  SetBit(kIsBayesian,false);
3444  break;
3445  case kFFC:
3447  SetBit(kIsBayesian,false);
3448  break;
3449  case kMidP:
3451  SetBit(kIsBayesian,false);
3452  break;
3453  case kBJeffrey:
3454  fBeta_alpha = 0.5;
3455  fBeta_beta = 0.5;
3456  SetBit(kIsBayesian,true);
3457  SetBit(kUseBinPrior,false);
3458  break;
3459  case kBUniform:
3460  fBeta_alpha = 1;
3461  fBeta_beta = 1;
3462  SetBit(kIsBayesian,true);
3463  SetBit(kUseBinPrior,false);
3464  break;
3465  case kBBayesian:
3466  SetBit(kIsBayesian,true);
3467  break;
3468  default:
3471  SetBit(kIsBayesian,false);
3472  }
3473 }
3474 
3475 ////////////////////////////////////////////////////////////////////////////////
3476 /// Sets the title
3477 ///
3478 /// Notes:
3479 /// - The titles of the internal histograms are set to "title + (total)"
3480 /// or "title + (passed)" respectively.
3481 /// - It is possible to label the axis of the histograms as usual (see
3482 /// TH1::SetTitle).
3483 ///
3484 /// Example: Setting the title to "My Efficiency" and label the axis
3485 /// pEff->SetTitle("My Efficiency;x label;eff");
3486 
3487 void TEfficiency::SetTitle(const char* title)
3488 {
3489 
3490  //setting the titles (looking for the first semicolon and insert the tokens there)
3491  TString title_passed = title;
3492  TString title_total = title;
3493  Ssiz_t pos = title_passed.First(";");
3494  if (pos != kNPOS) {
3495  title_passed.Insert(pos," (passed)");
3496  title_total.Insert(pos," (total)");
3497  }
3498  else {
3499  title_passed.Append(" (passed)");
3500  title_total.Append(" (total)");
3501  }
3502  fPassedHistogram->SetTitle(title_passed);
3503  fTotalHistogram->SetTitle(title_total);
3504 
3505  // strip (total) for the TEfficiency title
3506  // HIstogram SetTitle has already stripped the axis
3507  TString teffTitle = fTotalHistogram->GetTitle();
3508  teffTitle.ReplaceAll(" (total)","");
3509  TNamed::SetTitle(teffTitle);
3510 
3511 }
3512 
3513 ////////////////////////////////////////////////////////////////////////////////
3514 /// Sets the number of total events in the given global bin
3515 ///
3516 /// returns "true" if the number of total events has been updated
3517 /// otherwise "false" ist returned
3518 ///
3519 /// Note: - requires: fPassedHistogram->GetBinContent(bin) <= events
3520 
3522 {
3523  if(events >= fPassedHistogram->GetBinContent(bin)) {
3524  fTotalHistogram->SetBinContent(bin,events);
3525  return true;
3526  }
3527  else {
3528  Error("SetTotalEvents(Int_t,Int_t)","passed number of events (%.1lf) in bin %i is bigger than given number of total events %i",fPassedHistogram->GetBinContent(bin),bin,events);
3529  return false;
3530  }
3531 }
3532 
3533 ////////////////////////////////////////////////////////////////////////////////
3534 /// Sets the histogram containing all events
3535 ///
3536 /// The given histogram is cloned and stored internally as histogram containing
3537 /// all events. The given histogram has to be consistent with the current
3538 /// fPassedHistogram (see CheckConsistency(const TH1&,const TH1&)).
3539 /// The method returns whether the fTotalHistogram has been replaced (true) or
3540 /// not (false).
3541 ///
3542 /// Note: The list of associated functions fFunctions is cleared.
3543 ///
3544 /// Option:
3545 /// - "f": force the replacement without checking the consistency
3546 /// This can lead to inconsistent histograms and useless results
3547 /// or unexpected behaviour. But sometimes it might be the only
3548 /// way to change the histograms. If you use this option, you
3549 /// should ensure that the fPassedHistogram is replaced by a
3550 /// consistent one (with respect to rTotal) as well.
3551 
3553 {
3554  TString option = opt;
3555  option.ToLower();
3556 
3557  Bool_t bReplace = option.Contains("f");
3558 
3559  if(!bReplace)
3560  bReplace = CheckConsistency(*fPassedHistogram,rTotal);
3561 
3562  if(bReplace) {
3563  delete fTotalHistogram;
3564  Bool_t bStatus = TH1::AddDirectoryStatus();
3566  fTotalHistogram = (TH1*)(rTotal.Clone());
3568  TH1::AddDirectory(bStatus);
3569 
3570  if(fFunctions)
3571  fFunctions->Delete();
3572 
3573  //check whether both histograms are filled with weights
3574  bool useWeights = CheckWeights(*fPassedHistogram,rTotal);
3575  SetUseWeightedEvents(useWeights);
3576 
3577  return true;
3578  }
3579  else
3580  return false;
3581 }
3582 
3583 ////////////////////////////////////////////////////////////////////////////////
3584 
3586 {
3587  if (on && !TestBit(kUseWeights) )
3588  gROOT->Info("TEfficiency::SetUseWeightedEvents","Histograms are filled with weights");
3589 
3590  SetBit(kUseWeights,on);
3591 
3592  // no need to set sumw2 should be already there
3593  if (on) assert(fTotalHistogram->GetSumw2N() > 0 && fPassedHistogram->GetSumw2N() > 0 );
3594  //fTotalHistogram->Sumw2(on);
3595  //fPassedHistogram->Sumw2(on);
3596 }
3597 
3598 ////////////////////////////////////////////////////////////////////////////////
3599 /// Sets the global weight for this TEfficiency object
3600 ///
3601 /// Note: - weight has to be positive ( > 0)
3602 
3604 {
3605  if(weight > 0)
3606  fWeight = weight;
3607  else
3608  Warning("SetWeight","invalid weight %.2lf",weight);
3609 }
3610 
3611 ////////////////////////////////////////////////////////////////////////////////
3612 /**
3613 Calculates the boundaries for the frequentist Wilson interval
3614 
3615 \param[in] total number of total events
3616 \param[in] passed 0 <= number of passed events <= total
3617 \param[in] level confidence level
3618 \param[in] bUpper
3619  - true - upper boundary is returned
3620  - false - lower boundary is returned
3621 
3622 Calculation:
3623 \f{eqnarray*}{
3624  \alpha &=& 1 - \frac{level}{2}\\
3625  \kappa &=& \Phi^{-1}(1 - \alpha,1) ...\ normal\ quantile\ function\\
3626  mode &=& \frac{passed + \frac{\kappa^{2}}{2}}{total + \kappa^{2}}\\
3627  \Delta &=& \frac{\kappa}{total + \kappa^{2}} * \sqrt{passed (1 - \frac{passed}{total}) + \frac{\kappa^{2}}{4}}\\
3628  return &=& max(0,mode - \Delta)\ or\ min(1,mode + \Delta)
3629 \f}
3630 
3631 */
3632 
3634 {
3635  Double_t alpha = (1.0 - level)/2;
3636  if (total == 0) return (bUpper) ? 1 : 0;
3637  Double_t average = ((Double_t)passed) / total;
3638  Double_t kappa = ROOT::Math::normal_quantile(1 - alpha,1);
3639 
3640  Double_t mode = (passed + 0.5 * kappa * kappa) / (total + kappa * kappa);
3641  Double_t delta = kappa / (total + kappa*kappa) * std::sqrt(total * average
3642  * (1 - average) + kappa * kappa / 4);
3643  if(bUpper)
3644  return ((mode + delta) > 1) ? 1.0 : (mode + delta);
3645  else
3646  return ((mode - delta) < 0) ? 0.0 : (mode - delta);
3647 }
3648 
3649 ////////////////////////////////////////////////////////////////////////////////
3650 /// Addition operator
3651 ///
3652 /// adds the corresponding histograms:
3653 /// ~~~ {.cpp}
3654 /// lhs.GetTotalHistogram() + rhs.GetTotalHistogram()
3655 /// lhs.GetPassedHistogram() + rhs.GetPassedHistogram()
3656 /// ~~~
3657 /// the statistic option and the confidence level are taken from lhs
3658 
3659 const TEfficiency operator+(const TEfficiency& lhs,const TEfficiency& rhs)
3660 {
3661  TEfficiency tmp(lhs);
3662  tmp += rhs;
3663  return tmp;
3664 }
3665 
3666 #endif
TGraphAsymmErrors * CreateGraph(Option_t *opt="") const
Create the graph used be painted (for dim=1 TEfficiency) The return object is managed by the caller...
const int nx
Definition: kalman.C:16
Double_t * GetEXlow() const
double normal_quantile(double z, double sigma)
Inverse ( ) of the cumulative distribution function of the lower tail of the normal (Gaussian) distri...
void SetConfidenceLevel(Double_t level)
Sets the confidence level (0 < level < 1) The default value is 1-sigma :~ 0.683.
static Double_t Bayesian(Double_t total, Double_t passed, Double_t level, Double_t alpha, Double_t beta, Bool_t bUpper, Bool_t bShortest=false)
Calculates the boundaries for a Bayesian confidence interval (shortest or central interval depending ...
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Double_t fBeta_beta
Definition: TEfficiency.h:47
object has not been deleted
Definition: TObject.h:71
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3126
virtual void Paint(Option_t *option="")
Control routine to paint any kind of histograms.
Definition: TH1.cxx:5539
static Bool_t CheckEntries(const TH1 &pass, const TH1 &total, Option_t *opt="")
Checks whether bin contents are compatible with binomial statistics.
float xmin
Definition: THbookFile.cxx:93
virtual void Delete(Option_t *option="")
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:409
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
Definition: TH1.cxx:8253
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:847
long long Long64_t
Definition: RtypesCore.h:69
void Build(const char *name, const char *title)
Building standard data structure of a TEfficiency object.
Class to handle efficiency histograms.
Definition: TEfficiency.h:27
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
virtual void SetDirectory(TDirectory *dir)
By default when an histogram is created, it is added to the list of histogram objects in the current ...
Definition: TH1.cxx:8053
static Double_t BetaCentralInterval(Double_t level, Double_t alpha, Double_t beta, Bool_t bUpper)
Calculates the boundaries for a central confidence interval for a Beta distribution.
const char Option_t
Definition: RtypesCore.h:62
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition: TAxis.cxx:504
float ymin
Definition: THbookFile.cxx:93
double normal_quantile_c(double z, double sigma)
Inverse ( ) of the cumulative distribution function of the upper tail of the normal (Gaussian) distri...
const Ssiz_t kNPOS
Definition: RtypesCore.h:115
Binomial fitter for the division of two histograms.
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:640
const Double_t * GetArray() const
Definition: TArrayD.h:43
Bool_t TestBit(UInt_t f) const
Definition: TObject.h:159
virtual void SetBins(Int_t nx, Double_t xmin, Double_t xmax)
Redefine x axis parameters.
Definition: TH1.cxx:7883
static Bool_t CheckBinning(const TH1 &pass, const TH1 &total)
Checks binning for each axis.
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:131
THist< 1, float, THistStatContent, THistStatUncertainty > TH1F
Definition: THist.hxx:311
TEfficiency & operator+=(const TEfficiency &rhs)
Adds the histograms of another TEfficiency object to current histograms.
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4639
double beta_pdf(double x, double a, double b)
Probability density function of the beta distribution.
TAxis * GetYaxis() const
Get y axis of the graph.
Definition: TGraph.cxx:1602
Double_t GetBetaBeta(Int_t bin=-1) const
Definition: TEfficiency.h:107
virtual Int_t GetNbinsZ() const
Definition: TH1.h:279
const TEfficiency operator+(const TEfficiency &lhs, const TEfficiency &rhs)
Addition operator.
#define gROOT
Definition: TROOT.h:375
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition: TString.h:587
const Double_t kDefBetaBeta
Definition: TEfficiency.cxx:37
Basic string class.
Definition: TString.h:129
static Bool_t AddDirectoryStatus()
Static function: cannot be inlined on Windows/NT.
Definition: TH1.cxx:700
void FillGraph(TGraphAsymmErrors *graph, Option_t *opt) const
Fill the graph to be painted with information from TEfficiency Internal method called by TEfficiency:...
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1099
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
void Copy(TAttMarker &attmarker) const
Copy this marker attributes to a new TAttMarker.
Definition: TAttMarker.cxx:209
virtual void SetTitle(const char *title="")
Set graph title.
Definition: TGraph.cxx:2180
void SetBetaBinParameters(Int_t bin, Double_t alpha, Double_t beta)
Sets different shape parameter α and β for the prior distribution for each bin.
int nbins[3]
TFitResultPtr Fit(TF1 *f1, Option_t *opt="")
Fits the efficiency using the TBinomialEfficiencyFitter class.
~TEfficiency()
default destructor
virtual Double_t GetBinLowEdge(Int_t bin) const
Return bin lower edge for 1D histogram.
Definition: TH1.cxx:8264
EStatOption fStatisticOption
Definition: TEfficiency.h:57
Bool_t SetPassedEvents(Int_t bin, Int_t events)
Sets the number of passed events in the given global bin.
TString & Insert(Ssiz_t pos, const char *s)
Definition: TString.h:597
Double_t(* fBoundary)(Double_t, Double_t, Double_t, Bool_t)
Definition: TEfficiency.h:50
Float_t delta
TRObject operator()(const T1 &t1) const
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:687
virtual void GetStats(Double_t *stats) const
fill the array stats from the contents of this histogram The array stats must be correctly dimensione...
Definition: TH1.cxx:7021
virtual void Reset(Option_t *option="")
Reset this histogram: contents, errors, etc.
Definition: TH1.cxx:6419
double beta(double x, double y)
Calculates the beta function.
static void AddDirectory(Bool_t add=kTRUE)
Sets the flag controlling the automatic add of histograms in memory.
Definition: TH1.cxx:1218
Double_t GetWeight() const
Definition: TEfficiency.h:123
virtual void AppendPad(Option_t *option="")
Append graphics object to current pad.
Definition: TObject.cxx:112
Template class to wrap any C++ callable object which takes one argument i.e.
Marker Attributes class.
Definition: TAttMarker.h:19
TH1 * GetCopyPassedHisto() const
Returns a cloned version of fPassedHistogram.
Double_t fWeight
Definition: TEfficiency.h:59
virtual Int_t GetDimension() const
Definition: TH1.h:263
double sqrt(double)
static struct mg_connection * fc(struct mg_context *ctx)
Definition: civetweb.c:1956
TGraph with asymmetric error bars.
Double_t GetXmin() const
Definition: TAxis.h:133
EStatOption GetStatisticOption() const
Definition: TEfficiency.h:121
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:135
Double_t GetEfficiency(Int_t bin) const
Returns the efficiency in the given global bin.
static Double_t Combine(Double_t &up, Double_t &low, Int_t n, const Int_t *pass, const Int_t *total, Double_t alpha, Double_t beta, Double_t level=0.683, const Double_t *w=0, Option_t *opt="")
Calculates the combined efficiency and its uncertainties.
Fill Area Attributes class.
Definition: TAttFill.h:19
Double_t x[n]
Definition: legend1.C:17
virtual void Paint(Option_t *chopt="")
Draw this graph with its current attributes.
Definition: TGraph.cxx:1936
void Class()
Definition: Class.C:29
TGraphAsymmErrors * fPaintGraph
Definition: TEfficiency.h:54
static Double_t AgrestiCoull(Double_t total, Double_t passed, Double_t level, Bool_t bUpper)
Calculates the boundaries for the frequentist Agresti-Coull interval.
void Copy(TAttLine &attline) const
Copy this line attributes to a new TAttLine.
Definition: TAttLine.cxx:162
virtual void SaveLineAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1, Int_t widdef=1)
Save line attributes as C++ statement(s) on output stream out.
Definition: TAttLine.cxx:260
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
const int ny
Definition: kalman.C:17
static Bool_t CheckWeights(const TH1 &pass, const TH1 &total)
Check if both histogram are weighted.
TH1 * fPassedHistogram
temporary histogram for painting
Definition: TEfficiency.h:56
TH1 * fTotalHistogram
Definition: TEfficiency.h:58
TH2 * fPaintHisto
temporary graph for painting
Definition: TEfficiency.h:55
TList * GetListOfFunctions()
TString & Append(const char *cs)
Definition: TString.h:497
std::vector< std::vector< double > > Data
Double_t * fArray
Definition: TArrayD.h:30
TList * fFunctions
pointer to directory holding this TEfficiency object
Definition: TEfficiency.h:53
Ssiz_t First(char c) const
Find first occurrence of a character c.
Definition: TString.cxx:477
const Double_t sigma
Bool_t UsesBayesianStat() const
Definition: TEfficiency.h:156
virtual void SaveMarkerAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1, Int_t sizdef=1)
Save line attributes as C++ statement(s) on output stream out.
Definition: TAttMarker.cxx:244
virtual TArrayD * GetSumw2()
Definition: TH1.h:293
Bool_t SetPassedHistogram(const TH1 &rPassed, Option_t *opt)
Sets the histogram containing the passed events.
void SetTitle(const char *title)
Sets the title.
TH1F * h1
Definition: legend1.C:5
User class for performing function minimization.
static Double_t BetaMean(Double_t alpha, Double_t beta)
Compute the mean (average) of the beta distribution.
double beta_cdf(double x, double a, double b)
Cumulative distribution function of the beta distribution Upper tail of the integral of the beta_pdf...
The 3-D histogram classes derived from the 1-D histogram classes.
Definition: TH3.h:31
static Double_t MidPInterval(Double_t total, Double_t passed, Double_t level, Bool_t bUpper)
Calculates the boundaries using the mid-P binomial interval (Lancaster method) from B...
Bool_t SetBins(Int_t nx, Double_t xmin, Double_t xmax)
Set the bins for the underlined passed and total histograms If the class have been already filled the...
Double_t fConfLevel
pointer to a method calculating the boundaries of confidence intervals
Definition: TEfficiency.h:51
void SetStatisticOption(EStatOption option)
Sets the statistic option which affects the calculation of the confidence interval.
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TGraph.cxx:966
A doubly linked list.
Definition: TList.h:43
TH1 * GetCopyTotalHisto() const
Returns a cloned version of fTotalHistogram.
Bool_t SetTotalEvents(Int_t bin, Int_t events)
Sets the number of total events in the given global bin.
Bool_t AreEqualRel(Double_t af, Double_t bf, Double_t relPrec)
Definition: TMath.h:322
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a line.
Definition: TH1.cxx:2589
void SetWeight(Double_t weight)
Sets the global weight for this TEfficiency object.
THist< 3, double, THistStatContent, THistStatUncertainty > TH3D
Definition: THist.hxx:322
virtual void ExecuteEvent(Int_t event, Int_t px, Int_t py)
Execute action corresponding to one event.
Definition: TH1.cxx:3022
Int_t fN
Definition: TArray.h:38
float ymax
Definition: THbookFile.cxx:93
const double tol
virtual TObject * First() const
Return the first object in the list. Returns 0 when list is empty.
Definition: TList.cxx:561
static Bool_t FeldmanCousinsInterval(Double_t total, Double_t passed, Double_t level, Double_t &lower, Double_t &upper)
Calculates the interval boundaries using the frequentist methods of Feldman-Cousins.
const char * GetTitle() const
Returns title of object.
Definition: TAxis.h:129
Service class for 2-Dim histogram classes.
Definition: TH2.h:30
Class to manage histogram axis.
Definition: TAxis.h:30
void SetBetaAlpha(Double_t alpha)
Sets the shape parameter α.
SVector< double, 2 > v
Definition: Dict.h:5
if object ctor succeeded but object should not be used
Definition: TObject.h:65
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:679
Provides an indirection to the TFitResult class and with a semantics identical to a TFitResult pointe...
Definition: TFitResultPtr.h:31
static Double_t Normal(Double_t total, Double_t passed, Double_t level, Bool_t bUpper)
Returns the confidence limits for the efficiency supposing that the efficiency follows a normal distr...
void SetFunction(const ROOT::Math::IGenFunction &f, double xlow, double xup)
Sets function to be minimized.
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:436
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
Definition: TH1.cxx:8325
Collection abstract base class.
Definition: TCollection.h:42
unsigned int UInt_t
Definition: RtypesCore.h:42
virtual void SaveFillAttributes(std::ostream &out, const char *name, Int_t coldef=1, Int_t stydef=1001)
Save fill attributes as C++ statement(s) on output stream out.
Definition: TAttFill.cxx:232
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:873
Ssiz_t Length() const
Definition: TString.h:388
virtual void Append(TObject *obj, Bool_t replace=kFALSE)
Append object to this directory.
Definition: TDirectory.cxx:153
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a graph.
Int_t GetN() const
Definition: TGraph.h:122
const Double_t kDefWeight
Definition: TEfficiency.cxx:40
TAxis * GetYaxis()
Definition: TH1.h:301
Bin contents are average (used by Add)
Definition: TH1.h:154
void FillHistogram(TH2 *h2) const
Fill the 2d histogram to be painted with information from TEfficiency 2D Internal method called by TE...
float xmax
Definition: THbookFile.cxx:93
static double p1(double t, double a, double b)
void SavePrimitive(std::ostream &out, Option_t *opt="")
Have histograms fixed bins along each axis?
const Double_t kDefBetaAlpha
Definition: TEfficiency.cxx:36
virtual Int_t GetBin(Int_t binx, Int_t biny=0, Int_t binz=0) const
Return Global bin number corresponding to binx,y,z.
Definition: TH1.cxx:4541
Definition: graph.py:1
Int_t FindFixBin(Double_t x, Double_t y=0, Double_t z=0) const
Returns the global bin number containing the given values.
TGraphErrors * gr
Definition: legend1.C:25
constexpr Double_t E()
Definition: TMath.h:74
TAxis * GetXaxis() const
Get x axis of the graph.
Definition: TGraph.cxx:1592
const Bool_t kFALSE
Definition: RtypesCore.h:92
Double_t * GetX() const
Definition: TGraph.h:129
double beta_quantile_c(double x, double a, double b)
Inverse ( ) of the cumulative distribution function of the lower tail of the beta distribution (beta_...
virtual void SetName(const char *name)
Change the name of this histogram.
Definition: TH1.cxx:8073
Double_t At(Int_t i) const
Definition: TArrayD.h:79
static unsigned int total
int Ssiz_t
Definition: RtypesCore.h:63
virtual Int_t GetSumw2N() const
Definition: TH1.h:295
virtual Bool_t IsEmpty() const
Definition: TCollection.h:93
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width for 1D histogram.
Definition: TH1.cxx:8275
virtual TObject * Remove(TObject *)
Remove an object from the in-memory list.
#define ClassImp(name)
Definition: Rtypes.h:336
const TEfficiency::EStatOption kDefStatOpt
Definition: TEfficiency.cxx:39
void SetDirectory(TDirectory *dir)
Sets the directory holding this TEfficiency object.
TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases...
Definition: TGraph.cxx:1469
double Double_t
Definition: RtypesCore.h:55
void SetUseWeightedEvents(Bool_t on=kTRUE)
Bool_t SetTotalHistogram(const TH1 &rTotal, Option_t *opt)
Sets the histogram containing all events.
virtual void PaintStats(TF1 *fit)
Draw the stats.
Definition: TGraph.cxx:1963
Describe directory structure in memory.
Definition: TDirectory.h:34
Double_t GetEfficiencyErrorUp(Int_t bin) const
Returns the upper error on the efficiency in the given global bin.
Double_t fBeta_alpha
Definition: TEfficiency.h:46
virtual Int_t DistancetoPrimitive(Int_t px, Int_t py)
Compute distance from point px,py to a graph.
Definition: TGraph.cxx:789
Double_t GetConfidenceLevel() const
Definition: TEfficiency.h:108
Double_t y[n]
Definition: legend1.C:17
double func(double *x, double *p)
Definition: stressTF1.cxx:213
Double_t GetBetaAlpha(Int_t bin=-1) const
Definition: TEfficiency.h:106
virtual double FValMinimum() const
Return function value at current estimate of the minimum.
TEfficiency & operator=(const TEfficiency &rhs)
Assignment operator.
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:572
The TH1 histogram class.
Definition: TH1.h:56
void SetBetaBeta(Double_t beta)
Sets the shape parameter β.
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
virtual Double_t GetEntries() const
Return the current number of entries.
Definition: TH1.cxx:4054
THist< 2, double, THistStatContent, THistStatUncertainty > TH2D
Definition: THist.hxx:316
TDirectory * fDirectory
Definition: TEfficiency.h:52
Long64_t Merge(TCollection *list)
Merges the TEfficiency objects in the given list to the given TEfficiency object using the operator+=...
Bool_t IsNull() const
Definition: TString.h:385
Bool_t IsVariableBinSize() const
Definition: TAxis.h:136
virtual Bool_t Add(TF1 *h1, Double_t c1=1, Option_t *option="")
Performs the operation: this = this + c1*f1 if errors are defined (see TH1::Sumw2), errors are also recalculated.
Definition: TH1.cxx:770
TH2 * CreateHistogram(Option_t *opt="") const
Create the histogram used to be painted (for dim=2 TEfficiency) The return object is managed by the c...
static Double_t FeldmanCousins(Double_t total, Double_t passed, Double_t level, Bool_t bUpper)
Calculates the boundaries for the frequentist Feldman-Cousins interval.
TAxis * GetZaxis()
Definition: TH1.h:302
double beta_cdf_c(double x, double a, double b)
Complement of the cumulative distribution function of the beta distribution.
virtual double XMinimum() const
Return current estimate of the position of the minimum.
Mother of all ROOT objects.
Definition: TObject.h:37
virtual Int_t FindFixBin(Double_t x) const
Find bin number corresponding to abscissa x.
Definition: TAxis.cxx:405
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
Double_t GetEfficiencyErrorLow(Int_t bin) const
Returns the lower error on the efficiency in the given global bin.
Int_t GetGlobalBin(Int_t binx, Int_t biny=0, Int_t binz=0) const
Returns the global bin number which can be used as argument for the following functions: ...
void SetName(const char *name)
Sets the name.
TFitResultPtr Fit(TF1 *f1, Option_t *option="")
Carry out the fit of the given function to the given histograms.
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition: TGraph.cxx:2156
Double_t * GetY() const
Definition: TGraph.h:130
static Bool_t BetaShortestInterval(Double_t level, Double_t alpha, Double_t beta, Double_t &lower, Double_t &upper)
Calculates the boundaries for a shortest confidence interval for a Beta distribution.
void Draw(Option_t *opt="")
Draws the current TEfficiency object.
void Paint(Option_t *opt)
Paints this TEfficiency object.
virtual void Add(TObject *obj)
Definition: TList.h:77
Double_t * GetEYlow() const
Double_t * GetEYhigh() const
void SetNpx(int npx)
Set the number of point used to bracket root using a grid.
static Double_t ClopperPearson(Double_t total, Double_t passed, Double_t level, Bool_t bUpper)
Calculates the boundaries for the frequentist Clopper-Pearson interval.
1-Dim function class
Definition: TF1.h:150
void FillWeighted(Bool_t bPassed, Double_t weight, Double_t x, Double_t y=0, Double_t z=0)
This function is used for filling the two histograms with a weight.
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
Definition: TH1.cxx:2544
TF1 * f1
Definition: legend1.C:11
const Double_t kDefConfLevel
Definition: TEfficiency.cxx:38
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
THist< 1, double, THistStatContent, THistStatUncertainty > TH1D
Definition: THist.hxx:310
#define gPad
Definition: TVirtualPad.h:284
void Init(double alpha)
#define gDirectory
Definition: TDirectory.h:211
double result[121]
std::vector< std::pair< Double_t, Double_t > > fBeta_bin_params
Definition: TEfficiency.h:48
void ResetBit(UInt_t f)
Definition: TObject.h:158
Int_t GetDimension() const
returns the dimension of the current TEfficiency object
virtual void SetTitle(const char *title)
Change (i.e.
Definition: TH1.cxx:6028
virtual bool Minimize(int maxIter, double absTol=1.E-8, double relTol=1.E-10)
Find minimum position iterating until convergence specified by the absolute and relative tolerance or...
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content.
Definition: TH2.cxx:2471
virtual Int_t GetNbinsX() const
Definition: TH1.h:277
Double_t * GetEXhigh() const
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition: TObject.cxx:901
Int_t GetNbins() const
Definition: TAxis.h:121
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:364
virtual void Set(Int_t n)
Set number of points in the graph Existing coordinates are preserved New coordinates above fNpoints a...
Definition: TGraph.cxx:2105
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:155
THist< 2, float, THistStatContent, THistStatUncertainty > TH2F
Definition: THist.hxx:317
const Bool_t kTRUE
Definition: RtypesCore.h:91
double norm(double *x, double *p)
Definition: unuranDistr.cxx:40
virtual void SetNormFactor(Double_t factor=1)
Definition: TH1.h:381
Double_t GetXmax() const
Definition: TAxis.h:134
const Int_t n
Definition: legend1.C:16
Line Attributes class.
Definition: TAttLine.h:18
virtual void SetStats(Bool_t stats=kTRUE)
Set statistics option on/off.
Definition: TH1.cxx:8103
TEfficiency()
default constructor
const TArrayD * GetXbins() const
Definition: TAxis.h:130
static Double_t BetaMode(Double_t alpha, Double_t beta)
Compute the mode of the beta distribution.
void Fill(Bool_t bPassed, Double_t x, Double_t y=0, Double_t z=0)
This function is used for filling the two histograms.
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:859
TAxis * GetXaxis()
Definition: TH1.h:300
static Double_t Wilson(Double_t total, Double_t passed, Double_t level, Bool_t bUpper)
Calculates the boundaries for the frequentist Wilson interval.
static Bool_t CheckConsistency(const TH1 &pass, const TH1 &total, Option_t *opt="")
Checks the consistence of the given histograms.
virtual Int_t GetNbinsY() const
Definition: TH1.h:278
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
void Copy(TAttFill &attfill) const
Copy this fill attributes to a new TAttFill.
Definition: TAttFill.cxx:200
double beta_quantile(double x, double a, double b)
Inverse ( ) of the cumulative distribution function of the upper tail of the beta distribution (beta_...
const char * Data() const
Definition: TString.h:347
void Calculate(const double X, const double n)
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Definition: TObject.cxx:657
virtual void SetPointError(Double_t exl, Double_t exh, Double_t eyl, Double_t eyh)
Set ex and ey values for point pointed by the mouse.