Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
HybridCalculatorOriginal.cxx
Go to the documentation of this file.
1// @(#)root/roostats:$Id$
2
3/*************************************************************************
4 * Project: RooStats *
5 * Package: RooFit/RooStats *
6 * Authors: *
7 * Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke *
8 * Other author of this class: Danilo Piparo *
9 *************************************************************************
10 * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
11 * All rights reserved. *
12 * *
13 * For the licensing terms see $ROOTSYS/LICENSE. *
14 * For the list of contributors see $ROOTSYS/README/CREDITS. *
15 *************************************************************************/
16
17/** \class RooStats::HybridCalculatorOriginal
18 \ingroup Roostats
19
20
21HybridCalculatorOriginal class. This class is deprecated and it is replaced by
22the HybridCalculator.
23
24This is a fresh rewrite in RooStats of `RooStatsCms/LimitCalculator` developed
25by D. Piparo and G. Schott - Universitaet Karlsruhe
26
27The class is born from the need to have an implementation of the CLs
28method that could take advantage from the RooFit Package.
29The basic idea is the following:
30
31 - Instantiate an object specifying a signal+background model, a background model and a dataset.
32 - Perform toy MC experiments to know the distributions of -2lnQ
33 - Calculate the CLsb and CLs values as "integrals" of these distributions.
34
35The class allows the user to input models as RooAbsPdf ( TH1 object could be used
36by using the RooHistPdf class)
37
38The pdfs must be "extended": for more information please refer to
39http://roofit.sourceforge.net). The dataset can be entered as a
40RooAbsData objects.
41
42Unlike the TLimit Class a complete MC generation is performed at each step
43and not a simple Poisson fluctuation of the contents of the bins.
44Another innovation is the treatment of the nuisance parameters. The user
45can input in the constructor nuisance parameters.
46To include the information that we have about the nuisance parameters a prior
47PDF (RooAbsPdf) should be specified
48
49Different test statistic can be used (likelihood ratio, number of events or
50profile likelihood ratio. The default is the likelihood ratio.
51See the method SetTestStatistic.
52
53The number of toys to be generated is controlled by SetNumberOfToys(n).
54
55The result of the calculations is returned as a HybridResult object pointer.
56
57see also the following interesting references:
58
59 - Alex Read, "Presentation of search results: the CLs technique",
60 Journal of Physics G: Nucl. Part. Phys. 28 2693-2704 (2002).
61 see http://www.iop.org/EJ/abstract/0954-3899/28/10/313/
62
63 - Alex Read, "Modified Frequentist Analysis of Search Results (The CLs Method)" CERN 2000-005 (30 May 2000)
64
65 - V. Bartsch, G.Quast, "Expected signal observability at future experiments" CMS NOTE 2005/004
66
67 - TLimit
68*/
69
71
73
74#include "RooDataHist.h"
75#include "RooDataSet.h"
76#include "RooGlobalFunc.h"
77#include "RooNLLVar.h"
78#include "RooRealVar.h"
79#include "RooAbsData.h"
80#include "RooWorkspace.h"
81
82#include "TH1.h"
83
84using namespace std;
85
87
88using namespace RooStats;
89
90////////////////////////////////////////////////////////////////////////////////
91/// constructor with name and title
92
95 fSbModel(0),
96 fBModel(0),
97 fObservables(0),
98 fNuisanceParameters(0),
99 fPriorPdf(0),
100 fData(0),
101 fGenerateBinned(false),
102 fUsePriorPdf(false), fTmpDoExtended(true)
103{
104 // set default parameters
106 SetNumberOfToys(1000);
107}
108
109////////////////////////////////////////////////////////////////////////////////
110/// HybridCalculatorOriginal constructor without specifying a data set
111/// the user need to specify the models in the S+B case and B-only case,
112/// the list of observables of the model(s) (for MC-generation), the list of parameters
113/// that are marginalised and the prior distribution of those parameters
114
116 RooAbsPdf& bModel,
117 RooArgList& observables,
118 const RooArgSet* nuisance_parameters,
119 RooAbsPdf* priorPdf ,
120 bool GenerateBinned,
121 int testStatistics,
122 int numToys) :
123 fSbModel(&sbModel),
124 fBModel(&bModel),
125 fNuisanceParameters(nuisance_parameters),
126 fPriorPdf(priorPdf),
127 fData(0),
128 fGenerateBinned(GenerateBinned),
129 fUsePriorPdf(false),
130 fTmpDoExtended(true)
131{
132
133 // observables are managed by the class (they are copied in)
134 fObservables = new RooArgList(observables);
135 //Try to recover the information from the pdf's
136 //fObservables=new RooArgList("fObservables");
137 //fNuisanceParameters=new RooArgSet("fNuisanceParameters");
138 // if (priorPdf){
139
140
141 SetTestStatistic(testStatistics);
142 SetNumberOfToys(numToys);
143
144 if (priorPdf) UseNuisance(true);
145
146 // this->Print();
147 /* if ( _verbose ) */ //this->PrintMore("v"); /// TO DO: add the verbose mode
148}
149
150////////////////////////////////////////////////////////////////////////////////
151/// HybridCalculatorOriginal constructor for performing hypotesis test
152/// the user need to specify the data set, the models in the S+B case and B-only case.
153/// In case of treatment of nuisance parameter, the user need to specify the
154/// the list of parameters that are marginalised and the prior distribution of those parameters
155
157 RooAbsPdf& sbModel,
158 RooAbsPdf& bModel,
159 const RooArgSet* nuisance_parameters,
160 RooAbsPdf* priorPdf,
161 bool GenerateBinned,
162 int testStatistics,
163 int numToys) :
164 fSbModel(&sbModel),
165 fBModel(&bModel),
166 fObservables(0),
167 fNuisanceParameters(nuisance_parameters),
168 fPriorPdf(priorPdf),
169 fData(&data),
170 fGenerateBinned(GenerateBinned),
171 fUsePriorPdf(false),
172 fTmpDoExtended(true)
173{
174
175
176 SetTestStatistic(testStatistics);
177 SetNumberOfToys(numToys);
178
179 if (priorPdf) UseNuisance(true);
180}
181
182////////////////////////////////////////////////////////////////////////////////
183/// Constructor with a ModelConfig object representing the signal + background model and
184/// another model config representig the background only model
185/// a Prior pdf for the nuiscane parameter of the signal and background can be specified in
186/// the s+b model or the b model. If it is specified in the s+b model, the one of the s+b model will be used
187
189 const ModelConfig& sbModel,
190 const ModelConfig& bModel,
191 bool GenerateBinned,
192 int testStatistics,
193 int numToys) :
194 fSbModel(sbModel.GetPdf()),
195 fBModel(bModel.GetPdf()),
196 fObservables(0), // no need to set them - can be taken from the data
197 fNuisanceParameters((sbModel.GetNuisanceParameters()) ? sbModel.GetNuisanceParameters() : bModel.GetNuisanceParameters()),
198 fPriorPdf((sbModel.GetPriorPdf()) ? sbModel.GetPriorPdf() : bModel.GetPriorPdf()),
199 fData(&data),
200 fGenerateBinned(GenerateBinned),
201 fUsePriorPdf(false),
202 fTmpDoExtended(true)
203{
204
205 if (fPriorPdf) UseNuisance(true);
206
207 SetTestStatistic(testStatistics);
208 SetNumberOfToys(numToys);
209}
210
211////////////////////////////////////////////////////////////////////////////////
212/// HybridCalculatorOriginal destructor
213
215{
216 if (fObservables) delete fObservables;
217}
218
219////////////////////////////////////////////////////////////////////////////////
220/// Set the model describing the null hypothesis
221
223{
224 fBModel = model.GetPdf();
225 // only if it has not been set before
226 if (!fPriorPdf) fPriorPdf = model.GetPriorPdf();
228}
229
230////////////////////////////////////////////////////////////////////////////////
231/// Set the model describing the alternate hypothesis
232
234{
235 fSbModel = model.GetPdf();
236 fPriorPdf = model.GetPriorPdf();
238}
239
240////////////////////////////////////////////////////////////////////////////////
241/// set the desired test statistics:
242/// - index=1 : likelihood ratio: 2 * log( L_sb / L_b ) (DEFAULT)
243/// - index=2 : number of generated events
244/// - index=3 : profiled likelihood ratio
245/// if the index is different to any of those values, the default is used
246
248{
249 fTestStatisticsIdx = index;
250}
251
252////////////////////////////////////////////////////////////////////////////////
253/// first compute the test statistics for data and then prepare and run the toy-MC experiments
254
255HybridResult* HybridCalculatorOriginal::Calculate(TH1& data, unsigned int nToys, bool usePriors) const
256{
257
258 /// convert data TH1 histogram to a RooDataHist
259 auto dataHistName = std::string(GetName()) + "_roodatahist";
260 RooDataHist dataHist(dataHistName,"Data distribution as RooDataHist converted from TH1",*fObservables,&data);
261
262 HybridResult* result = Calculate(dataHist,nToys,usePriors);
263
264 return result;
265}
266
267////////////////////////////////////////////////////////////////////////////////
268/// first compute the test statistics for data and then prepare and run the toy-MC experiments
269
270HybridResult* HybridCalculatorOriginal::Calculate(RooAbsData& data, unsigned int nToys, bool usePriors) const
271{
272
273 double testStatData = 0;
274 if ( fTestStatisticsIdx==2 ) {
275 /// number of events used as test statistics
276 double nEvents = data.sumEntries();
277 testStatData = nEvents;
278 } else if ( fTestStatisticsIdx==3 ) {
279 /// profiled likelihood ratio used as test statistics
280 if ( fTmpDoExtended ) {
281 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,data,RooFit::CloneData(false),RooFit::Extended());
283 double sb_nll_val = sb_nll.getVal();
284 RooNLLVar b_nll("b_nll","b_nll",*fBModel,data,RooFit::CloneData(false),RooFit::Extended());
286 double b_nll_val = b_nll.getVal();
287 double m2lnQ = 2*(sb_nll_val-b_nll_val);
288 testStatData = m2lnQ;
289 } else {
290 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,data,RooFit::CloneData(false));
292 double sb_nll_val = sb_nll.getVal();
293 RooNLLVar b_nll("b_nll","b_nll",*fBModel,data,RooFit::CloneData(false));
295 double b_nll_val = b_nll.getVal();
296 double m2lnQ = 2*(sb_nll_val-b_nll_val);
297 testStatData = m2lnQ;
298 }
299 } else if ( fTestStatisticsIdx==1 ) {
300 /// likelihood ratio used as test statistics (default)
301 if ( fTmpDoExtended ) {
302 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,data,RooFit::Extended());
303 RooNLLVar b_nll("b_nll","b_nll",*fBModel,data,RooFit::Extended());
304 double m2lnQ = 2*(sb_nll.getVal()-b_nll.getVal());
305 testStatData = m2lnQ;
306 } else {
307 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,data);
308 RooNLLVar b_nll("b_nll","b_nll",*fBModel,data);
309 double m2lnQ = 2*(sb_nll.getVal()-b_nll.getVal());
310 testStatData = m2lnQ;
311 }
312 }
313
314 std::cout << "Test statistics has been evaluated for data\n";
315
316 HybridResult* result = Calculate(nToys,usePriors);
317 result->SetDataTestStatistics(testStatData);
318
319 return result;
320}
321
322////////////////////////////////////////////////////////////////////////////////
323
324HybridResult* HybridCalculatorOriginal::Calculate(unsigned int nToys, bool usePriors) const
325{
326 std::vector<double> bVals;
327 bVals.reserve(nToys);
328
329 std::vector<double> sbVals;
330 sbVals.reserve(nToys);
331
332 RunToys(bVals,sbVals,nToys,usePriors);
333
334 HybridResult* result;
335
336 TString name = "HybridResult_" + TString(GetName() );
337
338 if ( fTestStatisticsIdx==2 )
339 result = new HybridResult(name,sbVals,bVals,false);
340 else
341 result = new HybridResult(name,sbVals,bVals);
342
343 return result;
344}
345
346////////////////////////////////////////////////////////////////////////////////
347/// do the actual run-MC processing
348
349void HybridCalculatorOriginal::RunToys(std::vector<double>& bVals, std::vector<double>& sbVals, unsigned int nToys, bool usePriors) const
350{
351 std::cout << "HybridCalculatorOriginal: run " << nToys << " toy-MC experiments\n";
352 std::cout << "with test statistics index: " << fTestStatisticsIdx << "\n";
353 if (usePriors) std::cout << "marginalize nuisance parameters \n";
354
355 assert(nToys > 0);
356 assert(fBModel);
357 assert(fSbModel);
358 if (usePriors) {
359 assert(fPriorPdf);
360 assert(fNuisanceParameters);
361 }
362
363 std::vector<double> parameterValues; /// array to hold the initial parameter values
364 /// backup the initial values of the parameters that are varied by the prior MC-integration
365 int nParameters = (fNuisanceParameters) ? fNuisanceParameters->getSize() : 0;
366 RooArgList parametersList("parametersList"); /// transforms the RooArgSet in a RooArgList (needed for .at())
367 if (usePriors && nParameters>0) {
368 parametersList.add(*fNuisanceParameters);
369 parameterValues.resize(nParameters);
370 for (int iParameter=0; iParameter<nParameters; iParameter++) {
371 RooRealVar* oneParam = (RooRealVar*) parametersList.at(iParameter);
372 parameterValues[iParameter] = oneParam->getVal();
373 }
374 }
375
376 // create a cloned list of all parameters need in case of test statistics 3 where those
377 // changed by the best fit
378 RooArgSet originalSbParams;
379 RooArgSet originalBParams;
380 if (fTestStatisticsIdx == 3) {
383 if (sbparams) originalSbParams.addClone(*sbparams);
384 if (bparams) originalBParams.addClone(*bparams);
385 delete sbparams;
386 delete bparams;
387// originalSbParams.Print("V");
388// originalBParams.Print("V");
389 }
390
391
392 for (unsigned int iToy=0; iToy<nToys; iToy++) {
393
394 /// prints a progress report every 500 iterations
395 /// TO DO: add a global verbose flag
396 if ( /*verbose && */ iToy%500==0 ) {
397 std::cout << "....... toy number " << iToy << " / " << nToys << std::endl;
398 }
399
400 /// vary the value of the integrated parameters according to the prior pdf
401 if (usePriors && nParameters>0) {
402 /// generation from the prior pdf (TO DO: RooMCStudy could be used here)
404 for (int iParameter=0; iParameter<nParameters; iParameter++) {
405 RooRealVar* oneParam = (RooRealVar*) parametersList.at(iParameter);
406 oneParam->setVal(tmpValues->get()->getRealValue(oneParam->GetName()));
407 }
408 delete tmpValues;
409 }
410
411
412 /// generate the dataset in the B-only hypothesis
413 RooAbsData* bData;
414 if (fGenerateBinned)
416 else {
417 if ( fTmpDoExtended ) bData = static_cast<RooAbsData*> (fBModel->generate(*fObservables,RooFit::Extended()));
418 else bData = static_cast<RooAbsData*> (fBModel->generate(*fObservables,1));
419 }
420
421 /// work-around in case of an empty dataset (TO DO: need a debug in RooFit?)
422 bool bIsEmpty = false;
423 if (bData==NULL) {
424 bIsEmpty = true;
425 // if ( _verbose ) std::cout << "empty B-only dataset!\n";
426 RooDataSet* bDataDummy=new RooDataSet("bDataDummy","empty dataset",*fObservables);
427 bData = static_cast<RooAbsData*>(new RooDataHist ("bDataEmpty","",*fObservables,*bDataDummy));
428 delete bDataDummy;
429 }
430
431 /// generate the dataset in the S+B hypothesis
432 RooAbsData* sbData;
433 if (fGenerateBinned)
434 sbData = static_cast<RooAbsData*> (fSbModel->generateBinned(*fObservables,RooFit::Extended()));
435 else {
436 if ( fTmpDoExtended ) sbData = static_cast<RooAbsData*> (fSbModel->generate(*fObservables,RooFit::Extended()));
437 else sbData = static_cast<RooAbsData*> (fSbModel->generate(*fObservables,1));
438 }
439
440 /// work-around in case of an empty dataset (TO DO: need a debug in RooFit?)
441 bool sbIsEmpty = false;
442 if (sbData==NULL) {
443 sbIsEmpty = true;
444 // if ( _verbose ) std::cout << "empty S+B dataset!\n";
445 RooDataSet* sbDataDummy=new RooDataSet("sbDataDummy","empty dataset",*fObservables);
446 sbData = static_cast<RooAbsData*>(new RooDataHist ("sbDataEmpty","",*fObservables,*sbDataDummy));
447 delete sbDataDummy;
448 }
449
450 /// restore the parameters to their initial values
451 if (usePriors && nParameters>0) {
452 for (int iParameter=0; iParameter<nParameters; iParameter++) {
453 RooRealVar* oneParam = (RooRealVar*) parametersList.at(iParameter);
454 oneParam->setVal(parameterValues[iParameter]);
455 }
456 }
457
458 // test first the S+B hypothesis and the the B-only hypothesis
459 for (int hypoTested=0; hypoTested<=1; hypoTested++) {
460 RooAbsData* dataToTest = sbData;
461 bool dataIsEmpty = sbIsEmpty;
462 if ( hypoTested==1 ) { dataToTest = bData; dataIsEmpty = bIsEmpty; }
463 /// evaluate the test statistic in the tested hypothesis case
464 if ( fTestStatisticsIdx==2 ) { /// number of events used as test statistics
465 double nEvents = 0;
466 if ( !dataIsEmpty ) nEvents = dataToTest->numEntries();
467 if ( hypoTested==0 ) sbVals.push_back(nEvents);
468 else bVals.push_back(nEvents);
469 } else if ( fTestStatisticsIdx==3 ) { /// profiled likelihood ratio used as test statistics
470 if ( fTmpDoExtended ) {
471 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,*dataToTest,RooFit::CloneData(false),RooFit::Extended());
473 double sb_nll_val = sb_nll.getVal();
474 RooNLLVar b_nll("b_nll","b_nll",*fBModel,*dataToTest,RooFit::CloneData(false),RooFit::Extended());
476 double b_nll_val = b_nll.getVal();
477 double m2lnQ = -2*(b_nll_val-sb_nll_val);
478 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
479 else bVals.push_back(m2lnQ);
480 } else {
481 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,*dataToTest,RooFit::CloneData(false));
483 double sb_nll_val = sb_nll.getVal();
484 RooNLLVar b_nll("b_nll","b_nll",*fBModel,*dataToTest,RooFit::CloneData(false));
486 double b_nll_val = b_nll.getVal();
487 double m2lnQ = -2*(b_nll_val-sb_nll_val);
488 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
489 else bVals.push_back(m2lnQ);
490 }
491 } else if ( fTestStatisticsIdx==1 ) { /// likelihood ratio used as test statistics (default)
492 if ( fTmpDoExtended ) {
493 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,*dataToTest,RooFit::CloneData(false),RooFit::Extended());
494 RooNLLVar b_nll("b_nll","b_nll",*fBModel,*dataToTest,RooFit::CloneData(false),RooFit::Extended());
495 double m2lnQ = -2*(b_nll.getVal()-sb_nll.getVal());
496 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
497 else bVals.push_back(m2lnQ);
498 } else {
499 RooNLLVar sb_nll("sb_nll","sb_nll",*fSbModel,*dataToTest,RooFit::CloneData(false));
500 RooNLLVar b_nll("b_nll","b_nll",*fBModel,*dataToTest,RooFit::CloneData(false));
501 double m2lnQ = -2*(b_nll.getVal()-sb_nll.getVal());
502 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
503 else bVals.push_back(m2lnQ);
504 }
505 }
506 } // tested both hypotheses
507
508 /// delete the toy-MC datasets
509 delete sbData;
510 delete bData;
511
512 /// restore the parameters to their initial values in case fitting is done
513 if (fTestStatisticsIdx == 3) {
515 if (sbparams) {
516 assert(originalSbParams.getSize() == sbparams->getSize());
517 sbparams->assign(originalSbParams);
518 delete sbparams;
519 }
521 if (bparams) {
522 assert(originalBParams.getSize() == bparams->getSize());
523 bparams->assign(originalBParams);
524 delete bparams;
525 }
526 }
527
528
529
530 } /// end of loop over toy-MC experiments
531
532
533 /// restore the parameters to their initial values
534 if (usePriors && nParameters>0) {
535 for (int iParameter=0; iParameter<nParameters; iParameter++) {
536 RooRealVar* oneParam = (RooRealVar*) parametersList.at(iParameter);
537 oneParam->setVal(parameterValues[iParameter]);
538 }
539 }
540
541 return;
542}
543
544////////////////////////////////////////////////////////////////////////////////
545/// Print out some information about the input models
546
547void HybridCalculatorOriginal::PrintMore(const char* options) const
548{
549
550 if (fSbModel) {
551 std::cout << "Signal plus background model:\n";
552 fSbModel->Print(options);
553 }
554
555 if (fBModel) {
556 std::cout << "\nBackground model:\n";
557 fBModel->Print(options);
558 }
559
560 if (fObservables) {
561 std::cout << "\nObservables:\n";
562 fObservables->Print(options);
563 }
564
566 std::cout << "\nParameters being integrated:\n";
567 fNuisanceParameters->Print(options);
568 }
569
570 if (fPriorPdf) {
571 std::cout << "\nPrior PDF model for integration:\n";
572 fPriorPdf->Print(options);
573 }
574
575 return;
576}
577
578////////////////////////////////////////////////////////////////////////////////
579/// implementation of inherited methods from HypoTestCalculator
580
582 // perform the hypothesis test and return result of hypothesis test
583
584 // check first that everything needed is there
585 if (!DoCheckInputs()) return 0;
586 RooAbsData * treeData = dynamic_cast<RooAbsData *> (fData);
587 if (!treeData) {
588 std::cerr << "Error in HybridCalculatorOriginal::GetHypoTest - invalid data type - return NULL" << std::endl;
589 return 0;
590 }
591 bool usePrior = (fUsePriorPdf && fPriorPdf );
592 return Calculate( *treeData, fNToys, usePrior);
593}
594
595
597 if (!fData) {
598 std::cerr << "Error in HybridCalculatorOriginal - data have not been set" << std::endl;
599 return false;
600 }
601
602 // if observable have not been set take them from data
603 if (!fObservables && fData->get() ) fObservables = new RooArgList( *fData->get() );
604 if (!fObservables) {
605 std::cerr << "Error in HybridCalculatorOriginal - no observables" << std::endl;
606 return false;
607 }
608
609 if (!fSbModel) {
610 std::cerr << "Error in HybridCalculatorOriginal - S+B pdf has not been set " << std::endl;
611 return false;
612 }
613
614 if (!fBModel) {
615 std::cerr << "Error in HybridCalculatorOriginal - B pdf has not been set" << std::endl;
616 return false;
617 }
619 std::cerr << "Error in HybridCalculatorOriginal - nuisance parameters have not been set " << std::endl;
620 return false;
621 }
622 if (fUsePriorPdf && !fPriorPdf) {
623 std::cerr << "Error in HybridCalculatorOriginal - prior pdf has not been set " << std::endl;
624 return false;
625 }
626 return true;
627}
#define ClassImp(name)
Definition Rtypes.h:364
char name[80]
Definition TGX11.cxx:110
virtual void Print(Option_t *options=0) const
Print the object to the defaultPrintStream().
Definition RooAbsArg.h:337
RooArgSet * getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Int_t getSize() const
virtual RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE)
Add a clone of the specified argument to list.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
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.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:82
virtual const RooArgSet * get() const
Definition RooAbsData.h:128
virtual Double_t sumEntries() const =0
Return effective number of entries in dataset, i.e., sum all weights.
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
virtual RooDataHist * generateBinned(const RooArgSet &whatVars, Double_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()) const
As RooAbsPdf::generateBinned(const RooArgSet&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,...
Definition RooAbsPdf.h:107
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.
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())
See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,...
Definition RooAbsPdf.h:58
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:94
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:35
The RooDataHist is a container class to hold N-dimensional binned data.
Definition RooDataHist.h:45
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:36
virtual const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
Class RooNLLVar implements a -log(likelihood) calculation from a dataset and a PDF.
Definition RooNLLVar.h:30
RooRealVar represents a variable that can be changed from the outside.
Definition RooRealVar.h:39
virtual void setVal(Double_t value)
Set value of variable to 'value'.
HybridCalculatorOriginal class.
void PrintMore(const char *options) const
Print out some information about the input models.
virtual void SetAlternateModel(const ModelConfig &)
Set the model describing the alternate hypothesis.
virtual ~HybridCalculatorOriginal()
Destructor of HybridCalculator.
virtual HybridResult * GetHypoTest() const
inherited methods from HypoTestCalculator interface
void RunToys(std::vector< double > &bVals, std::vector< double > &sbVals, unsigned int nToys, bool usePriors) const
do the actual run-MC processing
virtual void SetNullModel(const ModelConfig &)
Set the model describing the null hypothesis.
HybridCalculatorOriginal(const char *name=0)
Dummy Constructor with only name.
void SetTestStatistic(int index)
set the desired test statistics: index=1 : 2 * log( L_sb / L_b ) (DEFAULT) index=2 : number of genera...
HybridResult * Calculate(TH1 &data, unsigned int nToys, bool usePriors) const
first compute the test statistics for data and then prepare and run the toy-MC experiments
Class encapsulating the result of the HybridCalculatorOriginal.
void SetDataTestStatistics(double testStat_data_val)
set the value of the test statistics on data
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition ModelConfig.h:30
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
RooAbsPdf * GetPriorPdf() const
get parameters prior pdf (return NULL if not existing)
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:58
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47
Basic string class.
Definition TString.h:136
RooCmdArg Strategy(Int_t code)
RooCmdArg Extended(Bool_t flag=kTRUE)
RooCmdArg Hesse(Bool_t flag=kTRUE)
RooCmdArg CloneData(Bool_t flag)
RooCmdArg PrintLevel(Int_t code)
Namespace for the RooStats classes.
Definition Asimov.h:19