126bool useProof =
false;
131void TwoSidedFrequentistUpperLimitWithBands(
const char *infile =
"",
const char *workspaceName =
"combined",
132 const char *modelConfigName =
"ModelConfig",
133 const char *dataName =
"obsData")
136 double confidenceLevel = 0.95;
139 double additionalToysFac = 0.5;
140 int nPointsToScan = 20;
146 const char *filename =
"";
147 if (!strcmp(infile,
"")) {
148 filename =
"results/example_combined_GaussExample_model.root";
153 cout <<
"HistFactory file cannot be generated on Windows - exit" << endl;
157 cout <<
"will run standard hist2workspace example" << endl;
158 gROOT->ProcessLine(
".! prepareHistFactory .");
159 gROOT->ProcessLine(
".! hist2workspace config/example.xml");
160 cout <<
"\n\n---------------------" << endl;
161 cout <<
"Done creating example input" << endl;
162 cout <<
"---------------------\n\n" << endl;
173 cout <<
"StandardRooStatsDemoMacro: Input file " << filename <<
" is not found" << endl;
183 cout <<
"workspace not found" << endl;
196 cout <<
"data or ModelConfig was not found" << endl;
200 cout <<
"Found data and ModelConfig:" << endl;
218 fc.SetConfidenceLevel(confidenceLevel);
219 fc.AdditionalNToysFactor(additionalToysFac);
221 fc.SetNBins(nPointsToScan);
222 fc.CreateConfBelt(
true);
244 fc.FluctuateNumDataEntries(
false);
246 cout <<
"Not sure what to do about this model" << endl;
261 cout <<
"will use global observables for unconditional ensemble" << endl;
271 cout <<
"\n95% interval on " << firstPOI->
GetName() <<
" is : [" << interval->
LowerLimit(*firstPOI) <<
", "
272 << interval->
UpperLimit(*firstPOI) <<
"] " << endl;
276 double observedUL = interval->
UpperLimit(*firstPOI);
277 firstPOI->
setVal(observedUL);
278 double obsTSatObsUL =
fc.GetTestStatSampler()->EvaluateTestStatistic(*data, tmpPOI);
285 TH1F *histOfThresholds =
298 histOfThresholds->
Fill(poiVal, arMax);
304 histOfThresholds->
Draw();
319 w->
saveSnapshot(
"paramsToGenerateData", *poiAndNuisance);
321 cout <<
"\nWill use these parameter points to generate pseudo data for bkg only" << endl;
322 paramsToGenerateData->
Print(
"v");
329 double CLbinclusive = 0;
332 TH1F *histOfUL =
new TH1F(
"histOfUL",
"", 100, 0, firstPOI->
getMax());
335 for (
int imc = 0; imc < nToyMC; ++imc) {
348 cout <<
"Not sure what to do about this model" << endl;
379 firstPOI->
setVal(observedUL);
380 double toyTSatObsUL =
fc.GetTestStatSampler()->EvaluateTestStatistic(*toyData, tmpPOI);
383 if (obsTSatObsUL < toyTSatObsUL)
384 CLb += (1.) / nToyMC;
385 if (obsTSatObsUL <= toyTSatObsUL)
386 CLbinclusive += (1.) / nToyMC;
395 double thisTS =
fc.GetTestStatSampler()->EvaluateTestStatistic(*toyData, tmpPOI);
400 if (thisTS <= arMax) {
401 thisUL = firstPOI->
getVal();
407 histOfUL->
Fill(thisUL);
415 c1->SaveAs(
"two-sided_upper_limit_output.pdf");
433 double band2sigDown = 0, band1sigDown = 0, bandMedian = 0, band1sigUp = 0, band2sigUp = 0;
434 for (
int i = 1; i <= cumulative->
GetNbinsX(); ++i) {
446 cout <<
"-2 sigma band " << band2sigDown << endl;
447 cout <<
"-1 sigma band " << band1sigDown <<
" [Power Constraint)]" << endl;
448 cout <<
"median of band " << bandMedian << endl;
449 cout <<
"+1 sigma band " << band1sigUp << endl;
450 cout <<
"+2 sigma band " << band2sigUp << endl;
453 cout <<
"\nobserved 95% upper-limit " << interval->
UpperLimit(*firstPOI) << endl;
454 cout <<
"CLb strict [P(toy>obs|0)] for observed 95% upper-limit " << CLb << endl;
455 cout <<
"CLb inclusive [P(toy>=obs|0)] for observed 95% upper-limit " << CLbinclusive << endl;
R__EXTERN TSystem * gSystem
static struct mg_connection * fc(struct mg_context *ctx)
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
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.
RooAbsArg * first() const
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.
virtual Int_t numEntries() const
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
Bool_t canBeExtended() const
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&,...
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
virtual RooAbsReal * createProfile(const RooArgSet ¶msOfInterest)
Create a RooProfileLL object that eliminates all nuisance parameters in the present function.
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
RooArgSet is a container object that can hold multiple RooAbsArg objects.
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
virtual TObject * clone(const char *newname) const
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
RooDataSet is a container class to hold unbinned data.
virtual const RooArgSet * get(Int_t index) const override
Return RooArgSet with coordinates of event 'index'.
RooRealVar represents a variable that can be changed from the outside.
virtual void setVal(Double_t value)
Set value of variable to 'value'.
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of 'global observables' – for RooStats tools.
ConfidenceBelt is a concrete implementation of the ConfInterval interface.
Double_t GetAcceptanceRegionMax(RooArgSet &, Double_t cl=-1., Double_t leftside=-1.)
The FeldmanCousins class (like the Feldman-Cousins technique) is essentially a specific configuration...
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
const RooArgSet * GetGlobalObservables() const
get RooArgSet for global observables (return NULL if not existing)
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return NULL if not existing)
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
virtual void Print(Option_t *option="") const override
overload the print method
const RooArgSet * GetObservables() const
get RooArgSet for observables (return NULL if not existing)
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
PointSetInterval is a concrete implementation of the ConfInterval interface.
Double_t UpperLimit(RooRealVar ¶m)
return upper limit on a given parameter
Double_t LowerLimit(RooRealVar ¶m)
return lower limit on a given parameter
ProfileLikelihoodTestStat is an implementation of the TestStatistic interface that calculates the pro...
Holds configuration options for proof and proof-lite.
ToyMCSampler is an implementation of the TestStatSampler interface.
void SetProofConfig(ProofConfig *pc=NULL)
virtual TestStatistic * GetTestStatistic(unsigned int i) const
virtual void SetGlobalObservables(const RooArgSet &o)
The RooWorkspace is a persistable container for RooFit projects.
RooAbsData * data(const char *name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found.
void Print(Option_t *opts=0) const
Print contents of the workspace.
Bool_t saveSnapshot(const char *name, const char *paramNames)
Save snapshot of values and attributes (including "Constant") of given parameters.
Bool_t loadSnapshot(const char *name)
Load the values and attributes of the parameters in the snapshot saved with the given name.
TObject * obj(const char *name) const
Return any type of object (RooAbsArg, RooAbsData or generic object) with given name)
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
1-D histogram with a float per channel (see TH1 documentation)}
virtual Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
TObject * Clone(const char *newname=0) const
Make a complete copy of the underlying object.
virtual Int_t GetNbinsX() const
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
virtual void SetContent(const Double_t *content)
Replace bin contents by the contents of array content.
virtual void SetMinimum(Double_t minimum=-1111)
virtual void Draw(Option_t *option="")
Draw this histogram with options.
virtual Double_t * GetIntegral()
Return a pointer to the array of bins integral.
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
virtual const char * GetName() const
Returns name of object.
virtual Bool_t AccessPathName(const char *path, EAccessMode mode=kFileExists)
Returns FALSE if one can access a file using the specified access mode.
RooCmdArg Extended(Bool_t flag=kTRUE)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Namespace for the RooStats classes.
Double_t SignificanceToPValue(Double_t Z)
returns p-value corresponding to a 1-sided significance
static constexpr double pc