Logo ROOT  
Reference Guide
RooGaussKronrodIntegrator1D.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 RooGaussKronrodIntegrator1D.cxx
19\class RooGaussKronrodIntegrator1D
20\ingroup Roofitcore
21
22RooGaussKronrodIntegrator1D implements the Gauss-Kronrod integration algorithm.
23
24An Gaussian quadrature method for numerical integration in which
25error is estimation based on evaluation at special points known as
26"Kronrod points." By suitably picking these points, abscissas from
27previous iterations can be reused as part of the new set of points,
28whereas usual Gaussian quadrature would require recomputation of
29all abscissas at each iteration.
30
31This class automatically handles (-inf,+inf) integrals by dividing
32the integration in three regions (-inf,-1), (-1,1), (1,inf) and
33calculating the 1st and 3rd term using a x -> 1/x coordinate
34transformation
35
36This class embeds the Gauss-Kronrod integrator from the GNU
37Scientific Library version 1.5 and applies the 10-, 21-, 43- and
3887-point rule in succession until the required target precision is
39reached
40**/
41
42#include <assert.h>
43#include <math.h>
44#include <float.h>
45#include "Riostream.h"
46#include "TMath.h"
48#include "RooArgSet.h"
49#include "RooRealVar.h"
50#include "RooNumber.h"
51#include "RooNumIntFactory.h"
53#include "RooMsgService.h"
54
55
56using namespace std;
57
59;
60
61
62// --- From GSL_MATH.h -------------------------------------------
64{
65 double (* function) (double x, void * params);
66 void * params;
67};
69#define GSL_FN_EVAL(F,x) (*((F)->function))(x,(F)->params)
70//----From GSL_INTEGRATION.h ---------------------------------------
72 double a, double b,
73 double epsabs, double epsrel,
74 double *result, double *abserr,
75 size_t * neval);
76//-------------------------------------------------------------------
77
78// register integrator class
79// create a derived class in order to call the protected method of the
80// RoodaptiveGaussKronrodIntegrator1D
81namespace RooFit_internal {
83
84 static void registerIntegrator()
85 {
86 auto &intFactory = RooNumIntFactory::instance();
88 }
89};
90// class used to register integrator at loafing time
93};
94
96}
97
98
99////////////////////////////////////////////////////////////////////////////////
100/// Register RooGaussKronrodIntegrator1D, its parameters and capabilities with RooNumIntConfig
101
103{
105 oocoutI(nullptr,Integration) << "RooGaussKronrodIntegrator1D has been registered" << std::endl;
106}
107
108
109
110////////////////////////////////////////////////////////////////////////////////
111/// coverity[UNINIT_CTOR]
112/// Default constructor
113
115{
116}
117
118
119
120////////////////////////////////////////////////////////////////////////////////
121/// Construct integral on 'function' using given configuration object. The integration
122/// range is taken from the definition in the function binding
123
126 _epsAbs(config.epsRel()),
127 _epsRel(config.epsAbs())
128{
131}
132
133
134
135////////////////////////////////////////////////////////////////////////////////
136/// Construct integral on 'function' using given configuration object in the given range
137
139 double xmin, double xmax, const RooNumIntConfig& config) :
141 _epsAbs(config.epsRel()),
142 _epsRel(config.epsAbs()),
143 _xmin(xmin),
144 _xmax(xmax)
145{
146 _useIntegrandLimits= false;
148}
149
150
151
152////////////////////////////////////////////////////////////////////////////////
153/// Clone integrator with given function and configuration. Needed for RooNumIntFactory
154
156{
157 return new RooGaussKronrodIntegrator1D(function,config) ;
158}
159
160
161
162////////////////////////////////////////////////////////////////////////////////
163/// Perform one-time initialization of integrator
164
166{
167 // Allocate coordinate buffer size after number of function dimensions
168 _x = new double[_function->getDimension()] ;
169
170 return checkLimits();
171}
172
173
174
175////////////////////////////////////////////////////////////////////////////////
176/// Destructor
177
179{
180 if (_x) {
181 delete[] _x ;
182 }
183}
184
185
186
187////////////////////////////////////////////////////////////////////////////////
188/// Change our integration limits. Return true if the new limits are
189/// ok, or otherwise false. Always returns false and does nothing
190/// if this object was constructed to always use our integrand's limits.
191
193{
195 oocoutE(nullptr,Eval) << "RooGaussKronrodIntegrator1D::setLimits: cannot override integrand's limits" << endl;
196 return false;
197 }
198 _xmin= *xmin;
199 _xmax= *xmax;
200 return checkLimits();
201}
202
203
204
205////////////////////////////////////////////////////////////////////////////////
206/// Check that our integration range is finite and otherwise return false.
207/// Update the limits from the integrand if requested.
208
210{
212 assert(0 != integrand() && integrand()->isValid());
215 }
216 return true ;
217}
218
219
220
222{
224 return instance->integrand(instance->xvec(x)) ;
225}
226
227
228
229////////////////////////////////////////////////////////////////////////////////
230/// Calculate and return integral
231
233{
234 assert(isValid());
235
236 // Copy yvec to xvec if provided
237 if (yvec) {
238 UInt_t i ; for (i=0 ; i<_function->getDimension()-1 ; i++) {
239 _x[i+1] = yvec[i] ;
240 }
241 }
242
243 // Setup glue function
246 F.params = this ;
247
248 // Return values
249 double result, error;
250 size_t neval = 0 ;
251
252 // Call GSL implementation of integeator
253 gsl_integration_qng (&F, _xmin, _xmax, _epsAbs, _epsRel, &result, &error, &neval);
254
255 return result;
256}
257
258
259// ----------------------------------------------------------------------------
260// ---------- Code below imported from GSL version 1.5 ------------------------
261// ----------------------------------------------------------------------------
262
263/*
264 *
265 * Copyright (C) 1996, 1997, 1998, 1999, 2000 Brian Gough
266 *
267 * This program is free software; you can redistribute it and/or modify
268 * it under the terms of the GNU General Public License as published by
269 * the Free Software Foundation; either version 2 of the License, or (at
270 * your option) any later version.
271 *
272 * This program is distributed in the hope that it will be useful, but
273 * WITHOUT ANY WARRANTY; without even the implied warranty of
274 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
275 * General Public License for more details.
276 *
277 * You should have received a copy of the GNU General Public License
278 * along with this program; if not, write to the Free Software
279 * Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
280 */
281
282#define GSL_SUCCESS 0
283#define GSL_EBADTOL 13 /* user specified an invalid tolerance */
284#define GSL_ETOL 14 /* failed to reach the specified tolerance */
285#define GSL_ERROR(a,b) oocoutE(nullptr,Eval) << "RooGaussKronrodIntegrator1D::integral() ERROR: " << a << endl ; return b ;
286#define GSL_DBL_MIN 2.2250738585072014e-308
287#define GSL_DBL_EPSILON 2.2204460492503131e-16
288
289
290// INCLUDED BELOW #include "qng.c"
291
293 double a, double b,
294 double epsabs, double epsrel,
295 double *result, double *abserr,
296 size_t * neval);
297
298
299
300// INCLUDED BELOW #include "err.c"
301static double rescale_error (double err, const double result_abs, const double result_asc) ;
302
303static double
304rescale_error (double err, const double result_abs, const double result_asc)
305{
306 err = fabs(err) ;
307
308 if (result_asc != 0 && err != 0)
309 {
310 double scale = TMath::Power((200 * err / result_asc), 1.5) ;
311
312 if (scale < 1)
313 {
314 err = result_asc * scale ;
315 }
316 else
317 {
318 err = result_asc ;
319 }
320 }
321 if (result_abs > GSL_DBL_MIN / (50 * GSL_DBL_EPSILON))
322 {
323 double min_err = 50 * GSL_DBL_EPSILON * result_abs ;
324
325 if (min_err > err)
326 {
327 err = min_err ;
328 }
329 }
330
331 return err ;
332}
333
334
335// INCLUDED BELOW #include "qng.h"
336/* Gauss-Kronrod-Patterson quadrature coefficients for use in
337 quadpack routine qng. These coefficients were calculated with
338 101 decimal digit arithmetic by L. W. Fullerton, Bell Labs, Nov
339 1981. */
340
341/* x1, abscissae common to the 10-, 21-, 43- and 87-point rule */
342static const double x1[5] = {
343 0.973906528517171720077964012084452,
344 0.865063366688984510732096688423493,
345 0.679409568299024406234327365114874,
346 0.433395394129247190799265943165784,
347 0.148874338981631210884826001129720
348} ;
349
350/* w10, weights of the 10-point formula */
351static const double w10[5] = {
352 0.066671344308688137593568809893332,
353 0.149451349150580593145776339657697,
354 0.219086362515982043995534934228163,
355 0.269266719309996355091226921569469,
356 0.295524224714752870173892994651338
357} ;
358
359/* x2, abscissae common to the 21-, 43- and 87-point rule */
360static const double x2[5] = {
361 0.995657163025808080735527280689003,
362 0.930157491355708226001207180059508,
363 0.780817726586416897063717578345042,
364 0.562757134668604683339000099272694,
365 0.294392862701460198131126603103866
366} ;
367
368/* w21a, weights of the 21-point formula for abscissae x1 */
369static const double w21a[5] = {
370 0.032558162307964727478818972459390,
371 0.075039674810919952767043140916190,
372 0.109387158802297641899210590325805,
373 0.134709217311473325928054001771707,
374 0.147739104901338491374841515972068
375} ;
376
377/* w21b, weights of the 21-point formula for abscissae x2 */
378static const double w21b[6] = {
379 0.011694638867371874278064396062192,
380 0.054755896574351996031381300244580,
381 0.093125454583697605535065465083366,
382 0.123491976262065851077958109831074,
383 0.142775938577060080797094273138717,
384 0.149445554002916905664936468389821
385} ;
386
387/* x3, abscissae common to the 43- and 87-point rule */
388static const double x3[11] = {
389 0.999333360901932081394099323919911,
390 0.987433402908088869795961478381209,
391 0.954807934814266299257919200290473,
392 0.900148695748328293625099494069092,
393 0.825198314983114150847066732588520,
394 0.732148388989304982612354848755461,
395 0.622847970537725238641159120344323,
396 0.499479574071056499952214885499755,
397 0.364901661346580768043989548502644,
398 0.222254919776601296498260928066212,
399 0.074650617461383322043914435796506
400} ;
401
402/* w43a, weights of the 43-point formula for abscissae x1, x3 */
403static const double w43a[10] = {
404 0.016296734289666564924281974617663,
405 0.037522876120869501461613795898115,
406 0.054694902058255442147212685465005,
407 0.067355414609478086075553166302174,
408 0.073870199632393953432140695251367,
409 0.005768556059769796184184327908655,
410 0.027371890593248842081276069289151,
411 0.046560826910428830743339154433824,
412 0.061744995201442564496240336030883,
413 0.071387267268693397768559114425516
414} ;
415
416/* w43b, weights of the 43-point formula for abscissae x3 */
417static const double w43b[12] = {
418 0.001844477640212414100389106552965,
419 0.010798689585891651740465406741293,
420 0.021895363867795428102523123075149,
421 0.032597463975345689443882222526137,
422 0.042163137935191811847627924327955,
423 0.050741939600184577780189020092084,
424 0.058379395542619248375475369330206,
425 0.064746404951445885544689259517511,
426 0.069566197912356484528633315038405,
427 0.072824441471833208150939535192842,
428 0.074507751014175118273571813842889,
429 0.074722147517403005594425168280423
430} ;
431
432/* x4, abscissae of the 87-point rule */
433static const double x4[22] = {
434 0.999902977262729234490529830591582,
435 0.997989895986678745427496322365960,
436 0.992175497860687222808523352251425,
437 0.981358163572712773571916941623894,
438 0.965057623858384619128284110607926,
439 0.943167613133670596816416634507426,
440 0.915806414685507209591826430720050,
441 0.883221657771316501372117548744163,
442 0.845710748462415666605902011504855,
443 0.803557658035230982788739474980964,
444 0.757005730685495558328942793432020,
445 0.706273209787321819824094274740840,
446 0.651589466501177922534422205016736,
447 0.593223374057961088875273770349144,
448 0.531493605970831932285268948562671,
449 0.466763623042022844871966781659270,
450 0.399424847859218804732101665817923,
451 0.329874877106188288265053371824597,
452 0.258503559202161551802280975429025,
453 0.185695396568346652015917141167606,
454 0.111842213179907468172398359241362,
455 0.037352123394619870814998165437704
456} ;
457
458/* w87a, weights of the 87-point formula for abscissae x1, x2, x3 */
459static const double w87a[21] = {
460 0.008148377384149172900002878448190,
461 0.018761438201562822243935059003794,
462 0.027347451050052286161582829741283,
463 0.033677707311637930046581056957588,
464 0.036935099820427907614589586742499,
465 0.002884872430211530501334156248695,
466 0.013685946022712701888950035273128,
467 0.023280413502888311123409291030404,
468 0.030872497611713358675466394126442,
469 0.035693633639418770719351355457044,
470 0.000915283345202241360843392549948,
471 0.005399280219300471367738743391053,
472 0.010947679601118931134327826856808,
473 0.016298731696787335262665703223280,
474 0.021081568889203835112433060188190,
475 0.025370969769253827243467999831710,
476 0.029189697756475752501446154084920,
477 0.032373202467202789685788194889595,
478 0.034783098950365142750781997949596,
479 0.036412220731351787562801163687577,
480 0.037253875503047708539592001191226
481} ;
482
483/* w87b, weights of the 87-point formula for abscissae x4 */
484static const double w87b[23] = {
485 0.000274145563762072350016527092881,
486 0.001807124155057942948341311753254,
487 0.004096869282759164864458070683480,
488 0.006758290051847378699816577897424,
489 0.009549957672201646536053581325377,
490 0.012329447652244853694626639963780,
491 0.015010447346388952376697286041943,
492 0.017548967986243191099665352925900,
493 0.019938037786440888202278192730714,
494 0.022194935961012286796332102959499,
495 0.024339147126000805470360647041454,
496 0.026374505414839207241503786552615,
497 0.028286910788771200659968002987960,
498 0.030052581128092695322521110347341,
499 0.031646751371439929404586051078883,
500 0.033050413419978503290785944862689,
501 0.034255099704226061787082821046821,
502 0.035262412660156681033782717998428,
503 0.036076989622888701185500318003895,
504 0.036698604498456094498018047441094,
505 0.037120549269832576114119958413599,
506 0.037334228751935040321235449094698,
507 0.037361073762679023410321241766599
508} ;
509
510
511int
513 double a, double b,
514 double epsabs, double epsrel,
515 double * result, double * abserr, size_t * neval)
516{
517 double fv1[5], fv2[5], fv3[5], fv4[5];
518 double savfun[21]; /* array of function values which have been computed */
519 double res10, res21, res43, res87; /* 10, 21, 43 and 87 point results */
520 double result_kronrod, err ;
521 double resabs; /* approximation to the integral of abs(f) */
522 double resasc; /* approximation to the integral of abs(f-i/(b-a)) */
523
524 const double half_length = 0.5 * (b - a);
525 const double abs_half_length = fabs (half_length);
526 const double center = 0.5 * (b + a);
527 const double f_center = GSL_FN_EVAL(f, center);
528
529 int k ;
530
531 if (epsabs <= 0 && (epsrel < 50 * GSL_DBL_EPSILON || epsrel < 0.5e-28))
532 {
533 * result = 0;
534 * abserr = 0;
535 * neval = 0;
536 GSL_ERROR ("tolerance cannot be acheived with given epsabs and epsrel",
538 };
539
540 /* Compute the integral using the 10- and 21-point formula. */
541
542 res10 = 0;
543 res21 = w21b[5] * f_center;
544 resabs = w21b[5] * fabs (f_center);
545
546 for (k = 0; k < 5; k++)
547 {
548 const double abscissa = half_length * x1[k];
549 const double fval1 = GSL_FN_EVAL(f, center + abscissa);
550 const double fval2 = GSL_FN_EVAL(f, center - abscissa);
551 const double fval = fval1 + fval2;
552 res10 += w10[k] * fval;
553 res21 += w21a[k] * fval;
554 resabs += w21a[k] * (fabs (fval1) + fabs (fval2));
555 savfun[k] = fval;
556 fv1[k] = fval1;
557 fv2[k] = fval2;
558 }
559
560 for (k = 0; k < 5; k++)
561 {
562 const double abscissa = half_length * x2[k];
563 const double fval1 = GSL_FN_EVAL(f, center + abscissa);
564 const double fval2 = GSL_FN_EVAL(f, center - abscissa);
565 const double fval = fval1 + fval2;
566 res21 += w21b[k] * fval;
567 resabs += w21b[k] * (fabs (fval1) + fabs (fval2));
568 savfun[k + 5] = fval;
569 fv3[k] = fval1;
570 fv4[k] = fval2;
571 }
572
573 resabs *= abs_half_length ;
574
575 {
576 const double mean = 0.5 * res21;
577
578 resasc = w21b[5] * fabs (f_center - mean);
579
580 for (k = 0; k < 5; k++)
581 {
582 resasc +=
583 (w21a[k] * (fabs (fv1[k] - mean) + fabs (fv2[k] - mean))
584 + w21b[k] * (fabs (fv3[k] - mean) + fabs (fv4[k] - mean)));
585 }
586 resasc *= abs_half_length ;
587 }
588
589 result_kronrod = res21 * half_length;
590
591 err = rescale_error ((res21 - res10) * half_length, resabs, resasc) ;
592
593 /* test for convergence. */
594
595 if (err < epsabs || err < epsrel * fabs (result_kronrod))
596 {
597 * result = result_kronrod ;
598 * abserr = err ;
599 * neval = 21;
600 return GSL_SUCCESS;
601 }
602
603 /* compute the integral using the 43-point formula. */
604
605 res43 = w43b[11] * f_center;
606
607 for (k = 0; k < 10; k++)
608 {
609 res43 += savfun[k] * w43a[k];
610 }
611
612 for (k = 0; k < 11; k++)
613 {
614 const double abscissa = half_length * x3[k];
615 const double fval = (GSL_FN_EVAL(f, center + abscissa)
616 + GSL_FN_EVAL(f, center - abscissa));
617 res43 += fval * w43b[k];
618 savfun[k + 10] = fval;
619 }
620
621 /* test for convergence */
622
623 result_kronrod = res43 * half_length;
624 err = rescale_error ((res43 - res21) * half_length, resabs, resasc);
625
626 if (err < epsabs || err < epsrel * fabs (result_kronrod))
627 {
628 * result = result_kronrod ;
629 * abserr = err ;
630 * neval = 43;
631 return GSL_SUCCESS;
632 }
633
634 /* compute the integral using the 87-point formula. */
635
636 res87 = w87b[22] * f_center;
637
638 for (k = 0; k < 21; k++)
639 {
640 res87 += savfun[k] * w87a[k];
641 }
642
643 for (k = 0; k < 22; k++)
644 {
645 const double abscissa = half_length * x4[k];
646 res87 += w87b[k] * (GSL_FN_EVAL(f, center + abscissa)
647 + GSL_FN_EVAL(f, center - abscissa));
648 }
649
650 /* test for convergence */
651
652 result_kronrod = res87 * half_length ;
653
654 err = rescale_error ((res87 - res43) * half_length, resabs, resasc);
655
656 if (err < epsabs || err < epsrel * fabs (result_kronrod))
657 {
658 * result = result_kronrod ;
659 * abserr = err ;
660 * neval = 87;
661 return GSL_SUCCESS;
662 }
663
664 /* failed to converge */
665
666 * result = result_kronrod ;
667 * abserr = err ;
668 * neval = 87;
669
670 // GSL_ERROR("failed to reach tolerance with highest-order rule", GSL_ETOL) ;
671 return GSL_ETOL ;
672}
#define f(i)
Definition: RSha256.hxx:104
#define GSL_ERROR(a, b)
static const double x2[5]
static double rescale_error(double err, const double result_abs, const double result_asc)
static const double w10[5]
#define GSL_EBADTOL
#define GSL_DBL_EPSILON
static const double w43a[10]
static const double x4[22]
static const double w21b[6]
static const double x1[5]
double RooGaussKronrodIntegrator1D_GSL_GlueFunction(double x, void *data)
static const double w87a[21]
static const double w21a[5]
static const double w43b[12]
static const double x3[11]
#define GSL_SUCCESS
#define GSL_FN_EVAL(F, x)
static const double w87b[23]
int gsl_integration_qng(const gsl_function *f, double a, double b, double epsabs, double epsrel, double *result, double *abserr, size_t *neval)
#define GSL_DBL_MIN
#define oocoutE(o, a)
Definition: RooMsgService.h:52
#define oocoutI(o, a)
Definition: RooMsgService.h:49
#define ClassImp(name)
Definition: Rtypes.h:375
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
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 b
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
float xmin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
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
RooAbsIntegrator is the abstract interface for integrators of real-valued functions that implement th...
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:57
RooGaussKronrodIntegrator1D implements the Gauss-Kronrod integration algorithm.
double integral(const double *yvec=0) override
Calculate and return integral.
double _xmax
Lower integration bound.
~RooGaussKronrodIntegrator1D() override
Destructor.
friend double RooGaussKronrodIntegrator1D_GSL_GlueFunction(double x, void *data)
bool initialize()
Perform one-time initialization of integrator.
RooGaussKronrodIntegrator1D()
coverity[UNINIT_CTOR] Default constructor
bool checkLimits() const override
Check that our integration range is finite and otherwise return false.
RooAbsIntegrator * clone(const RooAbsFunc &function, const RooNumIntConfig &config) const override
Clone integrator with given function and configuration. Needed for RooNumIntFactory.
bool setLimits(double *xmin, double *xmax) override
Change our integration limits.
static void registerIntegrator(RooNumIntFactory &fact)
Register RooGaussKronrodIntegrator1D, its parameters and capabilities with RooNumIntConfig.
RooNumIntConfig holds the configuration parameters of the various numeric integrators used by RooReal...
RooNumIntFactory is a factory to instantiate numeric integrators from a given function binding and a ...
bool storeProtoIntegrator(RooAbsIntegrator *proto, const RooArgSet &defConfig, const char *depName="")
Method accepting registration of a prototype numeric integrator along with a RooArgSet of its default...
static RooNumIntFactory & instance()
Static method returning reference to singleton instance of factory.
Double_t x[n]
Definition: legend1.C:17
#define F(x, y, z)
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition: RExports.h:167
static Roo_reg_GKInteg1D instance
static Roo_reg_AGKInteg1D instance
@ Integration
Definition: RooGlobalFunc.h:63
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:685
double(* function)(double x, void *params)
auto * a
Definition: textangle.C:12