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