Skip to content

Commit

Permalink
add ROOT cross-check
Browse files Browse the repository at this point in the history
  • Loading branch information
lukasheinrich committed Nov 4, 2020
1 parent 9e29825 commit a8c2f29
Show file tree
Hide file tree
Showing 2 changed files with 201 additions and 0 deletions.
55 changes: 55 additions & 0 deletions validation/StandardHypoTestDemo.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,55 @@
import ROOT

def StandardHypoTestDemo(
infile,
ntoys,
workspaceName = "combined",
modelSBName = "ModelConfig",
dataName = "obsData",
):

file = ROOT.TFile.Open(infile);
w = file.Get(workspaceName);
sbModel = w.obj(modelSBName);
data = w.data(dataName);

bModel = sbModel.Clone();
bModel.SetName(modelSBName+"B_only");
var = (bModel.GetParametersOfInterest().first());
oldval = var.getVal();
var.setVal(0);
bModel.SetSnapshot( ROOT.RooArgSet(var) );
var.setVal(oldval);

var = sbModel.GetParametersOfInterest().first();
oldval = var.getVal();
sbModel.SetSnapshot( ROOT.RooArgSet(var) );

profll = ROOT.RooStats.ProfileLikelihoodTestStat(bModel.GetPdf());
profll.SetOneSidedDiscovery(False);
profll.SetOneSided(True);
hypoCalc = ROOT.RooStats.FrequentistCalculator(data, bModel, sbModel);
# profll.SetPrintLevel(2)


print('by hand',profll.Evaluate(data,ROOT.RooArgSet(var) ))

hypoCalc.SetToys(ntoys, ntoys)

sampler = hypoCalc.GetTestStatSampler();
sampler.SetTestStatistic(profll);

paramPointNull = ROOT.RooArgSet(sbModel.GetParametersOfInterest())

htr = hypoCalc.GetHypoTest();
htr.SetPValueIsRightTail(True);
htr.SetBackgroundAsAlt(True);
htr.Print();
plot = ROOT.RooStats.HypoTestPlot(htr,100);
plot.SetLogYaxis(True);
plot.Draw();
ROOT.gPad.SaveAs("plot.png");

import sys
if __name__ == '__main__':
StandardHypoTestDemo(sys.argv[1],int(sys.argv[2]) if len(sys.argv) > 2 else 2000)
146 changes: 146 additions & 0 deletions validation/xmlimport_input_bkg.json
Original file line number Diff line number Diff line change
@@ -0,0 +1,146 @@
{
"channels": [
{
"name": "channel1",
"samples": [
{
"data": [
20.0,
10.0
],
"modifiers": [
{
"data": {
"hi": 1.05,
"lo": 0.95
},
"name": "syst1",
"type": "normsys"
},
{
"data": null,
"name": "SigXsecOverSM",
"type": "normfactor"
}
],
"name": "signal"
},
{
"data": [
100.0,
0.0
],
"modifiers": [
{
"data": null,
"name": "lumi",
"type": "lumi"
},
{
"data": [
5.000000074505806,
0.0
],
"name": "staterror_channel1",
"type": "staterror"
},
{
"data": {
"hi": 1.05,
"lo": 0.95
},
"name": "syst2",
"type": "normsys"
}
],
"name": "background1"
},
{
"data": [
0.0,
100.0
],
"modifiers": [
{
"data": null,
"name": "lumi",
"type": "lumi"
},
{
"data": [
0.0,
10.0
],
"name": "staterror_channel1",
"type": "staterror"
},
{
"data": {
"hi": 1.05,
"lo": 0.95
},
"name": "syst3",
"type": "normsys"
}
],
"name": "background2"
}
]
}
],
"measurements": [
{
"config": {
"parameters": [
{
"auxdata": [
1.0
],
"bounds": [
[
0.5,
1.5
]
],
"fixed": true,
"inits": [
1.0
],
"name": "lumi",
"sigmas": [
0.1
]
},
{
"bounds": [
[
0.0,
3.0
]
],
"inits": [
1.0
],
"name": "SigXsecOverSM"
},
{
"fixed": true,
"name": "syst1"
}
],
"poi": "SigXsecOverSM"
},
"name": "GaussExample"
}
],
"observations": [
{
"data": [
100.0,
100.0
],
"name": "channel1"
}
],
"version": "1.0.0"
}

0 comments on commit a8c2f29

Please sign in to comment.