Fitting and plotting in sub ranges.
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:17 AM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooPolynomial.h"
#include "RooAddPdf.h"
#include "RooFitResult.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH1.h"
using namespace RooFit;
Construct observables x
RooRealVar x("x", "x", -10, 10);
Construct gaussx(x,mx,1)
RooRealVar mx("mx", "mx", 0, -10, 10);
RooGaussian gx("gx", "gx", x, mx, 1.0);
Construct px = 1 (flat in x)
RooPolynomial px("px", "px", x);
Construct model = f*gx + (1-f)px
RooRealVar f("f", "f", 0., 1.);
RooAddPdf model("model", "model", RooArgList(gx, px), f);
Generated 10000 events in (x,y) from pdf model
std::unique_ptr<RooDataSet> modelData{model.generate(x, 10000)};
Fit pdf to all data
std::unique_ptr<RooFitResult> r_full{model.fitTo(*modelData, Save(true), PrintLevel(-1))};
[#1] INFO:Fitting -- RooAbsPdf::fitTo(model) fixing normalization set for coefficient determination to observables in data [#1] INFO:Fitting -- using CPU computation library compiled with -mavx2 [#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_model_modelData) Summation contains a RooNLLVar, using its error level [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Define "signal" range in x as [-3,3]
x.setRange("signal", -3, 3);
[#1] INFO:Eval -- RooRealVar::setRange(x) new range named 'signal' created with bounds [-3,3]
Fit pdf only to data in "signal" range
std::unique_ptr<RooFitResult> r_sig{model.fitTo(*modelData, Save(true), Range("signal"), PrintLevel(-1))};
[#1] INFO:Eval -- RooRealVar::setRange(x) new range named 'fit_nll_model_modelData' created with bounds [-3,3] [#1] INFO:Fitting -- RooAbsPdf::fitTo(model) fixing normalization set for coefficient determination to observables in data [#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_model_modelData) Summation contains a RooNLLVar, using its error level [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Make plot frame in x and add data and fitted model
RooPlot *frame = x.frame(Title("Fitting a sub range"));
modelData->plotOn(frame);
model.plotOn(frame, Range(""), LineStyle(kDashed), LineColor(kRed)); // Add shape in full ranged dashed
model.plotOn(frame); // By default only fitted range is shown
[#1] INFO:Plotting -- RooAbsPdf::plotOn(model) p.d.f was fitted in a subrange and no explicit NormRange() was specified. Plotting / normalising in fit range. To override, do one of the following - Clear the automatic fit range attribute: <pdf>.removeStringAttribute("fitrange"); - Explicitly specify the plotting range: Range("<rangeName>"). - Explicitly specify where to compute the normalisation: NormRange("<rangeName>"). The default (full) range can be denoted with Range("") / NormRange(""). [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) only plotting range '' [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) p.d.f. curve is normalized using explicit choice of ranges 'fit_nll_model_modelData' [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) p.d.f was fitted in a subrange and no explicit Range() and NormRange() was specified. Plotting / normalising in fit range. To override, do one of the following - Clear the automatic fit range attribute: <pdf>.removeStringAttribute("fitrange"); - Explicitly specify the plotting range: Range("<rangeName>"). - Explicitly specify where to compute the normalisation: NormRange("<rangeName>"). The default (full) range can be denoted with Range("") / NormRange(""). [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) only plotting range 'fit_nll_model_modelData' [#1] INFO:Plotting -- RooAbsPdf::plotOn(model) p.d.f. curve is normalized using explicit choice of ranges 'fit_nll_model_modelData'
Print fit results
cout << "result of fit on all data " << endl;
r_full->Print();
cout << "result of fit in in signal region (note increased error on signal fraction)" << endl;
r_sig->Print();
result of fit on all data RooFitResult: minimized FCN value: 25939.4, estimated distance to minimum: 3.77183e-06 covariance matrix quality: Full, accurate covariance matrix Status : MINIMIZE=0 HESSE=0 Floating Parameter FinalValue +/- Error -------------------- -------------------------- f 5.0441e-01 +/- 6.32e-03 mx -2.1605e-02 +/- 1.77e-02 result of fit in in signal region (note increased error on signal fraction) RooFitResult: minimized FCN value: 10339.5, estimated distance to minimum: 0.000279216 covariance matrix quality: Full, accurate covariance matrix Status : MINIMIZE=0 HESSE=0 Floating Parameter FinalValue +/- Error -------------------- -------------------------- f 4.8979e-01 +/- 1.62e-02 mx -2.1518e-02 +/- 1.79e-02
Draw frame on canvas
new TCanvas("rf203_ranges", "rf203_ranges", 600, 600);
gPad->SetLeftMargin(0.15);
frame->GetYaxis()->SetTitleOffset(1.4);
frame->Draw();
return;
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()