rf407_latextables

Data and categories: latex printing of lists and sets of RooArgSets

Author: Clemens Lange, Wouter Verkerke (C++ version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Thursday, December 08, 2022 at 12:02 PM.

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

Setup composite pdf

Declare observable x

In [2]:
x = ROOT.RooRealVar("x", "x", 0, 10)

Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters

In [3]:
mean = ROOT.RooRealVar("mean", "mean of gaussians", 5)
sigma1 = ROOT.RooRealVar("sigma1", "width of gaussians", 0.5)
sigma2 = ROOT.RooRealVar("sigma2", "width of gaussians", 1)
sig1 = ROOT.RooGaussian("sig1", "Signal component 1", x, mean, sigma1)
sig2 = ROOT.RooGaussian("sig2", "Signal component 2", x, mean, sigma2)
[#0] WARNING:InputArguments -- The parameter 'sigma1' with range [-1e+30, 1e+30] of the RooGaussian 'sig1' exceeds the safe range of (0, inf). Advise to limit its range.
[#0] WARNING:InputArguments -- The parameter 'sigma2' with range [-1e+30, 1e+30] of the RooGaussian 'sig2' exceeds the safe range of (0, inf). Advise to limit its range.

Sum the signal components into a composite signal pdf

In [4]:
sig1frac = ROOT.RooRealVar("sig1frac", "fraction of component 1 in signal", 0.8, 0.0, 1.0)
sig = ROOT.RooAddPdf("sig", "Signal", [sig1, sig2], [sig1frac])

Build Chebychev polynomial pdf

In [5]:
a0 = ROOT.RooRealVar("a0", "a0", 0.5, 0.0, 1.0)
a1 = ROOT.RooRealVar("a1", "a1", -0.2, 0.0, 1.0)
bkg1 = ROOT.RooChebychev("bkg1", "Background 1", x, [a0, a1])

Build expontential pdf

In [6]:
alpha = ROOT.RooRealVar("alpha", "alpha", -1)
bkg2 = ROOT.RooExponential("bkg2", "Background 2", x, alpha)

Sum the background components into a composite background pdf

In [7]:
bkg1frac = ROOT.RooRealVar("sig1frac", "fraction of component 1 in background", 0.2, 0.0, 1.0)
bkg = ROOT.RooAddPdf("bkg", "Signal", [bkg1, bkg2], [sig1frac])

Sum the composite signal and background

In [8]:
bkgfrac = ROOT.RooRealVar("bkgfrac", "fraction of background", 0.5, 0.0, 1.0)
model = ROOT.RooAddPdf("model", "g1+g2+a", [bkg, sig], [bkgfrac])

Make list of parameters before and after fit

Make list of model parameters

In [9]:
params = model.getParameters({x})

Save snapshot of prefit parameters

In [10]:
initParams = params.snapshot()

Do fit to data, obtain error estimates on parameters

In [11]:
data = model.generate({x}, 1000)
model.fitTo(data)
Out[11]:
<cppyy.gbl.RooFitResult object at 0x(nil)>
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#1] INFO:Minimization --  The following expressions have been identified as constant and will be precalculated and cached: (bkg2,sig1,sig2)
[#1] INFO:Minimization --  The following expressions will be evaluated in cache-and-track mode: (bkg1)
 **********
 **    1 **SET PRINT           1
 **********
 **********
 **    2 **SET NOGRAD
 **********
 PARAMETER DEFINITIONS:
    NO.   NAME         VALUE      STEP SIZE      LIMITS
     1 a0           5.00000e-01  1.00000e-01    0.00000e+00  1.00000e+00
     2 a1           0.00000e+00  1.00000e-01    0.00000e+00  1.00000e+00
 MINUIT WARNING IN PARAM DEF
 ============== STARTING VALUE IS AT LIMIT.
 MINUIT WARNING IN PARAMETR
 ============== VARIABLE2 IS AT ITS LOWER ALLOWED LIMIT.
 MINUIT WARNING IN PARAMETR
 ============== VARIABLE2 BROUGHT BACK INSIDE LIMITS.
     3 bkgfrac      5.00000e-01  1.00000e-01    0.00000e+00  1.00000e+00
     4 sig1frac     8.00000e-01  1.00000e-01    0.00000e+00  1.00000e+00
 **********
 **    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        2000           1
 **********
 FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
 MINUIT WARNING IN MIGrad    
 ============== VARIABLE2 IS AT ITS LOWER ALLOWED LIMIT.
 START MIGRAD MINIMIZATION.  STRATEGY  1.  CONVERGENCE WHEN EDM .LT. 1.00e-03
 FCN=1951.38 FROM MIGRAD    STATUS=INITIATE       36 CALLS          37 TOTAL
                     EDM= unknown      STRATEGY= 1      NO ERROR MATRIX       
  EXT PARAMETER               CURRENT GUESS       STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  a0           5.00000e-01   1.00000e-01   0.00000e+00   7.50773e-01
   2  a1           1.51202e-01   1.00000e-01   7.96807e-01  -2.24543e-02
   3  bkgfrac      5.00000e-01   1.00000e-01   0.00000e+00   3.17611e+01
   4  sig1frac     8.00000e-01   1.00000e-01   0.00000e+00  -1.99597e+00
                               ERR DEF= 0.5
 MIGRAD MINIMIZATION HAS CONVERGED.
 MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
 COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=1950.14 FROM MIGRAD    STATUS=CONVERGED     128 CALLS         129 TOTAL
                     EDM=7.1077e-08    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                   STEP         FIRST   
  NO.   NAME      VALUE            ERROR          SIZE      DERIVATIVE 
   1  a0           4.75556e-01   1.89453e-01   6.32546e-03  -6.37818e-04
   2  a1           2.74695e-01   1.35637e-01   7.26942e-03  -1.06943e-04
   3  bkgfrac      4.61270e-01   2.63715e-02   1.17680e-03   7.97710e-04
   4  sig1frac     7.92476e-01   5.50659e-02   2.13939e-03   7.46653e-05
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  4    ERR DEF=0.5
  3.776e-02 -8.088e-03 -1.409e-03 -8.885e-03 
 -8.088e-03  1.899e-02 -1.627e-03  1.195e-03 
 -1.409e-03 -1.627e-03  6.961e-04  6.425e-04 
 -8.885e-03  1.195e-03  6.425e-04  3.051e-03 
 PARAMETER  CORRELATION COEFFICIENTS  
       NO.  GLOBAL      1      2      3      4
        1  0.84589   1.000 -0.302 -0.275 -0.828
        2  0.63162  -0.302  1.000 -0.448  0.157
        3  0.68422  -0.275 -0.448  1.000  0.441
        4  0.85784  -0.828  0.157  0.441  1.000
 **********
 **    7 **SET ERR         0.5
 **********
 **********
 **    8 **SET PRINT           1
 **********
 **********
 **    9 **HESSE        2000
 **********
 COVARIANCE MATRIX CALCULATED SUCCESSFULLY
 FCN=1950.14 FROM HESSE     STATUS=OK             23 CALLS         152 TOTAL
                     EDM=7.04448e-08    STRATEGY= 1      ERROR MATRIX ACCURATE 
  EXT PARAMETER                                INTERNAL      INTERNAL  
  NO.   NAME      VALUE            ERROR       STEP SIZE       VALUE   
   1  a0           4.75556e-01   1.88347e-01   1.26509e-03  -4.89085e-02
   2  a1           2.74695e-01   1.35262e-01   2.90777e-04  -4.67449e-01
   3  bkgfrac      4.61270e-01   2.63173e-02   4.70719e-05  -7.75387e-02
   4  sig1frac     7.92476e-01   5.47555e-02   8.55758e-05   6.24820e-01
                               ERR DEF= 0.5
 EXTERNAL ERROR MATRIX.    NDIM=  25    NPAR=  4    ERR DEF=0.5
  3.730e-02 -7.950e-03 -1.395e-03 -8.757e-03 
 -7.950e-03  1.889e-02 -1.616e-03  1.165e-03 
 -1.395e-03 -1.616e-03  6.932e-04  6.378e-04 
 -8.757e-03  1.165e-03  6.378e-04  3.017e-03 
 PARAMETER  CORRELATION COEFFICIENTS  
       NO.  GLOBAL      1      2      3      4
        1  0.84381   1.000 -0.300 -0.274 -0.826
        2  0.62889  -0.300  1.000 -0.447  0.154
        3  0.68261  -0.274 -0.447  1.000  0.441
        4  0.85607  -0.826  0.154  0.441  1.000
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization

Print parameter list in LaTeX for (one column with names, column with values)

In [12]:
params.printLatex()
\begin{tabular}{lc}
$\verb+a0+ $ & $  0.5\pm 0.2$\\
$\verb+a1+ $ & $  0.3\pm 0.1$\\
$\verb+alpha+ $ & $ -1.00$\\
$\verb+bkgfrac+ $ & $  0.46\pm 0.03$\\
$\verb+mean+ $ & $  5$\\
$\verb+sig1frac+ $ & $  0.79\pm 0.05$\\
$\verb+sigma1+ $ & $  0.5$\\
$\verb+sigma2+ $ & $  1$\\
\end{tabular}

Print parameter list in LaTeX for (names values|names values)

In [13]:
params.printLatex(Columns=2)
\begin{tabular}{lc|lc}
$\verb+a0+ $ & $  0.5\pm 0.2$ & $\verb+mean+ $ & $  5$\\
$\verb+a1+ $ & $  0.3\pm 0.1$ & $\verb+sig1frac+ $ & $  0.79\pm 0.05$\\
$\verb+alpha+ $ & $ -1.00$ & $\verb+sigma1+ $ & $  0.5$\\
$\verb+bkgfrac+ $ & $  0.46\pm 0.03$ & $\verb+sigma2+ $ & $  1$\\
\end{tabular}

Print two parameter lists side by side (name values initvalues)

In [14]:
params.printLatex(Sibling=initParams)
\begin{tabular}{lcc}
$\verb+a0+ $ & $  0.5\pm 0.2$ & $ 0.5$\\
$\verb+a1+ $ & $  0.3\pm 0.1$ & $ 0$\\
$\verb+alpha+ $ & $ -1.00$ & $-1.00$\\
$\verb+bkgfrac+ $ & $  0.46\pm 0.03$ & $ 0.5$\\
$\verb+mean+ $ & $  5$ & $ 5$\\
$\verb+sig1frac+ $ & $  0.79\pm 0.05$ & $ 0.8$\\
$\verb+sigma1+ $ & $  0.5$ & $ 0.5$\\
$\verb+sigma2+ $ & $  1$ & $ 1$\\
\end{tabular}

Print two parameter lists side by side (name values initvalues|name values initvalues)

In [15]:
params.printLatex(Sibling=initParams, Columns=2)
\begin{tabular}{lcc|lcc}
$\verb+a0+ $ & $  0.5\pm 0.2$ & $ 0.5$ & $\verb+mean+ $ & $  5$ & $ 5$\\
$\verb+a1+ $ & $  0.3\pm 0.1$ & $ 0$ & $\verb+sig1frac+ $ & $  0.79\pm 0.05$ & $ 0.8$\\
$\verb+alpha+ $ & $ -1.00$ & $-1.00$ & $\verb+sigma1+ $ & $  0.5$ & $ 0.5$\\
$\verb+bkgfrac+ $ & $  0.46\pm 0.03$ & $ 0.5$ & $\verb+sigma2+ $ & $  1$ & $ 1$\\
\end{tabular}

Write LaTex table to file

In [16]:
params.printLatex(Sibling=initParams, OutputFile="rf407_latextables.tex")