Multidimensional models: working with parametrized ranges in a fit. This an example of a fit with an acceptance that changes per-event
pdf = exp(-t/tau)
with t[tmin,5]
where t
and tmin
are both observables in the dataset
Author: Wouter Verkerke
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, March 19, 2024 at 07:15 PM.
%%cpp -d
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooExponential.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "RooFitResult.h"
using namespace RooFit;
Declare observables
RooRealVar t("t", "t", 0, 5);
RooRealVar tmin("tmin", "tmin", 0, 0, 5);
Make parametrized range in t : [tmin,5]
t.setRange(tmin, RooConst(t.getMax()));
Make pdf
RooRealVar tau("tau", "tau", -1.54, -10, -0.1);
RooExponential model("model", "model", t, tau);
Generate complete dataset without acceptance cuts (for reference)
std::unique_ptr<RooDataSet> dall{model.generate(t, 10000)};
Generate a (fake) prototype dataset for acceptance limit values
std::unique_ptr<RooDataSet> tmp{RooGaussian("gmin", "gmin", tmin, 0.0, 0.5).generate(tmin, 5000)};
Generate dataset with t values that observe (t>tmin)
std::unique_ptr<RooDataSet> dacc{model.generate(t, ProtoData(*tmp))};
std::unique_ptr<RooFitResult> r{model.fitTo(*dacc, Save(), PrintLevel(-1))};
[#1] INFO:Fitting -- RooAbsPdf::fitTo(model_over_model_Int[t]) 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_over_model_Int[t]_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, add datasets and overlay model
RooPlot *frame = t.frame(Title("Fit to data with per-event acceptance"));
dall->plotOn(frame, MarkerColor(kRed), LineColor(kRed));
model.plotOn(frame);
dacc->plotOn(frame);
[#1] INFO:Plotting -- RooPlot::updateFitRangeNorm: New event count of 5000 will supersede previous event count of 10000 for normalization of PDF projections
Print fit results to demonstrate absence of bias
r->Print("v");
new TCanvas("rf314_paramranges", "rf314_paramranges", 600, 600);
gPad->SetLeftMargin(0.15);
frame->GetYaxis()->SetTitleOffset(1.6);
frame->Draw();
return;
RooFitResult: minimized FCN value: 2823.97, estimated distance to minimum: 3.17108e-08 covariance matrix quality: Full, accurate covariance matrix Status : MINIMIZE=0 HESSE=0 Floating Parameter InitialValue FinalValue +/- Error GblCorr. -------------------- ------------ -------------------------- -------- tau -1.5400e+00 -1.5335e+00 +/- 2.22e-02 <none>
Draw all canvases
%jsroot on
gROOT->GetListOfCanvases()->Draw()