37 using namespace RooStats;
41 HybridCalculatorOriginal::HybridCalculatorOriginal(const
char *name) :
46 fNuisanceParameters(0),
49 fGenerateBinned(
false),
50 fUsePriorPdf(false), fTmpDoExtended(true)
55 SetNumberOfToys(1000);
60 HybridCalculatorOriginal::HybridCalculatorOriginal(
RooAbsPdf& sbModel,
70 fNuisanceParameters(nuisance_parameters),
73 fGenerateBinned(GenerateBinned),
111 fNuisanceParameters(nuisance_parameters),
114 fGenerateBinned(GenerateBinned),
133 const ModelConfig& sbModel,
134 const ModelConfig& bModel,
138 fSbModel(sbModel.GetPdf()),
139 fBModel(bModel.GetPdf()),
141 fNuisanceParameters((sbModel.GetNuisanceParameters()) ? sbModel.GetNuisanceParameters() : bModel.GetNuisanceParameters()),
142 fPriorPdf((sbModel.GetPriorPdf()) ? sbModel.GetPriorPdf() : bModel.GetPriorPdf()),
144 fGenerateBinned(GenerateBinned),
217 double testStatData = 0;
227 double sb_nll_val = sb_nll.getVal();
230 double b_nll_val = b_nll.getVal();
231 double m2lnQ = 2*(sb_nll_val-b_nll_val);
232 testStatData = m2lnQ;
236 double sb_nll_val = sb_nll.getVal();
239 double b_nll_val = b_nll.getVal();
240 double m2lnQ = 2*(sb_nll_val-b_nll_val);
241 testStatData = m2lnQ;
248 double m2lnQ = 2*(sb_nll.getVal()-b_nll.getVal());
249 testStatData = m2lnQ;
253 double m2lnQ = 2*(sb_nll.getVal()-b_nll.getVal());
254 testStatData = m2lnQ;
258 std::cout <<
"Test statistics has been evaluated for data\n";
260 HybridResult* result =
Calculate(nToys,usePriors);
261 result->SetDataTestStatistics(testStatData);
270 std::vector<double> bVals;
271 bVals.reserve(nToys);
273 std::vector<double> sbVals;
274 sbVals.reserve(nToys);
276 RunToys(bVals,sbVals,nToys,usePriors);
283 result =
new HybridResult(name,sbVals,bVals,
false);
285 result =
new HybridResult(name,sbVals,bVals);
295 std::cout <<
"HybridCalculatorOriginal: run " << nToys <<
" toy-MC experiments\n";
297 if (usePriors) std::cout <<
"marginalize nuisance parameters \n";
307 std::vector<double> parameterValues;
311 if (usePriors && nParameters>0) {
313 parameterValues.resize(nParameters);
314 for (
int iParameter=0; iParameter<nParameters; iParameter++) {
316 parameterValues[iParameter] = oneParam->
getVal();
324 if (fTestStatisticsIdx == 3) {
327 if (sbparams) originalSbParams.
addClone(*sbparams);
328 if (bparams) originalBParams.
addClone(*bparams);
336 for (
unsigned int iToy=0; iToy<
nToys; iToy++) {
341 std::cout <<
"....... toy number " << iToy <<
" / " << nToys << std::endl;
345 if (usePriors && nParameters>0) {
348 for (
int iParameter=0; iParameter<nParameters; iParameter++) {
366 bool bIsEmpty =
false;
385 bool sbIsEmpty =
false;
395 if (usePriors && nParameters>0) {
396 for (
int iParameter=0; iParameter<nParameters; iParameter++) {
398 oneParam->
setVal(parameterValues[iParameter]);
403 for (
int hypoTested=0; hypoTested<=1; hypoTested++) {
405 bool dataIsEmpty = sbIsEmpty;
406 if ( hypoTested==1 ) { dataToTest = bData; dataIsEmpty = bIsEmpty; }
408 if ( fTestStatisticsIdx==2 ) {
410 if ( !dataIsEmpty ) nEvents = dataToTest->
numEntries();
411 if ( hypoTested==0 ) sbVals.push_back(nEvents);
412 else bVals.push_back(nEvents);
413 }
else if ( fTestStatisticsIdx==3 ) {
417 double sb_nll_val = sb_nll.
getVal();
420 double b_nll_val = b_nll.
getVal();
421 double m2lnQ = -2*(b_nll_val-sb_nll_val);
422 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
423 else bVals.push_back(m2lnQ);
427 double sb_nll_val = sb_nll.
getVal();
430 double b_nll_val = b_nll.
getVal();
431 double m2lnQ = -2*(b_nll_val-sb_nll_val);
432 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
433 else bVals.push_back(m2lnQ);
435 }
else if ( fTestStatisticsIdx==1 ) {
440 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
441 else bVals.push_back(m2lnQ);
446 if ( hypoTested==0 ) sbVals.push_back(m2lnQ);
447 else bVals.push_back(m2lnQ);
457 if (fTestStatisticsIdx == 3) {
461 *sbparams = originalSbParams;
467 *bparams = originalBParams;
478 if (usePriors && nParameters>0) {
479 for (
int iParameter=0; iParameter<nParameters; iParameter++) {
481 oneParam->
setVal(parameterValues[iParameter]);
495 std::cout <<
"Signal plus background model:\n";
500 std::cout <<
"\nBackground model:\n";
505 std::cout <<
"\nObservables:\n";
510 std::cout <<
"\nParameters being integrated:\n";
515 std::cout <<
"\nPrior PDF model for integration:\n";
531 std::cerr <<
"Error in HybridCalculatorOriginal::GetHypoTest - invalid data type - return NULL" << std::endl;
541 std::cerr <<
"Error in HybridCalculatorOriginal - data have not been set" << std::endl;
548 std::cerr <<
"Error in HybridCalculatorOriginal - no observables" << std::endl;
553 std::cerr <<
"Error in HybridCalculatorOriginal - S+B pdf has not been set " << std::endl;
558 std::cerr <<
"Error in HybridCalculatorOriginal - B pdf has not been set" << std::endl;
562 std::cerr <<
"Error in HybridCalculatorOriginal - nuisance parameters have not been set " << std::endl;
566 std::cerr <<
"Error in HybridCalculatorOriginal - prior pdf has not been set " << std::endl;
virtual Double_t sumEntries() const =0
HybridCalculatorOriginal(const char *name=0)
Dummy Constructor with only name.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
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())
Generate a new dataset containing the specified variables with events sampled from our distribution...
RooCmdArg CloneData(Bool_t flag)
RooCmdArg PrintLevel(Int_t code)
RooCmdArg Strategy(Int_t code)
const RooArgSet * fNuisanceParameters
unsigned int fTestStatisticsIdx
HybridCalculatorOriginal class.
ClassImp(TIterator) Bool_t TIterator return false
Compare two iterator objects.
RooCmdArg Extended(Bool_t flag=kTRUE)
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
void RunToys(std::vector< double > &bVals, std::vector< double > &sbVals, unsigned int nToys, bool usePriors) const
void SetNumberOfToys(unsigned int ntoys)
RooDataSet is a container class to hold N-dimensional binned data.
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
The TNamed class is the base class for all named ROOT classes.
virtual void SetAlternateModel(const ModelConfig &)
virtual void SetNullModel(const ModelConfig &)
Double_t getVal(const RooArgSet *set=0) const
virtual HybridResult * GetHypoTest() const
inherited methods from HypoTestCalculator interface
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
RooRealVar represents a fundamental (non-derived) real valued object.
virtual const RooArgSet * get() const
bool DoCheckInputs() const
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Class RooNLLVar implements a a -log(likelihood) calculation from a dataset and a PDF.
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
virtual Int_t numEntries() const
void UseNuisance(bool on=true)
HybridResult * Calculate(TH1 &data, unsigned int nToys, bool usePriors) const
RooArgList * fObservables
virtual const char * GetName() const
Returns name of object.
RooAbsData is the common abstract base class for binned and unbinned datasets.
RooDataSet is a container class to hold unbinned data.
RooCmdArg Hesse(Bool_t flag=kTRUE)
RooAbsArg * at(Int_t idx) const
void PrintMore(const char *options) const
ClassImp(RooStats::HybridCalculatorOriginal) using namespace RooStats
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
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...
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.
virtual RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE)
Add clone of specified element to an owning set.
RooAbsPdf * GetPriorPdf() const
get parameters prior pdf (return NULL if not existing)
virtual const RooArgSet * get(Int_t index) const
Return RooArgSet with coordinates of event 'index'.
void SetTestStatistic(int index)
set the desired test statistics: index=1 : 2 * log( L_sb / L_b ) (DEFAULT) index=2 : number of genera...
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.
virtual ~HybridCalculatorOriginal()
Destructor of HybridCalculator.