Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
HypoTestInverterPlot.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2// Author: Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke
3/*************************************************************************
4 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
5 * All rights reserved. *
6 * *
7 * For the licensing terms see $ROOTSYS/LICENSE. *
8 * For the list of contributors see $ROOTSYS/README/CREDITS. *
9 *************************************************************************/
10
11/** \class RooStats::HypoTestInverterPlot
12 \ingroup Roostats
13
14Class to plot a HypoTestInverterResult, the output of the HypoTestInverter calculator.
15
16It can be used to plot the obtained p-values ( CLb, CLs+b or CLs) for each scanned point, as well as
17the test statistic distributions (when a calculator based on pseudo-experiments is used) for the two
18hypotheses.
19
20*/
21
22// include header file of this class
24
25// include other header files
30
31#include "TGraphErrors.h"
32#include "TGraphAsymmErrors.h"
33#include "TMultiGraph.h"
34#include "TROOT.h"
35#include "TLine.h"
36#include "TAxis.h"
37#include "TLegend.h"
38#include "TVirtualPad.h"
39#include "TError.h"
41
42#include <cmath>
43
44using namespace std;
45
47
48using namespace RooStats;
49
50////////////////////////////////////////////////////////////////////////////////
51/// constructor from a HypoTestInverterResult class
52/// name and title are taken from the result class
53
55 TNamed( results->GetName(), results->GetTitle() ),
56 fResults(results)
57{
58}
59
60////////////////////////////////////////////////////////////////////////////////
61/// constructor with name and title from a HypoTestInverterResult class
62
64 const char* title,
65 HypoTestInverterResult* results ) :
66 TNamed( TString(name), TString(title) ),
67 fResults(results)
68{
69}
70
71////////////////////////////////////////////////////////////////////////////////
72/// Make the plot of the result of the scan using the observed data.
73/// By default plot CLs or CLsb depending if the flag UseCLs is set for the results
74/// that are passed to this instance.
75///
76/// \param opt Options according to following list:
77/// - Empty: Return CLs or CLs+b depending on the value of UseCLs.ƒ
78/// - "CLB": return CLb plot
79/// - "CLS+B" / "CLSPLUSB": return CLs+b plot independently of the flag
80/// - "CLS": return CLs plot independently of the flag
81
83{
84 TString option(opt);
85 option.ToUpper();
86 enum Mode_t { Default, CLb, CLsPlusb, CLs };
87 Mode_t type = Default;
88 if (option.Contains("CLB")) type = CLb;
89 else if (option.Contains("CLS+B") || option.Contains("CLSPLUSB")) type = CLsPlusb;
90 else if (option.Contains("CLS" )) type = CLs;
91
92 const int nEntries = fResults->ArraySize();
93
94 // sort the arrays based on the x values
95 std::vector<unsigned int> index(nEntries);
96 TMath::SortItr(fResults->fXValues.begin(), fResults->fXValues.end(), index.begin(), false);
97
98 // copy result in sorted arrays
99 std::vector<double> xArray;
100 std::vector<double> yArray;
101 std::vector<double> yErrArray;
102
103 for (int i=0; i<nEntries; i++) {
104 double CLVal = 0.;
105 double CLErr = 0.;
106 if (type == Default) {
107 CLVal = fResults->GetYValue(index[i]);
108 CLErr = fResults->GetYError(index[i]);
109 } else if (type == CLb) {
110 CLVal = fResults->CLb(index[i]);
111 CLErr = fResults->CLbError(index[i]);
112 } else if (type == CLsPlusb) {
113 CLVal = fResults->CLsplusb(index[i]);
114 CLErr = fResults->CLsplusbError(index[i]);
115 } else if (type == CLs) {
116 CLVal = fResults->CLs(index[i]);
117 CLErr = fResults->CLsError(index[i]);
118 }
119
120 if (CLVal < 0. || !std::isfinite(CLVal)) {
121 Warning("HypoTestInverterPlot::MakePlot", "Got a confidence level of %f at x=%f (failed fit?). Skipping this point.", CLVal, fResults->GetXValue(index[i]));
122 continue;
123 }
124
125 yArray.push_back(CLVal);
126 yErrArray.push_back(CLErr);
127 xArray.push_back(fResults->GetXValue(index[i]));
128 }
129
130 TGraphErrors* graph = new TGraphErrors(static_cast<Int_t>(xArray.size()),&xArray.front(),&yArray.front(),nullptr,&yErrArray.front());
131 TString pValueName = "CLs";
132 if (type == CLb ) pValueName = "CLb";
133 if (type == CLsPlusb || (type == Default && !fResults->fUseCLs) ) pValueName = "CLs+b";
134 TString name = pValueName + TString("_observed");
135 TString title = TString("Observed ") + pValueName;
136 graph->SetName(name);
137 graph->SetTitle(title);
138 graph->SetMarkerStyle(20);
139 graph->SetLineWidth(2);
140 return graph;
141}
142
143////////////////////////////////////////////////////////////////////////////////
144/// Make the expected plot and the bands
145/// nsig1 and nsig2 indicates the n-sigma value for the bands
146/// if nsig1 = 0 no band is drawn (only expected value)
147/// if nsig2 > nsig1 (default is nsig1=1 and nsig2=2) the second band is also drawn
148/// The first band is drawn in green while the second in yellow
149/// THe return result is a TMultiGraph object
150
152{
153 const int nEntries = fResults->ArraySize();
154 bool doFirstBand = (nsig1 > 0);
155 bool doSecondBand = (nsig2 > nsig1);
156
157 nsig1 = std::abs(nsig1);
158 nsig2 = std::abs(nsig2);
159
160 // sort the arrays based on the x values
161 std::vector<unsigned int> index(nEntries);
162 TMath::SortItr(fResults->fXValues.begin(), fResults->fXValues.end(), index.begin(), false);
163
164 // create the graphs
165 TGraph * g0 = new TGraph;
166 TString pValueName = "CLs";
167 if (!fResults->fUseCLs) pValueName = "CLs+b";
168 g0->SetTitle(TString::Format("Expected %s - Median",pValueName.Data()) );
169 TGraphAsymmErrors * g1 = nullptr;
170 TGraphAsymmErrors * g2 = nullptr;
171 if (doFirstBand) {
172 g1 = new TGraphAsymmErrors;
173 if (nsig1 - int(nsig1) < 0.01) {
174 g1->SetTitle(TString::Format("Expected %s #pm %d #sigma",pValueName.Data(),int(nsig1)) );
175 } else {
176 g1->SetTitle(TString::Format("Expected %s #pm %3.1f #sigma", pValueName.Data(), nsig1));
177 }
178 }
179 if (doSecondBand) {
180 g2 = new TGraphAsymmErrors;
181 if (nsig2 - int(nsig2) < 0.01) {
182 g2->SetTitle(TString::Format("Expected %s #pm %d #sigma",pValueName.Data(),int(nsig2)) );
183 } else {
184 g2->SetTitle(TString::Format("Expected %s #pm %3.1f #sigma", pValueName.Data(), nsig2));
185 }
186 }
187 double p[7];
188 double q[7];
189 p[0] = ROOT::Math::normal_cdf(-nsig2);
190 p[1] = ROOT::Math::normal_cdf(-nsig1);
191 p[2] = 0.5;
192 p[3] = ROOT::Math::normal_cdf(nsig1);
193 p[4] = ROOT::Math::normal_cdf(nsig2);
194
195 bool resultIsAsymptotic = ( !fResults->GetNullTestStatDist(0) && !fResults->GetAltTestStatDist(0) );
196 int np = 0;
197 for (int j=0; j<nEntries; ++j) {
198 int i = index[j]; // i is the order index
200 if ( !s) continue;
201 const std::vector<double> & values = s->GetSamplingDistribution();
202 // special case for asymptotic results (cannot use TMath::quantile in that case)
203 if (resultIsAsymptotic) {
204 double maxSigma = fResults->fgAsymptoticMaxSigma;
205 double dsig = 2* maxSigma/ (values.size() -1) ;
206 int i0 = (int) TMath::Floor ( ( -nsig2 + maxSigma )/dsig + 0.5);
207 int i1 = (int) TMath::Floor ( (-nsig1 + maxSigma )/dsig + 0.5);
208 int i2 = (int) TMath::Floor ( ( maxSigma)/dsig + 0.5);
209 int i3 = (int) TMath::Floor ( ( nsig1 + maxSigma)/dsig + 0.5);
210 int i4 = (int) TMath::Floor ( ( nsig2 + maxSigma)/dsig + 0.5);
211 q[0] = values[i0];
212 q[1] = values[i1];
213 q[2] = values[i2];
214 q[3] = values[i3];
215 q[4] = values[i4];
216 }
217 else {
218 double * x = const_cast<double *>(&values[0]); // need to change TMath::Quantiles
219 TMath::Quantiles(values.size(), 5, x,q,p,false);
220 }
221
222 g0->SetPoint(np, fResults->GetXValue(i), q[2]);
223 if (g1) {
224 g1->SetPoint(np, fResults->GetXValue(i), q[2]);
225 g1->SetPointEYlow(np, q[2] - q[1]); // -1 sigma error
226 g1->SetPointEYhigh(np, q[3] - q[2]);//+1 sigma error
227 }
228 if (g2) {
229 g2->SetPoint(np, fResults->GetXValue(i), q[2]);
230
231 g2->SetPointEYlow(np, q[2]-q[0]); // -2 -- -1 sigma error
232 g2->SetPointEYhigh(np, q[4]-q[2]);
233 }
234 if (s) delete s;
235 np++;
236 }
237
238 TString name = GetName() + TString("_expected");
239 TString title = TString("Expected ") + GetTitle();
240 TMultiGraph* graph = new TMultiGraph(name,title);
241
242 // set the graphics options and add in multi graph
243 // orderof adding is drawing order
244 if (g2) {
246 graph->Add(g2,"3");
247 }
248 if (g1) {
249 g1->SetFillColor(kGreen);
250 graph->Add(g1,"3");
251 }
252 g0->SetLineStyle(2);
253 g0->SetLineWidth(2);
254 graph->Add(g0,"L");
255
256 return graph;
257}
258
259////////////////////////////////////////////////////////////////////////////////
260/// destructor
261
263{
264}
265
266////////////////////////////////////////////////////////////////////////////////
267/// Draw the result in the current canvas
268/// Possible options:
269/// - SAME : draw in the current axis
270/// - OBS : draw only the observed plot
271/// - EXP : draw only the expected plot
272///
273/// - CLB : draw also the CLB
274/// - 2CL : draw both clsplusb and cls
275///
276/// default draw observed + expected with 1 and 2 sigma bands
277
279
280 TString option(opt);
281 option.ToUpper();
282 bool drawAxis = !option.Contains("SAME");
283 bool drawObs = option.Contains("OBS") || !option.Contains("EXP");
284 bool drawExp = option.Contains("EXP") || !option.Contains("OBS");
285 bool drawCLb = option.Contains("CLB");
286 bool draw2CL = option.Contains("2CL");
287
288 TGraphErrors * gobs = nullptr;
289 TGraph * gplot = nullptr;
290 if (drawObs) {
291 gobs = MakePlot();
292 // add object to top-level directory to avoid mem leak
293 if (gROOT) gROOT->Add(gobs);
294 if (drawAxis) {
295 gobs->Draw("APL");
296 gplot = gobs;
297 gplot->GetHistogram()->SetTitle( GetTitle() );
298 }
299 else gobs->Draw("PL");
300
301 }
302 TMultiGraph * gexp = nullptr;
303 if (drawExp) {
304 gexp = MakeExpectedPlot();
305 // add object to current directory to avoid mem leak
306 if (gROOT) gROOT->Add(gexp);
307 if (drawAxis && !drawObs) {
308 gexp->Draw("A");
309 if (gexp->GetHistogram()) gexp->GetHistogram()->SetTitle( GetTitle() );
310 gplot = static_cast<TGraph*>(gexp->GetListOfGraphs()->First());
311 }
312 else
313 gexp->Draw();
314
315 }
316
317 // draw also an horizontal line at the desired conf level
318 if (gplot) {
319 double alpha = 1.-fResults->ConfidenceLevel();
320 double x1 = gplot->GetXaxis()->GetXmin();
321 double x2 = gplot->GetXaxis()->GetXmax();
322 TLine * line = new TLine(x1, alpha, x2,alpha);
324 line->Draw();
325 // put axis labels
327 if (arg) gplot->GetXaxis()->SetTitle(arg->GetName());
328 gplot->GetYaxis()->SetTitle("p value");
329 }
330
331
332 TGraph *gclb = nullptr;
333 if (drawCLb) {
334 gclb = MakePlot("CLb");
335 if (gROOT) gROOT->Add(gclb);
336 gclb->SetMarkerColor(kBlue+4);
337 gclb->Draw("PL");
338 // draw in red observed cls or clsb
339 if (gobs) gobs->SetMarkerColor(kRed);
340 }
341 TGraph * gclsb = nullptr;
342 TGraph * gcls = nullptr;
343 if (draw2CL) {
344 if (fResults->fUseCLs) {
345 gclsb = MakePlot("CLs+b");
346 if (gROOT) gROOT->Add(gclsb);
347 gclsb->SetMarkerColor(kBlue);
348 gclsb->Draw("PL");
349 gclsb->SetLineStyle(3);
350 }
351 else {
352 gcls = MakePlot("CLs");
353 if (gROOT) gROOT->Add(gcls);
354 gcls->SetMarkerColor(kBlue);
355 gcls->Draw("PL");
356 gcls->SetLineStyle(3);
357 }
358 }
359 // draw again observed values otherwise will be covered by the bands
360 if (gobs) {
361 gobs->Draw("PL");
362 }
363
364
365 double y0 = 0.6;
366 double verticalSize = (gexp || draw2CL || drawCLb ) ? 0.3 : 0.15;
367 double y1 = y0 + verticalSize;
368 TLegend * l = new TLegend(0.6,y0,0.9,y1);
369 if (gobs) l->AddEntry(gobs,"","PEL");
370 if (gclsb) l->AddEntry(gclsb,"","PEL");
371 if (gcls) l->AddEntry(gcls,"","PEL");
372 if (gclb) l->AddEntry(gclb,"","PEL");
373 if (gexp) {
374 // loop in reverse order (opposite to drawing one)
375 int ngraphs = gexp->GetListOfGraphs()->GetSize();
376 for (int i = ngraphs-1; i>=0; --i) {
377 TObject * obj = gexp->GetListOfGraphs()->At(i);
378 TString lopt = "F";
379 if (i == ngraphs-1) lopt = "L";
380 if (obj) l->AddEntry(obj,"",lopt);
381 }
382 }
383 l->Draw();
384 // redraw the axis
385 if (gPad) gPad->RedrawAxis();
386
387}
388
389////////////////////////////////////////////////////////////////////////////////
390/// \param index Index of the result stored in HypoTestInverterResult
391/// \param type Type of the test (see below)
392/// \param nbins Number of bins
393/// - type =0 null and alt
394/// - type = 1 only null (S+B)
395/// - type = 2 only alt (B)
396
398 SamplingDistPlot * pl = nullptr;
399 if (type == 0) {
401 if (result)
402 pl = new HypoTestPlot(*result, nbins );
403 return pl;
404 }
405 if (type == 1) {
407 if (sbDist) {
408 pl = new SamplingDistPlot( nbins);
409 pl->AddSamplingDistribution(sbDist);
410 return pl;
411 }
412 }
413 if (type == 2) {
415 if (bDist) {
416 pl = new SamplingDistPlot( nbins);
417 pl->AddSamplingDistribution(bDist);
418 return pl;
419 }
420 }
421 return nullptr;
422}
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:377
@ kRed
Definition Rtypes.h:66
@ kGreen
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
@ kYellow
Definition Rtypes.h:66
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t np
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
Option_t Option_t TPoint TPoint const char x2
Option_t Option_t TPoint TPoint const char x1
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t type
Option_t Option_t TPoint TPoint const char y1
char name[80]
Definition TGX11.cxx:110
float * q
#define gROOT
Definition TROOT.h:406
#define gPad
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:77
RooAbsArg * first() const
Class to plot a HypoTestInverterResult, the output of the HypoTestInverter calculator.
HypoTestInverterResult * fResults
void Draw(Option_t *opt="") override
Draw the scan result in the current canvas Possible options: "" (default): draw observed + expected w...
SamplingDistPlot * MakeTestStatPlot(int index, int type=0, int nbins=100)
Plot the test statistic distributions.
HypoTestInverterPlot(HypoTestInverterResult *results)
constructor
TMultiGraph * MakeExpectedPlot(double sig1=1, double sig2=2)
Make the expected plot and the bands nsig1 and nsig2 indicates the n-sigma value for the bands if nsi...
TGraphErrors * MakePlot(Option_t *opt="")
return a TGraphErrors with the obtained observed p-values resultinf from the scan By default (Option ...
HypoTestInverterResult class holds the array of hypothesis test results and compute a confidence inte...
double GetYValue(int index) const
function to return the value of the confidence level for the i^th entry in the results
SamplingDistribution * GetSignalAndBackgroundTestStatDist(int index) const
get the signal and background test statistic distribution
int ArraySize() const
number of entries in the results array
double CLsplusbError(int index) const
return the observed CLsplusb value for the i-th entry
double CLsError(int index) const
return the observed CLb value for the i-th entry
double CLbError(int index) const
return the observed CLb value for the i-th entry
double GetXValue(int index) const
function to return the value of the parameter of interest for the i^th entry in the results
static double fgAsymptoticMaxSigma
max sigma value used to scan asymptotic expected p values
double CLs(int index) const
return the observed CLb value for the i-th entry
double GetYError(int index) const
function to return the estimated error on the value of the confidence level for the i^th entry in the...
double CLb(int index) const
return the observed CLb value for the i-th entry
TList fYObjects
list of HypoTestResult for each point
SamplingDistribution * GetExpectedPValueDist(int index) const
return expected distribution of p-values (Cls or Clsplusb)
double CLsplusb(int index) const
return the observed CLsplusb value for the i-th entry
SamplingDistribution * GetNullTestStatDist(int index) const
same in terms of alt and null
SamplingDistribution * GetAltTestStatDist(int index) const
SamplingDistribution * GetBackgroundTestStatDist(int index) const
get the background test statistic distribution
This class provides the plots for the result of a study performed with any of the HypoTestCalculatorG...
HypoTestResult is a base class for results from hypothesis tests.
This class provides simple and straightforward utilities to plot SamplingDistribution objects.
double AddSamplingDistribution(const SamplingDistribution *samplingDist, Option_t *drawOptions="NORMALIZE HIST")
adds the sampling distribution and returns the scale factor
This class simply holds a sampling distribution of some test statistic.
const std::vector< double > & GetSamplingDistribution() const
Get test statistics values.
RooArgSet fParameters
set containing the parameter of interest
double ConfidenceLevel() const override
return the confidence interval
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition TAttFill.h:37
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition TAttLine.h:42
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition TAttLine.h:43
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:40
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition TAttMarker.h:38
Double_t GetXmax() const
Definition TAxis.h:140
Double_t GetXmin() const
Definition TAxis.h:139
virtual Int_t GetSize() const
Return the capacity of the collection, i.e.
TGraph with asymmetric error bars.
virtual void SetPointEYlow(Int_t i, Double_t eyl)
Set EYlow for point i.
virtual void SetPointEYhigh(Int_t i, Double_t eyh)
Set EYhigh for point i.
A TGraphErrors is a TGraph with error bars.
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
virtual void SetPoint(Int_t i, Double_t x, Double_t y)
Set x and y values for point number i.
Definition TGraph.cxx:2319
void Draw(Option_t *chopt="") override
Draw this graph with its current attributes.
Definition TGraph.cxx:809
TAxis * GetXaxis() const
Get x axis of the graph.
Definition TGraph.cxx:1544
TAxis * GetYaxis() const
Get y axis of the graph.
Definition TGraph.cxx:1553
virtual TH1F * GetHistogram() const
Returns a pointer to the histogram used to draw the axis Takes into account the two following cases.
Definition TGraph.cxx:1406
void SetTitle(const char *title="") override
Change (i.e.
Definition TGraph.cxx:2374
void SetTitle(const char *title) override
Change/set the title.
Definition TH1.cxx:6682
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
Use the TLine constructor to create a simple line.
Definition TLine.h:22
TObject * First() const override
Return the first object in the list. Returns 0 when list is empty.
Definition TList.cxx:659
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:357
A TMultiGraph is a collection of TGraph (or derived) objects.
Definition TMultiGraph.h:34
TList * GetListOfGraphs() const
Definition TMultiGraph.h:68
TH1F * GetHistogram()
Returns a pointer to the histogram used to draw the axis.
void Draw(Option_t *chopt="") override
Draw this multigraph with its current attributes.
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
Mother of all ROOT objects.
Definition TObject.h:41
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:962
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:274
Basic string class.
Definition TString.h:139
const char * Data() const
Definition TString.h:378
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
Definition TString.cxx:2378
TLine * line
double normal_cdf(double x, double sigma=1, double x0=0)
Cumulative distribution function of the normal (Gaussian) distribution (lower tail).
Double_t x[n]
Definition legend1.C:17
Namespace for the RooStats classes.
Definition Asimov.h:19
void SortItr(Iterator first, Iterator last, IndexIterator index, Bool_t down=kTRUE)
Sort the n1 elements of the Short_t array defined by its iterators.
Definition TMathBase.h:406
Double_t Floor(Double_t x)
Rounds x downward, returning the largest integral value that is not greater than x.
Definition TMath.h:680
void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob, Bool_t isSorted=kTRUE, Int_t *index=nullptr, Int_t type=7)
Computes sample quantiles, corresponding to the given probabilities.
Definition TMath.cxx:1207
Definition graph.py:1
TLine l
Definition textangle.C:4