Loading [MathJax]/extensions/tex2jax.js
Logo ROOT   6.14/05
Reference Guide
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
rf703_effpdfprod.C File Reference

Detailed Description

View in nbviewer Open in SWAN 'SPECIAL PDFS' RooFit tutorial macro #703

Using a product of an (acceptance) efficiency and a p.d.f as p.d.f.

pict1_rf703_effpdfprod.C.png
Processing /mnt/build/workspace/root-makedoc-v614/rootspi/rdoc/src/v6-14-00-patches/tutorials/roofit/rf703_effpdfprod.C...
RooFit v3.60 -- Developed by Wouter Verkerke and David Kirkby
Copyright (C) 2000-2013 NIKHEF, University of California & Stanford University
All rights reserved, please read http://roofit.sourceforge.net/license.txt
[#1] INFO:NumericIntegration -- RooRealIntegral::init(modelEff_clone_Int[t]) using numeric integrator RooIntegrator1D to calculate Int(t)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(modelEff_clone_Int[t]) using numeric integrator RooIntegrator1D to calculate Int(t)
[#1] INFO:Minization -- RooMinimizer::optimizeConst: activating const optimization
[#1] INFO:Minization -- The following expressions have been identified as constant and will be precalculated and cached: (eff)
**********
** 1 **SET PRINT 1
**********
**********
** 2 **SET NOGRAD
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 tau -1.54000e+00 3.90000e-01 -4.00000e+00 -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=9696.74 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 3.90000e-01 2.09076e-01 1.80888e+02
ERR DEF= 0.5
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=9695.84 FROM MIGRAD STATUS=CONVERGED 12 CALLS 13 TOTAL
EDM=2.41592e-05 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 tau -1.55887e+00 1.40931e-02 5.05936e-04 6.58163e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 1 ERR DEF=0.5
1.986e-04
**********
** 7 **SET ERR 0.5
**********
**********
** 8 **SET PRINT 1
**********
**********
** 9 **HESSE 500
**********
COVARIANCE MATRIX CALCULATED SUCCESSFULLY
FCN=9695.84 FROM HESSE STATUS=OK 5 CALLS 18 TOTAL
EDM=2.41576e-05 STRATEGY= 1 ERROR MATRIX ACCURATE
EXT PARAMETER INTERNAL INTERNAL
NO. NAME VALUE ERROR STEP SIZE VALUE
1 tau -1.55887e+00 1.40931e-02 1.01187e-04 2.54601e-01
ERR DEF= 0.5
EXTERNAL ERROR MATRIX. NDIM= 25 NPAR= 1 ERR DEF=0.5
1.986e-04
[#1] INFO:Minization -- RooMinimizer::optimizeConst: deactivating const optimization
[#1] INFO:NumericIntegration -- RooRealIntegral::init(modelEff_clone_Int[t]) using numeric integrator RooIntegrator1D to calculate Int(t)
#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooGaussian.h"
#include "RooConstVar.h"
#include "RooExponential.h"
#include "RooEffProd.h"
#include "RooFormulaVar.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
using namespace RooFit ;
void rf703_effpdfprod()
{
// 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) ;
// Make pdf
RooRealVar tau("tau","tau",-1.54,-4,-0.1) ;
RooExponential model("model","model",t,tau) ;
// D e f i n e e f f i c i e n c y f u n c t i o n
// ---------------------------------------------------
// Use error function to simulate turn-on slope
RooFormulaVar eff("eff","0.5*(TMath::Erf((t-1)/0.5)+1)",t) ;
// D e f i n e d e c a y p d f w i t h e f f i c i e n c y
// ---------------------------------------------------------------
// Multiply pdf(t) with efficiency in t
RooEffProd modelEff("modelEff","model with efficiency",model,eff) ;
// P l o t e f f i c i e n c y , p d f
// ----------------------------------------
RooPlot* frame1 = t.frame(Title("Efficiency")) ;
eff.plotOn(frame1,LineColor(kRed)) ;
RooPlot* frame2 = t.frame(Title("Pdf with and without efficiency")) ;
model.plotOn(frame2,LineStyle(kDashed)) ;
modelEff.plotOn(frame2) ;
// G e n e r a t e t o y d a t a , f i t m o d e l E f f t o d a t a
// ------------------------------------------------------------------------------
// Generate events. If the input pdf has an internal generator, the internal generator
// is used and an accept/reject sampling on the efficiency is applied.
RooDataSet* data = modelEff.generate(t,10000) ;
// Fit pdf. The normalization integral is calculated numerically.
modelEff.fitTo(*data) ;
// Plot generated data and overlay fitted pdf
RooPlot* frame3 = t.frame(Title("Fitted pdf with efficiency")) ;
data->plotOn(frame3) ;
modelEff.plotOn(frame3) ;
TCanvas* c = new TCanvas("rf703_effpdfprod","rf703_effpdfprod",1200,400) ;
c->Divide(3) ;
c->cd(1) ; gPad->SetLeftMargin(0.15) ; frame1->GetYaxis()->SetTitleOffset(1.4) ; frame1->Draw() ;
c->cd(2) ; gPad->SetLeftMargin(0.15) ; frame2->GetYaxis()->SetTitleOffset(1.6) ; frame2->Draw() ;
c->cd(3) ; gPad->SetLeftMargin(0.15) ; frame3->GetYaxis()->SetTitleOffset(1.6) ; frame3->Draw() ;
}
Author
07/2008 - Wouter Verkerke

Definition in file rf703_effpdfprod.C.