Logo ROOT   6.12/07
Reference Guide
rs603_HLFactoryElaborateExample.C File Reference
Tutorials » RooStats Tutorials

Detailed Description

View in nbviewer Open in SWAN 'High Level Factory Example' RooStats tutorial macro #602

This tutorial shows an example of creating a combined model using the High Level model Factory.

pict1_rs603_HLFactoryElaborateExample.C.png
Processing /mnt/build/workspace/root-makedoc-v612/rootspi/rdoc/src/v6-12-00-patches/tutorials/roostats/rs603_HLFactoryElaborateExample.C...
RooFit v3.60 -- Developed by Wouter Verkerke and David Kirkby
Copyright (C) 2000-2013 NIKHEF, University of California & Stanford University
All rights reserved, please read http://roofit.sourceforge.net/license.txt
[#1] INFO:ObjectHandling -- RooWorkspace::exportToCint(HLFactoryComplexExample_ws) INFO: references to all objects in this workspace will be created in CINT in 'namespace HLFactoryComplexExample_ws'
[HLFactoryComplexExample] echo: In the middle!
[HLFactoryComplexExample] echo: At the end!
[#1] INFO:ObjectHandling -- RooWorkspace::import(rs603_ws) importing dataset data1
[#1] INFO:ObjectHandling -- RooWorkspace::import(rs603_ws) importing RooRealVar::x
[#1] INFO:ObjectHandling -- RooWorkspace::import(rs603_ws) importing dataset data2
-------------------------------------------------------------------
Rootfile and Workspace prepared
-------------------------------------------------------------------
[#1] INFO:ObjectHandling -- RooWorkspace::exportToCint(HLFactoryElaborateExample_ws) INFO: references to all objects in this workspace will be created in CINT in 'namespace HLFactoryElaborateExample_ws'
[HLFactoryElaborateExample] echo: In the middle!
[HLFactoryElaborateExample] echo: Now reading the included file!
[HLFactoryElaborateExample] echo: Including datasets in a Workspace in a Root file...
[#1] INFO:ObjectHandling -- RooWorkspace::import(HLFactoryElaborateExample_ws) importing dataset data1
[#1] INFO:ObjectHandling -- RooWorkspace::import(HLFactoryElaborateExample_ws) importing dataset data2
[HLFactoryElaborateExample] echo: At the end!
RooDataSet::data1[x] = 219 entries
RooCategory::HLFactoryElaborateExample_category = model2(idx = 1)
[#1] INFO:Minization -- p.d.f. provides expected number of events, including extended term in likelihood.
[#1] INFO:Minization -- RooMinimizer::optimizeConst: activating const optimization
RooAbsTestStatistic::initSimMode: creating slave calculator #0 for state model1 (219 dataset entries)
RooAbsTestStatistic::initSimMode: creating slave calculator #1 for state model2 (164 dataset entries)
[#1] INFO:Fitting -- RooAbsTestStatistic::initSimMode: created 2 slave calculators.
[#1] INFO:Minization -- The following expressions have been identified as constant and will be precalculated and cached: (flat1)
[#1] INFO:Minization -- The following expressions will be evaluated in cache-and-track mode: (gauss1)
[#1] INFO:Minization -- The following expressions have been identified as constant and will be precalculated and cached: (flat2)
[#1] INFO:Minization -- The following expressions will be evaluated in cache-and-track mode: (gauss2)
**********
** 1 **SET PRINT 1
**********
**********
** 2 **SET NOGRAD
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 mean1 5.00000e+01 1.00000e+01 0.00000e+00 1.00000e+02
2 mean2 8.00000e+01 1.00000e+01 0.00000e+00 1.00000e+02
3 nbkg1 1.00000e+02 5.00000e+01 0.00000e+00 1.00000e+03
4 nbkg2 8.00000e+01 4.00000e+01 0.00000e+00 1.00000e+03
5 nsig1 1.20000e+02 3.00000e+01 0.00000e+00 3.00000e+02
6 nsig2 9.00000e+01 4.00000e+01 0.00000e+00 4.00000e+02
**********
** 3 **SET ERR 0.5
**********
**********
** 4 **SET PRINT 1
**********
**********
** 5 **SET STR 1
**********
NOW USING STRATEGY 1: TRY TO BALANCE SPEED AGAINST RELIABILITY
**********
** 6 **MIGRAD 3000 1
**********
FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
START MIGRAD MINIMIZATION. STRATEGY 1. CONVERGENCE WHEN EDM .LT. 1.00e-03
FCN=168.565 FROM MIGRAD STATUS=INITIATE 24 CALLS 25 TOTAL
EDM= unknown STRATEGY= 1 NO ERROR MATRIX
EXT PARAMETER CURRENT GUESS STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 mean1 5.00000e+01 1.00000e+01 2.01358e-01 -1.82836e+02
2 mean2 8.00000e+01 1.00000e+01 2.57889e-01 1.04014e+01
3 nbkg1 1.00000e+02 5.00000e+01 1.72186e-01 2.51699e+01
4 nbkg2 8.00000e+01 4.00000e+01 1.52384e-01 2.98132e+01
5 nsig1 1.20000e+02 3.00000e+01 2.05758e-01 -9.05084e+00
6 nsig2 9.00000e+01 4.00000e+01 2.45245e-01 -5.18071e+00
ERR DEF= 0.5
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=165.577 FROM MIGRAD STATUS=CONVERGED 105 CALLS 106 TOTAL
EDM=1.82025e-06 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 mean1 5.06657e+01 4.11801e-01 7.32781e-05 -9.51127e-02
2 mean2 7.98756e+01 6.29188e-01 1.39674e-04 1.05318e-02
3 nbkg1 8.74181e+01 1.04372e+01 3.24289e-04 1.44260e-02
4 nbkg2 6.75739e+01 9.49233e+00 3.28828e-04 2.19919e-02
5 nsig1 1.31590e+02 1.23613e+01 7.28660e-04 3.75153e-03
6 nsig2 9.64299e+01 1.09011e+01 5.54180e-04 -3.77030e-03
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 6 ERR DEF=0.5
1.696e-01 0.000e+00 -9.815e-02 0.000e+00 9.813e-02 0.000e+00
0.000e+00 3.959e-01 0.000e+00 1.140e-01 0.000e+00 -1.141e-01
-9.815e-02 0.000e+00 1.090e+02 0.000e+00 -2.155e+01 0.000e+00
0.000e+00 1.140e-01 0.000e+00 9.015e+01 0.000e+00 -2.256e+01
9.813e-02 0.000e+00 -2.155e+01 0.000e+00 1.532e+02 0.000e+00
0.000e+00 -1.141e-01 0.000e+00 -2.256e+01 0.000e+00 1.190e+02
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2 3 4 5 6
1 0.02769 1.000 0.000 -0.023 0.000 0.019 0.000
2 0.02296 0.000 1.000 0.000 0.019 0.000 -0.017
3 0.16798 -0.023 0.000 1.000 0.000 -0.167 0.000
4 0.21836 0.000 0.019 0.000 1.000 0.000 -0.218
5 0.16755 0.019 0.000 -0.167 0.000 1.000 0.000
6 0.21817 0.000 -0.017 0.000 -0.218 0.000 1.000
**********
** 7 **SET ERR 0.5
**********
**********
** 8 **SET PRINT 1
**********
**********
** 9 **HESSE 3000
**********
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=165.577 FROM HESSE STATUS=OK 40 CALLS 146 TOTAL
EDM=1.82039e-06 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 mean1 5.06657e+01 4.11802e-01 1.46556e-05 1.33141e-02
2 mean2 7.98756e+01 6.29186e-01 2.79348e-05 6.40394e-01
3 nbkg1 8.74181e+01 1.04374e+01 6.48579e-05 -9.70492e-01
4 nbkg2 6.75739e+01 9.49253e+00 6.57656e-05 -1.04486e+00
5 nsig1 1.31590e+02 1.23615e+01 1.45732e-04 -1.23045e-01
6 nsig2 9.64299e+01 1.09013e+01 1.10836e-04 -5.44336e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 6 ERR DEF=0.5
1.696e-01 0.000e+00 -9.850e-02 0.000e+00 9.849e-02 0.000e+00
0.000e+00 3.959e-01 0.000e+00 1.133e-01 0.000e+00 -1.133e-01
-9.850e-02 0.000e+00 1.090e+02 0.000e+00 -2.157e+01 0.000e+00
0.000e+00 1.133e-01 0.000e+00 9.015e+01 0.000e+00 -2.257e+01
9.849e-02 0.000e+00 -2.157e+01 0.000e+00 1.532e+02 0.000e+00
0.000e+00 -1.133e-01 0.000e+00 -2.257e+01 0.000e+00 1.190e+02
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 1 2 3 4 5 6
1 0.02779 1.000 0.000 -0.023 0.000 0.019 0.000
2 0.02281 0.000 1.000 0.000 0.019 0.000 -0.017
3 0.16807 -0.023 0.000 1.000 0.000 -0.167 0.000
4 0.21846 0.000 0.019 0.000 1.000 0.000 -0.218
5 0.16763 0.019 0.000 -0.167 0.000 1.000 0.000
6 0.21827 0.000 -0.017 0.000 -0.218 0.000 1.000
[#1] INFO:Minization -- RooMinimizer::optimizeConst: deactivating const optimization
[#1] INFO:Plotting -- RooSimultaneous::plotOn(HLFactoryElaborateExample_sigbkg) plot on x represents a slice in the index category (HLFactoryElaborateExample_category)
[#1] INFO:Plotting -- RooAbsReal::plotOn(sb_model1) slice variable HLFactoryElaborateExample_category was not projected anyway
[#1] INFO:Plotting -- RooSimultaneous::plotOn(HLFactoryElaborateExample_sigbkg) plot on x represents a slice in the index category (HLFactoryElaborateExample_category)
[#1] INFO:Plotting -- RooAbsReal::plotOn(sb_model2) slice variable HLFactoryElaborateExample_category was not projected anyway
#include <fstream>
#include "TString.h"
#include "TFile.h"
#include "TROOT.h"
#include "RooGlobalFunc.h"
#include "RooWorkspace.h"
#include "RooRealVar.h"
#include "RooAbsPdf.h"
#include "RooDataSet.h"
#include "RooPlot.h"
#include "RooStats/HLFactory.h"
// use this order for safety on library loading
using namespace RooFit;
using namespace RooStats;
using namespace std;
void rs603_HLFactoryElaborateExample() {
// --- Prepare the 2 needed datacards for this example ---
TString card_name("rs603_card_WsMaker.rs");
ofstream ofile(card_name);
ofile << "// The simplest card for combination\n\n";
ofile << "gauss1 = Gaussian(x[0,100],mean1[50,0,100],4);\n";
ofile << "flat1 = Polynomial(x,0);\n";
ofile << "sb_model1 = SUM(nsig1[120,0,300]*gauss1 , nbkg1[100,0,1000]*flat1);\n\n";
ofile << "echo In the middle!;\n\n";
ofile << "gauss2 = Gaussian(x,mean2[80,0,100],5);\n";
ofile << "flat2 = Polynomial(x,0);\n";
ofile << "sb_model2 = SUM(nsig2[90,0,400]*gauss2 , nbkg2[80,0,1000]*flat2);\n\n";
ofile << "echo At the end!;\n";
ofile.close();
TString card_name2("rs603_card.rs");
ofstream ofile2(card_name2);
ofile2 << "// The simplest card for combination\n\n";
ofile2 << "gauss1 = Gaussian(x[0,100],mean1[50,0,100],4);\n";
ofile2 << "flat1 = Polynomial(x,0);\n";
ofile2 << "sb_model1 = SUM(nsig1[120,0,300]*gauss1 , nbkg1[100,0,1000]*flat1);\n\n";
ofile2 << "echo In the middle!;\n\n";
ofile2 << "gauss2 = Gaussian(x,mean2[80,0,100],5);\n";
ofile2 << "flat2 = Polynomial(x,0);\n";
ofile2 << "sb_model2 = SUM(nsig2[90,0,400]*gauss2 , nbkg2[80,0,1000]*flat2);\n\n";
ofile2 << "#include rs603_included_card.rs;\n\n";
ofile2 << "echo At the end!;\n";
ofile2.close();
TString card_name3("rs603_included_card.rs");
ofstream ofile3(card_name3);
ofile3 << "echo Now reading the included file!;\n\n";
ofile3 << "echo Including datasets in a Workspace in a Root file...;\n";
ofile3 << "data1 = import(rs603_infile.root,\n";
ofile3 << " rs603_ws,\n";
ofile3 << " data1);\n\n";
ofile3 << "data2 = import(rs603_infile.root,\n";
ofile3 << " rs603_ws,\n";
ofile3 << " data2);\n";
ofile3.close();
// --- Produce the two separate datasets into a WorkSpace ---
HLFactory hlf("HLFactoryComplexExample",
"rs603_card_WsMaker.rs",
false);
RooRealVar* x = static_cast<RooRealVar*>(hlf.GetWs()->arg("x"));
RooAbsPdf* pdf1 = hlf.GetWs()->pdf("sb_model1");
RooAbsPdf* pdf2 = hlf.GetWs()->pdf("sb_model2");
RooWorkspace w("rs603_ws");
RooDataSet* data1 = pdf1->generate(RooArgSet(*x),Extended());
data1->SetName("data1");
w.import(*data1);
RooDataSet* data2 = pdf2->generate(RooArgSet(*x),Extended());
data2->SetName("data2");
w.import(*data2);
// --- Write the WorkSpace into a rootfile ---
TFile outfile("rs603_infile.root","RECREATE");
w.Write();
outfile.Close();
cout << "-------------------------------------------------------------------\n"
<< " Rootfile and Workspace prepared \n"
<< "-------------------------------------------------------------------\n";
HLFactory hlf_2("HLFactoryElaborateExample",
"rs603_card.rs",
false);
x = hlf_2.GetWs()->var("x");
pdf1 = hlf_2.GetWs()->pdf("sb_model1");
pdf2 = hlf_2.GetWs()->pdf("sb_model2");
hlf_2.AddChannel("model1","sb_model1","flat1","data1");
hlf_2.AddChannel("model2","sb_model2","flat2","data2");
RooDataSet* data = hlf_2.GetTotDataSet();
RooAbsPdf* pdf = hlf_2.GetTotSigBkgPdf();
RooCategory* thecat = hlf_2.GetTotCategory();
// --- Perform extended ML fit of composite PDF to toy data ---
pdf->fitTo(*data) ;
// --- Plot toy data and composite PDF overlaid ---
RooPlot* xframe = x->frame() ;
data->plotOn(xframe);
thecat->setIndex(0);
pdf->plotOn(xframe,Slice(*thecat),ProjWData(*thecat,*data)) ;
thecat->setIndex(1);
pdf->plotOn(xframe,Slice(*thecat),ProjWData(*thecat,*data)) ;
gROOT->SetStyle("Plain");
xframe->Draw();
}
Author
Danilo Piparo

Definition in file rs603_HLFactoryElaborateExample.C.