Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooRombergIntegrator.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/**
18\file RooRombergIntegrator.cxx
19\class RooRombergIntegrator
20\ingroup Roofitcore
21
22Adaptive numerical integration algorithm.
23
24It uses Romberg's method with trapezoids or midpoints.
25The integrand is approximated by \f$ 1, 2, 4, 8, \ldots, 2^n \f$ trapezoids, and
26Richardson series acceleration is applied to the series. For refinement step \f$ n \f$, that is
27\f[
28 R(n,m) = \frac{1}{4^m - 1} \left( 4^m R(n, m-1) - R(n-1, m-1) \right)
29\f]
30The integrator will evaluate the first six refinements (i.e. 64 points) in one go,
31apply a five-point series acceleration, and successively add more steps until the
32desired precision is reached.
33**/
34
35#include "Riostream.h"
36
37#include "TClass.h"
39#include "RooArgSet.h"
40#include "RooRealVar.h"
41#include "RooNumber.h"
42#include "RooNumIntConfig.h"
43#include "RooNumIntFactory.h"
44#include "RooMsgService.h"
45
46#include <cassert>
47
48namespace {
49
50constexpr int nPoints = 5;
51
52////////////////////////////////////////////////////////////////////////////////
53/// Calculate the n-th stage of refinement of the Second Euler-Maclaurin
54/// summation rule which has the useful property of not evaluating the
55/// integrand at either of its endpoints but requires more function
56/// evaluations than the trapezoidal rule. This rule can be used with
57/// a suitable change of variables to estimate improper integrals.
58
59double addMidpoints(std::function<double(double)> const &func, double savedResult, int n, double xmin, double xmax)
60{
61 const double range = xmax - xmin;
62
63 if (n == 1) {
64 double xmid = 0.5 * (xmin + xmax);
65 return range * func(xmid);
66 }
67
68 int it = 1;
69 for (int j = 1; j < n - 1; j++) {
70 it *= 3;
71 }
72 double tnm = it;
73 double del = range / (3. * tnm);
74 double ddel = del + del;
75 double x = xmin + 0.5 * del;
76 double sum = 0;
77 for (int j = 1; j <= it; j++) {
78 sum += func(x);
79 x += ddel;
80 sum += func(x);
81 x += del;
82 }
83 return (savedResult + range * sum / tnm) / 3.;
84}
85
86////////////////////////////////////////////////////////////////////////////////
87/// Calculate the n-th stage of refinement of the extended trapezoidal
88/// summation rule. This is the most efficient rule for a well behaved
89/// integrands that can be evaluated over its entire range, including the
90/// endpoints.
91
92double addTrapezoids(std::function<double(double)> const &func, double savedResult, int n, double xmin, double xmax)
93{
94 const double range = xmax - xmin;
95
96 if (n == 1) {
97 // use a single trapezoid to cover the full range
98 return 0.5 * range * (func(xmin) + func(xmax));
99 }
100
101 // break the range down into several trapezoids using 2**(n-2)
102 // equally-spaced interior points
103 const int nInt = 1 << (n - 2);
104 const double del = range / nInt;
105
106 double sum = 0.;
107 // TODO Replace by batch computation
108 for (int j = 0; j < nInt; ++j) {
109 double x = xmin + (0.5 + j) * del;
110 sum += func(x);
111 }
112
113 return 0.5 * (savedResult + range * sum / nInt);
114}
115
116////////////////////////////////////////////////////////////////////////////////
117/// Extrapolate result to final value.
118
119std::pair<double, double> extrapolate(int n, double const *h, double const *s, double *c, double *d)
120{
121 double const *xa = &h[n - nPoints];
122 double const *ya = &s[n - nPoints];
123 int ns = 1;
124
125 double dif = std::abs(xa[1]);
126 for (int i = 1; i <= nPoints; i++) {
127 double dift = std::abs(xa[i]);
128 if (dift < dif) {
129 ns = i;
130 dif = dift;
131 }
132 c[i] = ya[i];
133 d[i] = ya[i];
134 }
135 double extrapError = 0.0;
136 double extrapValue = ya[ns--];
137 for (int m = 1; m < nPoints; m++) {
138 for (int i = 1; i <= nPoints - m; i++) {
139 double ho = xa[i];
140 double hp = xa[i + m];
141 double w = c[i + 1] - d[i];
142 double den = ho - hp;
143 if (den == 0.0) {
144 throw std::runtime_error("RooRombergIntegrator::extrapolate: internal error");
145 }
146 den = w / den;
147 d[i] = hp * den;
148 c[i] = ho * den;
149 }
150 extrapError = 2 * ns < (nPoints - m) ? c[ns + 1] : d[ns--];
151 extrapValue += extrapError;
152 }
153 return {extrapValue, extrapError};
154}
155
156} // namespace
157
158namespace RooFit {
159namespace Detail {
160
161std::pair<double, int> integrate1d(std::function<double(double)> func, bool doTrapezoid, int maxSteps, int minStepsZero,
162 int fixSteps, double epsAbs, double epsRel, bool doExtrap, double xmin, double xmax,
163 std::span<double> hArr, std::span<double> sArr)
164{
165 assert(int(hArr.size()) == maxSteps + 2);
166 assert(int(sArr.size()) == maxSteps + 2);
167
168 const double range = xmax - xmin;
169
170 // Small working arrays can be on the stack
171 std::array<double, nPoints + 1> cArr = {};
172 std::array<double, nPoints + 1> dArr = {};
173
174 hArr[1] = 1.0;
175 double zeroThresh = epsAbs / range;
176 for (int j = 1; j <= maxSteps; ++j) {
177 // refine our estimate using the appropriate summation rule
178 sArr[j] =
179 doTrapezoid ? addTrapezoids(func, sArr[j - 1], j, xmin, xmax) : addMidpoints(func, sArr[j - 1], j, xmin, xmax);
180
181 if (j >= minStepsZero) {
182 bool allZero(true);
183 for (int jj = 0; jj <= j; jj++) {
184 if (sArr[j] >= zeroThresh) {
185 allZero = false;
186 }
187 }
188 if (allZero) {
189 // cout << "Roo1DIntegrator(" << this << "): zero convergence at step " << j << ", value = " << 0 <<
190 // std::endl ;
191 return {0, j};
192 }
193 }
194
195 if (fixSteps > 0) {
196
197 // Fixed step mode, return result after fixed number of steps
198 if (j == fixSteps) {
199 // cout << "returning result at fixed step " << j << std::endl ;
200 return {sArr[j], j};
201 }
202
203 } else if (j >= nPoints) {
204
205 double extrapValue = 0.0;
206 double extrapError = 0.0;
207
208 // extrapolate the results of recent refinements and check for a stable result
209 if (doExtrap) {
210 std::tie(extrapValue, extrapError) = extrapolate(j, hArr.data(), sArr.data(), cArr.data(), dArr.data());
211 } else {
212 extrapValue = sArr[j];
213 extrapError = sArr[j] - sArr[j - 1];
214 }
215
216 if (std::abs(extrapError) <= epsRel * std::abs(extrapValue)) {
217 return {extrapValue, j};
218 }
219 if (std::abs(extrapError) <= epsAbs) {
220 return {extrapValue, j};
221 }
222 }
223 // update the step size for the next refinement of the summation
224 hArr[j + 1] = doTrapezoid ? hArr[j] / 4. : hArr[j] / 9.;
225 }
226
227 return {sArr[maxSteps], maxSteps};
228}
229
230} // namespace Detail
231} // namespace RooFit
232
233
234// Register this class with RooNumIntConfig
235
236////////////////////////////////////////////////////////////////////////////////
237/// Register integrator plugins, their parameters and capabilities with RooNumIntFactory.
238
240{
241 RooCategory sumRule("sumRule", "Summation Rule");
242 sumRule.defineType("Trapezoid", RooRombergIntegrator::Trapezoid);
243 sumRule.defineType("Midpoint", RooRombergIntegrator::Midpoint);
244 sumRule.setLabel("Trapezoid");
245 RooCategory extrap("extrapolation", "Extrapolation procedure");
246 extrap.defineType("None", 0);
247 extrap.defineType("Wynn-Epsilon", 1);
248 extrap.setLabel("Wynn-Epsilon");
249 RooRealVar maxSteps("maxSteps", "Maximum number of steps", 20);
250 RooRealVar minSteps("minSteps", "Minimum number of steps", 999);
251 RooRealVar fixSteps("fixSteps", "Fixed number of steps", 0);
252 RooRealVar numSeg("numSeg", "Number of segments", 3); // only for the segmented integrators
253
254 std::string name = "RooIntegrator1D";
255
256 auto creator = [](const RooAbsFunc &function, const RooNumIntConfig &config) {
257 return std::make_unique<RooRombergIntegrator>(function, config, 1, /*doSegmentation=*/false);
258 };
259
260 fact.registerPlugin(name, creator, {sumRule, extrap, maxSteps, minSteps, fixSteps},
261 /*canIntegrate1D=*/true,
262 /*canIntegrate2D=*/false,
263 /*canIntegrateND=*/false,
264 /*canIntegrateOpenEnded=*/false);
265
267
268 auto creator2d = [](const RooAbsFunc &function, const RooNumIntConfig &config) {
269 return std::make_unique<RooRombergIntegrator>(function, config, 2, /*doSegmentation=*/false);
270 };
271 std::string name2d = "RooIntegrator2D";
272 fact.registerPlugin(name2d, creator2d, {},
273 /*canIntegrate1D=*/false,
274 /*canIntegrate2D=*/true,
275 /*canIntegrateND=*/false,
276 /*canIntegrateOpenEnded=*/false,
277 /*depName=*/"RooIntegrator1D");
279
280 auto creatorSeg = [](const RooAbsFunc &function, const RooNumIntConfig &config) {
281 return std::make_unique<RooRombergIntegrator>(function, config, 1, /*doSegmentation=*/true);
282 };
283
284 fact.registerPlugin("RooSegmentedIntegrator1D", creatorSeg, {numSeg},
285 /*canIntegrate1D=*/true,
286 /*canIntegrate2D=*/false,
287 /*canIntegrateND=*/false,
288 /*canIntegrateOpenEnded=*/false,
289 /*depName=*/"RooIntegrator1D");
290
291 auto creatorSeg2d = [](const RooAbsFunc &function, const RooNumIntConfig &config) {
292 return std::make_unique<RooRombergIntegrator>(function, config, 2, /*doSegmentation=*/true);
293 };
294
295 fact.registerPlugin("RooSegmentedIntegrator2D", creatorSeg2d, {},
296 /*canIntegrate1D=*/false,
297 /*canIntegrate2D=*/true,
298 /*canIntegrateND=*/false,
299 /*canIntegrateOpenEnded=*/false,
300 /*depName=*/"RooSegmentedIntegrator1D");
301}
302
303////////////////////////////////////////////////////////////////////////////////
304/// Construct integrator on given function binding, using specified summation
305/// rule, maximum number of steps and conversion tolerance. The integration
306/// limits are taken from the function binding.
307
308RooRombergIntegrator::RooRombergIntegrator(const RooAbsFunc &function, SummationRule rule, int maxSteps, double eps)
309 : RooAbsIntegrator(function), _useIntegrandLimits(true), _rule(rule), _maxSteps(maxSteps), _epsAbs(eps), _epsRel(eps)
310{
311 _valid = initialize();
312}
313
314////////////////////////////////////////////////////////////////////////////////
315/// Construct integrator on given function binding for given range,
316/// using specified summation rule, maximum number of steps and
317/// conversion tolerance. The integration limits are taken from the
318/// function binding.
319
321 int maxSteps, double eps)
322 : RooAbsIntegrator(function),
323 _useIntegrandLimits(false),
324 _rule(rule),
325 _maxSteps(maxSteps),
326 _epsAbs(eps),
327 _epsRel(eps)
328{
329 _xmin.push_back(xmin);
330 _xmax.push_back(xmax);
331 _valid = initialize();
332}
333
334////////////////////////////////////////////////////////////////////////////////
335/// Construct integrator on given function binding, using specified
336/// configuration object. The integration limits are taken from the
337/// function binding
338
340 bool doSegmentation)
341 : RooAbsIntegrator(function, config.printEvalCounter()),
342 _nDim{nDim},
343 _epsAbs(config.epsAbs()),
344 _epsRel(config.epsRel())
345{
346 // Extract parameters from config object
347 const RooArgSet &configSet = config.getConfigSection("RooIntegrator1D");
348 _rule = (SummationRule)configSet.getCatIndex("sumRule", Trapezoid);
349 _maxSteps = (int)configSet.getRealValue("maxSteps", 20);
350 _minStepsZero = (int)configSet.getRealValue("minSteps", 999);
351 _fixSteps = (int)configSet.getRealValue("fixSteps", 0);
352 _doExtrap = (bool)configSet.getCatIndex("extrapolation", 1);
353 if (doSegmentation) {
354 _nSeg = (int)config.getConfigSection("RooSegmentedIntegrator1D").getRealValue("numSeg", 3);
355 _epsAbs /= std::sqrt(_nSeg);
356 _epsRel /= std::sqrt(_nSeg);
357 }
358
359 if (_fixSteps > _maxSteps) {
360 oocoutE(nullptr, Integration) << "RooRombergIntegrator::ctor() ERROR: fixSteps>maxSteps, fixSteps set to maxSteps"
361 << std::endl;
363 }
364
365 _useIntegrandLimits = true;
366 _valid = initialize();
367}
368
369////////////////////////////////////////////////////////////////////////////////
370/// Construct integrator on given function binding, using specified
371/// configuration object and integration range
372
374 const RooNumIntConfig &config, int nDim)
375 : RooAbsIntegrator(function, config.printEvalCounter()),
376 _useIntegrandLimits(false),
377 _nDim{nDim},
378 _epsAbs(config.epsAbs()),
379 _epsRel(config.epsRel())
380{
381 // Extract parameters from config object
382 const RooArgSet &configSet = config.getConfigSection("RooIntegrator1D");
383 _rule = (SummationRule)configSet.getCatIndex("sumRule", Trapezoid);
384 _maxSteps = (int)configSet.getRealValue("maxSteps", 20);
385 _minStepsZero = (int)configSet.getRealValue("minSteps", 999);
386 _fixSteps = (int)configSet.getRealValue("fixSteps", 0);
387 _doExtrap = (bool)configSet.getCatIndex("extrapolation", 1);
388
389 _xmin.push_back(xmin);
390 _xmax.push_back(xmax);
391 _valid = initialize();
392}
393
394////////////////////////////////////////////////////////////////////////////////
395/// Initialize the integrator
396
398{
399 // apply defaults if necessary
400 if (_maxSteps <= 0) {
401 _maxSteps = (_rule == Trapezoid) ? 20 : 14;
402 }
403
404 if (_epsRel <= 0)
405 _epsRel = 1e-6;
406 if (_epsAbs <= 0)
407 _epsAbs = 1e-6;
408
409 // check that the integrand is a valid function
410 if (!isValid()) {
411 oocoutE(nullptr, Integration) << "RooRombergIntegrator::initialize: cannot integrate invalid function"
412 << std::endl;
413 return false;
414 }
415
416 // Allocate coordinate buffer size after number of function dimensions
417 _x.resize(_function->getDimension());
418
419 // Allocate workspace for numerical integration engine
420 _wksp.resize(_nDim * 2 * _maxSteps + 4);
421
422 return checkLimits();
423}
424
425////////////////////////////////////////////////////////////////////////////////
426/// Change our integration limits. Return true if the new limits are
427/// ok, or otherwise false. Always returns false and does nothing
428/// if this object was constructed to always use our integrand's limits.
429
431{
433 oocoutE(nullptr, Integration) << "RooRombergIntegrator::setLimits: cannot override integrand's limits"
434 << std::endl;
435 return false;
436 }
437 _xmin.resize(_nDim);
438 _xmax.resize(_nDim);
439 for (int iDim = 0; iDim < _nDim; ++iDim) {
440 _xmin[iDim] = xmin[iDim];
441 _xmax[iDim] = xmax[iDim];
442 }
443 return checkLimits();
444}
445
446////////////////////////////////////////////////////////////////////////////////
447/// Check that our integration range is finite and otherwise return false.
448/// Update the limits from the integrand if requested.
449
451{
453 assert(nullptr != integrand() && integrand()->isValid());
454 const_cast<std::vector<double> &>(_xmin).resize(_nDim);
455 const_cast<std::vector<double> &>(_xmax).resize(_nDim);
456 for (int iDim = 0; iDim < _nDim; ++iDim) {
457 const_cast<double &>(_xmin[iDim]) = integrand()->getMinLimit(iDim);
458 const_cast<double &>(_xmax[iDim]) = integrand()->getMaxLimit(iDim);
459 }
460 }
461 for (int iDim = 0; iDim < _nDim; ++iDim) {
462 const double xmin = _xmin[iDim];
463 const double xmax = _xmax[iDim];
464 const double range = xmax - xmin;
465 if (range < 0.) {
466 oocoutE(nullptr, Integration) << "RooRombergIntegrator::checkLimits: bad range with min > max (_xmin[" << iDim
467 << "] = " << xmin << " _xmax[" << iDim << "] = " << xmax << ")" << std::endl;
468 return false;
469 }
471 return false;
472 }
473 }
474 return true;
475}
476
477double RooRombergIntegrator::integral(const double *yvec)
478{
479 // Copy yvec to xvec if provided
480 if (yvec) {
481 for (unsigned int i = 0; i < _function->getDimension() - 1; i++) {
482 _x[i + _nDim] = yvec[i];
483 }
484 }
485
486 return integral(_nDim - 1, _nSeg, _wksp);
487}
488
489////////////////////////////////////////////////////////////////////////////////
490/// Calculate numeric integral at given set of function binding parameters.
491
492double RooRombergIntegrator::integral(int iDim, int nSeg, std::span<double> wksp)
493{
494 assert(isValid());
495
496 const double xmin = _xmin[iDim];
497 const double xmax = _xmax[iDim];
498 const double range = xmax - xmin;
499
500 if (range == 0.)
501 return 0.;
502
503 // In case of segmentation, split this integral up in a loop.
504 if (nSeg > 1) {
505 const double segSize = (xmax - xmin) / nSeg;
506 double result = 0.0;
507 for (int iSeg = 0; iSeg < nSeg; iSeg++) {
508 _xmin[iDim] = xmin + iSeg * segSize;
509 _xmax[iDim] = xmin + (iSeg + 1) * segSize;
510 double part = integral(iDim, 1, wksp);
511 result += part;
512 }
513
514 // reset limits
515 _xmin[iDim] = xmin;
516 _xmax[iDim] = xmax;
517
518 return result;
519 }
520
521 // From the working array
522 std::size_t nWorkingArr = _maxSteps + 2;
523 double *hArr = wksp.data();
524 double *sArr = wksp.data() + nWorkingArr;
525
526 double output = 0.0;
527 int steps = 0;
528
529 std::span<double> nextWksp{wksp.data() + 2 * _maxSteps + 4, wksp.data() + wksp.size()};
530
531 auto func = [&](double x) {
532 _x[iDim] = x;
533
534 return iDim == 0 ? integrand(_x.data()) : integral(iDim - 1, _nSeg, nextWksp);
535 };
536
537 std::tie(output, steps) =
539 _doExtrap, xmin, xmax, {hArr, nWorkingArr}, {sArr, nWorkingArr});
540
541 if (steps == _maxSteps) {
542
543 oocoutW(nullptr, Integration) << "RooRombergIntegrator::integral: integral of " << _function->getName()
544 << " over range (" << xmin << "," << xmax << ") did not converge after "
545 << _maxSteps << " steps" << std::endl;
546 for (int j = 1; j <= _maxSteps; ++j) {
547 ooccoutW(nullptr, Integration) << " [" << j << "] h = " << hArr[j] << " , s = " << sArr[j] << std::endl;
548 }
549 }
550
551 return output;
552}
#define d(i)
Definition RSha256.hxx:102
#define c(i)
Definition RSha256.hxx:101
#define h(i)
Definition RSha256.hxx:106
#define e(i)
Definition RSha256.hxx:103
#define oocoutW(o, a)
#define oocoutE(o, a)
#define ooccoutW(o, a)
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t del
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
char name[80]
Definition TGX11.cxx:110
float xmin
float xmax
double getRealValue(const char *name, double defVal=0.0, bool verbose=false) const
Get value of a RooAbsReal stored in set with given name.
Int_t getCatIndex(const char *name, Int_t defVal=0, bool verbose=false) const
Get index value of a RooAbsCategory stored in set with given name.
Abstract interface for evaluating a real-valued function of one real variable and performing numerica...
Definition RooAbsFunc.h:27
virtual double getMaxLimit(UInt_t dimension) const =0
virtual double getMinLimit(UInt_t dimension) const =0
UInt_t getDimension() const
Definition RooAbsFunc.h:33
virtual const char * getName() const
Name of function binding.
Definition RooAbsFunc.h:65
Abstract interface for integrators of real-valued functions that implement the RooAbsFunc interface.
bool isValid() const
Is integrator in valid state.
const RooAbsFunc * _function
Pointer to function binding of integrand.
const RooAbsFunc * integrand() const
Return integrand function binding.
bool _valid
Is integrator in valid state?
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
Object to represent discrete states.
Definition RooCategory.h:28
bool defineType(const std::string &label)
Define a state with given name.
bool setLabel(const char *label, bool printError=true) override
Set value by specifying the name of the desired state.
Holds the configuration parameters of the various numeric integrators used by RooRealIntegral.
RooCategory & method2D()
const RooArgSet & getConfigSection(const char *name) const
Retrieve configuration information specific to integrator with given name.
RooCategory & method1D()
static RooNumIntConfig & defaultConfig()
Return reference to instance of default numeric integrator configuration object.
Factory to instantiate numeric integrators from a given function binding and a given configuration.
bool registerPlugin(std::string const &name, Creator const &creator, const RooArgSet &defConfig, bool canIntegrate1D, bool canIntegrate2D, bool canIntegrateND, bool canIntegrateOpenEnded, const char *depName="")
Method accepting registration of a prototype numeric integrator along with a RooArgSet of its default...
static constexpr int isInfinite(double x)
Return true if x is infinite by RooNumber internal specification.
Definition RooNumber.h:27
Variable that can be changed from the outside.
Definition RooRealVar.h:37
std::vector< double > _xmin
! Lower integration bounds
std::vector< double > _wksp
! Integrator workspace
RooRombergIntegrator(const RooAbsFunc &function, SummationRule rule=Trapezoid, int maxSteps=0, double eps=0)
Construct integrator on given function binding, using specified summation rule, maximum number of ste...
bool _useIntegrandLimits
If true limits of function binding are used.
static void registerIntegrator(RooNumIntFactory &fact)
Register integrator plugins, their parameters and capabilities with RooNumIntFactory.
double _epsRel
Relative convergence tolerance.
bool checkLimits() const override
Check that our integration range is finite and otherwise return false.
bool initialize()
Initialize the integrator.
double integral(const double *yvec=nullptr) override
bool _doExtrap
Apply conversion step?
int _nSeg
Number of segments.
std::vector< double > _xmax
! Upper integration bounds
bool setLimits(double *xmin, double *xmax) override
Change our integration limits.
int _fixSteps
Fixed number of steps.
double _epsAbs
Absolute convergence tolerance.
int _maxSteps
Maximum number of steps.
int _minStepsZero
Minimum number of steps to declare convergence to zero.
std::vector< double > _x
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
std::pair< double, int > integrate1d(std::function< double(double)> func, bool doTrapezoid, int maxSteps, int minStepsZero, int fixSteps, double epsAbs, double epsRel, bool doExtrap, double xmin, double xmax, std::span< double > hArr, std::span< double > sArr)
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition JSONIO.h:26
TMarker m
Definition textangle.C:8
static uint64_t sum(uint64_t i)
Definition Factory.cxx:2345
static void output()