rf211_paramconv.C File Reference

Detailed Description

Addition and convolution: working with a pdf with a convolution operator in terms of a parameter

This tutorial requires FFT3 to be enabled.

 
#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
   std::unique_ptr<RooDataHist> d{projModel->generateBinned(x, 1000)};
 
   // Fit pdf to toy data
   projModel->fitTo(*d, Verbose(), PrintLevel(-1));
 
   // Plot data and fitted pdf
   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");
}

[#1] INFO:Eval -- RooRealVar::setRange(mean) new range named 'refrange_fft_model' created with bounds [-3,3]
[#1] INFO:NumericIntegration -- RooRealIntegral::init(gx_Int[mean,x]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_mean_Int[mean]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#0] WARNING:Eval -- The FFT convolution 'model' will run with 50 bins. A decent accuracy for difficult convolutions is typically only reached with n >= 1000. Suggest to increase the number of bins of the observable 'mean'.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(model) creating new cache 0x37a80a50 with pdf gx_CONV_model_mean_CACHE_Obs[mean,x]_NORM_mean for nset (mean) with code 0
[#1] INFO:NumericIntegration -- RooRealIntegral::init(gx_Int[mean,x]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_mean_Int[mean]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#0] WARNING:Eval -- The FFT convolution 'model' will run with 50 bins. A decent accuracy for difficult convolutions is typically only reached with n >= 1000. Suggest to increase the number of bins of the observable 'mean'.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(model) creating new cache 0x37eecc00 with pdf gx_CONV_model_mean_CACHE_Obs[x,mean]_NORM_x_mean for nset (x,mean) with code 1
[#1] INFO:NumericIntegration -- RooRealIntegral::init(gx_Int[mean,x]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_mean_Int[mean]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#0] WARNING:Eval -- The FFT convolution 'model' will run with 50 bins. A decent accuracy for difficult convolutions is typically only reached with n >= 1000. Suggest to increase the number of bins of the observable 'mean'.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(model) creating new cache 0x37eecc00 with pdf gx_CONV_model_mean_CACHE_Obs[x,mean]_NORM_x_mean for nset (x,mean) with code 1
[#1] INFO:Fitting -- RooAbsPdf::fitTo(model_Int[mean]_Norm[mean,x]_wrapped_pdf) fixing normalization set for coefficient determination to observables in data
[#1] INFO:Fitting -- using generic CPU library compiled with no vectorizations
[#1] INFO:Fitting -- Creation of NLL object took 37.7932 ms
[#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_model_Int[mean]_Norm[mean,x]_wrapped_pdf_genData) Summation contains a RooNLLVar, using its error level
[#0] WARNING:Minimization -- RooAbsMinimizerFcn::synchronize: WARNING: no initial error estimate available for a: using 0.5
[#0] WARNING:Minimization -- RooAbsMinimizerFcn::synchronize: WARNING: no initial error estimate available for sigma: using 0.2
[#1] INFO:Minimization -- [fitFCN] No discrete parameters, performing continuous minimization only
sigma=0.5, [#1] INFO:NumericIntegration -- RooRealIntegral::init(model_mean_Int[mean]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#0] WARNING:Eval -- The FFT convolution 'model' will run with 50 bins. A decent accuracy for difficult convolutions is typically only reached with n >= 1000. Suggest to increase the number of bins of the observable 'mean'.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(model) creating new cache 0x38548d50 with pdf gx_CONV_model_mean_CACHE_Obs[mean] for nset () with code 2
[#1] INFO:FastEvaluations -- The class RooRealIntegral does not implement the faster batch evaluation interface. Consider requesting or implementing it to benefit from a speed up.
[#1] INFO:FastEvaluations -- The class RooWrapperPdf does not implement the faster batch evaluation interface. Consider requesting or implementing it to benefit from a speed up.
 
prevFCN = 2162.270958  a=2.017, 
prevFCN = 2162.270958  a=1.983, 
prevFCN = 2162.270958  a=2.172, 
prevFCN = 2163.393129  a=1.84, 
prevFCN = 2164.567188  a=2.017, 
prevFCN = 2162.270958  a=1.983, 
prevFCN = 2162.270958  a=2, sigma=0.5067, 
prevFCN = 2162.342065  sigma=0.4934, 
prevFCN = 2162.204536  sigma=0.5031, 
prevFCN = 2162.303734  sigma=0.4969, 
prevFCN = 2162.239214  sigma=0.4082, 
prevFCN = 2161.742001  sigma=0.4039, 
prevFCN = 2161.739445  a=2.172, 
prevFCN = 2162.012226  a=1.84, 
prevFCN = 2166.962347  a=2.017, 
prevFCN = 2161.739445  a=1.983, 
prevFCN = 2161.739445  a=2.172, 
prevFCN = 2162.012226  a=1.84, 
prevFCN = 2166.962347  a=2, sigma=0.4066, 
prevFCN = 2161.740792  sigma=0.4012, 
prevFCN = 2161.738979  a=2.067, sigma=0.3221, 
prevFCN = 2162.146113  a=2.024, sigma=0.3733, 
prevFCN = 2161.786096  a=2.01, sigma=0.3902, 
prevFCN = 2161.746127  a=2.005, sigma=0.3973, 
prevFCN = 2161.739845  a=2.002, sigma=0.4006, 
prevFCN = 2161.738993  a=2.002, sigma=0.4011, 
prevFCN = 2161.738978  a=2.019, 
prevFCN = 2161.738978  a=1.986, 
prevFCN = 2161.738978  a=2.174, 
prevFCN = 2161.983217  a=1.842, 
prevFCN = 2167.063346  a=2.019, 
prevFCN = 2161.738978  a=1.986, 
prevFCN = 2161.738978  a=2.002, sigma=0.404, 
prevFCN = 2161.739491  sigma=0.3981, 
prevFCN = 2161.739491  sigma=0.4011, 
prevFCN = 2161.738978  a=2.019, 
prevFCN = 2161.738978  a=1.986, 
prevFCN = 2161.738978  a=2.174, 
prevFCN = 2161.983217  a=1.842, 
prevFCN = 2167.063346  a=2.076, 
prevFCN = 2161.738978  a=1.931, 
prevFCN = 2163.330585  a=2.002, sigma=0.404, 
prevFCN = 2161.739491  sigma=0.3981, 
prevFCN = 2161.739491  a=2.003, sigma=0.4011, 
prevFCN = 2161.738978  a=2.001, 
prevFCN = 2161.738978  a=2.002, sigma=0.4017, 
prevFCN = 2161.738999  sigma=0.4005, 
prevFCN = 2161.738999  sigma=0.4012, 
prevFCN = 2161.738979  sigma=0.4009, 
prevFCN = 2161.738979  a=2.076, sigma=0.404, 
prevFCN = 2161.739491  a=2.039, sigma=0.4011, 
prevFCN = 2161.738978  a=2.02, sigma=0.4011, 
prevFCN = 2161.738978  a=2.229, sigma=0.4011, 
prevFCN = 2163.794407  a=2.132, sigma=0.4011, 
prevFCN = 2161.983321  a=2.085, sigma=0.4011, 
prevFCN = 2161.738978  a=2.062, sigma=0.4011, 
prevFCN = 2161.738978  a=2.046, sigma=0.4011, 
prevFCN = 2161.738978  a=2.031, 
prevFCN = 2161.738978  a=2.113, 
prevFCN = 2161.983247  a=1.966, 
prevFCN = 2163.330541  a=2.046, 
prevFCN = 2161.738978  a=2.031, 
prevFCN = 2161.738978  a=2.039, sigma=0.4022, 
prevFCN = 2161.73906  sigma=0.3999, 
prevFCN = 2161.73906  sigma=0.4011, 
prevFCN = 2161.738978  a=2.046, 
prevFCN = 2161.738978  a=2.031, 
prevFCN = 2161.738978  a=2.113, 
prevFCN = 2161.983247  a=1.966, 
prevFCN = 2163.330541  a=2.095, 
prevFCN = 2161.738978  a=1.984, 
prevFCN = 2161.738978  a=2.656, 
prevFCN = 2176.469586  a=1.55, 
prevFCN = 2182.341543  a=2.137, 
prevFCN = 2161.983247  a=1.944, 
prevFCN = 2163.330541  a=2.113, 
prevFCN = 2161.983247  a=1.967, 
prevFCN = 2163.330541  a=2.095, 
prevFCN = 2161.738978  a=1.984, 
prevFCN = 2161.738978  a=2.653, 
prevFCN = 2176.469586  a=1.552, 
prevFCN = 2182.341543  a=2.039, sigma=0.4022, 
prevFCN = 2161.73906  sigma=0.3999, 
prevFCN = 2161.73906  a=2.041, sigma=0.4011, 
prevFCN = 2161.738978  a=2.036, 
prevFCN = 2161.738978  a=2.039, sigma=0.4013, 
prevFCN = 2161.738982  sigma=0.4008, 
prevFCN = 2161.738982  sigma=0.4011, 
prevFCN = 2161.738978  sigma=0.401, 
prevFCN = 2161.738978  a=2.653, sigma=0.4022, 
prevFCN = 2176.498697  a=2.039, sigma=0.4011, [#1] INFO:NumericIntegration -- RooRealIntegral::init(gx_Int[mean,x]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_mean_Int[mean]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#0] WARNING:Eval -- The FFT convolution 'model' will run with 50 bins. A decent accuracy for difficult convolutions is typically only reached with n >= 1000. Suggest to increase the number of bins of the observable 'mean'.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(model) creating new cache 0x38541e10 with pdf gx_CONV_model_mean_CACHE_Obs[x,mean]_NORM_x_mean for nset (x,mean) with code 1
[#1] INFO:NumericIntegration -- RooRealIntegral::init(gx_Int[mean,x]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_mean_Int[mean]) using numeric integrator RooIntegrator1D to calculate Int(mean)
[#0] WARNING:Eval -- The FFT convolution 'model' will run with 50 bins. A decent accuracy for difficult convolutions is typically only reached with n >= 1000. Suggest to increase the number of bins of the observable 'mean'.
[#1] INFO:Caching -- RooAbsCachedPdf::getCache(model) creating new cache 0x38541e10 with pdf gx_CONV_model_mean_CACHE_Obs[x,mean]_NORM_x for nset (x) with code 3

Date

April 2009

Author

Wouter Verkerke

Definition in file rf211_paramconv.C.