Logo ROOT  
Reference Guide
RooRealIntegral.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 RooRealIntegral.cxx
19\class RooRealIntegral
20\ingroup Roofitcore
21
22RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
23The class performs none of the actual integration, but only manages the logic
24of what variables can be integrated analytically, accounts for eventual jacobian
25terms and defines what numerical integrations needs to be done to complement the
26analytical integral.
27The actual analytical integrations (if any) are done in the PDF themselves, the numerical
28integration is performed in the various implementations of the RooAbsIntegrator base class.
29**/
30
31#include "RooRealIntegral.h"
32
33#include "RooMsgService.h"
34#include "RooArgSet.h"
35#include "RooAbsRealLValue.h"
37#include "RooRealBinding.h"
38#include "RooRealAnalytic.h"
39#include "RooInvTransform.h"
40#include "RooSuperCategory.h"
41#include "RooNumIntFactory.h"
42#include "RooNumIntConfig.h"
43#include "RooNameReg.h"
45#include "RooConstVar.h"
46#include "RooDouble.h"
47#include "RooTrace.h"
48#include "RooHelpers.h"
49
50#include "TClass.h"
51
52#include <iostream>
53#include <memory>
54
55using namespace std;
56
58
59
61
62
63////////////////////////////////////////////////////////////////////////////////
64
66 _valid(false),
67 _respectCompSelect(true),
68 _funcNormSet(0),
69 _iconfig(0),
70 _mode(0),
71 _intOperMode(Hybrid),
72 _restartNumIntEngine(false),
73 _numIntEngine(0),
74 _numIntegrand(0),
75 _rangeName(0),
76 _params(0),
77 _cacheNum(false)
78{
80}
81
82
83
84////////////////////////////////////////////////////////////////////////////////
85/// Construct integral of 'function' over observables in 'depList'
86/// in range 'rangeName' with normalization observables 'funcNormSet'
87/// (for p.d.f.s). In the integral is performed to the maximum extent
88/// possible the internal (analytical) integrals advertised by function.
89/// The other integrations are performed numerically. The optional
90/// config object prescribes how these numeric integrations are configured.
91///
92
93RooRealIntegral::RooRealIntegral(const char *name, const char *title,
94 const RooAbsReal& function, const RooArgSet& depList,
95 const RooArgSet* funcNormSet, const RooNumIntConfig* config,
96 const char* rangeName) :
97 RooAbsReal(name,title),
98 _valid(true),
99 _respectCompSelect(true),
100 _sumList("!sumList","Categories to be summed numerically",this,false,false),
101 _intList("!intList","Variables to be integrated numerically",this,false,false),
102 _anaList("!anaList","Variables to be integrated analytically",this,false,false),
103 _jacList("!jacList","Jacobian product term",this,false,false),
104 _facList("!facList","Variables independent of function",this,false,true),
105 _function("!func","Function to be integrated",this,
106 const_cast<RooAbsReal&>(function),false,false),
107 _iconfig((RooNumIntConfig*)config),
108 _sumCat("!sumCat","SuperCategory for summation",this,false,false),
109 _mode(0),
110 _intOperMode(Hybrid),
111 _restartNumIntEngine(false),
112 _numIntEngine(0),
113 _numIntegrand(0),
114 _rangeName((TNamed*)RooNameReg::ptr(rangeName)),
115 _params(0),
116 _cacheNum(false)
117{
118 // A) Check that all dependents are lvalues
119 //
120 // B) Check if list of dependents can be re-expressed in
121 // lvalues that are higher in the expression tree
122 //
123 // C) Check for dependents that the PDF insists on integrating
124 // analytically iself
125 //
126 // D) Make list of servers that can be integrated analytically
127 // Add all parameters/dependents as value/shape servers
128 //
129 // E) Interact with function to make list of objects actually integrated analytically
130 //
131 // F) Make list of numerical integration variables consisting of:
132 // - Category dependents of RealLValues in analytical integration
133 // - Leaf nodes server lists of function server that are not analytically integrated
134 // - Make Jacobian list for analytically integrated RealLValues
135 //
136 // G) Split numeric list in integration list and summation list
137 //
138
139 oocxcoutI(&function,Integration) << "RooRealIntegral::ctor(" << GetName() << ") Constructing integral of function "
140 << function.GetName() << " over observables" << depList << " with normalization "
141 << (funcNormSet?*funcNormSet:RooArgSet()) << " with range identifier "
142 << (rangeName?rangeName:"<none>") << endl ;
143
144
145 // Choose same expensive object cache as integrand
146 setExpensiveObjectCache(function.expensiveObjectCache()) ;
147// cout << "RRI::ctor(" << GetName() << ") setting expensive object cache to " << &expensiveObjectCache() << " as taken from " << function.GetName() << endl ;
148
149 // Use objects integrator configuration if none is specified
150 if (!_iconfig) _iconfig = (RooNumIntConfig*) function.getIntegratorConfig() ;
151
152 // Save private copy of funcNormSet, if supplied, excluding factorizing terms
153 if (funcNormSet) {
154 _funcNormSet = new RooArgSet ;
155 for (const auto nArg : *funcNormSet) {
156 if (function.dependsOn(*nArg)) {
157 _funcNormSet->addClone(*nArg) ;
158 }
159 }
160 } else {
161 _funcNormSet = 0 ;
162 }
163
164 //_funcNormSet = funcNormSet ? (RooArgSet*)funcNormSet->snapshot(false) : 0 ;
165
166 // Make internal copy of dependent list
167 RooArgSet intDepList(depList) ;
168
169 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
170 // * A) Check that all dependents are lvalues and filter out any
171 // dependents that the PDF doesn't explicitly depend on
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
173
174 for (auto arg : intDepList) {
175 if(!arg->isLValue()) {
176 coutE(InputArguments) << ClassName() << "::" << GetName() << ": cannot integrate non-lvalue ";
177 arg->Print("1");
178 _valid= false;
179 }
180 if (!function.dependsOn(*arg)) {
181 RooAbsArg* argClone = (RooAbsArg*) arg->Clone() ;
182 _facListOwned.addOwned(*argClone) ;
183 _facList.add(*argClone) ;
184 addServer(*argClone,false,true) ;
185 }
186 }
187
188 if (_facList.getSize()>0) {
189 oocxcoutI(&function,Integration) << function.GetName() << ": Factorizing obserables are " << _facList << endl ;
190 }
191
192
193 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
194 // * B) Check if list of dependents can be re-expressed in *
195 // * lvalues that are higher in the expression tree *
196 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
197
198
199 // Initial fill of list of LValue branches
200 RooArgSet exclLVBranches("exclLVBranches") ;
201 RooArgSet branchList,branchListVD ;
202 function.branchNodeServerList(&branchList) ;
203
204 for (auto branch: branchList) {
205 RooAbsRealLValue *realArgLV = dynamic_cast<RooAbsRealLValue*>(branch) ;
206 RooAbsCategoryLValue *catArgLV = dynamic_cast<RooAbsCategoryLValue*>(branch) ;
207 if ((realArgLV && (realArgLV->isJacobianOK(intDepList)!=0)) || catArgLV) {
208 exclLVBranches.add(*branch) ;
209// cout << "exclv branch = " << endl ;
210// branch->printCompactTree() ;
211 }
212 if (dependsOnValue(*branch)) {
213 branchListVD.add(*branch) ;
214 } else {
215// cout << "value of self does not depend on branch " << branch->GetName() << endl ;
216 }
217 }
218 exclLVBranches.remove(depList,true,true) ;
219
220 // Initial fill of list of LValue leaf servers (put in intDepList)
221 RooArgSet exclLVServers("exclLVServers") ;
222 exclLVServers.add(intDepList) ;
223
224 // Obtain mutual exclusive dependence by iterative reduction
225 bool converged(false) ;
226 while(!converged) {
227 converged=true ;
228
229 // Reduce exclLVServers to only those serving exclusively exclLVBranches
230 std::vector<RooAbsArg*> toBeRemoved;
231 for (auto server : exclLVServers) {
232 if (!servesExclusively(server,exclLVBranches,branchListVD)) {
233 toBeRemoved.push_back(server);
234 converged=false ;
235 }
236 }
237 exclLVServers.remove(toBeRemoved.begin(), toBeRemoved.end());
238
239 // Reduce exclLVBranches to only those depending exclusively on exclLVservers
240 // Attention: counting loop, since erasing from container
241 for (std::size_t i=0; i < exclLVBranches.size(); ++i) {
242 const RooAbsArg* branch = exclLVBranches[i];
243 RooArgSet* brDepList = branch->getObservables(&intDepList) ;
244 RooArgSet bsList(*brDepList,"bsList") ;
245 delete brDepList ;
246 bsList.remove(exclLVServers,true,true) ;
247 if (bsList.getSize()>0) {
248 exclLVBranches.remove(*branch,true,true) ;
249 --i;
250 converged=false ;
251 }
252 }
253 }
254
255 // Eliminate exclLVBranches that do not depend on any LVServer
256 // Attention: Counting loop, since modifying container
257 for (std::size_t i=0; i < exclLVBranches.size(); ++i) {
258 const RooAbsArg* branch = exclLVBranches[i];
259 if (!branch->dependsOnValue(exclLVServers)) {
260 exclLVBranches.remove(*branch,true,true) ;
261 --i;
262 }
263 }
264
265 // Replace exclusive lvalue branch servers with lvalue branches
266 // WVE Don't do this for binned distributions - deal with this using numeric integration with transformed bin boundaroes
267 if (exclLVServers.getSize()>0 && !function.isBinnedDistribution(exclLVBranches)) {
268// cout << "activating LVservers " << exclLVServers << " for use in integration " << endl ;
269 intDepList.remove(exclLVServers) ;
270 intDepList.add(exclLVBranches) ;
271
272 //cout << "intDepList removing exclLVServers " << exclLVServers << endl ;
273 //cout << "intDepList adding exclLVBranches " << exclLVBranches << endl ;
274
275 }
276
277
278 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
279 // * C) Check for dependents that the PDF insists on integrating *
280 // analytically iself *
281 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
282
283 RooArgSet anIntOKDepList ;
284 for (auto arg : intDepList) {
285 if (function.forceAnalyticalInt(*arg)) {
286 anIntOKDepList.add(*arg) ;
287 }
288 }
289
290 if (anIntOKDepList.getSize()>0) {
291 oocxcoutI(&function,Integration) << function.GetName() << ": Observables that function forcibly requires to be integrated internally " << anIntOKDepList << endl ;
292 }
293
294
295 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
296 // * D) Make list of servers that can be integrated analytically *
297 // Add all parameters/dependents as value/shape servers *
298 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
299
300 for (const auto arg : function.servers()) {
301
302 //cout << "considering server" << arg->GetName() << endl ;
303
304 // Dependent or parameter?
305 if (!arg->dependsOnValue(intDepList)) {
306
307 //cout << " server does not depend on observables, adding server as value server to integral" << endl ;
308
309 if (function.dependsOnValue(*arg)) {
310 addServer(*arg,true,false) ;
311 }
312
313 continue ;
314
315 } else {
316
317 // Add final dependents of arg as shape servers
318 RooArgSet argLeafServers ;
319 arg->leafNodeServerList(&argLeafServers,0,false) ;
320
321 //arg->printCompactTree() ;
322 //cout << "leaf nodes of server are " << argLeafServers << " depList = " << depList << endl ;
323
324 // Skip arg if it is neither value or shape server
325 if (!arg->isValueServer(function) && !arg->isShapeServer(function)) {
326 //cout << " server is neither value not shape server of function, ignoring" << endl ;
327 continue ;
328 }
329
330 for (const auto leaf : argLeafServers) {
331
332 //cout << " considering leafnode " << leaf->GetName() << " of server " << arg->GetName() << endl ;
333
334 if (depList.find(leaf->GetName()) && function.dependsOnValue(*leaf)) {
335
336 RooAbsRealLValue* leaflv = dynamic_cast<RooAbsRealLValue*>(leaf) ;
337 if (leaflv && leaflv->getBinning(rangeName).isParameterized()) {
338 oocxcoutD(&function,Integration) << function.GetName() << " : Observable " << leaf->GetName() << " has parameterized binning, add value dependence of boundary objects rather than shape of leaf" << endl ;
339 if (leaflv->getBinning(rangeName).lowBoundFunc()) {
340 addServer(*leaflv->getBinning(rangeName).lowBoundFunc(),true,false) ;
341 }
342 if(leaflv->getBinning(rangeName).highBoundFunc()) {
343 addServer(*leaflv->getBinning(rangeName).highBoundFunc(),true,false) ;
344 }
345 } else {
346 oocxcoutD(&function,Integration) << function.GetName() << ": Adding observable " << leaf->GetName() << " of server "
347 << arg->GetName() << " as shape dependent" << endl ;
348 addServer(*leaf,false,true) ;
349 }
350 } else if (!depList.find(leaf->GetName())) {
351
352 if (function.dependsOnValue(*leaf)) {
353 oocxcoutD(&function,Integration) << function.GetName() << ": Adding parameter " << leaf->GetName() << " of server " << arg->GetName() << " as value dependent" << endl ;
354 addServer(*leaf,true,false) ;
355 } else {
356 oocxcoutD(&function,Integration) << function.GetName() << ": Adding parameter " << leaf->GetName() << " of server " << arg->GetName() << " as shape dependent" << endl ;
357 addServer(*leaf,false,true) ;
358 }
359 }
360 }
361 }
362
363 // If this dependent arg is self-normalized, stop here
364 //if (function.selfNormalized()) continue ;
365
366 bool depOK(false) ;
367 // Check for integratable AbsRealLValue
368
369 //cout << "checking server " << arg->IsA()->GetName() << "::" << arg->GetName() << endl ;
370
371 if (arg->isDerived()) {
372 RooAbsRealLValue *realArgLV = dynamic_cast<RooAbsRealLValue*>(arg) ;
373 RooAbsCategoryLValue *catArgLV = dynamic_cast<RooAbsCategoryLValue*>(arg) ;
374 //cout << "realArgLV = " << realArgLV << " intDepList = " << intDepList << endl ;
375 if ((realArgLV && intDepList.find(realArgLV->GetName()) && (realArgLV->isJacobianOK(intDepList)!=0)) || catArgLV) {
376
377 //cout << " arg " << arg->GetName() << " is derived LValue with valid jacobian" << endl ;
378
379 // Derived LValue with valid jacobian
380 depOK = true ;
381
382 // Now, check for overlaps
383 bool overlapOK = true ;
384 for (const auto otherArg : function.servers()) {
385 // skip comparison with self
386 if (arg==otherArg) continue ;
387 if (otherArg->IsA()==RooConstVar::Class()) continue ;
388 if (arg->overlaps(*otherArg,true)) {
389 //cout << "arg " << arg->GetName() << " overlaps with " << otherArg->GetName() << endl ;
390 //overlapOK=false ;
391 }
392 }
393 // coverity[DEADCODE]
394 if (!overlapOK) depOK=false ;
395
396 //cout << "overlap check returns OK=" << (depOK?"T":"F") << endl ;
397
398 }
399 } else {
400 // Fundamental types are always OK
401 depOK = true ;
402 }
403
404 // Add server to list of dependents that are OK for analytical integration
405 if (depOK) {
406 anIntOKDepList.add(*arg,true) ;
407 oocxcoutI(&function,Integration) << function.GetName() << ": Observable " << arg->GetName() << " is suitable for analytical integration (if supported by p.d.f)" << endl ;
408 }
409 }
410 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
411 // * E) interact with function to make list of objects actually integrated analytically *
412 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
413
414 RooArgSet anIntDepList ;
415
416 RooArgSet *anaSet = new RooArgSet( _anaList, Form("UniqueCloneOf_%s",_anaList.GetName()));
419 _anaList.add(*anaSet);
420 delete anaSet;
421
422 // Avoid confusion -- if mode is zero no analytical integral is defined regardless of contents of _anaListx
423 if (_mode==0) {
425 }
426
427 if (_mode!=0) {
428 oocxcoutI(&function,Integration) << function.GetName() << ": Function integrated observables " << _anaList << " internally with code " << _mode << endl ;
429 }
430
431 // when _funcNormSet is a nullptr a warning message appears for RooAddPdf functions
432 // This is not a problem since we do noty use the returned value from getVal()
433 // we then disable the produced warning message in the RooFit::Eval topic
434 std::unique_ptr<RooHelpers::LocalChangeMsgLevel> msgChanger;
435 if (_funcNormSet == nullptr) {
436 // remove only the RooFit::Eval message topic from current active streams
437 // passed level can be whatever if we provide a false as last argument
438 msgChanger = std::make_unique<RooHelpers::LocalChangeMsgLevel>(RooFit::WARNING, 0u, RooFit::Eval, false);
439 }
440
441 // WVE kludge: synchronize dset for use in analyticalIntegral
442 // LM : is this really needed ??
443 function.getVal(_funcNormSet) ;
444 // delete LocalChangeMsgLevel which will restore previous message level
445 msgChanger.reset(nullptr);
446
447 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
448 // * F) Make list of numerical integration variables consisting of: *
449 // * - Category dependents of RealLValues in analytical integration *
450 // * - Expanded server lists of server that are not analytically integrated *
451 // * Make Jacobian list with analytically integrated RealLValues *
452 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
453
454 RooArgSet numIntDepList ;
455
456
457 // Loop over actually analytically integrated dependents
458 for (const auto arg : _anaList) {
459
460 // Process only derived RealLValues
461 if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class()) && arg->isDerived() && !arg->isFundamental()) {
462
463 // Add to list of Jacobians to calculate
464 _jacList.add(*arg) ;
465
466 // Add category dependent of LValueReal used in integration
467 auto argDepList = std::unique_ptr<RooArgSet>(arg->getObservables(&intDepList));
468 for (const auto argDep : *argDepList) {
469 if (argDep->IsA()->InheritsFrom(RooAbsCategoryLValue::Class()) && intDepList.contains(*argDep)) {
470 numIntDepList.add(*argDep,true) ;
471 }
472 }
473 }
474 }
475
476
477 // If nothing was integrated analytically, swap back LVbranches for LVservers for subsequent numeric integration
478 if (_anaList.getSize()==0) {
479 if (exclLVServers.getSize()>0) {
480 //cout << "NUMINT phase analList is empty. exclLVServers = " << exclLVServers << endl ;
481 intDepList.remove(exclLVBranches) ;
482 intDepList.add(exclLVServers) ;
483 }
484 }
485 //cout << "NUMINT intDepList = " << intDepList << endl ;
486
487 // Loop again over function servers to add remaining numeric integrations
488 for (const auto arg : function.servers()) {
489
490 //cout << "processing server for numeric integration " << arg->IsA()->GetName() << "::" << arg->GetName() << endl ;
491
492 // Process only servers that are not treated analytically
493 if (!_anaList.find(arg->GetName()) && arg->dependsOn(intDepList)) {
494
495 // Process only derived RealLValues
496 if (dynamic_cast<RooAbsLValue*>(arg) && arg->isDerived() && intDepList.contains(*arg)) {
497 numIntDepList.add(*arg,true) ;
498 } else {
499
500 // WVE this will only get the observables, but not l-value transformations
501 // Expand server in final dependents
502 auto argDeps = std::unique_ptr<RooArgSet>(arg->getObservables(&intDepList));
503
504 if (argDeps->getSize()>0) {
505
506 // Add final dependents, that are not forcibly integrated analytically,
507 // to numerical integration list
508 for (const auto dep : *argDeps) {
509 if (!_anaList.find(dep->GetName())) {
510 numIntDepList.add(*dep,true) ;
511 }
512 }
513 }
514 }
515 }
516 }
517
518 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
519 // * G) Split numeric list in integration list and summation list *
520 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
521
522 // Split numeric integration list in summation and integration lists
523 for (const auto arg : numIntDepList) {
524 if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
525 _intList.add(*arg) ;
526 } else if (arg->IsA()->InheritsFrom(RooAbsCategoryLValue::Class())) {
527 _sumList.add(*arg) ;
528 }
529 }
530
531 if (_anaList.getSize()>0) {
532 oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _anaList << " are analytically integrated with code " << _mode << endl ;
533 }
534 if (_intList.getSize()>0) {
535 oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _intList << " are numerically integrated" << endl ;
536 }
537 if (_sumList.getSize()>0) {
538 oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _sumList << " are numerically summed" << endl ;
539 }
540
541
542 // Determine operating mode
543 if (numIntDepList.getSize()>0) {
544 // Numerical and optional Analytical integration
546 } else if (_anaList.getSize()>0) {
547 // Purely analytical integration
549 } else {
550 // No integration performed
552 }
553
554 // Determine auto-dirty status
556
557 // Create value caches for _intList and _sumList
560
561
562 if (_sumList.getSize()>0) {
563 RooSuperCategory *sumCat = new RooSuperCategory(Form("%s_sumCat",GetName()),"sumCat",_sumList) ;
564 _sumCat.addOwned(*sumCat) ;
565 }
566
568}
569
570
571
572////////////////////////////////////////////////////////////////////////////////
573/// Set appropriate cache operation mode for integral depending on cache operation
574/// mode of server objects
575
577{
578 // If any of our servers are is forcedDirty or a projectedDependent, then we need to be ADirty
579 for (const auto server : _serverList) {
580 if (server->isValueServer(*this)) {
581 RooArgSet leafSet ;
582 server->leafNodeServerList(&leafSet) ;
583 for (const auto leaf : leafSet) {
584 if (leaf->operMode()==ADirty && leaf->isValueServer(*this)) {
586 break ;
587 }
588 if (leaf->getAttribute("projectedDependent")) {
590 break ;
591 }
592 }
593 }
594 }
595}
596
597
598
599////////////////////////////////////////////////////////////////////////////////
600/// Utility function that returns true if 'object server' is a server
601/// to exactly one of the RooAbsArgs in 'exclLVBranches'
602
603bool RooRealIntegral::servesExclusively(const RooAbsArg* server,const RooArgSet& exclLVBranches, const RooArgSet& allBranches) const
604{
605 // Determine if given server serves exclusively exactly one of the given nodes in exclLVBranches
606
607 // Special case, no LV servers available
608 if (exclLVBranches.getSize()==0) return false ;
609
610 // If server has no clients and is not an LValue itself, return false
611 if (server->_clientList.empty() && exclLVBranches.find(server->GetName())) {
612 return false ;
613 }
614
615 // WVE must check for value relations only here!!!!
616
617 // Loop over all clients
618 Int_t numLVServ(0) ;
619 for (const auto client : server->valueClients()) {
620 // If client is not an LValue, recurse
621 if (!(exclLVBranches.find(client->GetName())==client)) {
622 if (allBranches.find(client->GetName())==client) {
623 if (!servesExclusively(client,exclLVBranches,allBranches)) {
624 // Client is a non-LValue that doesn't have an exclusive LValue server
625 return false ;
626 }
627 }
628 } else {
629 // Client is an LValue
630 numLVServ++ ;
631 }
632 }
633
634 return (numLVServ==1) ;
635}
636
637
638
639
640////////////////////////////////////////////////////////////////////////////////
641/// (Re)Initialize numerical integration engine if necessary. Return true if
642/// successful, or otherwise false.
643
645{
646 // if we already have an engine, check if it still works for the present limits.
647 if(0 != _numIntEngine) {
649 // otherwise, cleanup the old engine
650 delete _numIntEngine ;
651 _numIntEngine= 0;
652 if(0 != _numIntegrand) {
653 delete _numIntegrand;
654 _numIntegrand= 0;
655 }
656 }
657
658 // All done if there are no arguments to integrate numerically
659 if(0 == _intList.getSize()) return true;
660
661 // Bind the appropriate analytic integral (specified by _mode) of our RooRealVar object to
662 // those of its arguments that will be integrated out numerically.
663 if(_mode != 0) {
665 }
666 else {
668 }
669 if(0 == _numIntegrand || !_numIntegrand->isValid()) {
670 coutE(Integration) << ClassName() << "::" << GetName() << ": failed to create valid integrand." << endl;
671 return false;
672 }
673
674 // Create appropriate numeric integrator using factory
675 bool isBinned = _function.arg().isBinnedDistribution(_intList) ;
677
678 if(0 == _numIntEngine || !_numIntEngine->isValid()) {
679 coutE(Integration) << ClassName() << "::" << GetName() << ": failed to create valid integrator." << endl;
680 return false;
681 }
682
683 cxcoutI(NumIntegration) << "RooRealIntegral::init(" << GetName() << ") using numeric integrator "
684 << _numIntEngine->IsA()->GetName() << " to calculate Int" << _intList << endl ;
685
686 if (_intList.getSize()>3) {
687 cxcoutI(NumIntegration) << "RooRealIntegral::init(" << GetName() << ") evaluation requires " << _intList.getSize() << "-D numeric integration step. Evaluation may be slow, sufficient numeric precision for fitting & minimization is not guaranteed" << endl ;
688 }
689
690 _restartNumIntEngine = false ;
691 return true;
692}
693
694
695
696////////////////////////////////////////////////////////////////////////////////
697/// Copy constructor
698
700 RooAbsReal(other,name),
701 _valid(other._valid),
702 _respectCompSelect(other._respectCompSelect),
703 _sumList("!sumList",this,other._sumList),
704 _intList("!intList",this,other._intList),
705 _anaList("!anaList",this,other._anaList),
706 _jacList("!jacList",this,other._jacList),
707 _facList("!facList","Variables independent of function",this,false,true),
708 _function("!func",this,other._function),
709 _iconfig(other._iconfig),
710 _sumCat("!sumCat",this,other._sumCat),
711 _mode(other._mode),
712 _intOperMode(other._intOperMode),
713 _restartNumIntEngine(false),
714 _numIntEngine(0),
715 _numIntegrand(0),
716 _rangeName(other._rangeName),
717 _params(0),
718 _cacheNum(false)
719{
720 _funcNormSet = other._funcNormSet ? (RooArgSet*)other._funcNormSet->snapshot(false) : 0 ;
721
722 for (const auto arg : other._facList) {
723 RooAbsArg* argClone = (RooAbsArg*) arg->Clone() ;
724 _facListOwned.addOwned(*argClone) ;
725 _facList.add(*argClone) ;
726 addServer(*argClone,false,true) ;
727 }
728
729 other._intList.snapshot(_saveInt) ;
730 other._sumList.snapshot(_saveSum) ;
731
733}
734
735
736
737////////////////////////////////////////////////////////////////////////////////
738
740 // Destructor
741{
742 if (_numIntEngine) delete _numIntEngine ;
743 if (_numIntegrand) delete _numIntegrand ;
744 if (_funcNormSet) delete _funcNormSet ;
745 if (_params) delete _params ;
746
748}
749
750
751
752
753
754////////////////////////////////////////////////////////////////////////////////
755
756RooAbsReal* RooRealIntegral::createIntegral(const RooArgSet& iset, const RooArgSet* nset, const RooNumIntConfig* cfg, const char* rangeName) const
757{
758 // Handle special case of no integration with default algorithm
759 if (iset.getSize()==0) {
760 return RooAbsReal::createIntegral(iset,nset,cfg,rangeName) ;
761 }
762
763 // Special handling of integral of integral, return RooRealIntegral that represents integral over all dimensions in one pass
764 RooArgSet isetAll(iset) ;
765 isetAll.add(_sumList) ;
766 isetAll.add(_intList) ;
767 isetAll.add(_anaList) ;
768 isetAll.add(_facList) ;
769
770 const RooArgSet* newNormSet(0) ;
771 RooArgSet* tmp(0) ;
772 if (nset && !_funcNormSet) {
773 newNormSet = nset ;
774 } else if (!nset && _funcNormSet) {
775 newNormSet = _funcNormSet ;
776 } else if (nset && _funcNormSet) {
777 tmp = new RooArgSet ;
778 tmp->add(*nset) ;
779 tmp->add(*_funcNormSet,true) ;
780 newNormSet = tmp ;
781 }
782 RooAbsReal* ret = _function.arg().createIntegral(isetAll,newNormSet,cfg,rangeName) ;
783
784 if (tmp) {
785 delete tmp ;
786 }
787
788 return ret ;
789}
790
791
792
793////////////////////////////////////////////////////////////////////////////////
794/// Return value of object. If the cache is clean, return the
795/// cached value, otherwise recalculate on the fly and refill
796/// the cache
797
798double RooRealIntegral::getValV(const RooArgSet* nset) const
799{
800// // fast-track clean-cache processing
801// if (_operMode==AClean) {
802// return _value ;
803// }
804
805 if (nset && nset!=_lastNSet) {
806 ((RooAbsReal*) this)->setProxyNormSet(nset) ;
807 _lastNSet = (RooArgSet*) nset ;
808 }
809
811 _value = traceEval(nset) ;
812 }
813
814 return _value ;
815}
816
817
818
819
820
821////////////////////////////////////////////////////////////////////////////////
822/// Perform the integration and return the result
823
825{
827
828 double retVal(0) ;
829 switch (_intOperMode) {
830
831 case Hybrid:
832 {
833 // Cache numeric integrals in >1d expensive object cache
834 RooDouble* cacheVal(0) ;
837 }
838
839 if (cacheVal) {
840 retVal = *cacheVal ;
841 // cout << "using cached value of integral" << GetName() << endl ;
842 } else {
843
844
845 // Find any function dependents that are AClean
846 // and switch them temporarily to ADirty
847 bool origState = inhibitDirty() ;
848 setDirtyInhibit(true) ;
849
850 // try to initialize our numerical integration engine
851 if(!(_valid= initNumIntegrator())) {
852 coutE(Integration) << ClassName() << "::" << GetName()
853 << ":evaluate: cannot initialize numerical integrator" << endl;
854 return 0;
855 }
856
857 // Save current integral dependent values
860
861 // Evaluate sum/integral
862 retVal = sum() ;
863
864
865 // This must happen BEFORE restoring dependents, otherwise no dirty state propagation in restore step
866 setDirtyInhibit(origState) ;
867
868 // Restore integral dependent values
871
872 // Cache numeric integrals in >1d expensive object cache
874 RooDouble* val = new RooDouble(retVal) ;
876 // cout << "### caching value of integral" << GetName() << " in " << &expensiveObjectCache() << endl ;
877 }
878
879 }
880 break ;
881 }
882 case Analytic:
883 {
885 cxcoutD(Tracing) << "RooRealIntegral::evaluate_analytic(" << GetName()
886 << ")func = " << _function.arg().IsA()->GetName() << "::" << _function.arg().GetName()
887 << " raw = " << retVal << " _funcNormSet = " << (_funcNormSet?*_funcNormSet:RooArgSet()) << endl ;
888
889
890 break ;
891 }
892
893 case PassThrough:
894 {
895 //setDirtyInhibit(true) ;
896 retVal= _function.arg().getVal(_funcNormSet) ;
897 //setDirtyInhibit(false) ;
898 break ;
899 }
900 }
901
902
903 // Multiply answer with integration ranges of factorized variables
904 if (_facList.getSize()>0) {
905 for (const auto arg : _facList) {
906 // Multiply by fit range for 'real' dependents
907 if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
908 RooAbsRealLValue* argLV = (RooAbsRealLValue*)arg ;
909 retVal *= (argLV->getMax() - argLV->getMin()) ;
910 }
911 // Multiply by number of states for category dependents
912 if (arg->IsA()->InheritsFrom(RooAbsCategoryLValue::Class())) {
914 retVal *= argLV->numTypes() ;
915 }
916 }
917 }
918
919
920 if (dologD(Tracing)) {
921 cxcoutD(Tracing) << "RooRealIntegral::evaluate(" << GetName() << ") anaInt = " << _anaList << " numInt = " << _intList << _sumList << " mode = " ;
922 switch(_mode) {
923 case Hybrid: ccoutD(Tracing) << "Hybrid" ; break ;
924 case Analytic: ccoutD(Tracing) << "Analytic" ; break ;
925 case PassThrough: ccoutD(Tracing) << "PassThrough" ; break ;
926 }
927
928 ccxcoutD(Tracing) << "raw*fact = " << retVal << endl ;
929 }
930
931 return retVal ;
932}
933
934
935
936////////////////////////////////////////////////////////////////////////////////
937/// Return product of jacobian terms originating from analytical integration
938
940{
941 if (_jacList.getSize()==0) {
942 return 1 ;
943 }
944
945 double jacProd(1) ;
946 for (const auto elm : _jacList) {
947 auto arg = static_cast<const RooAbsRealLValue*>(elm);
948 jacProd *= arg->jacobian() ;
949 }
950
951 // Take fabs() here: if jacobian is negative, min and max are swapped and analytical integral
952 // will be positive, so must multiply with positive jacobian.
953 return fabs(jacProd) ;
954}
955
956
957
958////////////////////////////////////////////////////////////////////////////////
959/// Perform summation of list of category dependents to be integrated
960
962{
963 if (_sumList.getSize()!=0) {
964 // Add integrals for all permutations of categories summed over
965 double total(0) ;
966
968 for (const auto& nameIdx : *sumCat) {
969 sumCat->setIndex(nameIdx);
970 if (!_rangeName || sumCat->inRange(RooNameReg::str(_rangeName))) {
972 }
973 }
974
975 return total ;
976
977 } else {
978 // Simply return integral
979 double ret = integrate() / jacobianProduct() ;
980 return ret ;
981 }
982}
983
984
985////////////////////////////////////////////////////////////////////////////////
986/// Perform hybrid numerical/analytical integration over all real-valued dependents
987
989{
990 if (!_numIntEngine) {
991 // Trivial case, fully analytical integration
992 return ((RooAbsReal&)_function.arg()).analyticalIntegralWN(_mode,_funcNormSet,RooNameReg::str(_rangeName)) ;
993 } else {
994 return _numIntEngine->calculate() ;
995 }
996}
997
998
999
1000////////////////////////////////////////////////////////////////////////////////
1001/// Intercept server redirects and reconfigure internal object accordingly
1002
1004 bool /*mustReplaceAll*/, bool /*nameChange*/, bool /*isRecursive*/)
1005{
1006 _restartNumIntEngine = true ;
1007
1009
1010 // Update contents value caches for _intList and _sumList
1015
1016 // Delete parameters cache if we have one
1017 if (_params) {
1018 delete _params ;
1019 _params = 0 ;
1020 }
1021
1022 return false ;
1023}
1024
1025
1026
1027////////////////////////////////////////////////////////////////////////////////
1028
1030{
1031 if (!_params) {
1032 _params = new RooArgSet("params") ;
1033
1034 RooArgSet params ;
1035 for (const auto server : _serverList) {
1036 if (server->isValueServer(*this)) _params->add(*server) ;
1037 }
1038 }
1039
1040 return *_params ;
1041}
1042
1043
1044////////////////////////////////////////////////////////////////////////////////
1045/// Check if current value is valid
1046
1047bool RooRealIntegral::isValidReal(double /*value*/, bool /*printError*/) const
1048{
1049 return true ;
1050}
1051
1052////////////////////////////////////////////////////////////////////////////////
1053/// Check if component selection is allowed
1054
1056 return _respectCompSelect;
1057}
1058
1059////////////////////////////////////////////////////////////////////////////////
1060/// Set component selection to be allowed/forbidden
1061
1063 _respectCompSelect = allow;
1064}
1065
1066////////////////////////////////////////////////////////////////////////////////
1067/// Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the
1068/// integration operation
1069
1070void RooRealIntegral::printMetaArgs(ostream& os) const
1071{
1072
1073 if (intVars().getSize()!=0) {
1074 os << "Int " ;
1075 }
1076 os << _function.arg().GetName() ;
1077 if (_funcNormSet) {
1078 os << "_Norm" ;
1079 os << *_funcNormSet ;
1080 os << " " ;
1081 }
1082
1083 // List internally integrated observables and factorizing observables as analytically integrated
1084 RooArgSet tmp(_anaList) ;
1085 tmp.add(_facList) ;
1086 if (tmp.getSize()>0) {
1087 os << "d[Ana]" ;
1088 os << tmp ;
1089 os << " " ;
1090 }
1091
1092 // List numerically integrated and summed observables as numerically integrated
1093 RooArgSet tmp2(_intList) ;
1094 tmp2.add(_sumList) ;
1095 if (tmp2.getSize()>0) {
1096 os << " d[Num]" ;
1097 os << tmp2 ;
1098 os << " " ;
1099 }
1100}
1101
1102
1103
1104////////////////////////////////////////////////////////////////////////////////
1105/// Print the state of this object to the specified output stream.
1106
1107void RooRealIntegral::printMultiline(ostream& os, Int_t contents, bool verbose, TString indent) const
1108{
1110 os << indent << "--- RooRealIntegral ---" << endl;
1111 os << indent << " Integrates ";
1113 TString deeper(indent);
1114 deeper.Append(" ");
1115 os << indent << " operating mode is "
1116 << (_intOperMode==Hybrid?"Hybrid":(_intOperMode==Analytic?"Analytic":"PassThrough")) << endl ;
1117 os << indent << " Summed discrete args are " << _sumList << endl ;
1118 os << indent << " Numerically integrated args are " << _intList << endl;
1119 os << indent << " Analytically integrated args using mode " << _mode << " are " << _anaList << endl ;
1120 os << indent << " Arguments included in Jacobian are " << _jacList << endl ;
1121 os << indent << " Factorized arguments are " << _facList << endl ;
1122 os << indent << " Function normalization set " ;
1123 if (_funcNormSet)
1124 _funcNormSet->Print("1") ;
1125 else
1126 os << "<none>" ;
1127
1128 os << endl ;
1129}
1130
1131
1132
1133////////////////////////////////////////////////////////////////////////////////
1134/// Global switch to cache all integral values that integrate at least ndim dimensions numerically
1135
1137 _cacheAllNDim = ndim ;
1138}
1139
1140
1141////////////////////////////////////////////////////////////////////////////////
1142/// Return minimum dimensions of numeric integration for which values are cached.
1143
1145{
1146 return _cacheAllNDim ;
1147}
#define cxcoutI(a)
Definition: RooMsgService.h:89
#define cxcoutD(a)
Definition: RooMsgService.h:85
#define oocxcoutD(o, a)
Definition: RooMsgService.h:87
#define dologD(a)
Definition: RooMsgService.h:69
#define coutE(a)
Definition: RooMsgService.h:37
#define ccxcoutD(a)
Definition: RooMsgService.h:86
#define ccoutD(a)
Definition: RooMsgService.h:41
#define oocxcoutI(o, a)
Definition: RooMsgService.h:91
#define TRACE_DESTROY
Definition: RooTrace.h:24
#define TRACE_CREATE
Definition: RooTrace.h:23
int Int_t
Definition: RtypesCore.h:45
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
static unsigned int total
char name[80]
Definition: TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2447
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:77
RooExpensiveObjectCache & expensiveObjectCache() const
Definition: RooAbsArg.cxx:2354
TObject * Clone(const char *newname=0) const override
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:89
void setOperMode(OperMode mode, bool recurseADirty=true)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1876
bool dependsOnValue(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0) const
Check whether this object depends on values from an element in the serverList.
Definition: RooAbsArg.h:107
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:317
friend class RooArgSet
Definition: RooAbsArg.h:648
const RefCountList_t & valueClients() const
List of all value clients of this object. Value clients receive value updates.
Definition: RooAbsArg.h:194
virtual void setExpensiveObjectCache(RooExpensiveObjectCache &cache)
Definition: RooAbsArg.h:526
static void setDirtyInhibit(bool flag)
Control global dirty inhibit mode.
Definition: RooAbsArg.cxx:255
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:835
void addServer(RooAbsArg &server, bool valueProp=true, bool shapeProp=false, std::size_t refCount=1)
Register another RooAbsArg as a server to us, ie, declare that we depend on it.
Definition: RooAbsArg.cxx:379
bool isValueOrShapeDirtyAndClear() const
Definition: RooAbsArg.h:478
bool inhibitDirty() const
Delete watch flag.
Definition: RooAbsArg.cxx:114
RefCountList_t _clientList
Definition: RooAbsArg.h:657
RefCountList_t _serverList
Definition: RooAbsArg.h:656
virtual bool isParameterized() const
Interface function.
Definition: RooAbsBinning.h:80
virtual RooAbsReal * highBoundFunc() const
Return pointer to RooAbsReal parameterized upper bound, if any.
Definition: RooAbsBinning.h:88
virtual RooAbsReal * lowBoundFunc() const
Return pointer to RooAbsReal parameterized lower bound, if any.
Definition: RooAbsBinning.h:84
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
static TClass * Class()
Int_t numTypes(const char *=0) const
Return number of types defined (in range named rangeName if rangeName!=0)
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
bool contains(const RooAbsArg &var) const
Check if collection contains an argument with the same name as var.
void Print(Option_t *options=0) const override
This method must be overridden when a class wants to print itself.
Int_t getSize() const
Return the number of elements in the collection.
const char * GetName() const override
Returns name of object.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
Storage_t::size_type size() const
RooAbsArg * first() const
RooAbsArg * find(const char *name) const
Find object with given name in list.
bool isValid() const
Definition: RooAbsFunc.h:37
bool isValid() const
Is integrator in valid state.
double calculate(const double *yvec=0)
Calculate integral value with given array of parameter values.
virtual bool checkLimits() const
If true, finite limits are required on the observable range.
TClass * IsA() const override
Abstract base class for objects that are lvalues, i.e.
Definition: RooAbsLValue.h:26
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
virtual double getMax(const char *name=0) const
Get maximum of currently defined range.
virtual const RooAbsBinning & getBinning(const char *name=0, bool verbose=true, bool createOnTheFly=false) const =0
Retrive binning configuration with given name or default binning.
virtual bool isJacobianOK(const RooArgSet &depList) const
virtual double getMin(const char *name=0) const
Get minimum of currently defined range.
static TClass * Class()
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
virtual double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
Definition: RooAbsReal.cxx:420
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Variant of getAnalyticalIntegral that is also passed the normalization set that should be applied to ...
Definition: RooAbsReal.cxx:391
friend class RooRealBinding
Definition: RooAbsReal.h:480
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:553
RooArgSet * _lastNSet
!
Definition: RooAbsReal.h:585
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Structure printing.
Definition: RooAbsReal.cxx:495
double _value
Cache for current value of object.
Definition: RooAbsReal.h:484
double traceEval(const RooArgSet *set) const
Calculate current value of object, with error tracing wrapper.
Definition: RooAbsReal.cxx:359
virtual bool isBinnedDistribution(const RooArgSet &) const
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition: RooAbsReal.h:345
TClass * IsA() const override
Definition: RooAbsReal.h:588
static bool _globalSelectComp
Component selection flag for RooAbsPdf::plotCompOn.
Definition: RooAbsReal.h:566
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:180
void removeAll() override
Remove all argument inset using remove(const RooAbsArg&).
bool addOwned(RooAbsArg &var, bool silent=false) override
Overloaded RooCollection_t::addOwned() method insert object into owning set and registers object as s...
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
static TClass * Class()
RooDouble is a minimal implementation of a TObject holding a double value.
Definition: RooDouble.h:22
static TClass * Class()
const TObject * retrieveObject(const char *name, TClass *tclass, const RooArgSet &params)
Retrieve object from cache that was registered under given name with given parameters,...
bool registerObject(const char *ownerName, const char *objectName, TObject &cacheObject, TIterator *paramIter)
Register object associated with given name and given associated parameters with given values in cache...
RooNameReg is a registry for const char* names.
Definition: RooNameReg.h:25
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.cxx:102
RooNumIntConfig holds the configuration parameters of the various numeric integrators used by RooReal...
RooAbsIntegrator * createIntegrator(RooAbsFunc &func, const RooNumIntConfig &config, Int_t ndim=0, bool isBinned=false) const
Construct a numeric integrator instance that operates on function 'func' and is configured with 'conf...
static RooNumIntFactory & instance()
Static method returning reference to singleton instance of factory.
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer,...
Lightweight RooAbsFunc interface adaptor that binds an analytic integral of a RooAbsReal object (spec...
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
RooNumIntConfig * _iconfig
bool initNumIntegrator() const
(Re)Initialize numerical integration engine if necessary.
RooArgSet * _params
! cache for set of parameters
void setAllowComponentSelection(bool allow)
Set component selection to be allowed/forbidden.
RooRealProxy _function
Function being integration.
RooArgSet intVars() const
RooSetProxy _intList
Set of continuous observables over which is integrated numerically.
virtual double sum() const
Perform summation of list of category dependents to be integrated.
RooSetProxy _facList
Set of observables on which function does not depends, which are integrated nevertheless.
RooAbsFunc * _numIntegrand
! do not persist
static void setCacheAllNumeric(Int_t ndim)
Global switch to cache all integral values that integrate at least ndim dimensions numerically.
IntOperMode _intOperMode
integration operation mode
RooArgSet * _funcNormSet
Optional normalization set passed to function.
RooAbsIntegrator * _numIntEngine
! do not persist
bool _cacheNum
Cache integral if numeric.
double evaluate() const override
Perform the integration and return the result.
const RooArgSet & parameters() const
RooSetProxy _jacList
Set of lvalue observables over which is analytically integration that have a non-unit Jacobian.
bool isValidReal(double value, bool printError=false) const override
Check if current value is valid.
double getValV(const RooArgSet *set=0) const override
Return value of object.
RooArgSet _saveInt
! do not persist
RooSetProxy _anaList
Set of observables over which is integrated/summed analytically.
bool _restartNumIntEngine
! do not persist
bool servesExclusively(const RooAbsArg *server, const RooArgSet &exclLVBranches, const RooArgSet &allBranches) const
Utility function that returns true if 'object server' is a server to exactly one of the RooAbsArgs in...
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursive) override
Intercept server redirects and reconfigure internal object accordingly.
RooSetProxy _sumList
Set of discrete observable over which is summed numerically.
~RooRealIntegral() override
void printMetaArgs(std::ostream &os) const override
Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the...
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print the state of this object to the specified output stream.
virtual double integrate() const
Perform hybrid numerical/analytical integration over all real-valued dependents.
RooListProxy _sumCat
! do not persist
virtual double jacobianProduct() const
Return product of jacobian terms originating from analytical integration.
static Int_t getCacheAllNumeric()
Return minimum dimensions of numeric integration for which values are cached.
RooArgSet _saveSum
! do not persist
static Int_t _cacheAllNDim
! Cache all integrals with given numeric dimension
RooAbsReal * createIntegral(const RooArgSet &iset, const RooArgSet *nset=0, const RooNumIntConfig *cfg=0, const char *rangeName=0) const override
Create an object that represents the integral of the function over one or more observables listed in ...
void autoSelectDirtyMode()
Set appropriate cache operation mode for integral depending on cache operation mode of server objects...
RooArgSet _facListOwned
Owned components in _facList.
bool getAllowComponentSelection() const
Check if component selection is allowed.
bool empty() const
Check if empty.
The RooSuperCategory can join several RooAbsCategoryLValue objects into a single category.
bool setIndex(value_type index, bool printError=true) override
Set the value of the super category to the specified index.
bool inRange(const char *rangeName) const override
Check that all input category states are in the given range.
const T & arg() const
Return reference to object held in proxy.
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:130
Basic string class.
Definition: TString.h:136
TString & Append(const char *cs)
Definition: TString.h:564
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
@ NumIntegration
Definition: RooGlobalFunc.h:65
@ InputArguments
Definition: RooGlobalFunc.h:64
@ Integration
Definition: RooGlobalFunc.h:63