Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf903_numintcache.C File Reference

Detailed Description

View in nbviewer Open in SWAN
Numeric algorithm tuning: caching of slow numeric integrals and parameterization of slow numeric integrals

#include "RooRealVar.h"
#include "RooDataSet.h"
#include "RooDataHist.h"
#include "RooGaussian.h"
#include "TCanvas.h"
#include "TAxis.h"
#include "RooPlot.h"
#include "RooWorkspace.h"
#include "TFile.h"
#include "TH1.h"
using namespace RooFit;
RooWorkspace *getWorkspace(Int_t mode);
{
// Mode = 0 : Run plain fit (slow)
// Mode = 1 : Generate workspace with pre-calculated integral and store it on file (prepare for accelerated running)
// Mode = 2 : Run fit from previously stored workspace including cached integrals (fast, requires run in mode=1
// first)
// C r e a t e , s a v e o r l o a d w o r k s p a c e w i t h p . d . f .
// -----------------------------------------------------------------------------------
// Make/load workspace, exit here in mode 1
RooWorkspace *w1 = getWorkspace(mode);
if (mode == 1) {
// Show workspace that was created
w1->Print();
// Show plot of cached integral values
if (hhcache) {
new TCanvas("rf903_numintcache", "rf903_numintcache", 600, 600);
hhcache->createHistogram("a")->Draw();
} else {
Error("rf903_numintcache", "Cached histogram is not existing in workspace");
}
return;
}
// U s e p . d . f . f r o m w o r k s p a c e f o r g e n e r a t i o n a n d f i t t i n g
// -----------------------------------------------------------------------------------
// This is always slow (need to find maximum function value empirically in 3D space)
std::unique_ptr<RooDataSet> d{w1->pdf("model")->generate({*w1->var("x"), *w1->var("y"), *w1->var("z")}, 1000)};
// This is slow in mode 0, but fast in mode 1
w1->pdf("model")->fitTo(*d, Verbose(true), Timer(true), PrintLevel(-1));
// Projection on x (always slow as 2D integral over Y,Z at fitted value of a is not cached)
RooPlot *framex = w1->var("x")->frame(Title("Projection of 3D model on X"));
d->plotOn(framex);
w1->pdf("model")->plotOn(framex);
// Draw x projection on canvas
auto canv = new TCanvas("rf903_numintcache", "rf903_numintcache", 600, 600);
framex->Draw();
canv->Draw();
// Make workspace available on command line after macro finishes
gDirectory->Add(w1);
}
RooWorkspace *getWorkspace(Int_t mode)
{
// C r e a t e , s a v e o r l o a d w o r k s p a c e w i t h p . d . f .
// -----------------------------------------------------------------------------------
//
// Mode = 0 : Create workspace for plain running (no integral caching)
// Mode = 1 : Generate workspace with pre-calculated integral and store it on file
// Mode = 2 : Load previously stored workspace from file
if (mode != 2) {
// Create empty workspace workspace
w = new RooWorkspace("w", 1);
// Make a difficult to normalize pdf in 3 dimensions that is integrated numerically.
w->factory("EXPR::model('1/((x-a)*(x-a)+0.01)+1/((y-a)*(y-a)+0.01)+1/"
"((z-a)*(z-a)+0.01)',x[-1,1],y[-1,1],z[-1,1],a[-5,5])");
}
if (mode == 1) {
// Instruct model to pre-calculate normalization integral that integrate at least
// two dimensions numerically. In this specific case the integral value for
// all values of parameter 'a' are stored in a histogram and available for use
// in subsequent fitting and plotting operations (interpolation is applied)
// w->pdf("model")->setNormValueCaching(3) ;
w->pdf("model")->setStringAttribute("CACHEPARMINT", "x:y:z");
// Evaluate pdf once to trigger filling of cache
RooArgSet normSet(*w->var("x"), *w->var("y"), *w->var("z"));
w->pdf("model")->getVal(&normSet);
w->writeToFile("rf903_numintcache.root");
}
if (mode == 2) {
// Load preexisting workspace from file in mode==2
TFile *f = new TFile("rf903_numintcache.root");
w = (RooWorkspace *)f->Get("w");
}
// Return created or loaded workspace
return w;
}
#define d(i)
Definition RSha256.hxx:102
#define f(i)
Definition RSha256.hxx:104
int Int_t
Definition RtypesCore.h:45
#define gDirectory
Definition TDirectory.h:384
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Definition TError.cxx:185
Option_t Option_t TPoint TPoint const char mode
TH1 * createHistogram(const char *name, const RooAbsRealLValue &xvar, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Calls createHistogram(const char *name, const RooAbsRealLValue& xvar, const RooLinkedList& argList) c...
RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={}, const RooCmdArg &arg10={}) const override
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition RooAbsPdf.h:124
RooFit::OwningPtr< RooFitResult > fitTo(RooAbsData &data, CmdArgs_t const &... cmdArgs)
Fit PDF to given dataset.
Definition RooAbsPdf.h:157
RooFit::OwningPtr< RooDataSet > generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={})
See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,...
Definition RooAbsPdf.h:57
RooPlot * frame(const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Create a new RooPlot on the heap with a drawing frame initialized for this object,...
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
Container class to hold N-dimensional binned data.
Definition RooDataHist.h:39
const TObject * getObj(Int_t uniqueID)
Retrieve payload object of cache element with given unique ID.
Plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
void Draw(Option_t *options=nullptr) override
Draw this plot and all of the elements it contains.
Definition RooPlot.cxx:649
Persistable container for RooFit projects.
void Print(Option_t *opts=nullptr) const override
Print contents of the workspace.
RooAbsPdf * pdf(RooStringView name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
RooExpensiveObjectCache & expensiveObjectCache()
RooRealVar * var(RooStringView name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found.
The Canvas class.
Definition TCanvas.h:23
A ROOT file is composed of a header, followed by consecutive data records (TKey instances) with a wel...
Definition TFile.h:53
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3062
RooCmdArg PrintLevel(Int_t code)
RooCmdArg Verbose(bool flag=true)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition JSONIO.h:26
const char * Title
Definition TXMLSetup.cxx:68
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_Int[x,y,z]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(x,y,z)
[#1] INFO:Fitting -- RooAbsPdf::fitTo(model_over_model_Int[x,y,z]) 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[x,y,z]_modelData) Summation contains a RooNLLVar, using its error level
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#0] WARNING:Minimization -- RooAbsMinimizerFcn::synchronize: WARNING: no initial error estimate available for a: using 1
prevFCN = 1659.930708 a=0.02833, [#1] INFO:NumericIntegration -- RooRealIntegral::init(model_Int[x,y,z]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(x,y,z)
prevFCN = 1675.611563 a=-0.02833,
prevFCN = 1673.217894 a=0.002833,
prevFCN = 1660.205177 a=-0.002833,
prevFCN = 1659.94939 a=0.0002833,
prevFCN = 1659.944972 a=-0.0002833,
prevFCN = 1659.919376 a=-0.001237,
prevFCN = 1659.902781 a=-0.001089,
prevFCN = 1659.903175 a=-0.001384,
prevFCN = 1659.90318 a=-0.001237,
prevFCN = 1659.902781 a=-0.001089,
prevFCN = 1659.903175 a=-0.001384,
prevFCN = 1659.90318 a=-0.001207,
prevFCN = 1659.902797 a=-0.001266,
prevFCN = 1659.902798 [#1] INFO:Minimization -- Command timer: Real time 0:00:02, CP time 2.720
[#1] INFO:Minimization -- Session timer: Real time 0:00:02, CP time 2.720
a=-0.001237,
prevFCN = 1659.902781 a=-0.001207,
prevFCN = 1659.902797 a=-0.001266,
prevFCN = 1659.902798 a=-0.001231,
prevFCN = 1659.902782 a=-0.001243,
prevFCN = 1659.902782 [#1] INFO:Minimization -- Command timer: Real time 0:00:00, CP time 0.950
[#1] INFO:Minimization -- Session timer: Real time 0:00:03, CP time 3.670, 2 slices
a=-0.001237, [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
[#0] WARNING:NumericIntegration -- RooAdaptiveIntegratorND::dtor(model) WARNING: Number of suppressed warningings about integral evaluations where target precision was not reached is 17
[#1] INFO:Plotting -- RooAbsReal::plotOn(model) plot on x integrates over variables (y,z)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_Int[x,y,z]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(x,y,z)
[#1] INFO:NumericIntegration -- RooRealIntegral::init(model_Int[y,z]_Norm[x,y,z]) using numeric integrator RooAdaptiveIntegratorND to calculate Int(y,z)
[#0] WARNING:NumericIntegration -- RooAdaptiveIntegratorND::dtor(model) WARNING: Number of suppressed warningings about integral evaluations where target precision was not reached is 1
Date
July 2008
Author
Wouter Verkerke

Definition in file rf903_numintcache.C.