rf515_hfJSON

Code HistFactory Models in JSON.

With the HS3 standard, it is possible to code RooFit-Models of any kind as JSON files. In this tutorial, you can see how to code up a (simple) HistFactory-based model in JSON and import it into a RooWorkspace.

Author: Carsten Burgard
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, November 30, 2022 at 11:23 AM.

In [1]:
import ROOT
Welcome to JupyROOT 6.27/01

start by creating an empty workspace

In [2]:
ws = ROOT.RooWorkspace("workspace")

the RooJSONFactoryWSTool is responsible for importing and exporting things to and from your workspace

In [3]:
tool = ROOT.RooJSONFactoryWSTool(ws)

use it to import the information from your JSON file

In [4]:
tool.importJSON(ROOT.gROOT.GetTutorialDir().Data() + "/roofit/rf515_hfJSON.json")
ws.Print()
[#1] INFO:ObjectHandling -- RooWorkspace::import(workspace) importing RooRealVar::gamma_stat_channel1_bin_0
[#1] INFO:ObjectHandling -- RooWorkspace::import(workspace) importing RooRealVar::gamma_stat_channel1_bin_1
[#1] INFO:ObjectHandling -- RooWorkspace::import(workspace) importing dataset observed

RooWorkspace(workspace) workspace contents

variables
---------
(alpha_syst1,alpha_syst2,alpha_syst3,channelCat,gamma_stat_channel1_bin_0,gamma_stat_channel1_bin_1,mu,nom_alpha_syst1,nom_alpha_syst2,nom_alpha_syst3,nom_gamma_stat_channel1_bin_0,nom_gamma_stat_channel1_bin_1,obs_x_channel1,sigma_alpha_syst1,sigma_alpha_syst2,sigma_alpha_syst3,weight)

p.d.f.s
-------
RooProdPdf::channel1[ gamma_stat_channel1_bin_0_constraint * gamma_stat_channel1_bin_1_constraint * syst2_constraint * syst3_constraint * syst1_constraint * channel1_model(obs_x_channel1) ] = 0.190787
RooRealSumPdf::channel1_model[ background1_norm * background1 + background2_norm * background2 + signal_norm * signal ] = 240
RooPoisson::gamma_stat_channel1_bin_0_constraint[ x=nom_gamma_stat_channel1_bin_0 mean=gamma_stat_channel1_bin_0_poisMean ] = 0.019943
RooPoisson::gamma_stat_channel1_bin_1_constraint[ x=nom_gamma_stat_channel1_bin_1 mean=gamma_stat_channel1_bin_1_poisMean ] = 0.039861
RooSimultaneous::main[ indexCat=channelCat channel1=channel1 ] = 0.190787
RooGaussian::syst1_constraint[ x=alpha_syst1 mean=nom_alpha_syst1 sigma=sigma_alpha_syst1 ] = 1
RooGaussian::syst2_constraint[ x=alpha_syst2 mean=nom_alpha_syst2 sigma=sigma_alpha_syst2 ] = 1
RooGaussian::syst3_constraint[ x=alpha_syst3 mean=nom_alpha_syst3 sigma=sigma_alpha_syst3 ] = 1

functions
--------
RooProduct::background1[ background1_binWidth * mc_stat_channel1 * hist_background1 ] = 200
RooBinWidthFunction::background1_binWidth[ HistFuncForBinWidth=hist_background1 ] = 2
RooProduct::background1_norm[ overallSys_background1 ] = 1
RooProduct::background2[ background2_binWidth * mc_stat_channel1 * hist_background2 ] = 0
RooBinWidthFunction::background2_binWidth[ HistFuncForBinWidth=hist_background2 ] = 2
RooProduct::background2_norm[ overallSys_background2 ] = 1
RooProduct::gamma_stat_channel1_bin_0_poisMean[ gamma_stat_channel1_bin_0 * gamma_stat_channel1_bin_0_tau ] = 400
RooProduct::gamma_stat_channel1_bin_1_poisMean[ gamma_stat_channel1_bin_1 * gamma_stat_channel1_bin_1_tau ] = 100
RooHistFunc::hist_background1[ depList=(obs_x_channel1) ] = 100
RooHistFunc::hist_background2[ depList=(obs_x_channel1) ] = 0
RooHistFunc::hist_signal[ depList=(obs_x_channel1) ] = 20
ParamHistFunc::mc_stat_channel1[ ] = 1
RooStats::HistFactory::FlexibleInterpVar::overallSys_background1[ paramList=(alpha_syst2) ] = 1
RooStats::HistFactory::FlexibleInterpVar::overallSys_background2[ paramList=(alpha_syst3) ] = 1
RooStats::HistFactory::FlexibleInterpVar::overallSys_signal[ paramList=(alpha_syst1) ] = 1
RooProduct::signal[ signal_binWidth * hist_signal ] = 40
RooBinWidthFunction::signal_binWidth[ HistFuncForBinWidth=hist_signal ] = 2
RooProduct::signal_norm[ mu * overallSys_signal ] = 1

datasets
--------
RooDataSet::observed(channelCat,obs_x_channel1)

embedded datasets (in pdfs and functions)
-----------------------------------------
RooDataHist::dataHist_background1(obs_x_channel1)
RooDataHist::dataHist_background2(obs_x_channel1)
RooDataHist::dataHist_signal(obs_x_channel1)

parameter snapshots
-------------------
fromJSON = (mu=1,obs_x_channel1=1[C],gamma_stat_channel1_bin_0=1 +/- 0.05,gamma_stat_channel1_bin_1=1 +/- 0.1,alpha_syst2=0,nom_alpha_syst2=0[C],sigma_alpha_syst2=1[C],alpha_syst3=0,nom_alpha_syst3=0[C],sigma_alpha_syst3=1[C],alpha_syst1=0,nom_alpha_syst1=0[C],sigma_alpha_syst1=1[C],nom_gamma_stat_channel1_bin_0=400[C],nom_gamma_stat_channel1_bin_1=100[C],weight=0)

named sets
----------
ModelConfig_GlobalObservables:(nom_alpha_syst2,nom_alpha_syst3,nom_alpha_syst1,nom_gamma_stat_channel1_bin_0,nom_gamma_stat_channel1_bin_1)
ModelConfig_NuisParams:(gamma_stat_channel1_bin_0,gamma_stat_channel1_bin_1,alpha_syst2,alpha_syst3,alpha_syst1)
ModelConfig_Observables:(obs_x_channel1)
ModelConfig_POI:(mu)

generic objects
---------------
RooStats::ModelConfig::ModelConfig

now, you can easily use your workspace to run your fit (as you usually would) the model config is named after your pdf, i.e. _modelConfig

In [5]:
model = ws["ModelConfig"]
pdf = model.GetPdf()

we are fitting a clone of the model now, such that we are not double-fitting the model in the closure check

In [6]:
result = pdf.cloneTree().fitTo(ws["observed"], Save=True, GlobalObservables=model.GetGlobalObservables(), PrintLevel=-1)
result.Print()
[#1] INFO:Minimization -- p.d.f. provides expected number of events, including extended term in likelihood.
[#1] INFO:Minimization --  Including the following constraint terms in minimization: (gamma_stat_channel1_bin_0_constraint,gamma_stat_channel1_bin_1_constraint,syst2_constraint,syst3_constraint,syst1_constraint)
[#1] INFO:Minimization -- The following global observables have been defined and their values are taken from the model: (nom_alpha_syst2,nom_alpha_syst3,nom_alpha_syst1,nom_gamma_stat_channel1_bin_0,nom_gamma_stat_channel1_bin_1)
[#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_main_observed_with_constr) Summation contains a RooNLLVar, using its error level
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
RooAbsTestStatistic::initSimMode: creating slave calculator #0 for state channel1 (2 dataset entries)
[#1] INFO:Fitting -- RooAbsTestStatistic::initSimMode: created 1 slave calculators.
[#1] INFO:Minimization --  The following expressions have been identified as constant and will be precalculated and cached: (background1_binWidth,hist_background1,background2_binWidth,hist_background2,signal)
[#1] INFO:Minimization --  The following expressions will be evaluated in cache-and-track mode: (mc_stat_channel1)
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization

  RooFitResult: minimized FCN value: 16.4964, estimated distance to minimum: 3.65531e-05
                covariance matrix quality: Full, accurate covariance matrix
                Status : MINIMIZE=0 HESSE=0 

    Floating Parameter    FinalValue +/-  Error   
  --------------------  --------------------------
           alpha_syst1   -1.4417e-04 +/-  9.93e-01
           alpha_syst2   -7.5268e-03 +/-  9.83e-01
           alpha_syst3    1.2714e-02 +/-  9.48e-01
  gamma_stat_channel1_bin_0    9.9962e-01 +/-  4.94e-02
  gamma_stat_channel1_bin_1    1.0038e+00 +/-  8.03e-02
                    mu    1.1150e+00 +/-  6.04e-01

in the end, you can again write to json the result will be not completely identical to the JSON file you used as an input, but it will work just the same

In [7]:
tool.exportJSON("myWorkspace.json")
Out[7]:
True
[#0] WARNING:InputArguments -- RooAbsReal::createHistogram(hist_background1) WARNING extended mode requested for a non-pdf object, ignored
[#1] INFO:InputArguments -- RooAbsReal::createHistogram(hist_background1) INFO: Model has intrinsic binning definition, selecting that binning for the histogram
[#0] WARNING:InputArguments -- RooAbsReal::createHistogram(hist_background2) WARNING extended mode requested for a non-pdf object, ignored
[#1] INFO:InputArguments -- RooAbsReal::createHistogram(hist_background2) INFO: Model has intrinsic binning definition, selecting that binning for the histogram
[#0] WARNING:InputArguments -- RooAbsReal::createHistogram(hist_signal) WARNING extended mode requested for a non-pdf object, ignored
[#1] INFO:InputArguments -- RooAbsReal::createHistogram(hist_signal) INFO: Model has intrinsic binning definition, selecting that binning for the histogram

You can again import it if you want and check for closure

In [8]:
ws2 = ROOT.RooWorkspace("workspace")
tool2 = ROOT.RooJSONFactoryWSTool(ws2)
tool2.importJSON("myWorkspace.json")
model2 = ws2["main_modelConfig"]
result = model.GetPdf().fitTo(ws2["observed"], Save=True, GlobalObservables=model.GetGlobalObservables(), PrintLevel=-1)
result.Print()
[#1] INFO:ObjectHandling -- RooWorkspace::import(workspace) using existing copy of RooRealVar::gamma_stat_channel1_bin_0 for import of RooRealVar::gamma_stat_channel1_bin_0
[#1] INFO:ObjectHandling -- RooWorkspace::import(workspace) using existing copy of RooRealVar::gamma_stat_channel1_bin_1 for import of RooRealVar::gamma_stat_channel1_bin_1
[#1] INFO:ObjectHandling -- RooWorkspace::import(workspace) importing dataset observed
[#1] INFO:Minimization -- p.d.f. provides expected number of events, including extended term in likelihood.
[#1] INFO:Minimization --  Including the following constraint terms in minimization: (gamma_stat_channel1_bin_0_constraint,gamma_stat_channel1_bin_1_constraint,syst2_constraint,syst3_constraint,syst1_constraint)
[#1] INFO:Minimization -- The following global observables have been defined and their values are taken from the model: (nom_alpha_syst2,nom_alpha_syst3,nom_alpha_syst1,nom_gamma_stat_channel1_bin_0,nom_gamma_stat_channel1_bin_1)
[#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_main_observed_with_constr) Summation contains a RooNLLVar, using its error level
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
RooAbsTestStatistic::initSimMode: creating slave calculator #0 for state channel1 (2 dataset entries)
[#1] INFO:Fitting -- RooAbsTestStatistic::initSimMode: created 1 slave calculators.
[#1] INFO:Minimization --  The following expressions have been identified as constant and will be precalculated and cached: (background1_binWidth,hist_background1,background2_binWidth,hist_background2,signal)
[#1] INFO:Minimization --  The following expressions will be evaluated in cache-and-track mode: (mc_stat_channel1)
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization

  RooFitResult: minimized FCN value: 16.4964, estimated distance to minimum: 3.65532e-05
                covariance matrix quality: Full, accurate covariance matrix
                Status : MINIMIZE=0 HESSE=0 

    Floating Parameter    FinalValue +/-  Error   
  --------------------  --------------------------
           alpha_syst1   -1.4417e-04 +/-  9.93e-01
           alpha_syst2   -7.5268e-03 +/-  9.83e-01
           alpha_syst3    1.2714e-02 +/-  9.48e-01
  gamma_stat_channel1_bin_0    9.9962e-01 +/-  4.94e-02
  gamma_stat_channel1_bin_1    1.0038e+00 +/-  8.03e-02
                    mu    1.1150e+00 +/-  6.04e-01

Draw all canvases

In [9]:
from ROOT import gROOT 
gROOT.GetListOfCanvases().Draw()