rf604_constraints.py

Go to the documentation of this file.

## \file
## \ingroup tutorial_roofit
## \notebook -nodraw
## Likelihood and minimization: fitting with constraints
##
## \macro_code
## \macro_output
##
## \date February 2018
## \authors Clemens Lange, Wouter Verkerke (C++ version)
 
import ROOT
 
 
# Create model and dataset
# ----------------------------------------------
 
# Construct a Gaussian pdf
x = ROOT.RooRealVar("x", "x", -10, 10)
 
m = ROOT.RooRealVar("m", "m", 0, -10, 10)
s = ROOT.RooRealVar("s", "s", 2, 0.1, 10)
gauss = ROOT.RooGaussian("gauss", "gauss(x,m,s)", x, m, s)
 
# Construct a flat pdf (polynomial of 0th order)
poly = ROOT.RooPolynomial("poly", "poly(x)", x)
 
# model = f*gauss + (1-f)*poly
f = ROOT.RooRealVar("f", "f", 0.5, 0.0, 1.0)
model = ROOT.RooAddPdf("model", "model", [gauss, poly], [f])
 
# Generate small dataset for use in fitting below
d = model.generate({x}, 50)
 
# Create constraint pdf
# -----------------------------------------
 
# Construct Gaussian constraint pdf on parameter f at 0.8 with
# resolution of 0.1
fconstraint = ROOT.RooGaussian("fconstraint", "fconstraint", f, 0.8, 0.1)
 
# Method 1 - add internal constraint to model
# -------------------------------------------------------------------------------------
 
# Multiply constraint term with regular pdf using ROOT.RooProdPdf Specify in
# fitTo() that internal constraints on parameter f should be used
 
# Multiply constraint with pdf
modelc = ROOT.RooProdPdf("modelc", "model with constraint", [model, fconstraint])
 
# Fit model (without use of constraint term)
r1 = model.fitTo(d, Save=True, PrintLevel=-1)
 
# Fit modelc with constraint term on parameter f
r2 = modelc.fitTo(d, Constrain={f}, Save=True, PrintLevel=-1)
 
# Method 2 - specify external constraint when fitting
# ------------------------------------------------------------------------------------------
 
# Construct another Gaussian constraint pdf on parameter f at 0.8 with
# resolution of 0.1
fconstext = ROOT.RooGaussian("fconstext", "fconstext", f, 0.2, 0.1)
 
# Fit with external constraint
r3 = model.fitTo(d, ExternalConstraints={fconstext}, Save=True, PrintLevel=-1)
 
# Print the fit results
print("fit result without constraint (data generated at f=0.5)")
r1.Print("v")
print("fit result with internal constraint (data generated at f=0.5, is f=0.8+/-0.2)")
r2.Print("v")
print("fit result with (another) external constraint (data generated at f=0.5, is f=0.2+/-0.1)")
r3.Print("v")