ROOT   Reference Guide
Searching...
No Matches
combinedFit.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_fit
3## \notebook
4## Combined (simultaneous) fit of two histogram with separate functions
5## and some common parameters
6##
7## \macro_image
8## \macro_output
9## \macro_code
10##
11## \author Jonas Rembser, Lorenzo Moneta (C++ version)
12
13
14import ROOT
15import numpy as np
16
17
18# definition of shared parameter background function
19iparB = np.array([0, 2], dtype=np.int32) # exp amplitude in B histo and exp common parameter
20
21# signal + background function
22iparSB = np.array(
23 [
24 1, # exp amplitude in S+B histo
25 2, # exp common parameter
26 3, # Gaussian amplitude
27 4, # Gaussian mean
28 5, # Gaussian sigma
29 ],
30 dtype=np.int32,
31)
32
33# Create the GlobalCHi2 structure
34
35class GlobalChi2(object):
36 def __init__(self, f1, f2):
37 self._f1 = f1
38 self._f2 = f2
39
40 def __call__(self, par):
41 # parameter vector is first background (in common 1 and 2) and then is
42 # signal (only in 2)
43
44 # the zero-copy way to get a numpy array from a double *
45 par_arr = np.frombuffer(par, dtype=np.float64, count=6)
46
47 p1 = par_arr[iparB]
48 p2 = par_arr[iparSB]
49
50 return self._f1(p1) + self._f2(p2)
51
52
53hB = ROOT.TH1D("hB", "histo B", 100, 0, 100)
54hSB = ROOT.TH1D("hSB", "histo S+B", 100, 0, 100)
55
56fB = ROOT.TF1("fB", "expo", 0, 100)
57fB.SetParameters(1, -0.05)
58hB.FillRandom("fB")
59
60fS = ROOT.TF1("fS", "gaus", 0, 100)
61fS.SetParameters(1, 30, 5)
62
63hSB.FillRandom("fB", 2000)
64hSB.FillRandom("fS", 1000)
65
66# perform now global fit
67
68fSB = ROOT.TF1("fSB", "expo + gaus(2)", 0, 100)
69
71wfSB = ROOT.Math.WrappedMultiTF1(fSB, 1)
72
74rangeB = ROOT.Fit.DataRange()
75# set the data range
76rangeB.SetRange(10, 90)
77dataB = ROOT.Fit.BinData(opt, rangeB)
78ROOT.Fit.FillData(dataB, hB)
79
80rangeSB = ROOT.Fit.DataRange()
81rangeSB.SetRange(10, 50)
82dataSB = ROOT.Fit.BinData(opt, rangeSB)
83ROOT.Fit.FillData(dataSB, hSB)
84
85chi2_B = ROOT.Fit.Chi2Function(dataB, wfB)
86chi2_SB = ROOT.Fit.Chi2Function(dataSB, wfSB)
87
88globalChi2 = GlobalChi2(chi2_B, chi2_SB)
89
90fitter = ROOT.Fit.Fitter()
91
92Npar = 6
93par0 = np.array([5, 5, -0.1, 100, 30, 10])
94
95# create before the parameter settings in order to fix or set range on them
96fitter.Config().SetParamsSettings(6, par0)
97# fix 5-th parameter
98fitter.Config().ParSettings(4).Fix()
99# set limits on the third and 4-th parameter
100fitter.Config().ParSettings(2).SetLimits(-10, -1.0e-4)
101fitter.Config().ParSettings(3).SetLimits(0, 10000)
102fitter.Config().ParSettings(3).SetStepSize(5)
103
104fitter.Config().MinimizerOptions().SetPrintLevel(0)
106
107# we can't pass the Python object globalChi2 directly to FitFCN.
108# It needs to be wrapped in a ROOT::Math::Functor.
109globalChi2Functor = ROOT.Math.Functor(globalChi2, 6)
110
111# fit FCN function
112# (specify optionally data size and flag to indicate that is a chi2 fit)
113fitter.FitFCN(globalChi2Functor, 0, dataB.Size() + dataSB.Size(), True)
114result = fitter.Result()
115result.Print(ROOT.std.cout)
116
117c1 = ROOT.TCanvas("Simfit", "Simultaneous fit of two histograms", 10, 10, 700, 700)
118c1.Divide(1, 2)
119c1.cd(1)
120ROOT.gStyle.SetOptFit(1111)
121
122fB.SetFitResult(result, iparB)
123fB.SetRange(rangeB().first, rangeB().second)
124fB.SetLineColor(ROOT.kBlue)
126hB.Draw()
127
128c1.cd(2)
129fSB.SetFitResult(result, iparSB)
130fSB.SetRange(rangeSB().first, rangeSB().second)
131fSB.SetLineColor(ROOT.kRed)
133hSB.Draw()
134
135c1.SaveAs("combinedFit.png")
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition BinData.h:52
Chi2FCN class for binned fits using the least square methods.
Definition Chi2FCN.h:46
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition DataRange.h:35
Fitter class, entry point for performing all type of fits.
Definition Fitter.h:77
Documentation for class Functor class.
Definition Functor.h:47
Class to Wrap a ROOT Function class (like TF1) in a IParamMultiFunction interface of multi-dimensions...
void FillData(BinData &dv, const TH1 *hist, TF1 *func=nullptr)
fill the data vector from a TH1.
TMatrixT< Element > & Add(TMatrixT< Element > &target, Element scalar, const TMatrixT< Element > &source)
Modify addition: target += scalar * source.
DataOptions : simple structure holding the options on how the data are filled.
Definition DataOptions.h:28