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