Logo ROOT  
Reference Guide
rf314_paramfitrange.C File Reference

Detailed Description

View in nbviewer Open in SWAN 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

␛[1mRooFit v3.60 -- Developed by Wouter Verkerke and David Kirkby␛[0m
Copyright (C) 2000-2013 NIKHEF, University of California & Stanford University
All rights reserved, please read http://roofit.sourceforge.net/license.txt
[#1] INFO:Minization -- RooMinimizer::optimizeConst: activating const optimization
**********
** 1 **SET PRINT 1
**********
**********
** 2 **SET NOGRAD
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 tau -1.54000e+00 7.20000e-01 -1.00000e+01 -1.00000e-01
**********
** 3 **SET ERR 0.5
**********
**********
** 4 **SET PRINT 1
**********
**********
** 5 **SET STR 1
**********
NOW USING STRATEGY 1: TRY TO BALANCE SPEED AGAINST RELIABILITY
**********
** 6 **MIGRAD 500 1
**********
FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
START MIGRAD MINIMIZATION. STRATEGY 1. CONVERGENCE WHEN EDM .LT. 1.00e-03
FCN=2824.02 FROM MIGRAD STATUS=INITIATE 4 CALLS 5 TOTAL
EDM= unknown STRATEGY= 1 NO ERROR MATRIX
EXT PARAMETER CURRENT GUESS STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 tau -1.54000e+00 7.20000e-01 2.12947e-01 -4.56069e+01
ERR DEF= 0.5
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=2823.97 FROM MIGRAD STATUS=CONVERGED 12 CALLS 13 TOTAL
EDM=6.75021e-08 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 tau -1.53353e+00 2.21980e-02 2.34054e-04 4.07752e-02
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 1 ERR DEF=0.5
4.928e-04
**********
** 7 **SET ERR 0.5
**********
**********
** 8 **SET PRINT 1
**********
**********
** 9 **HESSE 500
**********
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=2823.97 FROM HESSE STATUS=OK 5 CALLS 18 TOTAL
EDM=6.74739e-08 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 tau -1.53353e+00 2.21980e-02 4.68108e-05 7.90063e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 1 ERR DEF=0.5
4.928e-04
[#1] INFO:Minization -- RooMinimizer::optimizeConst: deactivating const optimization
[#1] INFO:Plotting -- RooPlot::updateFitRangeNorm: New event count of 5000 will supercede previous event count of 10000 for normalization of PDF projections
RooFitResult: minimized FCN value: 2823.97, estimated distance to minimum: 6.74739e-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>
#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;
{
// D e f i n e o b s e r v a b l e s a n d d e c a y p d f
// ---------------------------------------------------------------
// 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);
// C r e a t e i n p u t d a t a
// ------------------------------------
// Generate complete dataset without acceptance cuts (for reference)
RooDataSet *dall = model.generate(t, 10000);
// Generate a (fake) prototype dataset for acceptance limit values
RooDataSet *tmp = RooGaussian("gmin", "gmin", tmin, RooConst(0), RooConst(0.5)).generate(tmin, 5000);
// Generate dataset with t values that observe (t>tmin)
RooDataSet *dacc = model.generate(t, ProtoData(*tmp));
// F i t p d f t o d a t a i n a c c e p t a n c e r e g i o n
// -----------------------------------------------------------------------
RooFitResult *r = model.fitTo(*dacc, Save());
// P l o t f i t t e d p d f o n f u l l a n d a c c e p t e d d a t a
// ---------------------------------------------------------------------------------
// 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);
// 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;
}
ROOT::R::TRInterface & r
Definition: Object.C:4
@ kRed
Definition: Rtypes.h:64
#define gPad
Definition: TVirtualPad.h:286
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Calls RooPlot* plotOn(RooPlot* frame, const RooLinkedList& cmdList) const ;.
Definition: RooAbsData.cxx:550
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,...
Definition: RooAbsPdf.h:55
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:31
Exponential PDF.
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Definition: RooFitResult.h:40
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:44
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1277
static RooPlot * frame(const RooAbsRealLValue &var, Double_t xmin, Double_t xmax, Int_t nBins)
Create a new frame for a given variable in x.
Definition: RooPlot.cxx:277
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:712
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:35
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:294
The Canvas class.
Definition: TCanvas.h:31
virtual void Print(Option_t *option="") const
This method must be overridden when a class wants to print itself.
Definition: TObject.cxx:550
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
RooCmdArg ProtoData(const RooDataSet &protoData, Bool_t randomizeOrder=kFALSE, Bool_t resample=kFALSE)
RooCmdArg MarkerColor(Color_t color)
RooConstVar & RooConst(Double_t val)
RooCmdArg Save(Bool_t flag=kTRUE)
RooCmdArg LineColor(Color_t color)
const char * Title
Definition: TXMLSetup.cxx:67
Author
07/2008 - Wouter Verkerke

Definition in file rf314_paramfitrange.C.