Logo ROOT   6.08/07
Reference Guide
rf211_paramconv.C File Reference

Detailed Description

View in nbviewer Open in SWAN 'ADDITION AND CONVOLUTION' RooFit tutorial macro #211

Working a with a p.d.f. with a convolution operator in terms of a parameter

(require ROOT to be compiled with –enable-fftw3)



Processing /mnt/build/workspace/root-makedoc-v608/rootspi/rdoc/src/v6-08-00-patches/tutorials/roofit/rf211_paramconv.C...
#include "RooRealVar.h"
#include "RooDataHist.h"
#include "RooGaussian.h"
#include "RooGenericPdf.h"
#include "RooFormulaVar.h"
#include "RooFFTConvPdf.h"
#include "RooPlot.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "TH2.h"
using namespace RooFit ;
void rf211_paramconv()
{
// S e t u p c o m p o n e n t p d f s
// ---------------------------------------
// Gaussian g(x ; mean,sigma)
RooRealVar x("x","x",-10,10) ;
RooRealVar mean("mean","mean",-3,3) ;
RooRealVar sigma("sigma","sigma",0.5,0.1,10) ;
RooGaussian modelx("gx","gx",x,mean,sigma) ;
// Block function in mean
RooRealVar a("a","a",2,1,10) ;
RooGenericPdf model_mean("model_mean","abs(mean)<a",RooArgList(mean,a)) ;
// Convolution in mean parameter model = g(x,mean,sigma) (x) block(mean)
x.setBins(1000,"cache") ;
mean.setBins(50,"cache") ;
RooFFTConvPdf model("model","model",mean,modelx,model_mean) ;
// Configure convolution to construct a 2-D cache in (x,mean)
// rather than a 1-d cache in mean that needs to be recalculated
// for each value of x
model.setCacheObservables(x) ;
model.setBufferFraction(1.0) ;
// Integrate model over mean projModel = Int model dmean
RooAbsPdf* projModel = model.createProjection(mean) ;
// Generate 1000 toy events
RooDataHist* d = projModel->generateBinned(x,1000) ;
// Fit p.d.f. to toy data
projModel->fitTo(*d,Verbose()) ;
// Plot data and fitted p.d.f.
RooPlot* frame = x.frame(Bins(25)) ;
d->plotOn(frame) ;
projModel->plotOn(frame) ;
// Make 2d histogram of model(x;mean)
TH1* hh = model.createHistogram("hh",x,Binning(50),YVar(mean,Binning(50)),ConditionalObservables(mean)) ;
hh->SetTitle("histogram of model(x|mean)") ;
hh->SetLineColor(kBlue) ;
// Draw frame on canvas
TCanvas* c = new TCanvas("rf211_paramconv","rf211_paramconv",800,400) ;
c->Divide(2) ;
c->cd(1) ; gPad->SetLeftMargin(0.15) ; frame->GetYaxis()->SetTitleOffset(1.4) ; frame->Draw() ;
c->cd(2) ; gPad->SetLeftMargin(0.20) ; hh->GetZaxis()->SetTitleOffset(2.5) ; hh->Draw("surf") ;
}
Author
04/2009 - Wouter Verkerke

Definition in file rf211_paramconv.C.