Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf901_numintconfig.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_roofit
3/// \notebook -nodraw
4/// Numeric algorithm tuning: configuration and customization of how numeric (partial) integrals are executed
5///
6/// \macro_code
7/// \macro_output
8///
9/// \date July 2008
10/// \author Wouter Verkerke
11
12#include "RooRealVar.h"
13#include "RooDataSet.h"
14#include "RooGaussian.h"
15#include "TCanvas.h"
16#include "TAxis.h"
17#include "RooPlot.h"
18#include "RooNumIntConfig.h"
19#include "RooLandau.h"
20#include "RooArgSet.h"
21#include <iomanip>
22using namespace RooFit;
23
25{
26
27 // A d j u s t g l o b a l 1 D i n t e g r a t i o n p r e c i s i o n
28 // ----------------------------------------------------------------------------
29
30 // Print current global default configuration for numeric integration strategies
32
33 // Example: Change global precision for 1D integrals from 1e-7 to 1e-6
34 //
35 // The relative epsilon (change as fraction of current best integral estimate) and
36 // absolute epsilon (absolute change w.r.t last best integral estimate) can be specified
37 // separately. For most pdf integrals the relative change criterium is the most important,
38 // however for certain non-pdf functions that integrate out to zero a separate absolute
39 // change criterium is necessary to declare convergence of the integral
40 //
41 // NB: This change is for illustration only. In general the precision should be at least 1e-7
42 // for normalization integrals for MINUIT to succeed.
43 //
46
47 // N u m e r i c i n t e g r a t i o n o f l a n d a u p d f
48 // ------------------------------------------------------------------
49
50 RooRealVar x("x", "x", -10, 10);
51 RooLandau landau("landau", "landau", x, 0.0, 0.1);
52
53 // Disable analytic integration from demonstration purposes
54 landau.forceNumInt(true);
55
56 // Activate debug-level messages for topic integration to be able to follow actions below
58
59 // Calculate integral over landau with default choice of numeric integrator
60 std::unique_ptr<RooAbsReal> intLandau{landau.createIntegral(x)};
61 double val = intLandau->getVal();
62 cout << " [1] int_dx landau(x) = " << setprecision(15) << val << endl;
63
64 // S a m e w i t h c u s t o m c o n f i g u r a t i o n
65 // -----------------------------------------------------------
66
67 // Construct a custom configuration which uses the adaptive Gauss-Kronrod technique
68 // for closed 1D integrals
70#ifdef R__HAS_MATHMORE
71 customConfig.method1D().setLabel("RooAdaptiveGaussKronrodIntegrator1D");
72#else
73 Warning("rf901_numintconfig","ROOT is built without Mathmore (GSL) support. Cannot use RooAdaptiveGaussKronrodIntegrator1D");
74#endif
75
76 // Calculate integral over landau with custom integral specification
77 std::unique_ptr<RooAbsReal> intLandau2{landau.createIntegral(x, NumIntConfig(customConfig))};
78 double val2 = intLandau2->getVal();
79 cout << " [2] int_dx landau(x) = " << val2 << endl;
80
81 // A d j u s t i n g d e f a u l t c o n f i g f o r a s p e c i f i c p d f
82 // -------------------------------------------------------------------------------------
83
84 // Another possibility: associate custom numeric integration configuration as default for object 'landau'
85 landau.setIntegratorConfig(customConfig);
86
87 // Calculate integral over landau custom numeric integrator specified as object default
88 std::unique_ptr<RooAbsReal> intLandau3{landau.createIntegral(x)};
89 double val3 = intLandau3->getVal();
90 cout << " [3] int_dx landau(x) = " << val3 << endl;
91
92 // Another possibility: Change global default for 1D numeric integration strategy on finite domains
93#ifdef R__HAS_MATHMORE
94 RooAbsReal::defaultIntegratorConfig()->method1D().setLabel("RooAdaptiveGaussKronrodIntegrator1D");
95
96 // A d j u s t i n g p a r a m e t e r s o f a s p e c i f i c t e c h n i q u e
97 // ---------------------------------------------------------------------------------------
98
99 // Adjust maximum number of steps of RooIntegrator1D in the global default configuration
100 RooAbsReal::defaultIntegratorConfig()->getConfigSection("RooIntegrator1D").setRealValue("maxSteps", 30);
101
102 // Example of how to change the parameters of a numeric integrator
103 // (Each config section is a RooArgSet with RooRealVars holding real-valued parameters
104 // and RooCategories holding parameters with a finite set of options)
105 customConfig.getConfigSection("RooAdaptiveGaussKronrodIntegrator1D").setRealValue("maxSeg", 50);
106 customConfig.getConfigSection("RooAdaptiveGaussKronrodIntegrator1D").setCatLabel("method", "15Points");
107
108 // Example of how to print set of possible values for "method" category
109 customConfig.getConfigSection("RooAdaptiveGaussKronrodIntegrator1D").find("method")->Print("v");
110#endif
111
112}
#define DEBUG
#define e(i)
Definition RSha256.hxx:103
void Warning(const char *location, const char *msgfmt,...)
Use this function in warning situations.
Definition TError.cxx:229
bool setRealValue(const char *name, double newVal=0.0, bool verbose=false)
Set value of a RooAbsRealLValue stored in set with given name to newVal No error messages are printed...
static RooNumIntConfig * defaultIntegratorConfig()
Returns the default numeric integration configuration for all RooAbsReals.
bool setLabel(const char *label, bool printError=true) override
Set value by specifying the name of the desired state.
Landau distribution p.d.f.
Definition RooLandau.h:24
static RooMsgService & instance()
Return reference to singleton instance.
Int_t addStream(RooFit::MsgLevel level, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={})
Add a message logging stream for message with given RooFit::MsgLevel or higher.
Holds the configuration parameters of the various numeric integrators used by RooRealIntegral.
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
void setEpsRel(double newEpsRel)
Set relative convergence criteria (convergence if std::abs(Err)/abs(Int)<newEpsRel)
const RooArgSet & getConfigSection(const char *name) const
Retrieve configuration information specific to integrator with given name.
RooCategory & method1D()
void setEpsAbs(double newEpsAbs)
Set absolute convergence criteria (convergence if std::abs(Err)<newEpsAbs)
Variable that can be changed from the outside.
Definition RooRealVar.h:37
RooCmdArg Topic(Int_t topic)
RooCmdArg NumIntConfig(const RooNumIntConfig &cfg)
Double_t x[n]
Definition legend1.C:17
double landau(double x, double mu, double sigma)
Definition MathFuncs.h:336
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition CodegenImpl.h:64