Likelihood and minimization: setting up a chi^2 fit to a binned dataset
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Wednesday, April 17, 2024 at 11:19 AM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooChebychev.h"
#include "RooAddPdf.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
using namespace RooFit;
Declare observable x
RooRealVar x("x", "x", 0, 10);
Create two Gaussian PDFs g1(x,mean1,sigma) anf g2(x,mean2,sigma) and their parameters
RooRealVar mean("mean", "mean of gaussians", 5);
RooRealVar sigma1("sigma1", "width of gaussians", 0.5);
RooRealVar sigma2("sigma2", "width of gaussians", 1);
RooGaussian sig1("sig1", "Signal component 1", x, mean, sigma1);
RooGaussian sig2("sig2", "Signal component 2", x, mean, sigma2);
[#0] WARNING:InputArguments -- The parameter 'sigma1' with range [-inf, inf] of the RooGaussian 'sig1' exceeds the safe range of (0, inf). Advise to limit its range. [#0] WARNING:InputArguments -- The parameter 'sigma2' with range [-inf, inf] of the RooGaussian 'sig2' exceeds the safe range of (0, inf). Advise to limit its range.
Build Chebychev polynomial pdf
RooRealVar a0("a0", "a0", 0.5, 0., 1.);
RooRealVar a1("a1", "a1", 0.2, 0., 1.);
RooChebychev bkg("bkg", "Background", x, RooArgSet(a0, a1));
Sum the signal components into a composite signal pdf
RooRealVar sig1frac("sig1frac", "fraction of component 1 in signal", 0.8, 0., 1.);
RooAddPdf sig("sig", "Signal", RooArgList(sig1, sig2), sig1frac);
Sum the composite signal and background
RooRealVar bkgfrac("bkgfrac", "fraction of background", 0.5, 0., 1.);
RooAddPdf model("model", "g1+g2+a", RooArgList(bkg, sig), bkgfrac);
std::unique_ptr<RooDataSet> d{model.generate(x, 10000)};
std::unique_ptr<RooDataHist> dh{d->binnedClone()};
Construct a chi^2 of the data and the model. When a pdf is used in a chi^2 fit, the probability density scaled by the number of events in the dataset to obtain the fit function If model is an extended pdf, the expected number events is used instead of the observed number of events.
model.chi2FitTo(*dh, {PrintLevel(-1)});
[#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: (sig1,sig2) [#1] INFO:Minimization -- The following expressions will be evaluated in cache-and-track mode: (bkg) [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
NB: It is also possible to fit a RooAbsReal function to a RooDataHist using chi2FitTo().
Note that entries with zero bins are not allowed for a proper chi^2 calculation and will give error messages
std::unique_ptr<RooAbsData> dsmall{d->reduce(EventRange(1, 100))};
std::unique_ptr<RooDataHist> dhsmall{static_cast<RooDataSet&>(*dsmall).binnedClone()};
std::unique_ptr<RooAbsReal> chi2_lowstat{model.createChi2(*dhsmall)};
cout << chi2_lowstat->getVal() << endl;
90.865