Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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
22Performs 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
34#include <RooAbsRealLValue.h>
35#include <RooArgSet.h>
36#include <RooConstVar.h>
37#include <RooDouble.h>
39#include <RooFuncWrapper.h>
40#include <RooHelpers.h>
41#include <RooInvTransform.h>
42#include <RooMsgService.h>
43#include <RooNameReg.h>
44#include <RooNumIntConfig.h>
45#include <RooNumIntFactory.h>
46#include <RooRealBinding.h>
47#include <RooSuperCategory.h>
48#include <RooTrace.h>
49
50#include "RooFitImplHelpers.h"
51
52#include <TClass.h>
53
54#include <iostream>
55#include <memory>
56
58
59namespace {
60
61struct ServerToAdd {
62 ServerToAdd(RooAbsArg *theArg, bool isShape) : arg{theArg}, isShapeServer{isShape} {}
63 RooAbsArg *arg = nullptr;
64 bool isShapeServer = false;
65};
66
67void addObservableToServers(RooAbsReal const &function, RooAbsArg &leaf, std::vector<ServerToAdd> &serversToAdd,
68 const char *rangeName)
69{
70 auto leaflv = dynamic_cast<RooAbsRealLValue *>(&leaf);
71 if (leaflv && leaflv->getBinning(rangeName).isParameterized()) {
72 oocxcoutD(&function, Integration)
73 << function.GetName() << " : Observable " << leaf.GetName()
74 << " has parameterized binning, add value dependence of boundary objects rather than shape of leaf"
75 << std::endl;
76 if (leaflv->getBinning(rangeName).lowBoundFunc()) {
77 serversToAdd.emplace_back(leaflv->getBinning(rangeName).lowBoundFunc(), false);
78 }
79 if (leaflv->getBinning(rangeName).highBoundFunc()) {
80 serversToAdd.emplace_back(leaflv->getBinning(rangeName).highBoundFunc(), false);
81 }
82 } else {
83 oocxcoutD(&function, Integration) << function.GetName() << ": Adding observable " << leaf.GetName()
84 << " as shape dependent" << std::endl;
85 serversToAdd.emplace_back(&leaf, true);
86 }
87}
88
89void addParameterToServers(RooAbsReal const &function, RooAbsArg &leaf, std::vector<ServerToAdd> &serversToAdd,
90 bool isShapeServer)
91{
92 if (!isShapeServer) {
93 oocxcoutD(&function, Integration) << function.GetName() << ": Adding parameter " << leaf.GetName()
94 << " as value dependent" << std::endl;
95 } else {
96 oocxcoutD(&function, Integration) << function.GetName() << ": Adding parameter " << leaf.GetName()
97 << " as shape dependent" << std::endl;
98 }
99 serversToAdd.emplace_back(&leaf, isShapeServer);
100}
101
102enum class MarkedState { Dependent, Independent, AlreadyAdded };
103
104/// Mark all args that recursively are value clients of "dep".
105void unmarkDepValueClients(RooAbsArg const &dep, RooArgSet const &args, std::vector<MarkedState> &marked)
106{
107 assert(args.size() == marked.size());
108 auto index = args.index(dep);
109 if (index >= 0) {
110 marked[index] = MarkedState::Dependent;
111 for (RooAbsArg *client : dep.valueClients()) {
112 unmarkDepValueClients(*client, args, marked);
113 }
114 }
115}
116
117std::vector<ServerToAdd>
118getValueAndShapeServers(RooAbsReal const &function, RooArgSet const &depList, const char *rangeName)
119{
120 std::vector<ServerToAdd> serversToAdd;
121
122 // Get the full computation graph and sort it topologically
123 RooArgList allArgsList;
124 function.treeNodeServerList(&allArgsList, nullptr, true, true, /*valueOnly=*/false, false);
125 RooArgSet allArgs{allArgsList};
126 allArgs.sortTopologically();
127
128 // Figure out what are all the value servers only
129 RooArgList allValueArgsList;
130 function.treeNodeServerList(&allValueArgsList, nullptr, true, true, /*valueOnly=*/true, false);
131 RooArgSet allValueArgs{allValueArgsList};
132
133 // All "marked" args will be added as value servers to the integral
134 std::vector<MarkedState> marked(allArgs.size(), MarkedState::Independent);
135 marked.back() = MarkedState::Dependent; // We don't want to consider the function itself
136
137 // Mark all args that are (indirect) value servers of the integration
138 // variable or the integration variable itself. If something was marked,
139 // it means the integration variable was in the compute graph and we will
140 // add it to the server list.
141 for (RooAbsArg *dep : depList) {
142 if (RooAbsArg * depInArgs = allArgs.find(dep->GetName())) {
143 unmarkDepValueClients(*depInArgs, allArgs, marked);
144 addObservableToServers(function, *depInArgs, serversToAdd, rangeName);
145 }
146 }
147
148 // We are adding all independent direct servers of the args depending on the
149 // integration variables
150 for (std::size_t i = 0; i < allArgs.size(); ++i) {
151 if(marked[i] == MarkedState::Dependent) {
152 for(RooAbsArg* server : allArgs[i]->servers()) {
153 int index = allArgs.index(server->GetName());
154 if(index >= 0 && marked[index] == MarkedState::Independent) {
155 addParameterToServers(function, *server, serversToAdd, !allValueArgs.find(*server));
156 marked[index] = MarkedState::AlreadyAdded;
157 }
158 }
159 }
160 }
161
162 return serversToAdd;
163}
164
165void fillAnIntOKDepList(RooAbsReal const &function, RooArgSet const &intDepList, RooArgSet &anIntOKDepList)
166{
167 for (const auto arg : function.servers()) {
168
169 // Dependent or parameter?
170 if (!arg->dependsOnValue(intDepList)) {
171 continue;
172 } else if (!arg->isValueServer(function) && !arg->isShapeServer(function)) {
173 // Skip arg if it is neither value or shape server
174 continue;
175 }
176
177 bool depOK(false);
178 // Check for integratable AbsRealLValue
179
180 if (arg->isDerived()) {
181 RooAbsRealLValue *realArgLV = dynamic_cast<RooAbsRealLValue *>(arg);
182 RooAbsCategoryLValue *catArgLV = dynamic_cast<RooAbsCategoryLValue *>(arg);
183 if ((realArgLV && intDepList.find(realArgLV->GetName()) && (realArgLV->isJacobianOK(intDepList) != 0)) ||
184 catArgLV) {
185
186 // Derived LValue with valid jacobian
187 depOK = true;
188
189 // Now, check for overlaps
190 bool overlapOK = true;
191 for (const auto otherArg : function.servers()) {
192 // skip comparison with self
193 if (arg == otherArg)
194 continue;
195 if (otherArg->IsA() == RooConstVar::Class())
196 continue;
197 if (arg->overlaps(*otherArg, true)) {
198 }
199 }
200 // coverity[DEADCODE]
201 if (!overlapOK)
202 depOK = false;
203 }
204 } else {
205 // Fundamental types are always OK
206 depOK = true;
207 }
208
209 // Add server to list of dependents that are OK for analytical integration
210 if (depOK) {
211 anIntOKDepList.add(*arg, true);
212 oocxcoutI(&function, Integration) << function.GetName() << ": Observable " << arg->GetName()
213 << " is suitable for analytical integration (if supported by p.d.f)"
214 << std::endl;
215 }
216 }
217}
218
219} // namespace
220
221
222using namespace std;
223
224
226
227
228////////////////////////////////////////////////////////////////////////////////
229
231{
233}
234
235
236
237////////////////////////////////////////////////////////////////////////////////
238/// Construct integral of 'function' over observables in 'depList'
239/// in range 'rangeName' with normalization observables 'funcNormSet'
240/// (for p.d.f.s). In the integral is performed to the maximum extent
241/// possible the internal (analytical) integrals advertised by function.
242/// The other integrations are performed numerically. The optional
243/// config object prescribes how these numeric integrations are configured.
244///
245/// \Note If pdf component selection was globally overridden to always include
246/// all components (either with RooAbsReal::globalSelectComp(bool) or a
247/// RooAbsReal::GlobalSelectComponentRAII), then any created integral will
248/// ignore component selections during its lifetime. This is especially useful
249/// when creating normalization or projection integrals.
250RooRealIntegral::RooRealIntegral(const char *name, const char *title,
251 const RooAbsReal& function, const RooArgSet& depList,
252 const RooArgSet* funcNormSet, const RooNumIntConfig* config,
253 const char* rangeName) :
254 RooAbsReal(name,title),
255 _valid(true),
256 _respectCompSelect{!_globalSelectComp},
257 _sumList("!sumList","Categories to be summed numerically",this,false,false),
258 _intList("!intList","Variables to be integrated numerically",this,false,false),
259 _anaList("!anaList","Variables to be integrated analytically",this,false,false),
260 _jacList("!jacList","Jacobian product term",this,false,false),
261 _facList("!facList","Variables independent of function",this,false,true),
262 _function("!func","Function to be integrated",this,false,false),
263 _iconfig((RooNumIntConfig*)config),
264 _sumCat("!sumCat","SuperCategory for summation",this,false,false),
265 _rangeName((TNamed*)RooNameReg::ptr(rangeName))
266{
267 // A) Check that all dependents are lvalues
268 //
269 // B) Check if list of dependents can be re-expressed in
270 // lvalues that are higher in the expression tree
271 //
272 // C) Check for dependents that the PDF insists on integrating
273 // analytically itself
274 //
275 // D) Make list of servers that can be integrated analytically
276 // Add all parameters/dependents as value/shape servers
277 //
278 // E) Interact with function to make list of objects actually integrated analytically
279 //
280 // F) Make list of numerical integration variables consisting of:
281 // - Category dependents of RealLValues in analytical integration
282 // - Leaf nodes server lists of function server that are not analytically integrated
283 // - Make Jacobian list for analytically integrated RealLValues
284 //
285 // G) Split numeric list in integration list and summation list
286 //
287
288 oocxcoutI(&function,Integration) << "RooRealIntegral::ctor(" << GetName() << ") Constructing integral of function "
289 << function.GetName() << " over observables" << depList << " with normalization "
290 << (funcNormSet?*funcNormSet:RooArgSet()) << " with range identifier "
291 << (rangeName?rangeName:"<none>") << std::endl ;
292
293
294 // Choose same expensive object cache as integrand
295 setExpensiveObjectCache(function.expensiveObjectCache()) ;
296// cout << "RRI::ctor(" << GetName() << ") setting expensive object cache to " << &expensiveObjectCache() << " as taken from " << function.GetName() << std::endl ;
297
298 // Use objects integrator configuration if none is specified
299 if (!_iconfig) _iconfig = (RooNumIntConfig*) function.getIntegratorConfig() ;
300
301 // Save private copy of funcNormSet, if supplied, excluding factorizing terms
302 if (funcNormSet) {
303 _funcNormSet = std::make_unique<RooArgSet>();
304 for (const auto nArg : *funcNormSet) {
305 if (function.dependsOn(*nArg)) {
306 _funcNormSet->addClone(*nArg) ;
307 }
308 }
309 }
310
311 //_funcNormSet = funcNormSet ? (RooArgSet*)funcNormSet->snapshot(false) : 0 ;
312
313 // Make internal copy of dependent list
314 RooArgSet intDepList(depList) ;
315
316 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
317 // * A) Check that all dependents are lvalues and filter out any
318 // dependents that the PDF doesn't explicitly depend on
319 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
320
321 for (auto arg : intDepList) {
322 if(!arg->isLValue()) {
323 coutE(InputArguments) << ClassName() << "::" << GetName() << ": cannot integrate non-lvalue ";
324 arg->Print("1");
325 _valid= false;
326 }
327 if (!function.dependsOn(*arg)) {
328 std::unique_ptr<RooAbsArg> argClone{static_cast<RooAbsArg*>(arg->Clone())};
329 _facList.add(*argClone) ;
330 addServer(*argClone,false,true) ;
331 _facListOwned.addOwned(std::move(argClone));
332 }
333 }
334
335 if (!_facList.empty()) {
336 oocxcoutI(&function,Integration) << function.GetName() << ": Factorizing obserables are " << _facList << std::endl ;
337 }
338
339
340 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
341 // * B) Check if list of dependents can be re-expressed in *
342 // * lvalues that are higher in the expression tree *
343 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
344
345
346 // Initial fill of list of LValue branches
347 RooArgSet exclLVBranches("exclLVBranches") ;
348 RooArgSet branchList,branchListVD ;
349 function.branchNodeServerList(&branchList) ;
350
351 for (auto branch: branchList) {
352 RooAbsRealLValue *realArgLV = dynamic_cast<RooAbsRealLValue*>(branch) ;
353 RooAbsCategoryLValue *catArgLV = dynamic_cast<RooAbsCategoryLValue*>(branch) ;
354 if ((realArgLV && (realArgLV->isJacobianOK(intDepList)!=0)) || catArgLV) {
355 exclLVBranches.add(*branch) ;
356// cout << "exclv branch = " << std::endl ;
357// branch->printCompactTree() ;
358 }
359 if (dependsOnValue(*branch)) {
360 branchListVD.add(*branch) ;
361 } else {
362// cout << "value of self does not depend on branch " << branch->GetName() << std::endl ;
363 }
364 }
365 exclLVBranches.remove(depList,true,true) ;
366
367 // Initial fill of list of LValue leaf servers (put in intDepList, but the
368 // instances that are in the actual computation graph of the function)
369 RooArgSet exclLVServers("exclLVServers") ;
370 function.getObservables(&intDepList, exclLVServers);
371
372 // Obtain mutual exclusive dependence by iterative reduction
373 bool converged(false) ;
374 while(!converged) {
375 converged=true ;
376
377 // Reduce exclLVServers to only those serving exclusively exclLVBranches
378 std::vector<RooAbsArg*> toBeRemoved;
379 for (auto server : exclLVServers) {
380 if (!servesExclusively(server,exclLVBranches,branchListVD)) {
381 toBeRemoved.push_back(server);
382 converged=false ;
383 }
384 }
385 exclLVServers.remove(toBeRemoved.begin(), toBeRemoved.end());
386
387 // Reduce exclLVBranches to only those depending exclusively on exclLVservers
388 // Attention: counting loop, since erasing from container
389 for (std::size_t i=0; i < exclLVBranches.size(); ++i) {
390 const RooAbsArg* branch = exclLVBranches[i];
391 RooArgSet brDepList;
392 branch->getObservables(&intDepList, brDepList);
393 RooArgSet bsList(brDepList,"bsList") ;
394 bsList.remove(exclLVServers,true,true) ;
395 if (!bsList.empty()) {
396 exclLVBranches.remove(*branch,true,true) ;
397 --i;
398 converged=false ;
399 }
400 }
401 }
402
403 // Eliminate exclLVBranches that do not depend on any LVServer
404 // Attention: Counting loop, since modifying container
405 for (std::size_t i=0; i < exclLVBranches.size(); ++i) {
406 const RooAbsArg* branch = exclLVBranches[i];
407 if (!branch->dependsOnValue(exclLVServers)) {
408 exclLVBranches.remove(*branch,true,true) ;
409 --i;
410 }
411 }
412
413 // Replace exclusive lvalue branch servers with lvalue branches
414 // WVE Don't do this for binned distributions - deal with this using numeric integration with transformed bin boundaroes
415 if (!exclLVServers.empty() && !function.isBinnedDistribution(exclLVBranches)) {
416// cout << "activating LVservers " << exclLVServers << " for use in integration " << std::endl ;
417 intDepList.remove(exclLVServers) ;
418 intDepList.add(exclLVBranches) ;
419
420 //cout << "intDepList removing exclLVServers " << exclLVServers << std::endl ;
421 //cout << "intDepList adding exclLVBranches " << exclLVBranches << std::endl ;
422
423 }
424
425
426 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
427 // * C) Check for dependents that the PDF insists on integrating *
428 // analytically itself *
429 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
430
431 RooArgSet anIntOKDepList ;
432 for (auto arg : intDepList) {
433 if (function.forceAnalyticalInt(*arg)) {
434 anIntOKDepList.add(*arg) ;
435 }
436 }
437
438 if (!anIntOKDepList.empty()) {
439 oocxcoutI(&function,Integration) << function.GetName() << ": Observables that function forcibly requires to be integrated internally " << anIntOKDepList << std::endl ;
440 }
441
442
443 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
444 // * D) Make list of servers that can be integrated analytically *
445 // Add all parameters/dependents as value/shape servers *
446 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
447
448 auto serversToAdd = getValueAndShapeServers(function, depList, rangeName);
449 fillAnIntOKDepList(function, intDepList, anIntOKDepList);
450 // We will not add the servers just now, because it makes only sense to add
451 // them once we have made sure that this integral is not operating in
452 // pass-through mode. It will be done at the end of this constructor.
453
454 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
455 // * E) interact with function to make list of objects actually integrated analytically *
456 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
457
458 _mode = function.getAnalyticalIntegralWN(anIntOKDepList,_anaList,_funcNormSet.get(),RooNameReg::str(_rangeName)) ;
459
460 // Avoid confusion -- if mode is zero no analytical integral is defined regardless of contents of _anaList
461 if (_mode==0) {
463 }
464
465 if (_mode!=0) {
466 oocxcoutI(&function,Integration) << function.GetName() << ": Function integrated observables " << _anaList << " internally with code " << _mode << std::endl ;
467 }
468
469 // WVE kludge: synchronize dset for use in analyticalIntegral
470 // LM : I think this is needed only if _funcNormSet is not an empty set
471 if (_funcNormSet && !_funcNormSet->empty()) {
472 function.getVal(_funcNormSet.get()) ;
473 }
474
475 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
476 // * F) Make list of numerical integration variables consisting of: *
477 // * - Category dependents of RealLValues in analytical integration *
478 // * - Expanded server lists of server that are not analytically integrated *
479 // * Make Jacobian list with analytically integrated RealLValues *
480 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
481
482 RooArgSet numIntDepList ;
483
484
485 // Loop over actually analytically integrated dependents
486 for (const auto arg : _anaList) {
487
488 // Process only derived RealLValues
489 if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class()) && arg->isDerived() && !arg->isFundamental()) {
490
491 // Add to list of Jacobians to calculate
492 _jacList.add(*arg) ;
493
494 // Add category dependent of LValueReal used in integration
495 auto argDepList = std::unique_ptr<RooArgSet>(arg->getObservables(&intDepList));
496 for (const auto argDep : *argDepList) {
497 if (argDep->IsA()->InheritsFrom(RooAbsCategoryLValue::Class()) && intDepList.contains(*argDep)) {
498 numIntDepList.add(*argDep,true) ;
499 }
500 }
501 }
502 }
503
504
505 // If nothing was integrated analytically, swap back LVbranches for LVservers for subsequent numeric integration
506 if (_anaList.empty()) {
507 if (!exclLVServers.empty()) {
508 //cout << "NUMINT phase analList is empty. exclLVServers = " << exclLVServers << std::endl ;
509 intDepList.remove(exclLVBranches) ;
510 intDepList.add(exclLVServers) ;
511 }
512 }
513 //cout << "NUMINT intDepList = " << intDepList << std::endl ;
514
515 // Loop again over function servers to add remaining numeric integrations
516 for (const auto arg : function.servers()) {
517
518 //cout << "processing server for numeric integration " << arg->ClassName() << "::" << arg->GetName() << std::endl ;
519
520 // Process only servers that are not treated analytically
521 if (!_anaList.find(arg->GetName()) && arg->dependsOn(intDepList)) {
522
523 // Process only derived RealLValues
524 if (dynamic_cast<RooAbsLValue*>(arg) && arg->isDerived() && intDepList.contains(*arg)) {
525 numIntDepList.add(*arg,true) ;
526 } else {
527
528 // WVE this will only get the observables, but not l-value transformations
529 // Expand server in final dependents
530 auto argDeps = std::unique_ptr<RooArgSet>(arg->getObservables(&intDepList));
531
532 if (!argDeps->empty()) {
533
534 // Add final dependents, that are not forcibly integrated analytically,
535 // to numerical integration list
536 for (const auto dep : *argDeps) {
537 if (!_anaList.find(dep->GetName())) {
538 numIntDepList.add(*dep,true) ;
539 }
540 }
541 }
542 }
543 }
544 }
545
546 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
547 // * G) Split numeric list in integration list and summation list *
548 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
549
550 // Split numeric integration list in summation and integration lists
551 for (const auto arg : numIntDepList) {
552 if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
553 _intList.add(*arg) ;
554 } else if (arg->IsA()->InheritsFrom(RooAbsCategoryLValue::Class())) {
555 _sumList.add(*arg) ;
556 }
557 }
558
559 if (!_anaList.empty()) {
560 oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _anaList << " are analytically integrated with code " << _mode << std::endl ;
561 }
562 if (!_intList.empty()) {
563 oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _intList << " are numerically integrated" << std::endl ;
564 }
565 if (!_sumList.empty()) {
566 oocxcoutI(&function,Integration) << function.GetName() << ": Observables " << _sumList << " are numerically summed" << std::endl ;
567 }
568
569
570 // Determine operating mode
571 if (!numIntDepList.empty()) {
572 // Numerical and optional Analytical integration
574 } else if (!_anaList.empty()) {
575 // Purely analytical integration
577 } else {
578 // No integration performed, where the function is a direct value server
580 _function._valueServer = true;
581 }
582 // We are only setting the function proxy now that it's clear if it's a value
583 // server or not.
584 _function.setArg(const_cast<RooAbsReal&>(function));
585
586 // Determine auto-dirty status
588
589 // Create value caches for _intList and _sumList
592
593
594 if (!_sumList.empty()) {
595 _sumCat.addOwned(std::make_unique<RooSuperCategory>(Form("%s_sumCat",GetName()),"sumCat",_sumList));
596 }
597
598 // Only if we are not in pass-through mode we need to add the shape and value
599 // servers separately.
601 for(auto const& toAdd : serversToAdd) {
602 addServer(*toAdd.arg, !toAdd.isShapeServer, toAdd.isShapeServer);
603 }
604 }
605
607}
608
609
610
611////////////////////////////////////////////////////////////////////////////////
612/// Set appropriate cache operation mode for integral depending on cache operation
613/// mode of server objects
614
616{
617 // If any of our servers are is forcedDirty or a projectedDependent, then we need to be ADirty
618 for (const auto server : _serverList) {
619 if (server->isValueServer(*this)) {
620 RooArgSet leafSet ;
621 server->leafNodeServerList(&leafSet) ;
622 for (const auto leaf : leafSet) {
623 if (leaf->operMode()==ADirty && leaf->isValueServer(*this)) {
625 break ;
626 }
627 if (leaf->getAttribute("projectedDependent")) {
629 break ;
630 }
631 }
632 }
633 }
634}
635
636
637
638////////////////////////////////////////////////////////////////////////////////
639/// Utility function that returns true if 'object server' is a server
640/// to exactly one of the RooAbsArgs in 'exclLVBranches'
641
642bool RooRealIntegral::servesExclusively(const RooAbsArg* server,const RooArgSet& exclLVBranches, const RooArgSet& allBranches) const
643{
644 // Determine if given server serves exclusively exactly one of the given nodes in exclLVBranches
645
646 // Special case, no LV servers available
647 if (exclLVBranches.empty()) return false ;
648
649 // If server has no clients and is not an LValue itself, return false
650 if (server->_clientList.empty() && exclLVBranches.find(server->GetName())) {
651 return false ;
652 }
653
654 // WVE must check for value relations only here!!!!
655
656 // Loop over all clients
657 Int_t numLVServ(0) ;
658 for (const auto client : server->valueClients()) {
659 // If client is not an LValue, recurse
660 if (!(exclLVBranches.find(client->GetName())==client)) {
661 if (allBranches.find(client->GetName())==client) {
662 if (!servesExclusively(client,exclLVBranches,allBranches)) {
663 // Client is a non-LValue that doesn't have an exclusive LValue server
664 return false ;
665 }
666 }
667 } else {
668 // Client is an LValue
669 numLVServ++ ;
670 }
671 }
672
673 return (numLVServ==1) ;
674}
675
676
677
678
679////////////////////////////////////////////////////////////////////////////////
680/// (Re)Initialize numerical integration engine if necessary. Return true if
681/// successful, or otherwise false.
682
684{
685 // if we already have an engine, check if it still works for the present limits.
686 if(_numIntEngine) {
687 if(_numIntEngine->isValid() && _numIntEngine->checkLimits() && !_restartNumIntEngine ) return true;
688 // otherwise, cleanup the old engine
689 _numIntEngine.reset();
690 _numIntegrand.reset();
691 }
692
693 // All done if there are no arguments to integrate numerically
694 if(_intList.empty()) return true;
695
696 // Bind the appropriate analytic integral of our RooRealVar object to
697 // those of its arguments that will be integrated out numerically.
698 if(_mode != 0) {
699 std::unique_ptr<RooAbsReal> analyticalPart{_function->createIntegral(_anaList, *funcNormSet(), RooNameReg::str(_rangeName))};
700 _numIntegrand = std::make_unique<RooRealBinding>(*analyticalPart,_intList,nullptr,false,_rangeName);
701 const_cast<RooRealIntegral*>(this)->addOwnedComponents(std::move(analyticalPart));
702 }
703 else {
704 _numIntegrand = std::make_unique<RooRealBinding>(*_function,_intList,funcNormSet(),false,_rangeName);
705 }
706 if(nullptr == _numIntegrand || !_numIntegrand->isValid()) {
707 coutE(Integration) << ClassName() << "::" << GetName() << ": failed to create valid integrand." << std::endl;
708 return false;
709 }
710
711 // Create appropriate numeric integrator using factory
712 bool isBinned = _function->isBinnedDistribution(_intList) ;
714
715 if(_numIntEngine == nullptr || !_numIntEngine->isValid()) {
716 coutE(Integration) << ClassName() << "::" << GetName() << ": failed to create valid integrator." << std::endl;
717 return false;
718 }
719
720 cxcoutI(NumIntegration) << "RooRealIntegral::init(" << GetName() << ") using numeric integrator "
721 << _numIntEngine->ClassName() << " to calculate Int" << _intList << std::endl ;
722
723 if (_intList.size()>3) {
724 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 ;
725 }
726
727 _restartNumIntEngine = false ;
728 return true;
729}
730
731
732
733////////////////////////////////////////////////////////////////////////////////
734/// Copy constructor
735
737 RooAbsReal(other,name),
738 _valid(other._valid),
739 _respectCompSelect(other._respectCompSelect),
740 _sumList("!sumList",this,other._sumList),
741 _intList("!intList",this,other._intList),
742 _anaList("!anaList",this,other._anaList),
743 _jacList("!jacList",this,other._jacList),
744 _facList("!facList","Variables independent of function",this,false,true),
745 _function("!func",this,other._function),
746 _iconfig(other._iconfig),
747 _sumCat("!sumCat",this,other._sumCat),
748 _mode(other._mode),
749 _intOperMode(other._intOperMode),
750 _restartNumIntEngine(false),
751 _rangeName(other._rangeName),
752 _cacheNum(false)
753{
754 if(other._funcNormSet) {
755 _funcNormSet = std::make_unique<RooArgSet>();
756 other._funcNormSet->snapshot(*_funcNormSet, false);
757 }
758
759 for (const auto arg : other._facList) {
760 std::unique_ptr<RooAbsArg> argClone{static_cast<RooAbsArg*>(arg->Clone())};
761 _facList.add(*argClone) ;
762 addServer(*argClone,false,true) ;
763 _facListOwned.addOwned(std::move(argClone));
764 }
765
766 other._intList.snapshot(_saveInt) ;
767 other._sumList.snapshot(_saveSum) ;
768
770}
771
772
773////////////////////////////////////////////////////////////////////////////////
774
776{
778}
779
780
781////////////////////////////////////////////////////////////////////////////////
782
783RooFit::OwningPtr<RooAbsReal> RooRealIntegral::createIntegral(const RooArgSet& iset, const RooArgSet* nset, const RooNumIntConfig* cfg, const char* rangeName) const
784{
785 // Handle special case of no integration with default algorithm
786 if (iset.empty()) {
787 return RooAbsReal::createIntegral(iset,nset,cfg,rangeName) ;
788 }
789
790 // Special handling of integral of integral, return RooRealIntegral that represents integral over all dimensions in one pass
791 RooArgSet isetAll(iset) ;
792 isetAll.add(_sumList) ;
793 isetAll.add(_intList) ;
794 isetAll.add(_anaList) ;
795 isetAll.add(_facList) ;
796
797 const RooArgSet* newNormSet(nullptr) ;
798 std::unique_ptr<RooArgSet> tmp;
799 if (nset && !_funcNormSet) {
800 newNormSet = nset ;
801 } else if (!nset && _funcNormSet) {
802 newNormSet = _funcNormSet.get();
803 } else if (nset && _funcNormSet) {
804 tmp = std::make_unique<RooArgSet>();
805 tmp->add(*nset) ;
806 tmp->add(*_funcNormSet,true) ;
807 newNormSet = tmp.get();
808 }
809 return _function->createIntegral(isetAll,newNormSet,cfg,rangeName);
810}
811
812
813
814////////////////////////////////////////////////////////////////////////////////
815/// Return value of object. If the cache is clean, return the
816/// cached value, otherwise recalculate on the fly and refill
817/// the cache
818
819double RooRealIntegral::getValV(const RooArgSet* nset) const
820{
821// // fast-track clean-cache processing
822// if (_operMode==AClean) {
823// return _value ;
824// }
825
826 if (nset && nset->uniqueId().value() != _lastNormSetId) {
827 const_cast<RooRealIntegral*>(this)->setProxyNormSet(nset);
828 _lastNormSetId = nset->uniqueId().value();
829 }
830
832 _value = traceEval(nset) ;
833 }
834
835 return _value ;
836}
837
838
839
840
841
842////////////////////////////////////////////////////////////////////////////////
843/// Perform the integration and return the result
844
846{
848
849 double retVal(0) ;
850 switch (_intOperMode) {
851
852 case Hybrid:
853 {
854 // Cache numeric integrals in >1d expensive object cache
855 RooDouble* cacheVal(nullptr) ;
858 }
859
860 if (cacheVal) {
861 retVal = *cacheVal ;
862 // cout << "using cached value of integral" << GetName() << std::endl ;
863 } else {
864
865
866 // Find any function dependents that are AClean
867 // and switch them temporarily to ADirty
868 bool origState = inhibitDirty() ;
869 setDirtyInhibit(true) ;
870
871 // try to initialize our numerical integration engine
872 if(!(_valid= initNumIntegrator())) {
873 coutE(Integration) << ClassName() << "::" << GetName()
874 << ":evaluate: cannot initialize numerical integrator" << std::endl;
875 return 0;
876 }
877
878 // Save current integral dependent values
881
882 // Evaluate sum/integral
883 retVal = sum() ;
884
885
886 // This must happen BEFORE restoring dependents, otherwise no dirty state propagation in restore step
887 setDirtyInhibit(origState) ;
888
889 // Restore integral dependent values
892
893 // Cache numeric integrals in >1d expensive object cache
895 RooDouble* val = new RooDouble(retVal) ;
897 // cout << "### caching value of integral" << GetName() << " in " << &expensiveObjectCache() << std::endl ;
898 }
899
900 }
901 break ;
902 }
903 case Analytic:
904 {
906 cxcoutD(Tracing) << "RooRealIntegral::evaluate_analytic(" << GetName()
907 << ")func = " << _function->ClassName() << "::" << _function->GetName()
908 << " raw = " << retVal << " _funcNormSet = " << (_funcNormSet?*_funcNormSet:RooArgSet()) << std::endl ;
909
910
911 break ;
912 }
913
914 case PassThrough:
915 {
916 // In pass through mode, the RooRealIntegral should have registered the
917 // function as a value server, because we directly depend on its value.
918 assert(_function.isValueServer());
919 // There should be no other servers besides the actual function and the
920 // factorized observables that the function doesn't depend on but are
921 // integrated over later.
922 assert(servers().size() == _facList.size() + 1);
923
924 //setDirtyInhibit(true) ;
925 retVal= _function->getVal(funcNormSet());
926 //setDirtyInhibit(false) ;
927 break ;
928 }
929 }
930
931
932 // Multiply answer with integration ranges of factorized variables
933 if (!_facList.empty()) {
934 for (const auto arg : _facList) {
935 // Multiply by fit range for 'real' dependents
936 if (arg->IsA()->InheritsFrom(RooAbsRealLValue::Class())) {
937 RooAbsRealLValue* argLV = (RooAbsRealLValue*)arg ;
938 retVal *= (argLV->getMax(intRange()) - argLV->getMin(intRange())) ;
939 }
940 // Multiply by number of states for category dependents
941 if (arg->IsA()->InheritsFrom(RooAbsCategoryLValue::Class())) {
943 retVal *= argLV->numTypes() ;
944 }
945 }
946 }
947
948
949 if (dologD(Tracing)) {
950 cxcoutD(Tracing) << "RooRealIntegral::evaluate(" << GetName() << ") anaInt = " << _anaList << " numInt = " << _intList << _sumList << " mode = " ;
951 switch(_intOperMode) {
952 case Hybrid: ccoutD(Tracing) << "Hybrid" ; break ;
953 case Analytic: ccoutD(Tracing) << "Analytic" ; break ;
954 case PassThrough: ccoutD(Tracing) << "PassThrough" ; break ;
955 }
956
957 ccxcoutD(Tracing) << "raw*fact = " << retVal << std::endl ;
958 }
959
960 return retVal ;
961}
962
963
964
965////////////////////////////////////////////////////////////////////////////////
966/// Return product of jacobian terms originating from analytical integration
967
969{
970 if (_jacList.empty()) {
971 return 1 ;
972 }
973
974 double jacProd(1) ;
975 for (const auto elm : _jacList) {
976 auto arg = static_cast<const RooAbsRealLValue*>(elm);
977 jacProd *= arg->jacobian() ;
978 }
979
980 // Take std::abs() here: if jacobian is negative, min and max are swapped and analytical integral
981 // will be positive, so must multiply with positive jacobian.
982 return std::abs(jacProd) ;
983}
984
985
986
987////////////////////////////////////////////////////////////////////////////////
988/// Perform summation of list of category dependents to be integrated
989
991{
992 if (!_sumList.empty()) {
993 // Add integrals for all permutations of categories summed over
994 double total(0) ;
995
997 for (const auto& nameIdx : *sumCat) {
998 sumCat->setIndex(nameIdx);
999 if (!_rangeName || sumCat->inRange(RooNameReg::str(_rangeName))) {
1001 }
1002 }
1003
1004 return total ;
1005
1006 } else {
1007 // Simply return integral
1008 double ret = integrate() / jacobianProduct() ;
1009 return ret ;
1010 }
1011}
1012
1013
1014////////////////////////////////////////////////////////////////////////////////
1015/// Perform hybrid numerical/analytical integration over all real-valued dependents
1016
1018{
1019 if (!_numIntEngine) {
1020 // Trivial case, fully analytical integration
1022 } else {
1023 return _numIntEngine->calculate() ;
1024 }
1025}
1026
1027
1028
1029////////////////////////////////////////////////////////////////////////////////
1030/// Intercept server redirects and reconfigure internal object accordingly
1031
1033 bool mustReplaceAll, bool nameChange, bool isRecursive)
1034{
1035 _restartNumIntEngine = true ;
1036
1038
1039 // Update contents value caches for _intList and _sumList
1044
1045 // Delete parameters cache if we have one
1046 _params.reset();
1047
1048 return RooAbsReal::redirectServersHook(newServerList, mustReplaceAll, nameChange, isRecursive);
1049}
1050
1051
1052
1053////////////////////////////////////////////////////////////////////////////////
1054
1056{
1057 if (!_params) {
1058 _params = std::make_unique<RooArgSet>("params") ;
1059
1060 RooArgSet params ;
1061 for (const auto server : _serverList) {
1062 if (server->isValueServer(*this)) _params->add(*server) ;
1063 }
1064 }
1065
1066 return *_params ;
1067}
1068
1069
1070////////////////////////////////////////////////////////////////////////////////
1071/// Check if current value is valid
1072
1073bool RooRealIntegral::isValidReal(double /*value*/, bool /*printError*/) const
1074{
1075 return true ;
1076}
1077
1078////////////////////////////////////////////////////////////////////////////////
1079/// Check if component selection is allowed
1080
1082 return _respectCompSelect;
1083}
1084
1085////////////////////////////////////////////////////////////////////////////////
1086/// Set component selection to be allowed/forbidden
1087
1089 _respectCompSelect = allow;
1090}
1091
1093{
1094 if (_sumList.empty() && _intList.empty()) {
1096 return;
1097 }
1098
1099 if (intVars().size() != 1 || _intList.size() != 1) {
1100 std::stringstream errorMsg;
1101 errorMsg << "Only analytical integrals and 1D numeric integrals are supported for AD for class"
1102 << _function.GetName();
1103 coutE(Minimization) << errorMsg.str() << std::endl;
1104 throw std::runtime_error(errorMsg.str().c_str());
1105 }
1106
1107 auto &intVar = static_cast<RooAbsRealLValue &>(*_intList[0]);
1108
1109 RooFuncWrapper wrapper{GetName(), GetTitle(), *_function, {}, nullptr, nullptr, false};
1110
1111 RooArgSet params;
1112 _function->getParameters(nullptr, params);
1113
1114 std::string paramsName = ctx.getTmpVarName();
1115
1116 std::stringstream ss;
1117
1118 std::string resName = RooFit::Detail::makeValidVarName(GetName()) + "Result";
1119 ctx.addResult(this, resName);
1120 ctx.addToGlobalScope("double " + resName + " = 0.0;\n");
1121
1122 ss << "double " << paramsName << "[] = {";
1123 std::string args;
1124 int intVarIdx = 0;
1125 for (RooAbsArg *param : params) {
1126 // Fill the integration variable with dummy value for now. This will then
1127 // be reset in the sampling loop.
1128 if (param->namePtr() == intVar.namePtr()) {
1129 args += "0.0,";
1130 } else if (!param->isConstant()) {
1131 args += ctx.getResult(*param) + ",";
1132 intVarIdx++;
1133 }
1134 }
1135 if (!args.empty()) {
1136 args.pop_back();
1137 }
1138 ss << args << "};\n";
1139
1140 // TODO: once Clad has support for higher-order functions (follow also the
1141 // Clad issue #637), we could refactor this code into an actual function
1142 // instead of hardcoding it here as a string.
1143 ss << "{\n"
1144 << " const int n = 1000; // number of sampling points\n"
1145 << " double d = " << intVar.getMax(intRange()) << " - " << intVar.getMin(intRange()) << ";\n"
1146 << " double eps = d / n;\n"
1147 << " for (int i = 0; i < n; ++i) {\n"
1148 << " " << paramsName << "[" << intVarIdx << "] = " << intVar.getMin(intRange()) << " + eps * i;\n"
1149 // TODO: the second "paramsName" should be nullptr, but until Clad issue
1150 // 636 is fixed, we have to pass a non-nullptr dummy.
1151 << " " << resName << " += " << " + " << ctx.buildCall(wrapper.funcName(), paramsName, paramsName) << ";\n"
1152 << " }\n"
1153 << " " << resName << " *= " << " d / n;\n"
1154 << "}\n";
1155
1156 ctx.addToGlobalScope(ss.str());
1157}
1158
1159////////////////////////////////////////////////////////////////////////////////
1160/// Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the
1161/// integration operation
1162
1163void RooRealIntegral::printMetaArgs(ostream& os) const
1164{
1165
1166 if (!intVars().empty()) {
1167 os << "Int " ;
1168 }
1169 os << _function->GetName() ;
1170 if (_funcNormSet) {
1171 os << "_Norm" << *_funcNormSet << " " ;
1172 }
1173
1174 // List internally integrated observables and factorizing observables as analytically integrated
1175 RooArgSet tmp(_anaList) ;
1176 tmp.add(_facList) ;
1177 if (!tmp.empty()) {
1178 os << "d[Ana]" << tmp << " ";
1179 }
1180
1181 // List numerically integrated and summed observables as numerically integrated
1182 RooArgSet tmp2(_intList) ;
1183 tmp2.add(_sumList) ;
1184 if (!tmp2.empty()) {
1185 os << " d[Num]" << tmp2 << " ";
1186 }
1187}
1188
1189
1190
1191////////////////////////////////////////////////////////////////////////////////
1192/// Print the state of this object to the specified output stream.
1193
1194void RooRealIntegral::printMultiline(ostream& os, Int_t contents, bool verbose, TString indent) const
1195{
1196 RooAbsReal::printMultiline(os,contents,verbose,indent) ;
1197 os << indent << "--- RooRealIntegral ---" << std::endl;
1198 os << indent << " Integrates ";
1200 TString deeper(indent);
1201 deeper.Append(" ");
1202 os << indent << " operating mode is "
1203 << (_intOperMode==Hybrid?"Hybrid":(_intOperMode==Analytic?"Analytic":"PassThrough")) << std::endl ;
1204 os << indent << " Summed discrete args are " << _sumList << std::endl ;
1205 os << indent << " Numerically integrated args are " << _intList << std::endl;
1206 os << indent << " Analytically integrated args using mode " << _mode << " are " << _anaList << std::endl ;
1207 os << indent << " Arguments included in Jacobian are " << _jacList << std::endl ;
1208 os << indent << " Factorized arguments are " << _facList << std::endl ;
1209 os << indent << " Function normalization set " ;
1210 if (_funcNormSet)
1211 _funcNormSet->Print("1") ;
1212 else
1213 os << "<none>" ;
1214
1215 os << std::endl ;
1216}
1217
1218
1219
1220////////////////////////////////////////////////////////////////////////////////
1221/// Global switch to cache all integral values that integrate at least ndim dimensions numerically
1222
1224 _cacheAllNDim = ndim ;
1225}
1226
1227
1228////////////////////////////////////////////////////////////////////////////////
1229/// Return minimum dimensions of numeric integration for which values are cached.
1230
1232{
1233 return _cacheAllNDim ;
1234}
1235
1236
1237std::unique_ptr<RooAbsArg>
1239{
1241}
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
#define cxcoutI(a)
#define cxcoutD(a)
#define oocxcoutD(o, a)
#define dologD(a)
#define coutE(a)
#define ccxcoutD(a)
#define ccoutD(a)
#define oocxcoutI(o, a)
#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:377
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:2467
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:79
RooExpensiveObjectCache & expensiveObjectCache() const
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.
void setOperMode(OperMode mode, bool recurseADirty=true)
Set the operation mode of this node.
RooFit::OwningPtr< RooArgSet > getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
const RefCountList_t & valueClients() const
List of all value clients of this object. Value clients receive value updates.
Definition RooAbsArg.h:198
virtual void setExpensiveObjectCache(RooExpensiveObjectCache &cache)
Definition RooAbsArg.h:503
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
static void setDirtyInhibit(bool flag)
Control global dirty inhibit mode.
virtual std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const
const RefCountList_t & servers() const
List of all servers of this object.
Definition RooAbsArg.h:208
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:108
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.
bool getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
bool isValueOrShapeDirtyAndClear() const
Definition RooAbsArg.h:455
bool inhibitDirty() const
Delete watch flag.
void setProxyNormSet(const RooArgSet *nset)
Forward a change in the cached normalization argset to all the registered proxies.
RefCountList_t _clientList
Definition RooAbsArg.h:635
RefCountList_t _serverList
Definition RooAbsArg.h:634
bool isValueServer(const RooAbsArg &arg) const
Check if this is serving values to arg.
Definition RooAbsArg.h:225
OperMode operMode() const
Query the operation mode of this node.
Definition RooAbsArg.h:484
Abstract base class for objects that represent a discrete value that can be set from the outside,...
static TClass * Class()
Int_t numTypes(const char *=nullptr) const
Return number of types defined (in range named rangeName if rangeName!=nullptr)
Abstract container object that can hold multiple RooAbsArg objects.
RooFit::UniqueId< RooAbsCollection > const & uniqueId() const
Returns a unique ID that is different for every instantiated RooAbsCollection.
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.
Int_t getSize() const
Return the number of elements in the collection.
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.
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
Abstract base class for objects that are lvalues, i.e.
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.
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Create an object that represents the integral of the function over one or more observables listed in ...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Structure printing.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursiveStep) override
Function that is called at the end of redirectServers().
double _value
Cache for current value of object.
Definition RooAbsReal.h:543
double traceEval(const RooArgSet *set) const
Calculate current value of object, with error tracing wrapper.
virtual std::string buildCallToAnalyticIntegral(Int_t code, const char *rangeName, RooFit::Detail::CodeSquashContext &ctx) const
This function defines the analytical integral translation for the class.
RooFit::UniqueId< RooArgSet >::Value_t _lastNormSetId
Component selection flag for RooAbsPdf::plotCompOn.
Definition RooAbsReal.h:550
virtual double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
virtual bool isBinnedDistribution(const RooArgSet &) const
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition RooAbsReal.h:353
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:60
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition RooArgSet.h:178
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,...
A class to maintain the context for squashing of RooFit models into code.
std::string buildCall(std::string const &funcname, Args_t const &...args)
Build the code to call the function with name funcname, passing some arguments.
void addResult(RooAbsArg const *key, std::string const &value)
A function to save an expression that includes/depends on the result of the input node.
std::string const & getResult(RooAbsArg const &arg)
Gets the result for the given node using the node name.
std::string getTmpVarName()
Get a unique variable name to be used in the generated code.
void addToGlobalScope(std::string const &str)
Adds the given string to the string block that will be emitted at the top of the squashed function.
A wrapper class to store a C++ function of type 'double (*)(double*, double*)'.
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...
std::unique_ptr< RooAbsIntegrator > createIntegrator(RooAbsFunc &func, const RooNumIntConfig &config, Int_t ndim=0, bool isBinned=false) const
Construct a numeric integrator instance that operates on function 'func' and is configured with 'conf...
static RooNumIntFactory & instance()
Static method returning reference to singleton instance of factory.
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer,...
Performs hybrid numerical/analytical integrals of RooAbsReal objects.
RooNumIntConfig * _iconfig
bool initNumIntegrator() const
(Re)Initialize numerical integration engine if necessary.
RooArgSet const * funcNormSet() const
RooFit::OwningPtr< 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 listed in ...
void translate(RooFit::Detail::CodeSquashContext &ctx) const override
This function defines a translation for each RooAbsReal based object that can be used to express the ...
void setAllowComponentSelection(bool allow)
Set component selection to be allowed/forbidden.
RooRealProxy _function
Function being integration.
RooArgSet intVars() const
RooSetProxy _intList
Set of continuous observables over which is integrated numerically.
virtual double sum() const
Perform summation of list of category dependents to be integrated.
RooSetProxy _facList
Set of observables on which function does not depends, which are integrated nevertheless.
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.
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
void autoSelectDirtyMode()
Set appropriate cache operation mode for integral depending on cache operation mode of server objects...
RooArgSet _facListOwned
Owned components in _facList.
const char * intRange() const
bool getAllowComponentSelection() const
Check if component selection is allowed.
Joins 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.
const T & arg() const
Return reference to object held in proxy.
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:207
Basic string class.
Definition TString.h:139
TString & Append(const char *cs)
Definition TString.h:576
void function(const Char_t *name_, T fun, const Char_t *docstring=0)
Definition RExports.h:167
std::string makeValidVarName(std::string const &in)
T * OwningPtr
An alias for raw pointers for indicating that the return type of a RooFit function is an owning point...
Definition Config.h:43
constexpr Value_t value() const
Return numerical value of ID.
Definition UniqueId.h:59