ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
rs401d_FeldmanCousins.C
Go to the documentation of this file.
1 /////////////////////////////////////////////////////////////////////////
2 //
3 // 'Neutrino Oscillation Example from Feldman & Cousins'
4 // author: Kyle Cranmer
5 // date March 2009
6 //
7 // This tutorial shows a more complex example using the FeldmanCousins utility
8 // to create a confidence interval for a toy neutrino oscillation experiment.
9 // The example attempts to faithfully reproduce the toy example described in Feldman & Cousins'
10 // original paper, Phys.Rev.D57:3873-3889,1998.
11 //
12 // to run it:
13 // .x tutorials/roostats/rs401d_FeldmanCousins.C+
14 /////////////////////////////////////////////////////////////////////////
15 
16 #include "RooGlobalFunc.h"
17 #include "RooStats/ConfInterval.h"
24 #include "RooStats/MCMCInterval.h"
25 
26 #include "RooDataSet.h"
27 #include "RooDataHist.h"
28 #include "RooRealVar.h"
29 #include "RooConstVar.h"
30 #include "RooAddition.h"
31 #include "RooProduct.h"
32 #include "RooProdPdf.h"
33 #include "RooAddPdf.h"
34 
35 #include "TROOT.h"
36 #include "RooPolynomial.h"
37 #include "RooRandom.h"
38 
39 #include "RooNLLVar.h"
40 #include "RooProfileLL.h"
41 
42 #include "RooPlot.h"
43 
44 #include "TCanvas.h"
45 #include "TH1F.h"
46 #include "TH2F.h"
47 #include "TTree.h"
48 #include "TMarker.h"
49 #include "TStopwatch.h"
50 
51 #include <iostream>
52 
53 // PDF class created for this macro
54 #if !defined(__CINT__) || defined(__MAKECINT__)
55 #include "../tutorials/roostats/NuMuToNuE_Oscillation.h"
56 #include "../tutorials/roostats/NuMuToNuE_Oscillation.cxx" // so that it can be executed directly
57 #else
58 #include "../tutorials/roostats/NuMuToNuE_Oscillation.cxx+" // so that it can be executed directly
59 #endif
60 
61 // use this order for safety on library loading
62 using namespace RooFit ;
63 using namespace RooStats ;
64 
65 
66 void rs401d_FeldmanCousins(bool doFeldmanCousins=false, bool doMCMC = true)
67 {
68 
69  // to time the macro
70  TStopwatch t;
71  t.Start();
72 
73 
74  /*
75  Taken from Feldman & Cousins paper, Phys.Rev.D57:3873-3889,1998.
76  e-Print: physics/9711021 (see page 13.)
77 
78  Quantum mechanics dictates that the probability of such a transformation is given by the formula
79  P (νµ → ν e ) = sin^2 (2θ) sin^2 (1.27 ∆m^2 L /E )
80  where P is the probability for a νµ to transform into a νe , L is the distance in km between
81  the creation of the neutrino from meson decay and its interaction in the detector, E is the
82  neutrino energy in GeV, and ∆m^2 = |m^2− m^2 | in (eV/c^2 )^2 .
83 
84  To demonstrate how this works in practice, and how it compares to alternative approaches
85  that have been used, we consider a toy model of a typical neutrino oscillation experiment.
86  The toy model is defined by the following parameters: Mesons are assumed to decay to
87  neutrinos uniformly in a region 600 m to 1000 m from the detector. The expected background
88  from conventional νe interactions and misidentified νµ interactions is assumed to be 100
89  events in each of 5 energy bins which span the region from 10 to 60 GeV. We assume that
90  the νµ flux is such that if P (νµ → ν e ) = 0.01 averaged over any bin, then that bin would
91  have an expected additional contribution of 100 events due to νµ → ν e oscillations.
92  */
93 
94  // Make signal model model
95  RooRealVar E("E","", 15,10,60,"GeV");
96  RooRealVar L("L","", .800,.600, 1.0,"km"); // need these units in formula
97  RooRealVar deltaMSq("deltaMSq","#Delta m^{2}",40,1,300,"eV/c^{2}");
98  RooRealVar sinSq2theta("sinSq2theta","sin^{2}(2#theta)", .006,.0,.02);
99  //RooRealVar deltaMSq("deltaMSq","#Delta m^{2}",40,20,70,"eV/c^{2}");
100  // RooRealVar sinSq2theta("sinSq2theta","sin^{2}(2#theta)", .006,.001,.01);
101  // PDF for oscillation only describes deltaMSq dependence, sinSq2theta goes into sigNorm
102  // 1) The code for this PDF was created by issuing these commands
103  // root [0] RooClassFactory x
104  // root [1] x.makePdf("NuMuToNuE_Oscillation","L,E,deltaMSq","","pow(sin(1.27*deltaMSq*L/E),2)")
105  NuMuToNuE_Oscillation PnmuTone("PnmuTone","P(#nu_{#mu} #rightarrow #nu_{e}",L,E,deltaMSq);
106 
107  // only E is observable, so create the signal model by integrating out L
108  RooAbsPdf* sigModel = PnmuTone.createProjection(L);
109 
110  // create \int dE' dL' P(E',L' | \Delta m^2).
111  // Given RooFit will renormalize the PDF in the range of the observables,
112  // the average probability to oscillate in the experiment's acceptance
113  // needs to be incorporated into the extended term in the likelihood.
114  // Do this by creating a RooAbsReal representing the integral and divide by
115  // the area in the E-L plane.
116  // The integral should be over "primed" observables, so we need
117  // an independent copy of PnmuTone not to interfere with the original.
118 
119  // Independent copy for Integral
120  RooRealVar EPrime("EPrime","", 15,10,60,"GeV");
121  RooRealVar LPrime("LPrime","", .800,.600, 1.0,"km"); // need these units in formula
122  NuMuToNuE_Oscillation PnmuTonePrime("PnmuTonePrime","P(#nu_{#mu} #rightarrow #nu_{e}",
123  LPrime,EPrime,deltaMSq);
124  RooAbsReal* intProbToOscInExp = PnmuTonePrime.createIntegral(RooArgSet(EPrime,LPrime));
125 
126  // Getting the flux is a bit tricky. It is more celear to include a cross section term that is not
127  // explicitly refered to in the text, eg.
128  // # events in bin = flux * cross-section for nu_e interaction in E bin * average prob nu_mu osc. to nu_e in bin
129  // let maxEventsInBin = flux * cross-section for nu_e interaction in E bin
130  // maxEventsInBin * 1% chance per bin = 100 events / bin
131  // therefore maxEventsInBin = 10,000.
132  // for 5 bins, this means maxEventsTot = 50,000
133  RooConstVar maxEventsTot("maxEventsTot","maximum number of sinal events",50000);
134  RooConstVar inverseArea("inverseArea","1/(#Delta E #Delta L)",
135  1./(EPrime.getMax()-EPrime.getMin())/(LPrime.getMax()-LPrime.getMin()));
136 
137  // sigNorm = maxEventsTot * (\int dE dL prob to oscillate in experiment / Area) * sin^2(2\theta)
138  RooProduct sigNorm("sigNorm", "", RooArgSet(maxEventsTot, *intProbToOscInExp, inverseArea, sinSq2theta));
139  // bkg = 5 bins * 100 events / bin
140  RooConstVar bkgNorm("bkgNorm","normalization for background",500);
141 
142  // flat background (0th order polynomial, so no arguments for coefficients)
143  RooPolynomial bkgEShape("bkgEShape","flat bkg shape", E);
144 
145  // total model
146  RooAddPdf model("model","",RooArgList(*sigModel,bkgEShape),
147  RooArgList(sigNorm,bkgNorm));
148 
149  // for debugging, check model tree
150  // model.printCompactTree();
151  // model.graphVizTree("model.dot");
152 
153 
154  // turn off some messages
158 
159 
160  //////////////////////////////////////////////
161  // n events in data to data, simply sum of sig+bkg
162  Int_t nEventsData = bkgNorm.getVal()+sigNorm.getVal();
163  cout << "generate toy data with nEvents = " << nEventsData << endl;
164  // adjust random seed to get a toy dataset similar to one in paper.
165  // Found by trial and error (3 trials, so not very "fine tuned")
167  // create a toy dataset
168  RooDataSet* data = model.generate(RooArgSet(E), nEventsData);
169 
170  /////////////////////////////////////////////
171  // make some plots
172  TCanvas* dataCanvas = new TCanvas("dataCanvas");
173  dataCanvas->Divide(2,2);
174 
175  // plot the PDF
176  dataCanvas->cd(1);
177  TH1* hh = PnmuTone.createHistogram("hh",E,Binning(40),YVar(L,Binning(40)),Scaling(kFALSE)) ;
178  hh->SetLineColor(kBlue) ;
179  hh->SetTitle("True Signal Model");
180  hh->Draw("surf");
181 
182  // plot the data with the best fit
183  dataCanvas->cd(2);
184  RooPlot* Eframe = E.frame();
185  data->plotOn(Eframe);
186  model.fitTo(*data, Extended());
187  model.plotOn(Eframe);
188  model.plotOn(Eframe,Components(*sigModel),LineColor(kRed));
189  model.plotOn(Eframe,Components(bkgEShape),LineColor(kGreen));
190  model.plotOn(Eframe);
191  Eframe->SetTitle("toy data with best fit model (and sig+bkg components)");
192  Eframe->Draw();
193 
194  // plot the likelihood function
195  dataCanvas->cd(3);
196  RooNLLVar nll("nll", "nll", model, *data, Extended());
197  RooProfileLL pll("pll", "", nll, RooArgSet(deltaMSq, sinSq2theta));
198  // TH1* hhh = nll.createHistogram("hhh",sinSq2theta,Binning(40),YVar(deltaMSq,Binning(40))) ;
199  TH1* hhh = pll.createHistogram("hhh",sinSq2theta,Binning(40),YVar(deltaMSq,Binning(40)),Scaling(kFALSE)) ;
200  hhh->SetLineColor(kBlue) ;
201  hhh->SetTitle("Likelihood Function");
202  hhh->Draw("surf");
203 
204  dataCanvas->Update();
205 
206 
207 
208  //////////////////////////////////////////////////////////
209  //////// show use of Feldman-Cousins utility in RooStats
210  // set the distribution creator, which encodes the test statistic
211  RooArgSet parameters(deltaMSq, sinSq2theta);
212  RooWorkspace* w = new RooWorkspace();
213 
214  ModelConfig modelConfig;
215  modelConfig.SetWorkspace(*w);
216  modelConfig.SetPdf(model);
217  modelConfig.SetParametersOfInterest(parameters);
218 
219  RooStats::FeldmanCousins fc(*data, modelConfig);
220  fc.SetTestSize(.1); // set size of test
221  fc.UseAdaptiveSampling(true);
222  fc.SetNBins(10); // number of points to test per parameter
223 
224  // use the Feldman-Cousins tool
225  ConfInterval* interval = 0;
226  if(doFeldmanCousins)
227  interval = fc.GetInterval();
228 
229 
230  ///////////////////////////////////////////////////////////////////
231  ///////// show use of ProfileLikeihoodCalculator utility in RooStats
232  RooStats::ProfileLikelihoodCalculator plc(*data, modelConfig);
233  plc.SetTestSize(.1);
234 
235  ConfInterval* plcInterval = plc.GetInterval();
236 
237  ///////////////////////////////////////////////////////////////////
238  ///////// show use of MCMCCalculator utility in RooStats
239  MCMCInterval* mcInt = NULL;
240 
241  if (doMCMC) {
242  // turn some messages back on
245 
246  TStopwatch mcmcWatch;
247  mcmcWatch.Start();
248 
249  RooArgList axisList(deltaMSq, sinSq2theta);
250  MCMCCalculator mc(*data, modelConfig);
251  mc.SetNumIters(5000);
252  mc.SetNumBurnInSteps(100);
253  mc.SetUseKeys(true);
254  mc.SetTestSize(.1);
255  mc.SetAxes(axisList); // set which is x and y axis in posterior histogram
256  //mc.SetNumBins(50);
257  mcInt = (MCMCInterval*)mc.GetInterval();
258 
259  mcmcWatch.Stop();
260  mcmcWatch.Print();
261  }
262  ////////////////////////////////////////////
263  // make plot of resulting interval
264 
265  dataCanvas->cd(4);
266 
267  // first plot a small dot for every point tested
268  if (doFeldmanCousins) {
269  RooDataHist* parameterScan = (RooDataHist*) fc.GetPointsToScan();
270  TH2F* hist = (TH2F*) parameterScan->createHistogram("sinSq2theta:deltaMSq",30,30);
271  // hist->Draw();
272  TH2F* forContour = (TH2F*)hist->Clone();
273 
274  // now loop through the points and put a marker if it's in the interval
275  RooArgSet* tmpPoint;
276  // loop over points to test
277  for(Int_t i=0; i<parameterScan->numEntries(); ++i){
278  // get a parameter point from the list of points to test.
279  tmpPoint = (RooArgSet*) parameterScan->get(i)->clone("temp");
280 
281  if (interval){
282  if (interval->IsInInterval( *tmpPoint ) ) {
283  forContour->SetBinContent( hist->FindBin(tmpPoint->getRealValue("sinSq2theta"),
284  tmpPoint->getRealValue("deltaMSq")), 1);
285  }else{
286  forContour->SetBinContent( hist->FindBin(tmpPoint->getRealValue("sinSq2theta"),
287  tmpPoint->getRealValue("deltaMSq")), 0);
288  }
289  }
290 
291 
292  delete tmpPoint;
293  }
294 
295  if (interval){
296  Double_t level=0.5;
297  forContour->SetContour(1,&level);
298  forContour->SetLineWidth(2);
299  forContour->SetLineColor(kRed);
300  forContour->Draw("cont2,same");
301  }
302  }
303 
304  MCMCIntervalPlot* mcPlot = NULL;
305  if (mcInt) {
306  cout << "MCMC actual confidence level: "
307  << mcInt->GetActualConfidenceLevel() << endl;
308  mcPlot = new MCMCIntervalPlot(*mcInt);
309  mcPlot->SetLineColor(kMagenta);
310  mcPlot->Draw();
311  }
312  dataCanvas->Update();
313 
314  LikelihoodIntervalPlot plotInt((LikelihoodInterval*)plcInterval);
315  plotInt.SetTitle("90% Confidence Intervals");
316  if (mcInt)
317  plotInt.Draw("same");
318  else
319  plotInt.Draw();
320  dataCanvas->Update();
321 
322  /// print timing info
323  t.Stop();
324  t.Print();
325 
326 
327 }
328 
void Print(Option_t *option="") const
Print the real and cpu time passed between the start and stop events.
Definition: TStopwatch.cxx:217
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:52
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Definition: TStopwatch.cxx:56
void Draw(const Option_t *options=0)
draw the likelihood interval or contour for the 1D case a RooPlot is drawn by default of the profiled...
This class provides simple and straightforward utilities to plot a MCMCInterval object.
Class RooProfileLL implements the profile likelihood estimator for a given likelihood and set of para...
Definition: RooProfileLL.h:26
RooCmdArg LineColor(Color_t color)
LikelihoodInterval is a concrete implementation of the RooStats::ConfInterval interface.
virtual void SetWorkspace(RooWorkspace &ws)
Definition: ModelConfig.h:88
RooArgList L(const RooAbsArg &v1)
Definition: Rtypes.h:61
virtual void SetTestSize(Double_t size)
set the size of the test (rate of Type I error) ( Eg. 0.05 for a 95% Confidence Interval) ...
virtual void SetTestSize(Double_t size)
set the size of the test (rate of Type I error) ( Eg. 0.05 for a 95% Confidence Interval) ...
This class provides simple and straightforward utilities to plot a LikelihoodInterval object...
ProfileLikelihoodCalculator is a concrete implementation of CombinedCalculator (the interface class f...
virtual void SetNumBurnInSteps(Int_t numBurnInSteps)
set the number of steps in the chain to discard as burn-in, starting from the first ...
virtual Bool_t IsInInterval(const RooArgSet &) const =0
check if given point is in the interval
virtual Double_t getMin(const char *name=0) const
int Int_t
Definition: RtypesCore.h:41
void Draw(const Option_t *options=NULL)
const Bool_t kFALSE
Definition: Rtypes.h:92
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition: TObject.cxx:254
void SetTitle(const char *name)
Set the title of the RooPlot to 'title'.
Definition: RooPlot.cxx:1099
Definition: Rtypes.h:61
static RooMsgService & instance()
Return reference to singleton instance.
virtual void SetSeed(UInt_t seed=0)
Set the random generator seed.
Definition: TRandom.cxx:568
RooCmdArg Extended(Bool_t flag=kTRUE)
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
Definition: TH1.cxx:2565
virtual const RooArgSet * get() const
Definition: RooDataHist.h:77
virtual LikelihoodInterval * GetInterval() const
Return a likelihood interval.
RooAbsData * GetPointsToScan()
void setStreamStatus(Int_t id, Bool_t active)
(De)Activate stream with given unique ID
static struct mg_connection * fc(struct mg_context *ctx)
Definition: civetweb.c:839
RooDataSet is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:40
void Stop()
Stop the stopwatch.
Definition: TStopwatch.cxx:75
virtual void SetPdf(const RooAbsPdf &pdf)
Set the Pdf, add to the the workspace if not already there.
Definition: ModelConfig.h:97
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition: RooRandom.cxx:54
RooConstVar represent a constant real-valued object.
Definition: RooConstVar.h:25
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Plot dataset on specified frame.
Definition: RooAbsData.cxx:626
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Return a p.d.f that represent a projection of this p.d.f integrated over given observables.
Definition: RooAbsPdf.cxx:2989
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
virtual void SetLineColor(Color_t lcolor)
Definition: TAttLine.h:54
TThread * t[5]
Definition: threadsh1.C:13
Class RooNLLVar implements a a -log(likelihood) calculation from a dataset and a PDF.
Definition: RooNLLVar.h:26
virtual void SetNumIters(Int_t numIters)
set the number of iterations to run the metropolis algorithm
virtual MCMCInterval * GetInterval() const
Main interface to get a ConfInterval.
virtual void SetTestSize(Double_t size)
set the size of the test (rate of Type I error) ( Eg. 0.05 for a 95% Confidence Interval) ...
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2878
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none(), const RooCmdArg &arg9=RooCmdArg::none(), const RooCmdArg &arg10=RooCmdArg::none()) const
Plot (project) PDF on specified frame.
Definition: RooAbsPdf.h:105
virtual void SetAxes(RooArgList &axes)
set which variables to put on each axis
void SetNBins(Int_t bins)
void SetLineColor(Color_t color)
Double_t E()
Definition: TMath.h:54
tuple w
Definition: qtexample.py:51
TH1 * createHistogram(const char *name, const RooAbsRealLValue &xvar, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create and fill a ROOT histogram TH1,TH2 or TH3 with the values of this dataset.
Definition: RooAbsData.cxx:735
RooPlot * frame(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
RooProduct a RooAbsReal implementation that represent the product of a given set of other RooAbsReal ...
Definition: RooProduct.h:32
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
Definition: RooPlot.h:41
The FeldmanCousins class (like the Feldman-Cousins technique) is essentially a specific configuration...
ConfInterval is an interface class for a generic interval in the RooStats framework.
Definition: ConfInterval.h:44
virtual TObject * clone(const char *newname) const
Definition: RooArgSet.h:82
virtual void SetUseKeys(Bool_t useKeys)
set whether to use kernel estimation to determine the interval
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooCmdArg Components(const RooArgSet &compSet)
RooCmdArg YVar(const RooAbsRealLValue &var, const RooCmdArg &arg=RooCmdArg::none())
The TH1 histogram class.
Definition: TH1.h:80
virtual PointSetInterval * GetInterval() const
Main interface to get a ConfInterval (will be a PointSetInterval)
TH1 * createHistogram(const char *varNameList, Int_t xbins=0, Int_t ybins=0, Int_t zbins=0) const
Create and fill a ROOT histogram TH1,TH2 or TH3 with the values of this function for the variables wi...
virtual Int_t numEntries() const
Return the number of bins.
RooCmdArg Scaling(Bool_t flag)
virtual Double_t getMax(const char *name=0) const
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
Generate a new dataset containing the specified variables with events sampled from our distribution...
Definition: RooAbsPdf.cxx:1702
#define NULL
Definition: Rtypes.h:82
virtual RooFitResult * fitTo(RooAbsData &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Fit PDF to given dataset.
Definition: RooAbsPdf.cxx:1056
RooPolynomial implements a polynomial p.d.f of the form By default coefficient a_0 is chosen to be 1...
Definition: RooPolynomial.h:28
virtual void SetParametersOfInterest(const RooArgSet &set)
Definition: ModelConfig.h:115
Definition: Rtypes.h:61
virtual void SetTitle(const char *title)
Change (i.e.
Definition: TH1.cxx:6268
const Bool_t kTRUE
Definition: Rtypes.h:91
Bayesian Calculator estimating an interval or a credible region using the Markov-Chain Monte Carlo me...
virtual void SetTitle(const char *title="")
Change (i.e. set) the title of the TNamed.
Definition: TNamed.cxx:152
void rs401d_FeldmanCousins(bool doFeldmanCousins=false, bool doMCMC=true)
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:42
void UseAdaptiveSampling(bool flag=true)
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:503
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:559
RooCmdArg Binning(const RooAbsBinning &binning)
Double_t getRealValue(const char *name, Double_t defVal=0, Bool_t verbose=kFALSE) const
Get value of a RooAbsReal stored in set with given name.
Definition: RooArgSet.cxx:527
Stopwatch class.
Definition: TStopwatch.h:30