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