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