Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooImproperIntegrator1D.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 RooImproperIntegrator1D.cxx
19\class RooImproperIntegrator1D
20\ingroup Roofitcore
21
22Special numeric integrator that can handle integrals over open domains.
23To this end the range is cut in up three pieces: [-inf,-1],[-1,+1] and [+1,inf]
24and the outer two pieces, if required are calculated using a 1/x transform
25**/
26
27
28#include "RooFit.h"
29
31#include "RooIntegrator1D.h"
32#include "RooInvTransform.h"
33#include "RooNumber.h"
34#include "RooNumIntFactory.h"
35#include "RooArgSet.h"
36#include "RooMsgService.h"
37
38#include "Riostream.h"
39#include <math.h>
40#include "TClass.h"
41
42
43
44using namespace std;
45
47;
48
49// Register this class with RooNumIntConfig
50
51////////////////////////////////////////////////////////////////////////////////
52/// Register RooImproperIntegrator1D, its parameters and capabilities with RooNumIntFactory
53
55{
57 fact.storeProtoIntegrator(proto,RooArgSet(),RooIntegrator1D::Class()->GetName()) ;
58}
59
60
61
62////////////////////////////////////////////////////////////////////////////////
63/// Default constructor
64
66 _case(ClosedBothEnds), _xmin(-10), _xmax(10), _useIntegrandLimits(kTRUE),
67 _origFunc(0), _function(0), _integrator1(0), _integrator2(0), _integrator3(0)
68{
69}
70
71
72////////////////////////////////////////////////////////////////////////////////
73/// Constructor with function binding. The integration range is taken from the
74/// definition in the function binding
75
77 RooAbsIntegrator(function),
78 _useIntegrandLimits(kTRUE),
79 _origFunc((RooAbsFunc*)&function),
80 _function(0),
81 _integrator1(0),
82 _integrator2(0),
83 _integrator3(0)
84{
85 initialize(&function) ;
86}
87
88
89
90////////////////////////////////////////////////////////////////////////////////
91/// Constructor with function binding and configuration object. The integration range is taken
92/// from the definition in the function binding
93
95 RooAbsIntegrator(function),
96 _useIntegrandLimits(kTRUE),
97 _origFunc((RooAbsFunc*)&function),
98 _function(0),
99 _config(config),
100 _integrator1(0),
101 _integrator2(0),
102 _integrator3(0)
103{
104 initialize(&function) ;
105}
106
107
108
109////////////////////////////////////////////////////////////////////////////////
110/// Constructor with function binding, definition of integration range and configuration object
111
113 RooAbsIntegrator(function),
114 _xmin(xmin),
115 _xmax(xmax),
116 _useIntegrandLimits(kFALSE),
117 _origFunc((RooAbsFunc*)&function),
118 _function(0),
119 _config(config),
120 _integrator1(0),
121 _integrator2(0),
122 _integrator3(0)
123{
124 initialize(&function) ;
125}
126
127
128
129////////////////////////////////////////////////////////////////////////////////
130/// Return clone of integrator with given function and configuration. Needed by RooNumIntFactory.
131
133{
134 return new RooImproperIntegrator1D(function,config) ;
135}
136
137
138
139////////////////////////////////////////////////////////////////////////////////
140/// Initialize the integrator, construct and initialize subintegrators
141
143{
144 if(!isValid()) {
145 oocoutE((TObject*)0,Integration) << "RooImproperIntegrator: cannot integrate invalid function" << endl;
146 return;
147 }
148 // Create a new function object that uses the change of vars: x -> 1/x
149 if (function) {
150 _function= new RooInvTransform(*function);
151 } else {
152 function = _origFunc ;
153 if (_integrator1) {
154 delete _integrator1 ;
155 _integrator1 = 0 ;
156 }
157 if (_integrator2) {
158 delete _integrator2 ;
159 _integrator2 = 0 ;
160 }
161 if (_integrator3) {
162 delete _integrator3 ;
163 _integrator3 = 0 ;
164 }
165 }
166
167 // partition the integration range into subranges that can each be
168 // handled by RooIntegrator1D
169 switch(_case= limitsCase()) {
170 case ClosedBothEnds:
171 // both limits are finite: use the plain trapezoid integrator
173 break;
174 case OpenBothEnds:
175 // both limits are infinite: integrate over (-1,+1) using
176 // the plain trapezoid integrator...
178 // ...and integrate the infinite tails using the midpoint integrator
181 break;
183 // xmax >= 0 so integrate from (-inf,-1) and (-1,xmax)
186 break;
187 case OpenBelow:
188 // xmax < 0 so integrate from (-inf,xmax)
190 break;
192 // xmin <= 0 so integrate from (xmin,+1) and (+1,+inf)
195 break;
196 case OpenAbove:
197 // xmin > 0 so integrate from (xmin,+inf)
199 break;
200 case Invalid:
201 default:
202 _valid= kFALSE;
203 }
204}
205
206
207////////////////////////////////////////////////////////////////////////////////
208/// Destructor
209
211{
212 if(0 != _integrator1) delete _integrator1;
213 if(0 != _integrator2) delete _integrator2;
214 if(0 != _integrator3) delete _integrator3;
215 if(0 != _function) delete _function;
216}
217
218
219////////////////////////////////////////////////////////////////////////////////
220/// Change our integration limits. Return kTRUE if the new limits are
221/// ok, or otherwise kFALSE. Always returns kFALSE and does nothing
222/// if this object was constructed to always use our integrand's limits.
223
225{
227 oocoutE((TObject*)0,Integration) << "RooIntegrator1D::setLimits: cannot override integrand's limits" << endl;
228 return kFALSE;
229 }
230
231 _xmin= *xmin;
232 _xmax= *xmax;
233 return checkLimits();
234}
235
236
237////////////////////////////////////////////////////////////////////////////////
238/// Check if the limits are valid. For this integrator all limit configurations
239/// are valid, but if the limits change between two calculate() calls it
240/// may be necessary to reconfigure (e.g. if an open ended range becomes
241/// a closed range
242
244{
245 // Has either limit changed?
247 if(_xmin == integrand()->getMinLimit(0) &&
248 _xmax == integrand()->getMaxLimit(0)) return kTRUE;
249 }
250
251 // The limits have changed: can we use the same strategy?
252 if(limitsCase() != _case) {
253 // Reinitialize embedded integrators, will automatically propagate new limits
254 const_cast<RooImproperIntegrator1D*>(this)->initialize() ;
255 return kTRUE ;
256 }
257
258 // Reuse our existing integrators by updating their limits
259 switch(_case) {
260 case ClosedBothEnds:
262 break;
263 case OpenBothEnds:
264 // nothing has changed
265 break;
268 break;
269 case OpenBelow:
271 break;
274 break;
275 case OpenAbove:
277 break;
278 case Invalid:
279 default:
280 return kFALSE;
281 }
282 return kTRUE;
283}
284
285
286////////////////////////////////////////////////////////////////////////////////
287/// Classify the type of limits we have: OpenBothEnds,ClosedBothEnds,OpenBelow or OpenAbove.
288
290{
291 // Analyze the specified limits to determine which case applies.
292 if(0 == integrand() || !integrand()->isValid()) return Invalid;
293
297 }
298
301 if(!inf1 && !inf2) {
302 // both limits are finite
303 return ClosedBothEnds;
304 }
305 else if(inf1 && inf2) {
306 // both limits are infinite
307 return OpenBothEnds;
308 }
309 else if(inf1) { // inf2==false
310 if(_xmax >= 0) {
311 return OpenBelowSpansZero;
312 }
313 else {
314 return OpenBelow;
315 }
316 }
317 else { // inf1==false && inf2==true
318 if(_xmin <= 0) {
319 return OpenAboveSpansZero;
320 }
321 else {
322 return OpenAbove;
323 }
324 }
325 // return Invalid; OSF-CC: Statement unreachable
326}
327
328
329////////////////////////////////////////////////////////////////////////////////
330/// Calculate the integral at the given parameter values of the function binding
331
333{
334 Double_t result(0);
335 if(0 != _integrator1) result+= _integrator1->integral(yvec);
336 if(0 != _integrator2) result+= _integrator2->integral(yvec);
337 if(0 != _integrator3) result+= _integrator3->integral(yvec);
338 return result;
339}
#define oocoutE(o, a)
const Bool_t kFALSE
Definition RtypesCore.h:101
const Bool_t kTRUE
Definition RtypesCore.h:100
#define ClassImp(name)
Definition Rtypes.h:364
float xmin
float xmax
const char * proto
Definition civetweb.c:16613
Abstract interface for evaluating a real-valued function of one real variable and performing numerica...
Definition RooAbsFunc.h:27
virtual Double_t getMinLimit(UInt_t dimension) const =0
virtual Double_t getMaxLimit(UInt_t dimension) const =0
RooAbsIntegrator is the abstract interface for integrators of real-valued functions that implement th...
const RooAbsFunc * integrand() const
Bool_t isValid() const
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:35
Special numeric integrator that can handle integrals over open domains.
LimitsCase limitsCase() const
Classify the type of limits we have: OpenBothEnds,ClosedBothEnds,OpenBelow or OpenAbove.
RooImproperIntegrator1D()
Default constructor.
virtual Bool_t checkLimits() const
Check if the limits are valid.
virtual RooAbsIntegrator * clone(const RooAbsFunc &function, const RooNumIntConfig &config) const
Return clone of integrator with given function and configuration. Needed by RooNumIntFactory.
void initialize(const RooAbsFunc *function=0)
Initialize the integrator, construct and initialize subintegrators.
virtual Double_t integral(const Double_t *yvec=0)
Calculate the integral at the given parameter values of the function binding.
virtual ~RooImproperIntegrator1D()
Destructor.
static void registerIntegrator(RooNumIntFactory &fact)
Register RooImproperIntegrator1D, its parameters and capabilities with RooNumIntFactory.
Bool_t setLimits(Double_t *xmin, Double_t *xmax)
Change our integration limits.
RooIntegrator1D implements an adaptive one-dimensional numerical integration algorithm.
Bool_t setLimits(Double_t *xmin, Double_t *xmax)
Change our integration limits.
virtual Double_t integral(const Double_t *yvec=0)
Calculate numeric integral at given set of function binding parameters.
Lightweight function binding that returns the inverse of an input function binding.
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_t 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 Int_t isInfinite(Double_t x)
Return true if x is infinite by RooNumBer internal specification.
Definition RooNumber.cxx:58
Mother of all ROOT objects.
Definition TObject.h:41
virtual const char * GetName() const
Returns name of object.
Definition TObject.cxx:429