Logo ROOT  
Reference Guide
RooAbsArg.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/** \class RooAbsArg
19 \ingroup Roofitcore
20
21RooAbsArg is the common abstract base class for objects that
22represent a value and a "shape" in RooFit. Values or shapes usually depend on values
23or shapes of other RooAbsArg instances. Connecting several RooAbsArg in
24a computation graph models an expression tree that can be evaluated.
25
26### Building a computation graph of RooFit objects
27Therefore, RooAbsArg provides functionality to connect objects of type RooAbsArg into
28a computation graph to pass values between those objects.
29A value can e.g. be a real-valued number, (instances of RooAbsReal), or an integer, that is,
30catgory index (instances of RooAbsCategory). The third subclass of RooAbsArg is RooStringVar,
31but it is rarely used.
32
33The "shapes" that a RooAbsArg can possess can e.g. be the definition
34range of an observable, or how many states a category object has. In computations,
35values are expected to change often, while shapes remain mostly constant
36(unless e.g. a new range is set for an observable).
37
38Nodes of a computation graph are connected using instances of RooAbsProxy.
39If Node B declares a member `RooTemplateProxy<TypeOfNodeA>`, Node A will be
40registered as a server of values to Node B, and Node B will know that it is
41a client of node A. Using functions like dependsOn(), or getObservables()
42/ getParameters(), the relation of `A --> B` can be queried. Using graphVizTree(),
43one can create a visualisation of the expression tree.
44
45
46An instance of RooAbsArg can have named attributes. It also has flags
47to indicate that either its value or its shape were changed (= it is dirty).
48RooAbsArg provides functionality to manage client/server relations in
49a computation graph (\ref clientServerInterface), and helps propagating
50value/shape changes through the graph. RooAbsArg implements interfaces
51for inspecting client/server relationships (\ref clientServerInterface) and
52setting/clearing/querying named attributes.
53
54### Caching of values
55The values of nodes in the computation graph are cached in RooFit. If
56a value is used in two nodes of a graph, it doesn't need to be recomputed. If
57a node acquires a new value, it notifies its consumers ("clients") that
58their cached values are dirty. See the functions in \ref optimisationInterface
59for details.
60A node uses its isValueDirty() and isShapeDirty() functions to decide if a
61computation is necessary. Caching can be vetoed globally by setting a
62bit using setDirtyInhibit(). This will make computations slower, but all the
63nodes of the computation graph will be evaluated irrespective of whether their
64state is clean or dirty. Using setOperMode(), caching can also be enabled/disabled
65for single nodes.
66
67*/
68
69#include "RooFit.h"
70
71#include "TBuffer.h"
72#include "TClass.h"
74#include "strlcpy.h"
75
76#include "RooSecondMoment.h"
77#include "RooWorkspace.h"
78
79#include "RooMsgService.h"
80#include "RooAbsArg.h"
81#include "RooArgSet.h"
82#include "RooArgProxy.h"
83#include "RooSetProxy.h"
84#include "RooListProxy.h"
85#include "RooAbsData.h"
87#include "RooAbsRealLValue.h"
88#include "RooTrace.h"
89#include "RooRealIntegral.h"
90#include "RooConstVar.h"
92#include "RooAbsDataStore.h"
93#include "RooResolutionModel.h"
94#include "RooVectorDataStore.h"
95#include "RooTreeDataStore.h"
96#include "RooHelpers.h"
97
98#include <algorithm>
99#include <cstring>
100#include <fstream>
101#include <iostream>
102#include <memory>
103#include <sstream>
104
105using namespace std;
106
107#if (__GNUC__==3&&__GNUC_MINOR__==2&&__GNUC_PATCHLEVEL__==3)
108char* operator+( streampos&, char* );
109#endif
110
112;
113
117
118std::map<RooAbsArg*,std::unique_ptr<TRefArray>> RooAbsArg::_ioEvoList;
119std::stack<RooAbsArg*> RooAbsArg::_ioReadStack ;
120
121
122////////////////////////////////////////////////////////////////////////////////
123/// Default constructor
124
126 : TNamed(), _deleteWatch(kFALSE), _valueDirty(kTRUE), _shapeDirty(kTRUE), _operMode(Auto), _fast(kFALSE), _ownedComponents(nullptr),
127 _prohibitServerRedirect(kFALSE), _namePtr(0), _isConstant(kFALSE), _localNoInhibitDirty(kFALSE),
128 _myws(0)
129{
131
132}
133
134////////////////////////////////////////////////////////////////////////////////
135/// Create an object with the specified name and descriptive title.
136/// The newly created object has no clients or servers and has its
137/// dirty flags set.
138
139RooAbsArg::RooAbsArg(const char *name, const char *title)
140 : TNamed(name, title), _deleteWatch(kFALSE), _valueDirty(kTRUE), _shapeDirty(kTRUE), _operMode(Auto), _fast(kFALSE),
141 _ownedComponents(0), _prohibitServerRedirect(kFALSE), _namePtr(0), _isConstant(kFALSE),
142 _localNoInhibitDirty(kFALSE), _myws(0)
143{
144 if (name == nullptr || strlen(name) == 0) {
145 throw std::logic_error("Each RooFit object needs a name. "
146 "Objects representing the same entity (e.g. an observable 'x') are identified using their name.");
147 }
149}
150
151////////////////////////////////////////////////////////////////////////////////
152/// Copy constructor transfers all boolean and string properties of the original
153/// object. Transient properties and client-server links are not copied
154
155RooAbsArg::RooAbsArg(const RooAbsArg &other, const char *name)
156 : TNamed(other.GetName(), other.GetTitle()), RooPrintable(other), _boolAttrib(other._boolAttrib),
157 _stringAttrib(other._stringAttrib), _deleteWatch(other._deleteWatch), _operMode(Auto), _fast(kFALSE),
158 _ownedComponents(0), _prohibitServerRedirect(kFALSE), _namePtr(other._namePtr),
159 _isConstant(other._isConstant), _localNoInhibitDirty(other._localNoInhibitDirty), _myws(0)
160{
161 // Use name in argument, if supplied
162 if (name) {
165 } else {
166 // Same name, don't recalculate name pointer (expensive)
167 TNamed::SetName(other.GetName()) ;
168 _namePtr = other._namePtr ;
169 }
170
171 // Copy server list by hand
172 Bool_t valueProp, shapeProp ;
173 for (const auto server : other._serverList) {
174 valueProp = server->_clientListValue.containsByNamePtr(&other);
175 shapeProp = server->_clientListShape.containsByNamePtr(&other);
176 addServer(*server,valueProp,shapeProp) ;
177 }
178
179 setValueDirty() ;
180 setShapeDirty() ;
181
182 //setAttribute(Form("CloneOf(%08x)",&other)) ;
183 //cout << "RooAbsArg::cctor(" << this << ") #bools = " << _boolAttrib.size() << " #strings = " << _stringAttrib.size() << endl ;
184
185}
186
187////////////////////////////////////////////////////////////////////////////////
188/// Assign all boolean and string properties of the original
189/// object. Transient properties and client-server links are not assigned.
191 TNamed::operator=(other);
193 _boolAttrib = other._boolAttrib;
196 _operMode = other._operMode;
197 _fast = other._fast;
198 _ownedComponents = nullptr;
200 _namePtr = other._namePtr;
201 _isConstant = other._isConstant;
203 _myws = nullptr;
204
205 bool valueProp, shapeProp;
206 for (const auto server : other._serverList) {
207 valueProp = server->_clientListValue.containsByNamePtr(&other);
208 shapeProp = server->_clientListShape.containsByNamePtr(&other);
209 addServer(*server,valueProp,shapeProp) ;
210 }
211
214
215 return *this;
216}
217
218
219////////////////////////////////////////////////////////////////////////////////
220/// Destructor.
221
223{
224 // Notify all servers that they no longer need to serve us
225 while (!_serverList.empty()) {
227 }
228
229 // Notify all clients that they are in limbo
230 std::vector<RooAbsArg*> clientListTmp(_clientList.begin(), _clientList.end()); // have to copy, as we invalidate iterators
232 for (auto client : clientListTmp) {
233 client->setAttribute("ServerDied") ;
234 TString attr("ServerDied:");
235 attr.Append(GetName());
236 attr.Append(Form("(%zx)",(size_t)this)) ;
237 client->setAttribute(attr.Data());
238 client->removeServer(*this,kTRUE);
239
240 if (_verboseDirty) {
241
242 if (first) {
243 cxcoutD(Tracing) << "RooAbsArg::dtor(" << GetName() << "," << this << ") DeleteWatch: object is being destroyed" << endl ;
244 first = kFALSE ;
245 }
246
247 cxcoutD(Tracing) << fName << "::" << ClassName() << ":~RooAbsArg: dependent \""
248 << client->GetName() << "\" should have been deleted first" << endl ;
249 }
250 }
251
252 if (_ownedComponents) {
253 delete _ownedComponents ;
254 _ownedComponents = 0 ;
255 }
256
257}
258
259
260////////////////////////////////////////////////////////////////////////////////
261/// Control global dirty inhibit mode. When set to true no value or shape dirty
262/// flags are propagated and cache is always considered to be dirty.
263
265{
266 _inhibitDirty = flag ;
267}
268
269
270////////////////////////////////////////////////////////////////////////////////
271/// Activate verbose messaging related to dirty flag propagation
272
274{
275 _verboseDirty = flag ;
276}
277
278////////////////////////////////////////////////////////////////////////////////
279/// Check if this object was created as a clone of 'other'
280
282{
283 return (getAttribute(Form("CloneOf(%zx)",(size_t)&other)) ||
284 other.getAttribute(Form("CloneOf(%zx)",(size_t)this))) ;
285}
286
287
288////////////////////////////////////////////////////////////////////////////////
289/// Set (default) or clear a named boolean attribute of this object.
290
292{
293 // Preserve backward compatibility - any strong
294 if(string("Constant")==name) {
295 _isConstant = value ;
296 }
297
298 if (value) {
299 _boolAttrib.insert(name) ;
300 } else {
301 set<string>::iterator iter = _boolAttrib.find(name) ;
302 if (iter != _boolAttrib.end()) {
303 _boolAttrib.erase(iter) ;
304 }
305
306 }
307
308}
309
310
311////////////////////////////////////////////////////////////////////////////////
312/// Check if a named attribute is set. By default, all attributes are unset.
313
315{
316 return (_boolAttrib.find(name) != _boolAttrib.end()) ;
317}
318
319
320////////////////////////////////////////////////////////////////////////////////
321/// Associate string 'value' to this object under key 'key'
322
323void RooAbsArg::setStringAttribute(const Text_t* key, const Text_t* value)
324{
325 if (value) {
326 _stringAttrib[key] = value ;
327 } else {
328 _stringAttrib.erase(key) ;
329 }
330}
331
332////////////////////////////////////////////////////////////////////////////////
333/// Get string attribute mapped under key 'key'. Returns null pointer
334/// if no attribute exists under that key
335
337{
338 map<string,string>::const_iterator iter = _stringAttrib.find(key) ;
339 if (iter!=_stringAttrib.end()) {
340 return iter->second.c_str() ;
341 } else {
342 return 0 ;
343 }
344}
345
346
347////////////////////////////////////////////////////////////////////////////////
348/// Set (default) or clear a named boolean attribute of this object.
349
351{
352 if (value) {
353
354 _boolAttribTransient.insert(name) ;
355
356 } else {
357
358 set<string>::iterator iter = _boolAttribTransient.find(name) ;
359 if (iter != _boolAttribTransient.end()) {
360 _boolAttribTransient.erase(iter) ;
361 }
362
363 }
364
365}
366
367
368////////////////////////////////////////////////////////////////////////////////
369/// Check if a named attribute is set. By default, all attributes
370/// are unset.
371
373{
374 return (_boolAttribTransient.find(name) != _boolAttribTransient.end()) ;
375}
376
377
378
379
380////////////////////////////////////////////////////////////////////////////////
381/// Register another RooAbsArg as a server to us, ie, declare that
382/// we depend on it.
383/// \param server The server to be registered.
384/// \param valueProp In addition to the basic client-server relationship, declare dependence on the server's value.
385/// \param valueProp In addition to the basic client-server relationship, declare dependence on the server's shape.
386/// \param refCount Optionally add with higher reference count (if multiple components depend on it)
387
388void RooAbsArg::addServer(RooAbsArg& server, Bool_t valueProp, Bool_t shapeProp, std::size_t refCount)
389{
391 cxcoutF(LinkStateMgmt) << "RooAbsArg::addServer(" << this << "," << GetName()
392 << "): PROHIBITED SERVER ADDITION REQUESTED: adding server " << server.GetName()
393 << "(" << &server << ") for " << (valueProp?"value ":"") << (shapeProp?"shape":"") << endl ;
394 throw std::logic_error("PROHIBITED SERVER ADDITION REQUESTED in RooAbsArg::addServer");
395 }
396
397 cxcoutD(LinkStateMgmt) << "RooAbsArg::addServer(" << this << "," << GetName() << "): adding server " << server.GetName()
398 << "(" << &server << ") for " << (valueProp?"value ":"") << (shapeProp?"shape":"") << endl ;
399
400 if (server.operMode()==ADirty && operMode()!=ADirty && valueProp) {
402 }
403
404
405 // LM: use hash tables for larger lists
406// if (_serverList.GetSize() > 999 && _serverList.getHashTableSize() == 0) _serverList.setHashTableSize(1000);
407// if (server._clientList.GetSize() > 999 && server._clientList.getHashTableSize() == 0) server._clientList.setHashTableSize(1000);
408// if (server._clientListValue.GetSize() > 999 && server._clientListValue.getHashTableSize() == 0) server._clientListValue.setHashTableSize(1000);
409
410 // Add server link to given server
411 _serverList.Add(&server, refCount) ;
412
413 server._clientList.Add(this, refCount);
414 if (valueProp) server._clientListValue.Add(this, refCount);
415 if (shapeProp) server._clientListShape.Add(this, refCount);
416}
417
418
419
420////////////////////////////////////////////////////////////////////////////////
421/// Register a list of RooAbsArg as servers to us by calling
422/// addServer() for each arg in the list
423
424void RooAbsArg::addServerList(RooAbsCollection& serverList, Bool_t valueProp, Bool_t shapeProp)
425{
426 _serverList.reserve(_serverList.size() + serverList.size());
427
428 for (const auto arg : serverList) {
429 addServer(*arg,valueProp,shapeProp) ;
430 }
431}
432
433
434
435////////////////////////////////////////////////////////////////////////////////
436/// Unregister another RooAbsArg as a server to us, ie, declare that
437/// we no longer depend on its value and shape.
438
440{
442 cxcoutF(LinkStateMgmt) << "RooAbsArg::addServer(" << this << "," << GetName() << "): PROHIBITED SERVER REMOVAL REQUESTED: removing server "
443 << server.GetName() << "(" << &server << ")" << endl ;
444 assert(0) ;
445 }
446
447 if (_verboseDirty) {
448 cxcoutD(LinkStateMgmt) << "RooAbsArg::removeServer(" << GetName() << "): removing server "
449 << server.GetName() << "(" << &server << ")" << endl ;
450 }
451
452 // Remove server link to given server
453 _serverList.Remove(&server, force) ;
454
455 server._clientList.Remove(this, force) ;
456 server._clientListValue.Remove(this, force) ;
457 server._clientListShape.Remove(this, force) ;
458}
459
460
461////////////////////////////////////////////////////////////////////////////////
462/// Replace 'oldServer' with 'newServer'
463
464void RooAbsArg::replaceServer(RooAbsArg& oldServer, RooAbsArg& newServer, Bool_t propValue, Bool_t propShape)
465{
466 Int_t count = _serverList.refCount(&oldServer);
467 removeServer(oldServer, kTRUE);
468 addServer(newServer, propValue, propShape, count);
469}
470
471
472////////////////////////////////////////////////////////////////////////////////
473/// Change dirty flag propagation mask for specified server
474
475void RooAbsArg::changeServer(RooAbsArg& server, Bool_t valueProp, Bool_t shapeProp)
476{
477 if (!_serverList.containsByNamePtr(&server)) {
478 coutE(LinkStateMgmt) << "RooAbsArg::changeServer(" << GetName() << "): Server "
479 << server.GetName() << " not registered" << endl ;
480 return ;
481 }
482
483 // This condition should not happen, but check anyway
484 if (!server._clientList.containsByNamePtr(this)) {
485 coutE(LinkStateMgmt) << "RooAbsArg::changeServer(" << GetName() << "): Server "
486 << server.GetName() << " doesn't have us registered as client" << endl ;
487 return ;
488 }
489
490 // Remove all propagation links, then reinstall requested ones ;
491 Int_t vcount = server._clientListValue.refCount(this) ;
492 Int_t scount = server._clientListShape.refCount(this) ;
493 server._clientListValue.RemoveAll(this) ;
494 server._clientListShape.RemoveAll(this) ;
495 if (valueProp) {
496 server._clientListValue.Add(this, vcount) ;
497 }
498 if (shapeProp) {
499 server._clientListShape.Add(this, scount) ;
500 }
501}
502
503
504
505////////////////////////////////////////////////////////////////////////////////
506/// Fill supplied list with all leaf nodes of the arg tree, starting with
507/// ourself as top node. A leaf node is node that has no servers declared.
508
509void RooAbsArg::leafNodeServerList(RooAbsCollection* list, const RooAbsArg* arg, Bool_t recurseNonDerived) const
510{
511 treeNodeServerList(list,arg,kFALSE,kTRUE,kFALSE,recurseNonDerived) ;
512}
513
514
515
516////////////////////////////////////////////////////////////////////////////////
517/// Fill supplied list with all branch nodes of the arg tree starting with
518/// ourself as top node. A branch node is node that has one or more servers declared.
519
520void RooAbsArg::branchNodeServerList(RooAbsCollection* list, const RooAbsArg* arg, Bool_t recurseNonDerived) const
521{
522 treeNodeServerList(list,arg,kTRUE,kFALSE,kFALSE,recurseNonDerived) ;
523}
524
525
526////////////////////////////////////////////////////////////////////////////////
527/// Fill supplied list with nodes of the arg tree, following all server links,
528/// starting with ourself as top node.
529/// \param[in] list Output list
530/// \param[in] arg Start searching at this element of the tree.
531/// \param[in] doBranch Add branch nodes to the list.
532/// \param[in] doLeaf Add leaf nodes to the list.
533/// \param[in] valueOnly Only check if an element is a value server (no shape server).
534/// \param[in] recurseFundamental
535
536void RooAbsArg::treeNodeServerList(RooAbsCollection* list, const RooAbsArg* arg, Bool_t doBranch, Bool_t doLeaf, Bool_t valueOnly, Bool_t recurseFundamental) const
537{
538// if (arg==0) {
539// cout << "treeNodeServerList(" << GetName() << ") doBranch=" << (doBranch?"T":"F") << " doLeaf = " << (doLeaf?"T":"F") << " valueOnly=" << (valueOnly?"T":"F") << endl ;
540// }
541
542 if (!arg) {
543 list->reserve(10);
544 arg=this ;
545 }
546
547 // Decide if to add current node
548 if ((doBranch&&doLeaf) ||
549 (doBranch&&arg->isDerived()) ||
550 (doLeaf&&arg->isFundamental()&&(!(recurseFundamental&&arg->isDerived()))) ||
551 (doLeaf && !arg->isFundamental() && !arg->isDerived())) {
552
553 list->add(*arg,kTRUE) ;
554 }
555
556 // Recurse if current node is derived
557 if (arg->isDerived() && (!arg->isFundamental() || recurseFundamental)) {
558 for (const auto server : arg->_serverList) {
559
560 // Skip non-value server nodes if requested
561 Bool_t isValueSrv = server->_clientListValue.containsByNamePtr(arg);
562 if (valueOnly && !isValueSrv) {
563 continue ;
564 }
565 treeNodeServerList(list,server,doBranch,doLeaf,valueOnly,recurseFundamental) ;
566 }
567 }
568}
569
570
571////////////////////////////////////////////////////////////////////////////////
572/// Create a list of leaf nodes in the arg tree starting with
573/// ourself as top node that don't match any of the names of the variable list
574/// of the supplied data set (the dependents). The caller of this
575/// function is responsible for deleting the returned argset.
576/// The complement of this function is getObservables()
577
578RooArgSet* RooAbsArg::getParameters(const RooAbsData* set, Bool_t stripDisconnected) const
579{
580 return getParameters(set?set->get():0,stripDisconnected) ;
581}
582
583
584////////////////////////////////////////////////////////////////////////////////
585/// INTERNAL helper function for getParameters()
586
587void RooAbsArg::addParameters(RooArgSet& params, const RooArgSet* nset,Bool_t stripDisconnected) const
588{
589 // INTERNAL helper function for getParameters()
590
591 RooArgSet nodeParamServers ;
592 RooArgSet nodeBranchServers ;
593 for (const auto server : _serverList) {
594 if (server->isValueServer(*this)) {
595 if (server->isFundamental()) {
596 if (!nset || !server->dependsOn(*nset)) {
597 nodeParamServers.add(*server) ;
598 }
599 } else {
600 nodeBranchServers.add(*server) ;
601 }
602 }
603 }
604
605 // Allow pdf to strip parameters from list before adding it
606 getParametersHook(nset,&nodeParamServers,stripDisconnected) ;
607
608 // Add parameters of this node to the combined list
609 params.add(nodeParamServers,kTRUE) ;
610
611 // Now recurse into branch servers
612 for (const auto server : nodeBranchServers) {
613 server->addParameters(params,nset) ;
614 }
615}
616
617
618////////////////////////////////////////////////////////////////////////////////
619/// Create a list of leaf nodes in the arg tree starting with
620/// ourself as top node that don't match any of the names the args in the
621/// supplied argset. The caller of this function is responsible
622/// for deleting the returned argset. The complement of this function
623/// is getObservables().
624
625RooArgSet* RooAbsArg::getParameters(const RooArgSet* observables, bool stripDisconnected) const {
626 auto * outputSet = new RooArgSet;
627 getParameters(observables, *outputSet, stripDisconnected);
628 return outputSet;
629}
630
631
632////////////////////////////////////////////////////////////////////////////////
633/// Fills a list with leaf nodes in the arg tree starting with
634/// ourself as top node that don't match any of the names the args in the
635/// supplied argset. Returns `true` only if something went wrong.
636/// The complement of this function is getObservables().
637/// \param[in] observables Set of leafs to ignore because they are observables and not parameters.
638/// \param[out] outputSet Output set.
639/// \param[in] stripDisconnected Allow pdf to strip parameters from list before adding it.
640
641bool RooAbsArg::getParameters(const RooArgSet* observables, RooArgSet& outputSet, bool stripDisconnected) const
642{
644
645 // Check for cached parameter set
646 if (_myws) {
647 auto nsetObs = getColonSeparatedNameString(observables ? *observables : RooArgSet());
648 const RooArgSet *paramSet = _myws->set(Form("CACHE_PARAMS_OF_PDF_%s_FOR_OBS_%s", GetName(), nsetObs.c_str()));
649 if (paramSet) {
650 outputSet.add(*paramSet);
651 return false;
652 }
653 }
654
655 outputSet.clear();
656 outputSet.setName("parameters");
657
658 addParameters(outputSet, observables, stripDisconnected);
659
660 outputSet.sort();
661
662 // Cache parameter set
663 if (_myws && outputSet.size() > 10) {
664 auto nsetObs = getColonSeparatedNameString(observables ? *observables : RooArgSet());
665 _myws->defineSetInternal(Form("CACHE_PARAMS_OF_PDF_%s_FOR_OBS_%s", GetName(), nsetObs.c_str()), outputSet);
666 // cout << " caching parameters in workspace for pdf " << IsA()->GetName() << "::" << GetName() << endl ;
667 }
668
669 return false;
670}
671
672
673////////////////////////////////////////////////////////////////////////////////
674/// Create a list of leaf nodes in the arg tree starting with
675/// ourself as top node that match any of the names of the variable list
676/// of the supplied data set (the dependents). The caller of this
677/// function is responsible for deleting the returned argset.
678/// The complement of this function is getParameters().
679
681{
682 if (!set) return new RooArgSet ;
683
684 return getObservables(set->get()) ;
685}
686
687
688////////////////////////////////////////////////////////////////////////////////
689/// Create a list of leaf nodes in the arg tree starting with
690/// ourself as top node that match any of the names the args in the
691/// supplied argset. The caller of this function is responsible
692/// for deleting the returned argset. The complement of this function
693/// is getParameters().
694
695RooArgSet* RooAbsArg::getObservables(const RooArgSet* dataList, bool valueOnly) const
696{
697 auto depList = new RooArgSet;
698 getObservables(dataList, *depList, valueOnly);
699 return depList;
700}
701
702
703////////////////////////////////////////////////////////////////////////////////
704/// Create a list of leaf nodes in the arg tree starting with
705/// ourself as top node that match any of the names the args in the
706/// supplied argset.
707/// Returns `true` only if something went wrong.
708/// The complement of this function is getParameters().
709/// \param[in] dataList Set of leaf nodes to match.
710/// \param[out] outputSet Output set.
711/// \param[in] valueOnly If this parameter is true, we only match leafs that
712/// depend on the value of any arg in `dataList`.
713
714bool RooAbsArg::getObservables(const RooAbsCollection* dataList, RooArgSet& outputSet, bool valueOnly) const
715{
716 outputSet.clear();
717 outputSet.setName("dependents");
718
719 if (!dataList) return false;
720
721 // Make iterator over tree leaf node list
722 RooArgSet leafList("leafNodeServerList") ;
723 treeNodeServerList(&leafList,0,kFALSE,kTRUE,valueOnly) ;
724
725 if (valueOnly) {
726 for (const auto arg : leafList) {
727 if (arg->dependsOnValue(*dataList) && arg->isLValue()) {
728 outputSet.add(*arg) ;
729 }
730 }
731 } else {
732 for (const auto arg : leafList) {
733 if (arg->dependsOn(*dataList) && arg->isLValue()) {
734 outputSet.add(*arg) ;
735 }
736 }
737 }
738
739 return false;
740}
741
742
743////////////////////////////////////////////////////////////////////////////////
744/// Create a RooArgSet with all components (branch nodes) of the
745/// expression tree headed by this object.
747{
748 TString name(GetName()) ;
749 name.Append("_components") ;
750
751 RooArgSet* set = new RooArgSet(name) ;
753
754 return set ;
755}
756
757
758
759////////////////////////////////////////////////////////////////////////////////
760/// Overloadable function in which derived classes can implement
761/// consistency checks of the variables. If this function returns
762/// true, indicating an error, the fitter or generator will abort.
763
765{
766 return kFALSE ;
767}
768
769
770////////////////////////////////////////////////////////////////////////////////
771/// Recursively call checkObservables on all nodes in the expression tree
772
774{
775 RooArgSet nodeList ;
776 treeNodeServerList(&nodeList) ;
777 RooFIter iter = nodeList.fwdIterator() ;
778
779 RooAbsArg* arg ;
780 Bool_t ret(kFALSE) ;
781 while((arg=iter.next())) {
782 if (arg->getAttribute("ServerDied")) {
783 coutE(LinkStateMgmt) << "RooAbsArg::recursiveCheckObservables(" << GetName() << "): ERROR: one or more servers of node "
784 << arg->GetName() << " no longer exists!" << endl ;
785 arg->Print("v") ;
786 ret = kTRUE ;
787 }
788 ret |= arg->checkObservables(nset) ;
789 }
790
791 return ret ;
792}
793
794
795////////////////////////////////////////////////////////////////////////////////
796/// Test whether we depend on (ie, are served by) any object in the
797/// specified collection. Uses the dependsOn(RooAbsArg&) member function.
798
799Bool_t RooAbsArg::dependsOn(const RooAbsCollection& serverList, const RooAbsArg* ignoreArg, Bool_t valueOnly) const
800{
801 // Test whether we depend on (ie, are served by) any object in the
802 // specified collection. Uses the dependsOn(RooAbsArg&) member function.
803
804 for (auto server : serverList) {
805 if (dependsOn(*server,ignoreArg,valueOnly)) {
806 return kTRUE;
807 }
808 }
809 return kFALSE;
810}
811
812
813////////////////////////////////////////////////////////////////////////////////
814/// Test whether we depend on (ie, are served by) the specified object.
815/// Note that RooAbsArg objects are considered equivalent if they have
816/// the same name.
817
818Bool_t RooAbsArg::dependsOn(const RooAbsArg& testArg, const RooAbsArg* ignoreArg, Bool_t valueOnly) const
819{
820 if (this==ignoreArg) return kFALSE ;
821
822 // First check if testArg is self
823 //if (!TString(testArg.GetName()).CompareTo(GetName())) return kTRUE ;
824 if (testArg.namePtr()==namePtr()) return kTRUE ;
825
826
827 // Next test direct dependence
828 RooAbsArg* foundServer = findServer(testArg) ;
829 if (foundServer) {
830
831 // Return true if valueOnly is FALSE or if server is value server, otherwise keep looking
832 if ( !valueOnly || foundServer->isValueServer(*this)) {
833 return kTRUE ;
834 }
835 }
836
837 // If not, recurse
838 for (const auto server : _serverList) {
839 if ( !valueOnly || server->isValueServer(*this)) {
840 if (server->dependsOn(testArg,ignoreArg,valueOnly)) {
841 return kTRUE ;
842 }
843 }
844 }
845
846 return kFALSE ;
847}
848
849
850
851////////////////////////////////////////////////////////////////////////////////
852/// Test if any of the nodes of tree are shared with that of the given tree
853
854Bool_t RooAbsArg::overlaps(const RooAbsArg& testArg, Bool_t valueOnly) const
855{
856 RooArgSet list("treeNodeList") ;
857 treeNodeServerList(&list) ;
858
859 return valueOnly ? testArg.dependsOnValue(list) : testArg.dependsOn(list) ;
860}
861
862
863
864////////////////////////////////////////////////////////////////////////////////
865/// Test if any of the dependents of the arg tree (as determined by getObservables)
866/// overlaps with those of the testArg.
867
869{
870 return observableOverlaps(dset->get(),testArg) ;
871}
872
873
874////////////////////////////////////////////////////////////////////////////////
875/// Test if any of the dependents of the arg tree (as determined by getObservables)
876/// overlaps with those of the testArg.
877
878Bool_t RooAbsArg::observableOverlaps(const RooArgSet* nset, const RooAbsArg& testArg) const
879{
880 RooArgSet* depList = getObservables(nset) ;
881 Bool_t ret = testArg.dependsOn(*depList) ;
882 delete depList ;
883 return ret ;
884}
885
886
887
888////////////////////////////////////////////////////////////////////////////////
889/// Mark this object as having changed its value, and propagate this status
890/// change to all of our clients. If the object is not in automatic dirty
891/// state propagation mode, this call has no effect.
892
894{
895 _allBatchesDirty = true;
896
897 if (_operMode!=Auto || _inhibitDirty) return ;
898
899 // Handle no-propagation scenarios first
900 if (_clientListValue.size() == 0) {
902 return ;
903 }
904
905 // Cyclical dependency interception
906 if (source==0) {
907 source=this ;
908 } else if (source==this) {
909 // Cyclical dependency, abort
910 coutE(LinkStateMgmt) << "RooAbsArg::setValueDirty(" << GetName()
911 << "): cyclical dependency detected, source = " << source->GetName() << endl ;
912 //assert(0) ;
913 return ;
914 }
915
916 // Propagate dirty flag to all clients if this is a down->up transition
917 if (_verboseDirty) {
918 cxcoutD(LinkStateMgmt) << "RooAbsArg::setValueDirty(" << (source?source->GetName():"self") << "->" << GetName() << "," << this
919 << "): dirty flag " << (_valueDirty?"already ":"") << "raised" << endl ;
920 }
921
923
924
925 for (auto client : _clientListValue) {
926 client->setValueDirty(source) ;
927 }
928
929
930}
931
932
933////////////////////////////////////////////////////////////////////////////////
934/// Mark this object as having changed its shape, and propagate this status
935/// change to all of our clients.
936
938{
939 if (_verboseDirty) {
940 cxcoutD(LinkStateMgmt) << "RooAbsArg::setShapeDirty(" << GetName()
941 << "): dirty flag " << (_shapeDirty?"already ":"") << "raised" << endl ;
942 }
943
944 if (_clientListShape.empty()) {
946 return ;
947 }
948
949 // Set 'dirty' shape state for this object and propagate flag to all its clients
950 if (source==0) {
951 source=this ;
952 } else if (source==this) {
953 // Cyclical dependency, abort
954 coutE(LinkStateMgmt) << "RooAbsArg::setShapeDirty(" << GetName()
955 << "): cyclical dependency detected" << endl ;
956 return ;
957 }
958
959 // Propagate dirty flag to all clients if this is a down->up transition
961
962 for (auto client : _clientListShape) {
963 client->setShapeDirty(source) ;
964 client->setValueDirty(source) ;
965 }
966
967}
968
969
970
971////////////////////////////////////////////////////////////////////////////////
972/// Replace all direct servers of this object with the new servers in `newServerList`.
973/// This substitutes objects that we receive values from with new objects that have the same name.
974/// \see recursiveRedirectServers() Use recursive version if servers that are only indirectly serving this object should be replaced as well.
975/// \see redirectServers() If only the direct servers of an object need to be replaced.
976///
977/// Note that changing the types of objects is generally allowed, but can be wrong if the interface of an object changes.
978/// For example, one can reparametrise a model by substituting a variable with a function:
979/// \f[
980/// f(x\, |\, a) = a \cdot x \rightarrow f(x\, |\, b) = (2.1 \cdot b) \cdot x
981/// \f]
982/// If an object, however, expects a PDF, and this is substituted with a function that isn't normalised, wrong results might be obtained
983/// or it might even crash the program. The types of the objects being substituted are not checked.
984///
985/// \param[in] newSetOrig Set of new servers that should be used instead of the current servers.
986/// \param[in] mustReplaceAll A warning is printed and error status is returned if not all servers could be
987/// substituted successfully.
988/// \param[in] nameChange If false, an object named "x" is only replaced with an object also named "x" in `newSetOrig`.
989/// If the object in `newSet` is called differently, set `nameChange` to true and use setAttribute() on the x object:
990/// ```
991/// objectToReplaceX.setAttribute("ORIGNAME:x")
992/// ```
993/// Now, the renamed object will be selected based on the attribute "ORIGNAME:<name>".
994/// \param[in] isRecursionStep Internal switch used when called from recursiveRedirectServers().
995Bool_t RooAbsArg::redirectServers(const RooAbsCollection& newSetOrig, Bool_t mustReplaceAll, Bool_t nameChange, Bool_t isRecursionStep)
996{
997 // Trivial case, no servers
998 if (_serverList.empty()) return kFALSE ;
999 if (newSetOrig.getSize()==0) return kFALSE ;
1000
1001 // Strip any non-matching removal nodes from newSetOrig
1002 RooAbsCollection* newSet ;
1003
1004 if (nameChange) {
1005 newSet = new RooArgSet ;
1006 for (auto arg : newSetOrig) {
1007
1008 if (string("REMOVAL_DUMMY")==arg->GetName()) {
1009
1010 if (arg->getAttribute("REMOVE_ALL")) {
1011 newSet->add(*arg) ;
1012 } else if (arg->getAttribute(Form("REMOVE_FROM_%s",getStringAttribute("ORIGNAME")))) {
1013 newSet->add(*arg) ;
1014 }
1015 } else {
1016 newSet->add(*arg) ;
1017 }
1018 }
1019 } else {
1020 newSet = (RooAbsCollection*) &newSetOrig ;
1021 }
1022
1023 // Replace current servers with new servers with the same name from the given list
1024 Bool_t ret(kFALSE) ;
1025
1026 //Copy original server list to not confuse the iterator while deleting
1027 std::vector<RooAbsArg*> origServerList, origServerValue, origServerShape;
1028 auto origSize = _serverList.size();
1029 origServerList.reserve(origSize);
1030 origServerValue.reserve(origSize);
1031
1032 for (const auto oldServer : _serverList) {
1033 origServerList.push_back(oldServer) ;
1034
1035 // Retrieve server side link state information
1036 if (oldServer->_clientListValue.containsByNamePtr(this)) {
1037 origServerValue.push_back(oldServer) ;
1038 }
1039 if (oldServer->_clientListShape.containsByNamePtr(this)) {
1040 origServerShape.push_back(oldServer) ;
1041 }
1042 }
1043
1044 // Delete all previously registered servers
1045 for (auto oldServer : origServerList) {
1046
1047 RooAbsArg * newServer= oldServer->findNewServer(*newSet, nameChange);
1048
1049 if (newServer && _verboseDirty) {
1050 cxcoutD(LinkStateMgmt) << "RooAbsArg::redirectServers(" << (void*)this << "," << GetName() << "): server " << oldServer->GetName()
1051 << " redirected from " << oldServer << " to " << newServer << endl ;
1052 }
1053
1054 if (!newServer) {
1055 if (mustReplaceAll) {
1056 coutE(LinkStateMgmt) << "RooAbsArg::redirectServers(" << (void*)this << "," << GetName() << "): server " << oldServer->GetName()
1057 << " (" << (void*)oldServer << ") not redirected" << (nameChange?"[nameChange]":"") << endl ;
1058 ret = kTRUE ;
1059 }
1060 continue ;
1061 }
1062
1063 auto findByNamePtr = [&oldServer](const RooAbsArg * item) {
1064 return oldServer->namePtr() == item->namePtr();
1065 };
1066 bool propValue = std::any_of(origServerValue.begin(), origServerValue.end(), findByNamePtr);
1067 bool propShape = std::any_of(origServerShape.begin(), origServerShape.end(), findByNamePtr);
1068
1069 if (newServer != this) {
1070 replaceServer(*oldServer,*newServer,propValue,propShape) ;
1071 }
1072 }
1073
1074
1075 setValueDirty() ;
1076 setShapeDirty() ;
1077
1078 // Process the proxies
1079 for (int i=0 ; i<numProxies() ; i++) {
1080 RooAbsProxy* p = getProxy(i) ;
1081 if (!p) continue ;
1082 Bool_t ret2 = p->changePointer(*newSet,nameChange,kFALSE) ;
1083
1084 if (mustReplaceAll && !ret2) {
1085 auto ap = dynamic_cast<const RooArgProxy*>(p);
1086 coutE(LinkStateMgmt) << "RooAbsArg::redirectServers(" << GetName()
1087 << "): ERROR, proxy '" << p->name()
1088 << "' with arg '" << (ap ? ap->absArg()->GetName() : "<could not cast>") << "' could not be adjusted" << endl;
1089 ret = kTRUE ;
1090 }
1091 }
1092
1093
1094 // Optional subclass post-processing
1095 for (Int_t i=0 ;i<numCaches() ; i++) {
1096 ret |= getCache(i)->redirectServersHook(*newSet,mustReplaceAll,nameChange,isRecursionStep) ;
1097 }
1098 ret |= redirectServersHook(*newSet,mustReplaceAll,nameChange,isRecursionStep) ;
1099
1100 if (nameChange) {
1101 delete newSet ;
1102 }
1103
1104 return ret ;
1105}
1106
1107////////////////////////////////////////////////////////////////////////////////
1108/// Find the new server in the specified set that matches the old server.
1109///
1110/// \param[in] newSet Search this set by name for a new server.
1111/// \param[in] If true, search for an item with the bool attribute "ORIGNAME:<oldName>" set.
1112/// Use `<object>.setAttribute("ORIGNAME:<oldName>")` to set this attribute.
1113/// \return Pointer to the new server or `nullptr` if there's no unique match.
1115{
1116 RooAbsArg *newServer = 0;
1117 if (!nameChange) {
1118 newServer = newSet.find(*this) ;
1119 }
1120 else {
1121 // Name changing server redirect:
1122 // use 'ORIGNAME:<oldName>' attribute instead of name of new server
1123 TString nameAttrib("ORIGNAME:") ;
1124 nameAttrib.Append(GetName()) ;
1125
1126 RooArgSet* tmp = (RooArgSet*) newSet.selectByAttrib(nameAttrib,kTRUE) ;
1127 if(0 != tmp) {
1128
1129 // Check if any match was found
1130 if (tmp->getSize()==0) {
1131 delete tmp ;
1132 return 0 ;
1133 }
1134
1135 // Check if match is unique
1136 if(tmp->getSize()>1) {
1137 coutF(LinkStateMgmt) << "RooAbsArg::redirectServers(" << GetName() << "): FATAL Error, " << tmp->getSize() << " servers with "
1138 << nameAttrib << " attribute" << endl ;
1139 tmp->Print("v") ;
1140 assert(0) ;
1141 }
1142
1143 // use the unique element in the set
1144 newServer= tmp->first();
1145 delete tmp ;
1146 }
1147 }
1148 return newServer;
1149}
1150
1151
1152////////////////////////////////////////////////////////////////////////////////
1153/// Recursively replace all servers with the new servers in `newSet`.
1154/// This substitutes objects that we receive values from (also indirectly through other objects) with new objects that have the same name.
1155///
1156/// *Copied from redirectServers:*
1157///
1158/// \copydetails RooAbsArg::redirectServers
1159Bool_t RooAbsArg::recursiveRedirectServers(const RooAbsCollection& newSet, Bool_t mustReplaceAll, Bool_t nameChange, Bool_t recurseInNewSet)
1160{
1161 // Cyclic recursion protection
1162 static std::set<const RooAbsArg*> callStack;
1163 {
1164 std::set<const RooAbsArg*>::iterator it = callStack.lower_bound(this);
1165 if (it != callStack.end() && this == *it) {
1166 return kFALSE;
1167 } else {
1168 callStack.insert(it, this);
1169 }
1170 }
1171
1172 // Do not recurse into newset if not so specified
1173// if (!recurseInNewSet && newSet.contains(*this)) {
1174// return kFALSE ;
1175// }
1176
1177
1178 // Apply the redirectServers function recursively on all branch nodes in this argument tree.
1179 Bool_t ret(kFALSE) ;
1180
1181 cxcoutD(LinkStateMgmt) << "RooAbsArg::recursiveRedirectServers(" << this << "," << GetName() << ") newSet = " << newSet << " mustReplaceAll = "
1182 << (mustReplaceAll?"T":"F") << " nameChange = " << (nameChange?"T":"F") << " recurseInNewSet = " << (recurseInNewSet?"T":"F") << endl ;
1183
1184 // Do redirect on self (identify operation as recursion step)
1185 ret |= redirectServers(newSet,mustReplaceAll,nameChange,kTRUE) ;
1186
1187 // Do redirect on servers
1188 for (const auto server : _serverList){
1189 ret |= server->recursiveRedirectServers(newSet,mustReplaceAll,nameChange,recurseInNewSet) ;
1190 }
1191
1192 callStack.erase(this);
1193 return ret ;
1194}
1195
1196
1197
1198////////////////////////////////////////////////////////////////////////////////
1199/// Register an RooArgProxy in the proxy list. This function is called by owned
1200/// proxies upon creation. After registration, this arg wil forward pointer
1201/// changes from serverRedirects and updates in cached normalization sets
1202/// to the proxies immediately after they occur. The proxied argument is
1203/// also added as value and/or shape server
1204
1206{
1207 // Every proxy can be registered only once
1208 if (_proxyList.FindObject(&proxy)) {
1209 coutE(LinkStateMgmt) << "RooAbsArg::registerProxy(" << GetName() << "): proxy named "
1210 << proxy.GetName() << " for arg " << proxy.absArg()->GetName()
1211 << " already registered" << endl ;
1212 return ;
1213 }
1214
1215// cout << (void*)this << " " << GetName() << ": registering proxy "
1216// << (void*)&proxy << " with name " << proxy.name() << " in mode "
1217// << (proxy.isValueServer()?"V":"-") << (proxy.isShapeServer()?"S":"-") << endl ;
1218
1219 // Register proxied object as server
1220 if (proxy.absArg()) {
1221 addServer(*proxy.absArg(),proxy.isValueServer(),proxy.isShapeServer()) ;
1222 }
1223
1224 // Register proxy itself
1225 _proxyList.Add(&proxy) ;
1226}
1227
1228
1229////////////////////////////////////////////////////////////////////////////////
1230/// Remove proxy from proxy list. This functions is called by owned proxies
1231/// upon their destruction.
1232
1234{
1235 _proxyList.Remove(&proxy) ;
1237}
1238
1239
1240
1241////////////////////////////////////////////////////////////////////////////////
1242/// Register an RooSetProxy in the proxy list. This function is called by owned
1243/// proxies upon creation. After registration, this arg wil forward pointer
1244/// changes from serverRedirects and updates in cached normalization sets
1245/// to the proxies immediately after they occur.
1246
1248{
1249 // Every proxy can be registered only once
1250 if (_proxyList.FindObject(&proxy)) {
1251 coutE(LinkStateMgmt) << "RooAbsArg::registerProxy(" << GetName() << "): proxy named "
1252 << proxy.GetName() << " already registered" << endl ;
1253 return ;
1254 }
1255
1256 // Register proxy itself
1257 _proxyList.Add(&proxy) ;
1258}
1259
1260
1261
1262////////////////////////////////////////////////////////////////////////////////
1263/// Remove proxy from proxy list. This functions is called by owned proxies
1264/// upon their destruction.
1265
1267{
1268 _proxyList.Remove(&proxy) ;
1270}
1271
1272
1273
1274////////////////////////////////////////////////////////////////////////////////
1275/// Register an RooListProxy in the proxy list. This function is called by owned
1276/// proxies upon creation. After registration, this arg wil forward pointer
1277/// changes from serverRedirects and updates in cached normalization sets
1278/// to the proxies immediately after they occur.
1279
1281{
1282 // Every proxy can be registered only once
1283 if (_proxyList.FindObject(&proxy)) {
1284 coutE(LinkStateMgmt) << "RooAbsArg::registerProxy(" << GetName() << "): proxy named "
1285 << proxy.GetName() << " already registered" << endl ;
1286 return ;
1287 }
1288
1289 // Register proxy itself
1290 Int_t nProxyOld = _proxyList.GetEntries() ;
1291 _proxyList.Add(&proxy) ;
1292 if (_proxyList.GetEntries()!=nProxyOld+1) {
1293 cout << "RooAbsArg::registerProxy(" << GetName() << ") proxy registration failure! nold=" << nProxyOld << " nnew=" << _proxyList.GetEntries() << endl ;
1294 }
1295}
1296
1297
1298
1299////////////////////////////////////////////////////////////////////////////////
1300/// Remove proxy from proxy list. This functions is called by owned proxies
1301/// upon their destruction.
1302
1304{
1305 _proxyList.Remove(&proxy) ;
1307}
1308
1309
1310
1311////////////////////////////////////////////////////////////////////////////////
1312/// Return the nth proxy from the proxy list.
1313
1315{
1316 // Cross cast: proxy list returns TObject base pointer, we need
1317 // a RooAbsProxy base pointer. C++ standard requires
1318 // a dynamic_cast for this.
1319 return dynamic_cast<RooAbsProxy*> (_proxyList.At(index)) ;
1320}
1321
1322
1323
1324////////////////////////////////////////////////////////////////////////////////
1325/// Return the number of registered proxies.
1326
1328{
1329 return _proxyList.GetEntriesFast();
1330}
1331
1332
1333
1334////////////////////////////////////////////////////////////////////////////////
1335/// Forward a change in the cached normalization argset
1336/// to all the registered proxies.
1337
1339{
1340 for (int i=0 ; i<numProxies() ; i++) {
1341 RooAbsProxy* p = getProxy(i) ;
1342 if (!p) continue ;
1343 getProxy(i)->changeNormSet(nset) ;
1344 }
1345}
1346
1347
1348
1349////////////////////////////////////////////////////////////////////////////////
1350/// Overloadable function for derived classes to implement
1351/// attachment as branch to a TTree
1352
1354{
1355 coutE(Contents) << "RooAbsArg::attachToTree(" << GetName()
1356 << "): Cannot be attached to a TTree" << endl ;
1357}
1358
1359
1360
1361////////////////////////////////////////////////////////////////////////////////
1362/// WVE (08/21/01) Probably obsolete now
1363
1365{
1366 return kTRUE ;
1367}
1368
1369
1370
1371
1372////////////////////////////////////////////////////////////////////////////////
1373/// Print object name
1374
1375void RooAbsArg::printName(ostream& os) const
1376{
1377 os << GetName() ;
1378}
1379
1380
1381
1382////////////////////////////////////////////////////////////////////////////////
1383/// Print object title
1384
1385void RooAbsArg::printTitle(ostream& os) const
1386{
1387 os << GetTitle() ;
1388}
1389
1390
1391
1392////////////////////////////////////////////////////////////////////////////////
1393/// Print object class name
1394
1395void RooAbsArg::printClassName(ostream& os) const
1396{
1397 os << IsA()->GetName() ;
1398}
1399
1400
1401void RooAbsArg::printAddress(ostream& os) const
1402{
1403 // Print addrss of this RooAbsArg
1404 os << this ;
1405}
1406
1407
1408
1409////////////////////////////////////////////////////////////////////////////////
1410/// Print object arguments, ie its proxies
1411
1412void RooAbsArg::printArgs(ostream& os) const
1413{
1414 // Print nothing if there are no dependencies
1415 if (numProxies()==0) return ;
1416
1417 os << "[ " ;
1418 for (Int_t i=0 ; i<numProxies() ; i++) {
1419 RooAbsProxy* p = getProxy(i) ;
1420 if (p==0) continue ;
1421 if (!TString(p->name()).BeginsWith("!")) {
1422 p->print(os) ;
1423 os << " " ;
1424 }
1425 }
1426 printMetaArgs(os) ;
1427 os << "]" ;
1428}
1429
1430
1431
1432////////////////////////////////////////////////////////////////////////////////
1433/// Define default contents to print
1434
1436{
1437 return kName|kClassName|kValue|kArgs ;
1438}
1439
1440
1441
1442////////////////////////////////////////////////////////////////////////////////
1443/// Implement multi-line detailed printing
1444
1445void RooAbsArg::printMultiline(ostream& os, Int_t /*contents*/, Bool_t /*verbose*/, TString indent) const
1446{
1447 os << indent << "--- RooAbsArg ---" << endl;
1448 // dirty state flags
1449 os << indent << " Value State: " ;
1450 switch(_operMode) {
1451 case ADirty: os << "FORCED DIRTY" ; break ;
1452 case AClean: os << "FORCED clean" ; break ;
1453 case Auto: os << (isValueDirty() ? "DIRTY":"clean") ; break ;
1454 }
1455 os << endl
1456 << indent << " Shape State: " << (isShapeDirty() ? "DIRTY":"clean") << endl;
1457 // attribute list
1458 os << indent << " Attributes: " ;
1459 printAttribList(os) ;
1460 os << endl ;
1461 // our memory address (for x-referencing with client addresses of other args)
1462 os << indent << " Address: " << (void*)this << endl;
1463 // client list
1464 os << indent << " Clients: " << endl;
1465 for (const auto client : _clientList) {
1466 os << indent << " (" << (void*)client << ","
1467 << (_clientListValue.containsByNamePtr(client)?"V":"-")
1468 << (_clientListShape.containsByNamePtr(client)?"S":"-")
1469 << ") " ;
1470 client->printStream(os,kClassName|kTitle|kName,kSingleLine);
1471 }
1472
1473 // server list
1474 os << indent << " Servers: " << endl;
1475 for (const auto server : _serverList) {
1476 os << indent << " (" << (void*)server << ","
1477 << (server->_clientListValue.containsByNamePtr(this)?"V":"-")
1478 << (server->_clientListShape.containsByNamePtr(this)?"S":"-")
1479 << ") " ;
1480 server->printStream(os,kClassName|kName|kTitle,kSingleLine);
1481 }
1482
1483 // proxy list
1484 os << indent << " Proxies: " << endl ;
1485 for (int i=0 ; i<numProxies() ; i++) {
1486 RooAbsProxy* proxy=getProxy(i) ;
1487 if (!proxy) continue ;
1488 if (proxy->IsA()->InheritsFrom(RooArgProxy::Class())) {
1489 os << indent << " " << proxy->name() << " -> " ;
1490 RooAbsArg* parg = ((RooArgProxy*)proxy)->absArg() ;
1491 if (parg) {
1492 parg->printStream(os,kName,kSingleLine) ;
1493 } else {
1494 os << " (empty)" << endl ; ;
1495 }
1496 } else {
1497 os << indent << " " << proxy->name() << " -> " ;
1498 os << endl ;
1499 TString moreIndent(indent) ;
1500 moreIndent.Append(" ") ;
1501 ((RooSetProxy*)proxy)->printStream(os,kName,kStandard,moreIndent.Data()) ;
1502 }
1503 }
1504}
1505
1506
1507////////////////////////////////////////////////////////////////////////////////
1508/// Print object tree structure
1509
1510void RooAbsArg::printTree(ostream& os, TString /*indent*/) const
1511{
1512 const_cast<RooAbsArg*>(this)->printCompactTree(os) ;
1513}
1514
1515
1516////////////////////////////////////////////////////////////////////////////////
1517/// Ostream operator
1518
1519ostream& operator<<(ostream& os, RooAbsArg const& arg)
1520{
1521 arg.writeToStream(os,kTRUE) ;
1522 return os ;
1523}
1524
1525////////////////////////////////////////////////////////////////////////////////
1526/// Istream operator
1527
1528istream& operator>>(istream& is, RooAbsArg &arg)
1529{
1530 arg.readFromStream(is,kTRUE,kFALSE) ;
1531 return is ;
1532}
1533
1534////////////////////////////////////////////////////////////////////////////////
1535/// Print the attribute list
1536
1537void RooAbsArg::printAttribList(ostream& os) const
1538{
1539 set<string>::const_iterator iter = _boolAttrib.begin() ;
1540 Bool_t first(kTRUE) ;
1541 while (iter != _boolAttrib.end()) {
1542 os << (first?" [":",") << *iter ;
1543 first=kFALSE ;
1544 ++iter ;
1545 }
1546 if (!first) os << "] " ;
1547}
1548
1549
1550////////////////////////////////////////////////////////////////////////////////
1551/// Bind this node to objects in `set`.
1552/// Search the set for objects that have the same name as our servers, and
1553/// attach ourselves to those. After this operation, this node is computing its
1554/// values based on the new servers. This can be used to e.g. read values from
1555// a dataset.
1556
1557
1559{
1560 RooArgSet branches;
1561 branchNodeServerList(&branches,0,true);
1562
1563 for(auto const& branch : branches) {
1564 branch->redirectServers(set,false,false);
1565 }
1566}
1567
1568
1569////////////////////////////////////////////////////////////////////////////////
1570/// Replace server nodes with names matching the dataset variable names
1571/// with those data set variables, making this PDF directly dependent on the dataset.
1572
1574{
1575 attachArgs(*data.get());
1576}
1577
1578
1579////////////////////////////////////////////////////////////////////////////////
1580/// Replace server nodes with names matching the dataset variable names
1581/// with those data set variables, making this PDF directly dependent on the dataset
1582
1584{
1585 attachArgs(*dstore.get());
1586}
1587
1588
1589////////////////////////////////////////////////////////////////////////////////
1590/// Utility function used by TCollection::Sort to compare contained TObjects
1591/// We implement comparison by name, resulting in alphabetical sorting by object name.
1592
1594{
1595 return strcmp(GetName(),other->GetName()) ;
1596}
1597
1598
1599
1600////////////////////////////////////////////////////////////////////////////////
1601/// Print information about current value dirty state information.
1602/// If depth flag is true, information is recursively printed for
1603/// all nodes in this arg tree.
1604
1606{
1607 if (depth) {
1608
1609 RooArgSet branchList ;
1610 branchNodeServerList(&branchList) ;
1611 RooFIter bIter = branchList.fwdIterator() ;
1612 RooAbsArg* branch ;
1613 while((branch=bIter.next())) {
1614 branch->printDirty(kFALSE) ;
1615 }
1616
1617 } else {
1618 cout << GetName() << " : " ;
1619 switch (_operMode) {
1620 case AClean: cout << "FORCED clean" ; break ;
1621 case ADirty: cout << "FORCED DIRTY" ; break ;
1622 case Auto: cout << "Auto " << (isValueDirty()?"DIRTY":"clean") ;
1623 }
1624 cout << endl ;
1625 }
1626}
1627
1628
1629////////////////////////////////////////////////////////////////////////////////
1630/// Activate cache mode optimization with given definition of observables.
1631/// The cache operation mode of all objects in the expression tree will
1632/// modified such that all nodes that depend directly or indirectly on
1633/// any of the listed observables will be set to ADirty, as they are
1634/// expected to change every time. This save change tracking overhead for
1635/// nodes that are a priori known to change every time
1636
1638{
1639 RooLinkedList proc;
1640 RooArgSet opt ;
1641 optimizeCacheMode(observables,opt,proc) ;
1642
1643 coutI(Optimization) << "RooAbsArg::optimizeCacheMode(" << GetName() << ") nodes " << opt << " depend on observables, "
1644 << "changing cache operation mode from change tracking to unconditional evaluation" << endl ;
1645}
1646
1647
1648////////////////////////////////////////////////////////////////////////////////
1649/// Activate cache mode optimization with given definition of observables.
1650/// The cache operation mode of all objects in the expression tree will
1651/// modified such that all nodes that depend directly or indirectly on
1652/// any of the listed observables will be set to ADirty, as they are
1653/// expected to change every time. This save change tracking overhead for
1654/// nodes that are a priori known to change every time
1655
1656void RooAbsArg::optimizeCacheMode(const RooArgSet& observables, RooArgSet& optimizedNodes, RooLinkedList& processedNodes)
1657{
1658 // Optimization applies only to branch nodes, not to leaf nodes
1659 if (!isDerived()) {
1660 return ;
1661 }
1662
1663
1664 // Terminate call if this node was already processed (tree structure may be cyclical)
1665 // LM : RooLinkedList::findArg looks by name and not but by object pointer,
1666 // should one use RooLinkedList::FindObject (look byt pointer) instead of findArg when
1667 // tree contains nodes with the same name ?
1668 // Add an info message if the require node does not exist but a different node already exists with same name
1669
1670 if (processedNodes.FindObject(this))
1671 return;
1672
1673 // check if findArgs returns something different (i.e. a different node with same name) when
1674 // this node has not been processed (FindObject returns a null pointer)
1675 auto obj = processedNodes.findArg(this);
1676 assert(obj != this); // obj == this cannot happen
1677 if (obj)
1678 // here for nodes with duplicate names
1679 cxcoutI(Optimization) << "RooAbsArg::optimizeCacheMode(" << GetName()
1680 << " node " << this << " exists already as " << obj << " but with the SAME name !" << endl;
1681
1682 processedNodes.Add(this);
1683
1684 // Set cache mode operator to 'AlwaysDirty' if we depend on any of the given observables
1685 if (dependsOnValue(observables)) {
1686
1687 if (dynamic_cast<RooRealIntegral*>(this)) {
1688 cxcoutI(Integration) << "RooAbsArg::optimizeCacheMode(" << GetName() << ") integral depends on value of one or more observables and will be evaluated for every event" << endl ;
1689 }
1690 optimizedNodes.add(*this,kTRUE) ;
1691 if (operMode()==AClean) {
1692 } else {
1693 setOperMode(ADirty,kTRUE) ; // WVE propagate flag recursively to top of tree
1694 }
1695 } else {
1696 }
1697 // Process any RooAbsArgs contained in any of the caches of this object
1698 for (Int_t i=0 ;i<numCaches() ; i++) {
1699 getCache(i)->optimizeCacheMode(observables,optimizedNodes,processedNodes) ;
1700 }
1701
1702 // Forward calls to all servers
1703 for (const auto server : _serverList) {
1704 server->optimizeCacheMode(observables,optimizedNodes,processedNodes) ;
1705 }
1706
1707}
1708
1709////////////////////////////////////////////////////////////////////////////////
1710/// Find branch nodes with all-constant parameters, and add them to the list of
1711/// nodes that can be cached with a dataset in a test statistic calculation
1712
1714{
1715 RooLinkedList proc ;
1716 Bool_t ret = findConstantNodes(observables,cacheList,proc) ;
1717
1718 // If node can be optimized and hasn't been identified yet, add it to the list
1719 coutI(Optimization) << "RooAbsArg::findConstantNodes(" << GetName() << "): components "
1720 << cacheList << " depend exclusively on constant parameters and will be precalculated and cached" << endl ;
1721
1722 return ret ;
1723}
1724
1725
1726
1727////////////////////////////////////////////////////////////////////////////////
1728/// Find branch nodes with all-constant parameters, and add them to the list of
1729/// nodes that can be cached with a dataset in a test statistic calculation
1730
1731Bool_t RooAbsArg::findConstantNodes(const RooArgSet& observables, RooArgSet& cacheList, RooLinkedList& processedNodes)
1732{
1733 // Caching only applies to branch nodes
1734 if (!isDerived()) {
1735 return kFALSE;
1736 }
1737
1738 // Terminate call if this node was already processed (tree structure may be cyclical)
1739 if (processedNodes.findArg(this)) {
1740 return kFALSE ;
1741 } else {
1742 processedNodes.Add(this) ;
1743 }
1744
1745 // Check if node depends on any non-constant parameter
1746 Bool_t canOpt(kTRUE) ;
1747 RooArgSet* paramSet = getParameters(observables) ;
1748 RooFIter iter = paramSet->fwdIterator() ;
1749 RooAbsArg* param ;
1750 while((param = iter.next())) {
1751 if (!param->isConstant()) {
1752 canOpt=kFALSE ;
1753 break ;
1754 }
1755 }
1756 delete paramSet ;
1757
1758
1759 if (getAttribute("NeverConstant")) {
1760 canOpt = kFALSE ;
1761 }
1762
1763 if (canOpt) {
1764 setAttribute("ConstantExpression") ;
1765 }
1766
1767 // If yes, list node eligible for caching, if not test nodes one level down
1768 if (canOpt||getAttribute("CacheAndTrack")) {
1769
1770 if (!cacheList.find(*this) && dependsOnValue(observables) && !observables.find(*this) ) {
1771
1772 // Add to cache list
1773 cxcoutD(Optimization) << "RooAbsArg::findConstantNodes(" << GetName() << ") adding self to list of constant nodes" << endl ;
1774
1775 if (canOpt) setAttribute("ConstantExpressionCached") ;
1776 cacheList.add(*this,kFALSE) ;
1777 }
1778 }
1779
1780 if (!canOpt) {
1781
1782 // If not, see if next level down can be cached
1783 for (const auto server : _serverList) {
1784 if (server->isDerived()) {
1785 server->findConstantNodes(observables,cacheList,processedNodes) ;
1786 }
1787 }
1788 }
1789
1790 // Forward call to all cached contained in current object
1791 for (Int_t i=0 ;i<numCaches() ; i++) {
1792 getCache(i)->findConstantNodes(observables,cacheList,processedNodes) ;
1793 }
1794
1795 return kFALSE ;
1796}
1797
1798
1799
1800
1801////////////////////////////////////////////////////////////////////////////////
1802/// Interface function signaling a request to perform constant term
1803/// optimization. This default implementation takes no action other than to
1804/// forward the calls to all servers
1805
1807{
1808 for (const auto server : _serverList) {
1809 server->constOptimizeTestStatistic(opcode,doAlsoTrackingOpt) ;
1810 }
1811}
1812
1813
1814////////////////////////////////////////////////////////////////////////////////
1815/// Change cache operation mode to given mode. If recurseAdirty
1816/// is true, then a mode change to AlwaysDirty will automatically
1817/// be propagated recursively to all client nodes
1818
1819void RooAbsArg::setOperMode(OperMode mode, Bool_t recurseADirty)
1820{
1821 // Prevent recursion loops
1822 if (mode==_operMode) return ;
1823
1824 _operMode = mode ;
1825 _fast = ((mode==AClean) || dynamic_cast<RooRealVar*>(this)!=0 || dynamic_cast<RooConstVar*>(this)!=0 ) ;
1826 for (Int_t i=0 ;i<numCaches() ; i++) {
1827 getCache(i)->operModeHook() ;
1828 }
1829 operModeHook() ;
1830
1831 // Propagate to all clients
1832 if (mode==ADirty && recurseADirty) {
1833 for (auto clientV : _clientListValue) {
1834 clientV->setOperMode(mode) ;
1835 }
1836 }
1837}
1838
1839
1840////////////////////////////////////////////////////////////////////////////////
1841/// Print tree structure of expression tree on stdout, or to file if filename is specified.
1842/// If namePat is not "*", only nodes with names matching the pattern will be printed.
1843/// The client argument is used in recursive calls to properly display the value or shape nature
1844/// of the client-server links. It should be zero in calls initiated by users.
1845
1846void RooAbsArg::printCompactTree(const char* indent, const char* filename, const char* namePat, RooAbsArg* client)
1847{
1848 if (filename) {
1849 ofstream ofs(filename) ;
1850 printCompactTree(ofs,indent,namePat,client) ;
1851 } else {
1852 printCompactTree(cout,indent,namePat,client) ;
1853 }
1854}
1855
1856
1857////////////////////////////////////////////////////////////////////////////////
1858/// Print tree structure of expression tree on given ostream.
1859/// If namePat is not "*", only nodes with names matching the pattern will be printed.
1860/// The client argument is used in recursive calls to properly display the value or shape nature
1861/// of the client-server links. It should be zero in calls initiated by users.
1862
1863void RooAbsArg::printCompactTree(ostream& os, const char* indent, const char* namePat, RooAbsArg* client)
1864{
1865 if ( !namePat || TString(GetName()).Contains(namePat)) {
1866 os << indent << this ;
1867 if (client) {
1868 os << "/" ;
1869 if (isValueServer(*client)) os << "V" ; else os << "-" ;
1870 if (isShapeServer(*client)) os << "S" ; else os << "-" ;
1871 }
1872 os << " " ;
1873
1874 os << IsA()->GetName() << "::" << GetName() << " = " ;
1875 printValue(os) ;
1876
1877 if (!_serverList.empty()) {
1878 switch(operMode()) {
1879 case Auto: os << " [Auto," << (isValueDirty()?"Dirty":"Clean") << "] " ; break ;
1880 case AClean: os << " [ACLEAN] " ; break ;
1881 case ADirty: os << " [ADIRTY] " ; break ;
1882 }
1883 }
1884 os << endl ;
1885
1886 for (Int_t i=0 ;i<numCaches() ; i++) {
1888 }
1890 }
1891
1892 TString indent2(indent) ;
1893 indent2 += " " ;
1894 for (const auto arg : _serverList) {
1895 arg->printCompactTree(os,indent2,namePat,this) ;
1896 }
1897}
1898
1899
1900////////////////////////////////////////////////////////////////////////////////
1901/// Print tree structure of expression tree on given ostream, only branch nodes are printed.
1902/// Lead nodes (variables) will not be shown
1903///
1904/// If namePat is not "*", only nodes with names matching the pattern will be printed.
1905
1906void RooAbsArg::printComponentTree(const char* indent, const char* namePat, Int_t nLevel)
1907{
1908 if (nLevel==0) return ;
1909 if (isFundamental()) return ;
1910 RooResolutionModel* rmodel = dynamic_cast<RooResolutionModel*>(this) ;
1911 if (rmodel && rmodel->isConvolved()) return ;
1912 if (InheritsFrom("RooConstVar")) return ;
1913
1914 if ( !namePat || TString(GetName()).Contains(namePat)) {
1915 cout << indent ;
1916 Print() ;
1917 }
1918
1919 TString indent2(indent) ;
1920 indent2 += " " ;
1921 for (const auto arg : _serverList) {
1922 arg->printComponentTree(indent2.Data(),namePat,nLevel-1) ;
1923 }
1924}
1925
1926
1927////////////////////////////////////////////////////////////////////////////////
1928/// Construct a mangled name from the actual name that
1929/// is free of any math symbols that might be interpreted by TTree
1930
1932{
1933 // Check for optional alternate name of branch for this argument
1934 TString rawBranchName = GetName() ;
1935 if (getStringAttribute("BranchName")) {
1936 rawBranchName = getStringAttribute("BranchName") ;
1937 }
1938
1939 TString cleanName(rawBranchName) ;
1940 cleanName.ReplaceAll("/","D") ;
1941 cleanName.ReplaceAll("-","M") ;
1942 cleanName.ReplaceAll("+","P") ;
1943 cleanName.ReplaceAll("*","X") ;
1944 cleanName.ReplaceAll("[","L") ;
1945 cleanName.ReplaceAll("]","R") ;
1946 cleanName.ReplaceAll("(","L") ;
1947 cleanName.ReplaceAll(")","R") ;
1948 cleanName.ReplaceAll("{","L") ;
1949 cleanName.ReplaceAll("}","R") ;
1950
1951 return cleanName;
1952}
1953
1954
1955////////////////////////////////////////////////////////////////////////////////
1956/// Hook function interface for object to insert additional information
1957/// when printed in the context of a tree structure. This default
1958/// implementation prints nothing
1959
1960void RooAbsArg::printCompactTreeHook(ostream&, const char *)
1961{
1962}
1963
1964
1965////////////////////////////////////////////////////////////////////////////////
1966/// Register RooAbsCache with this object. This function is called
1967/// by RooAbsCache constructors for objects that are a datamember
1968/// of this RooAbsArg. By registering itself the RooAbsArg is aware
1969/// of all its cache data members and will forward server change
1970/// and cache mode change calls to the cache objects, which in turn
1971/// can forward them their contents
1972
1974{
1975 _cacheList.push_back(&cache) ;
1976}
1977
1978
1979////////////////////////////////////////////////////////////////////////////////
1980/// Unregister a RooAbsCache. Called from the RooAbsCache destructor
1981
1983{
1984 _cacheList.erase(std::remove(_cacheList.begin(), _cacheList.end(), &cache),
1985 _cacheList.end());
1986}
1987
1988
1989////////////////////////////////////////////////////////////////////////////////
1990/// Return number of registered caches
1991
1993{
1994 return _cacheList.size() ;
1995}
1996
1997
1998////////////////////////////////////////////////////////////////////////////////
1999/// Return registered cache object by index
2000
2002{
2003 return _cacheList[index] ;
2004}
2005
2006
2007////////////////////////////////////////////////////////////////////////////////
2008/// Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
2009
2011{
2012 return getParameters(RooArgSet(),stripDisconnected) ;
2013}
2014
2015
2016////////////////////////////////////////////////////////////////////////////////
2017/// Return ancestors in cloning chain of this RooAbsArg. NOTE: Returned pointers
2018/// are not guaranteed to be 'live', so do not dereference without proper caution
2019
2021{
2022 RooLinkedList retVal ;
2023
2024 set<string>::const_iterator iter= _boolAttrib.begin() ;
2025 while(iter != _boolAttrib.end()) {
2026 if (TString(*iter).BeginsWith("CloneOf(")) {
2027 char buf[128] ;
2028 strlcpy(buf,iter->c_str(),128) ;
2029 strtok(buf,"(") ;
2030 char* ptrToken = strtok(0,")") ;
2031 RooAbsArg* ptr = (RooAbsArg*) strtoll(ptrToken,0,16) ;
2032 retVal.Add(ptr) ;
2033 }
2034 }
2035
2036 return retVal ;
2037}
2038
2039
2040////////////////////////////////////////////////////////////////////////////////
2041/// Create a GraphViz .dot file visualizing the expression tree headed by
2042/// this RooAbsArg object. Use the GraphViz tool suite to make e.g. a gif
2043/// or ps file from the .dot file.
2044/// If a node derives from RooAbsReal, its current (unnormalised) value is
2045/// printed as well.
2046///
2047/// Based on concept developed by Kyle Cranmer.
2048
2049void RooAbsArg::graphVizTree(const char* fileName, const char* delimiter, bool useTitle, bool useLatex)
2050{
2051 ofstream ofs(fileName) ;
2052 if (!ofs) {
2053 coutE(InputArguments) << "RooAbsArg::graphVizTree() ERROR: Cannot open graphViz output file with name " << fileName << endl ;
2054 return ;
2055 }
2056 graphVizTree(ofs, delimiter, useTitle, useLatex) ;
2057}
2058
2059////////////////////////////////////////////////////////////////////////////////
2060/// Write the GraphViz representation of the expression tree headed by
2061/// this RooAbsArg object to the given ostream.
2062/// If a node derives from RooAbsReal, its current (unnormalised) value is
2063/// printed as well.
2064///
2065/// Based on concept developed by Kyle Cranmer.
2066
2067void RooAbsArg::graphVizTree(ostream& os, const char* delimiter, bool useTitle, bool useLatex)
2068{
2069 if (!os) {
2070 coutE(InputArguments) << "RooAbsArg::graphVizTree() ERROR: output stream provided as input argument is in invalid state" << endl ;
2071 }
2072
2073 // silent warning messages coming when evaluating a RooAddPdf without a normalization set
2075
2076 // Write header
2077 os << "digraph \"" << GetName() << "\"{" << endl ;
2078
2079 // First list all the tree nodes
2080 RooArgSet nodeSet ;
2081 treeNodeServerList(&nodeSet) ;
2082 RooFIter iter = nodeSet.fwdIterator() ;
2083 RooAbsArg* node ;
2084
2085 // iterate over nodes
2086 while((node=iter.next())) {
2087 string nodeName = node->GetName();
2088 string nodeTitle = node->GetTitle();
2089 string nodeLabel = (useTitle && !nodeTitle.empty()) ? nodeTitle : nodeName;
2090
2091 // if using latex, replace ROOT's # with normal latex backslash
2092 string::size_type position = nodeLabel.find("#") ;
2093 while(useLatex && position!=nodeLabel.npos){
2094 nodeLabel.replace(position, 1, "\\");
2095 }
2096
2097 string typeFormat = "\\texttt{";
2098 string nodeType = (useLatex) ? typeFormat+node->IsA()->GetName()+"}" : node->IsA()->GetName();
2099
2100 if (auto realNode = dynamic_cast<RooAbsReal*>(node)) {
2101 nodeLabel += delimiter + std::to_string(realNode->getVal());
2102 }
2103
2104 os << "\"" << nodeName << "\" [ color=" << (node->isFundamental()?"blue":"red")
2105 << ", label=\"" << nodeType << delimiter << nodeLabel << "\"];" << endl ;
2106
2107 }
2108
2109 // Get set of all server links
2110 set<pair<RooAbsArg*,RooAbsArg*> > links ;
2111 graphVizAddConnections(links) ;
2112
2113 // And write them out
2114 for(auto const& link : links) {
2115 os << "\"" << link.first->GetName() << "\" -> \"" << link.second->GetName() << "\";" << endl ;
2116 }
2117
2118 // Write trailer
2119 os << "}" << endl ;
2120
2121}
2122
2123////////////////////////////////////////////////////////////////////////////////
2124/// Utility function that inserts all point-to-point client-server connections
2125/// between any two RooAbsArgs in the expression tree headed by this object
2126/// in the linkSet argument.
2127
2128void RooAbsArg::graphVizAddConnections(set<pair<RooAbsArg*,RooAbsArg*> >& linkSet)
2129{
2130 for (const auto server : _serverList) {
2131 linkSet.insert(make_pair(this,server)) ;
2132 server->graphVizAddConnections(linkSet) ;
2133 }
2134}
2135
2136
2137
2138// //_____________________________________________________________________________
2139// TGraphStruct* RooAbsArg::graph(Bool_t useFactoryTag, Double_t textSize)
2140// {
2141// // Return a TGraphStruct object filled with the tree structure of the pdf object
2142
2143// TGraphStruct* theGraph = new TGraphStruct() ;
2144
2145// // First list all the tree nodes
2146// RooArgSet nodeSet ;
2147// treeNodeServerList(&nodeSet) ;
2148// TIterator* iter = nodeSet.createIterator() ;
2149// RooAbsArg* node ;
2150
2151
2152// // iterate over nodes
2153// while((node=(RooAbsArg*)iter->Next())) {
2154
2155// // Skip node that represent numeric constants
2156// if (node->IsA()->InheritsFrom(RooConstVar::Class())) continue ;
2157
2158// string nodeName ;
2159// if (useFactoryTag && node->getStringAttribute("factory_tag")) {
2160// nodeName = node->getStringAttribute("factory_tag") ;
2161// } else {
2162// if (node->isFundamental()) {
2163// nodeName = node->GetName();
2164// } else {
2165// ostringstream oss ;
2166// node->printStream(oss,(node->defaultPrintContents(0)&(~kValue)),node->defaultPrintStyle(0)) ;
2167// nodeName= oss.str() ;
2168// // nodeName = Form("%s::%s",node->IsA()->GetName(),node->GetName());
2169
2170// }
2171// }
2172// if (strncmp(nodeName.c_str(),"Roo",3)==0) {
2173// nodeName = string(nodeName.c_str()+3) ;
2174// }
2175// node->setStringAttribute("graph_name",nodeName.c_str()) ;
2176
2177// TGraphNode* gnode = theGraph->AddNode(nodeName.c_str(),nodeName.c_str()) ;
2178// gnode->SetLineWidth(2) ;
2179// gnode->SetTextColor(node->isFundamental()?kBlue:kRed) ;
2180// gnode->SetTextSize(textSize) ;
2181// }
2182// delete iter ;
2183
2184// // Get set of all server links
2185// set<pair<RooAbsArg*,RooAbsArg*> > links ;
2186// graphVizAddConnections(links) ;
2187
2188// // And insert them into the graph
2189// set<pair<RooAbsArg*,RooAbsArg*> >::iterator liter = links.begin() ;
2190// for( ; liter != links.end() ; ++liter ) {
2191
2192// TGraphNode* n1 = (TGraphNode*)theGraph->GetListOfNodes()->FindObject(liter->first->getStringAttribute("graph_name")) ;
2193// TGraphNode* n2 = (TGraphNode*)theGraph->GetListOfNodes()->FindObject(liter->second->getStringAttribute("graph_name")) ;
2194// if (n1 && n2) {
2195// TGraphEdge* edge = theGraph->AddEdge(n1,n2) ;
2196// edge->SetLineWidth(2) ;
2197// }
2198// }
2199
2200// return theGraph ;
2201// }
2202
2203
2204
2205// //_____________________________________________________________________________
2206// Bool_t RooAbsArg::inhibitDirty()
2207// {
2208// // Return current status of the inhibitDirty global flag. If true
2209// // no dirty state change tracking occurs and all caches are considered
2210// // to be always dirty.
2211// return _inhibitDirty ;
2212// }
2213
2214
2215////////////////////////////////////////////////////////////////////////////////
2216/// Take ownership of the contents of 'comps'
2217
2219{
2220 if (!_ownedComponents) {
2221 _ownedComponents = new RooArgSet("owned components") ;
2222 }
2223 return _ownedComponents->addOwned(comps) ;
2224}
2225
2226
2227
2228////////////////////////////////////////////////////////////////////////////////
2229/// Clone tree expression of objects. All tree nodes will be owned by
2230/// the head node return by cloneTree()
2231
2232RooAbsArg* RooAbsArg::cloneTree(const char* newname) const
2233{
2234 // Clone tree using snapshot
2235 RooArgSet* clonedNodes = (RooArgSet*) RooArgSet(*this).snapshot(kTRUE) ;
2236
2237 // Find the head node in the cloneSet
2238 RooAbsArg* head = clonedNodes->find(*this) ;
2239 assert(head);
2240
2241 // Remove the head node from the cloneSet
2242 // To release it from the set ownership
2243 clonedNodes->remove(*head) ;
2244
2245 // Add the set as owned component of the head
2246 head->addOwnedComponents(*clonedNodes) ;
2247
2248 // Delete intermediate container
2249 clonedNodes->releaseOwnership() ;
2250 delete clonedNodes ;
2251
2252 // Adjust name of head node if requested
2253 if (newname) {
2254 head->TNamed::SetName(newname) ;
2255 head->_namePtr = (TNamed*) RooNameReg::instance().constPtr(newname) ;
2256 }
2257
2258 // Return the head
2259 return head ;
2260}
2261
2262
2263
2264////////////////////////////////////////////////////////////////////////////////
2265
2267{
2268 if (dynamic_cast<RooTreeDataStore*>(&store)) {
2269 attachToTree(((RooTreeDataStore&)store).tree()) ;
2270 } else if (dynamic_cast<RooVectorDataStore*>(&store)) {
2272 }
2273}
2274
2275
2276
2277////////////////////////////////////////////////////////////////////////////////
2278
2280{
2281 if (_eocache) {
2282 return *_eocache ;
2283 } else {
2285 }
2286}
2287
2288
2289////////////////////////////////////////////////////////////////////////////////
2290
2292{
2293 string suffix ;
2294
2295 RooArgSet branches ;
2296 branchNodeServerList(&branches) ;
2297 RooFIter iter = branches.fwdIterator();
2298 RooAbsArg* arg ;
2299 while((arg=iter.next())) {
2300 const char* tmp = arg->cacheUniqueSuffix() ;
2301 if (tmp) suffix += tmp ;
2302 }
2303 return Form("%s",suffix.c_str()) ;
2304}
2305
2306
2307////////////////////////////////////////////////////////////////////////////////
2308
2310{
2311 RooArgSet branches ;
2312 branchNodeServerList(&branches) ;
2313 for(auto const& arg : branches) {
2314 for (auto const& arg2 : arg->_cacheList) {
2315 arg2->wireCache() ;
2316 }
2317 }
2318}
2319
2320
2321
2322////////////////////////////////////////////////////////////////////////////////
2323
2324void RooAbsArg::SetName(const char* name)
2325{
2328 if (newPtr != _namePtr) {
2329 //cout << "Rename '" << _namePtr->GetName() << "' to '" << name << "' (set flag in new name)" << endl;
2330 _namePtr = newPtr;
2332 }
2333}
2334
2335
2336
2337
2338////////////////////////////////////////////////////////////////////////////////
2339
2340void RooAbsArg::SetNameTitle(const char *name, const char *title)
2341{
2342 TNamed::SetNameTitle(name,title) ;
2344 if (newPtr != _namePtr) {
2345 //cout << "Rename '" << _namePtr->GetName() << "' to '" << name << "' (set flag in new name)" << endl;
2346 _namePtr = newPtr;
2348 }
2349}
2350
2351
2352////////////////////////////////////////////////////////////////////////////////
2353/// Stream an object of class RooAbsArg.
2354
2355void RooAbsArg::Streamer(TBuffer &R__b)
2356{
2357 if (R__b.IsReading()) {
2358 _ioReadStack.push(this) ;
2359 R__b.ReadClassBuffer(RooAbsArg::Class(),this);
2360 _ioReadStack.pop() ;
2362 _isConstant = getAttribute("Constant") ;
2363 } else {
2365 }
2366}
2367
2368////////////////////////////////////////////////////////////////////////////////
2369/// Method called by workspace container to finalize schema evolution issues
2370/// that cannot be handled in a single ioStreamer pass.
2371///
2372/// A second pass is typically needed when evolving data member of RooAbsArg-derived
2373/// classes that are container classes with references to other members, which may
2374/// not yet be 'live' in the first ioStreamer() evolution pass.
2375///
2376/// Classes may overload this function, but must call the base method in the
2377/// overloaded call to ensure base evolution is handled properly
2378
2380{
2381 // Handling of v5-v6 migration (TRefArray _proxyList --> RooRefArray _proxyList)
2382 auto iter = _ioEvoList.find(this);
2383 if (iter != _ioEvoList.end()) {
2384
2385 // Transfer contents of saved TRefArray to RooRefArray now
2387 _proxyList.Expand(iter->second->GetEntriesFast());
2388 for (int i = 0; i < iter->second->GetEntriesFast(); i++) {
2389 _proxyList.Add(iter->second->At(i));
2390 }
2391 // Delete TRefArray and remove from list
2392 _ioEvoList.erase(iter);
2393 }
2394}
2395
2396
2397
2398
2399////////////////////////////////////////////////////////////////////////////////
2400/// Method called by workspace container to finalize schema evolution issues
2401/// that cannot be handled in a single ioStreamer pass. This static finalize method
2402/// is called after ioStreamerPass2() is called on each directly listed object
2403/// in the workspace. It's purpose is to complete schema evolution of any
2404/// objects in the workspace that are not directly listed as content elements
2405/// (e.g. analytical convolution tokens )
2406
2408{
2409 // Handling of v5-v6 migration (TRefArray _proxyList --> RooRefArray _proxyList)
2410 for (const auto& iter : _ioEvoList) {
2411
2412 // Transfer contents of saved TRefArray to RooRefArray now
2413 if (!iter.first->_proxyList.GetEntriesFast())
2414 iter.first->_proxyList.Expand(iter.second->GetEntriesFast());
2415 for (int i = 0; i < iter.second->GetEntriesFast(); i++) {
2416 iter.first->_proxyList.Add(iter.second->At(i));
2417 }
2418 }
2419
2420 _ioEvoList.clear();
2421}
2422
2423
2427}
2428
2429////////////////////////////////////////////////////////////////////////////////
2430/// Stream an object of class RooRefArray.
2431
2432void RooRefArray::Streamer(TBuffer &R__b)
2433{
2434 UInt_t R__s, R__c;
2435 if (R__b.IsReading()) {
2436
2437 Version_t R__v = R__b.ReadVersion(&R__s, &R__c); if (R__v) { }
2438
2439 // Make temporary refArray and read that from the streamer
2440 auto refArray = std::make_unique<TRefArray>();
2441 refArray->Streamer(R__b) ;
2442 R__b.CheckByteCount(R__s, R__c, refArray->IsA());
2443
2444 // Schedule deferred processing of TRefArray into proxy list
2445 RooAbsArg::_ioEvoList[RooAbsArg::_ioReadStack.top()] = std::move(refArray);
2446
2447 } else {
2448
2449 R__c = R__b.WriteVersion(RooRefArray::IsA(), kTRUE);
2450
2451 // Make a temporary refArray and write that to the streamer
2452 TRefArray refArray(GetEntriesFast());
2453 TIterator* iter = MakeIterator() ;
2454 TObject* tmpObj ; while ((tmpObj = iter->Next())) {
2455 refArray.Add(tmpObj) ;
2456 }
2457 delete iter ;
2458
2459 refArray.Streamer(R__b) ;
2460 R__b.SetByteCount(R__c, kTRUE) ;
2461
2462 }
2463}
2464
2465/// Print at the prompt
2466namespace cling {
2467std::string printValue(RooAbsArg *raa)
2468{
2469 std::stringstream s;
2470 if (0 == *raa->GetName() && 0 == *raa->GetTitle()) {
2471 s << "An instance of " << raa->ClassName() << ".";
2472 return s.str();
2473 }
2474 raa->printStream(s, raa->defaultPrintContents(""), raa->defaultPrintStyle(""));
2475 return s.str();
2476}
2477} // namespace cling
void Class()
Definition: Class.C:29
istream & operator>>(istream &is, RooAbsArg &arg)
Istream operator.
Definition: RooAbsArg.cxx:1528
ostream & operator<<(ostream &os, RooAbsArg const &arg)
Ostream operator.
Definition: RooAbsArg.cxx:1519
#define coutI(a)
Definition: RooMsgService.h:30
#define cxcoutI(a)
Definition: RooMsgService.h:85
#define cxcoutD(a)
Definition: RooMsgService.h:81
#define coutF(a)
Definition: RooMsgService.h:34
#define coutE(a)
Definition: RooMsgService.h:33
#define cxcoutF(a)
int Int_t
Definition: RtypesCore.h:45
char Text_t
Definition: RtypesCore.h:62
short Version_t
Definition: RtypesCore.h:65
unsigned int UInt_t
Definition: RtypesCore.h:46
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:364
static void indent(ostringstream &buf, int indent_level)
char name[80]
Definition: TGX11.cxx:110
Binding & operator=(OUT(*fun)(void))
char * Form(const char *fmt,...)
TString operator+(const TString &s1, const TString &s2)
Use the special concatenation constructor.
Definition: TString.cxx:1499
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooExpensiveObjectCache & expensiveObjectCache() const
Definition: RooAbsArg.cxx:2279
virtual Bool_t checkObservables(const RooArgSet *nset) const
Overloadable function in which derived classes can implement consistency checks of the variables.
Definition: RooAbsArg.cxx:764
virtual RooAbsArg * cloneTree(const char *newname=0) const
Clone tree expression of objects.
Definition: RooAbsArg.cxx:2232
RooRefArray _proxyList
Definition: RooAbsArg.h:613
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:295
Bool_t _fast
Definition: RooAbsArg.h:682
virtual Bool_t readFromStream(std::istream &is, Bool_t compact, Bool_t verbose=kFALSE)=0
virtual Bool_t isValid() const
WVE (08/21/01) Probably obsolete now.
Definition: RooAbsArg.cxx:1364
TNamed * _namePtr
Definition: RooAbsArg.h:691
static void verboseDirty(Bool_t flag)
Activate verbose messaging related to dirty flag propagation.
Definition: RooAbsArg.cxx:273
void attachToStore(RooAbsDataStore &store)
Attach this argument to the data store such that it reads data from there.
Definition: RooAbsArg.cxx:2266
virtual Int_t defaultPrintContents(Option_t *opt) const
Define default contents to print.
Definition: RooAbsArg.cxx:1435
const TNamed * namePtr() const
De-duplicated pointer to this object's name.
Definition: RooAbsArg.h:544
const char * aggregateCacheUniqueSuffix() const
Definition: RooAbsArg.cxx:2291
Bool_t redirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t isRecursionStep=kFALSE)
Replace all direct servers of this object with the new servers in newServerList.
Definition: RooAbsArg.cxx:995
void changeServer(RooAbsArg &server, Bool_t valueProp, Bool_t shapeProp)
Change dirty flag propagation mask for specified server.
Definition: RooAbsArg.cxx:475
Bool_t isValueServer(const RooAbsArg &arg) const
Check if this is serving values to arg.
Definition: RooAbsArg.h:217
virtual Bool_t redirectServersHook(const RooAbsCollection &, Bool_t, Bool_t, Bool_t)
Definition: RooAbsArg.h:263
RooWorkspace * _myws
Prevent 'AlwaysDirty' mode for this node.
Definition: RooAbsArg.h:699
void printCompactTree(const char *indent="", const char *fileName=0, const char *namePat=0, RooAbsArg *client=0)
Print tree structure of expression tree on stdout, or to file if filename is specified.
Definition: RooAbsArg.cxx:1846
void attachDataStore(const RooAbsDataStore &set)
Replace server nodes with names matching the dataset variable names with those data set variables,...
Definition: RooAbsArg.cxx:1583
virtual void printClassName(std::ostream &os) const
Print object class name.
Definition: RooAbsArg.cxx:1395
RooArgSet * _ownedComponents
Definition: RooAbsArg.h:685
void setShapeDirty()
Notify that a shape-like property (e.g. binning) has changed.
Definition: RooAbsArg.h:496
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:799
void registerProxy(RooArgProxy &proxy)
Register an RooArgProxy in the proxy list.
Definition: RooAbsArg.cxx:1205
void attachArgs(const RooAbsCollection &set)
Bind this node to objects in set.
Definition: RooAbsArg.cxx:1558
Bool_t isShapeServer(const RooAbsArg &arg) const
Check if this is serving shape to arg.
Definition: RooAbsArg.h:225
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string 'value' to this object under key 'key'.
Definition: RooAbsArg.cxx:323
static void ioStreamerPass2Finalize()
Method called by workspace container to finalize schema evolution issues that cannot be handled in a ...
Definition: RooAbsArg.cxx:2407
void leafNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with all leaf nodes of the arg tree, starting with ourself as top node.
Definition: RooAbsArg.cxx:509
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Implement multi-line detailed printing.
Definition: RooAbsArg.cxx:1445
Bool_t isCloneOf(const RooAbsArg &other) const
Check if this object was created as a clone of 'other'.
Definition: RooAbsArg.cxx:281
Bool_t isShapeDirty() const
Definition: RooAbsArg.h:417
static void setDirtyInhibit(Bool_t flag)
Control global dirty inhibit mode.
Definition: RooAbsArg.cxx:264
Bool_t findConstantNodes(const RooArgSet &observables, RooArgSet &cacheList)
Find branch nodes with all-constant parameters, and add them to the list of nodes that can be cached ...
Definition: RooAbsArg.cxx:1713
void graphVizAddConnections(std::set< std::pair< RooAbsArg *, RooAbsArg * > > &)
Utility function that inserts all point-to-point client-server connections between any two RooAbsArgs...
Definition: RooAbsArg.cxx:2128
void unRegisterProxy(RooArgProxy &proxy)
Remove proxy from proxy list.
Definition: RooAbsArg.cxx:1233
void printComponentTree(const char *indent="", const char *namePat=0, Int_t nLevel=999)
Print tree structure of expression tree on given ostream, only branch nodes are printed.
Definition: RooAbsArg.cxx:1906
friend class RooArgSet
Definition: RooAbsArg.h:600
Bool_t inhibitDirty() const
Delete watch flag.
Definition: RooAbsArg.cxx:116
std::set< std::string > _boolAttrib
Definition: RooAbsArg.h:635
void addServer(RooAbsArg &server, Bool_t valueProp=kTRUE, Bool_t shapeProp=kFALSE, std::size_t refCount=1)
Register another RooAbsArg as a server to us, ie, declare that we depend on it.
Definition: RooAbsArg.cxx:388
virtual void printTree(std::ostream &os, TString indent="") const
Print object tree structure.
Definition: RooAbsArg.cxx:1510
void unRegisterCache(RooAbsCache &cache)
Unregister a RooAbsCache. Called from the RooAbsCache destructor.
Definition: RooAbsArg.cxx:1982
virtual Bool_t isFundamental() const
Is this object a fundamental type that can be added to a dataset? Fundamental-type subclasses overrid...
Definition: RooAbsArg.h:243
virtual void printTitle(std::ostream &os) const
Print object title.
Definition: RooAbsArg.cxx:1385
RefCountList_t _clientListValue
Definition: RooAbsArg.h:611
void printAttribList(std::ostream &os) const
Transient boolean attributes (not copied in ctor)
Definition: RooAbsArg.cxx:1537
virtual void Print(Option_t *options=0) const
Print the object to the defaultPrintStream().
Definition: RooAbsArg.h:323
const Text_t * getStringAttribute(const Text_t *key) const
Get string attribute mapped under key 'key'.
Definition: RooAbsArg.cxx:336
void treeNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t doBranch=kTRUE, Bool_t doLeaf=kTRUE, Bool_t valueOnly=kFALSE, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with nodes of the arg tree, following all server links, starting with ourself as t...
Definition: RooAbsArg.cxx:536
void setAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:291
virtual void writeToStream(std::ostream &os, Bool_t compact) const =0
virtual void printAddress(std::ostream &os) const
Print class name of object.
Definition: RooAbsArg.cxx:1401
void setTransientAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:350
virtual void printCompactTreeHook(std::ostream &os, const char *ind="")
Hook function interface for object to insert additional information when printed in the context of a ...
Definition: RooAbsArg.cxx:1960
Bool_t _isConstant
De-duplicated name pointer. This will be equal for all objects with the same name.
Definition: RooAbsArg.h:692
Bool_t isValueDirty() const
Definition: RooAbsArg.h:422
Bool_t overlaps(const RooAbsArg &testArg, Bool_t valueOnly=kFALSE) const
Test if any of the nodes of tree are shared with that of the given tree.
Definition: RooAbsArg.cxx:854
Int_t Compare(const TObject *other) const
Utility function used by TCollection::Sort to compare contained TObjects We implement comparison by n...
Definition: RooAbsArg.cxx:1593
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2010
virtual void ioStreamerPass2()
Method called by workspace container to finalize schema evolution issues that cannot be handled in a ...
Definition: RooAbsArg.cxx:2379
void wireAllCaches()
Definition: RooAbsArg.cxx:2309
Bool_t _valueDirty
Definition: RooAbsArg.h:677
virtual ~RooAbsArg()
Destructor.
Definition: RooAbsArg.cxx:222
virtual Bool_t isDerived() const
Does value or shape of this arg depend on any other arg?
Definition: RooAbsArg.h:92
static Bool_t _inhibitDirty
Definition: RooAbsArg.h:666
virtual const char * cacheUniqueSuffix() const
Definition: RooAbsArg.h:499
RefCountListLegacyIterator_t * makeLegacyIterator(const RefCountList_t &list) const
Definition: RooAbsArg.cxx:2425
Bool_t _localNoInhibitDirty
Cached isConstant status.
Definition: RooAbsArg.h:694
virtual void constOptimizeTestStatistic(ConstOpCode opcode, Bool_t doAlsoTrackingOpt=kTRUE)
Interface function signaling a request to perform constant term optimization.
Definition: RooAbsArg.cxx:1806
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1819
virtual void attachToTree(TTree &t, Int_t bufSize=32000)=0
Overloadable function for derived classes to implement attachment as branch to a TTree.
Definition: RooAbsArg.cxx:1353
OperMode _operMode
Mark batches as dirty (only meaningful for RooAbsReal).
Definition: RooAbsArg.h:681
RooLinkedList getCloningAncestors() const
Return ancestors in cloning chain of this RooAbsArg.
Definition: RooAbsArg.cxx:2020
void setValueDirty()
Mark the element dirty. This forces a re-evaluation when a value is requested.
Definition: RooAbsArg.h:491
Bool_t addOwnedComponents(const RooArgSet &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2218
Bool_t recursiveCheckObservables(const RooArgSet *nset) const
Recursively call checkObservables on all nodes in the expression tree.
Definition: RooAbsArg.cxx:773
RooAbsCache * getCache(Int_t index) const
Return registered cache object by index.
Definition: RooAbsArg.cxx:2001
void printDirty(Bool_t depth=kTRUE) const
Print information about current value dirty state information.
Definition: RooAbsArg.cxx:1605
bool _allBatchesDirty
Definition: RooAbsArg.h:679
static Bool_t _verboseDirty
Definition: RooAbsArg.h:665
void registerCache(RooAbsCache &cache)
Register RooAbsCache with this object.
Definition: RooAbsArg.cxx:1973
virtual void optimizeCacheMode(const RooArgSet &observables)
Activate cache mode optimization with given definition of observables.
Definition: RooAbsArg.cxx:1637
RefCountList_t _clientListShape
Definition: RooAbsArg.h:610
virtual void attachToVStore(RooVectorDataStore &vstore)=0
TString cleanBranchName() const
Construct a mangled name from the actual name that is free of any math symbols that might be interpre...
Definition: RooAbsArg.cxx:1931
void SetNameTitle(const char *name, const char *title)
Set all the TNamed parameters (name and title).
Definition: RooAbsArg.cxx:2340
Bool_t _prohibitServerRedirect
Set of owned component.
Definition: RooAbsArg.h:687
void removeServer(RooAbsArg &server, Bool_t force=kFALSE)
Unregister another RooAbsArg as a server to us, ie, declare that we no longer depend on its value and...
Definition: RooAbsArg.cxx:439
Int_t numProxies() const
Return the number of registered proxies.
Definition: RooAbsArg.cxx:1327
Bool_t getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
Definition: RooAbsArg.cxx:314
virtual void printName(std::ostream &os) const
Print object name.
Definition: RooAbsArg.cxx:1375
void replaceServer(RooAbsArg &oldServer, RooAbsArg &newServer, Bool_t valueProp, Bool_t shapeProp)
Replace 'oldServer' with 'newServer'.
Definition: RooAbsArg.cxx:464
virtual void getParametersHook(const RooArgSet *, RooArgSet *, Bool_t) const
Definition: RooAbsArg.h:567
virtual void printMetaArgs(std::ostream &) const
Definition: RooAbsArg.h:333
void addServerList(RooAbsCollection &serverList, Bool_t valueProp=kTRUE, Bool_t shapeProp=kFALSE)
Register a list of RooAbsArg as servers to us by calling addServer() for each arg in the list.
Definition: RooAbsArg.cxx:424
Bool_t dependsOnValue(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0) const
Check whether this object depends on values from an element in the serverList.
Definition: RooAbsArg.h:102
void setProxyNormSet(const RooArgSet *nset)
Forward a change in the cached normalization argset to all the registered proxies.
Definition: RooAbsArg.cxx:1338
RefCountList_t _clientList
Definition: RooAbsArg.h:609
RooAbsProxy * getProxy(Int_t index) const
Return the nth proxy from the proxy list.
Definition: RooAbsArg.cxx:1314
RooAbsArg & operator=(const RooAbsArg &other)
Assign all boolean and string properties of the original object.
Definition: RooAbsArg.cxx:190
RefCountList_t _serverList
Definition: RooAbsArg.h:608
RooExpensiveObjectCache * _eocache
Prohibit server redirects – Debugging tool.
Definition: RooAbsArg.h:689
RooArgSet * getComponents() const
Create a RooArgSet with all components (branch nodes) of the expression tree headed by this object.
Definition: RooAbsArg.cxx:746
std::set< std::string > _boolAttribTransient
Definition: RooAbsArg.h:637
Bool_t _shapeDirty
Definition: RooAbsArg.h:678
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...
Definition: RooAbsArg.cxx:578
RooAbsArg * findNewServer(const RooAbsCollection &newSet, Bool_t nameChange) const
Find the new server in the specified set that matches the old server.
Definition: RooAbsArg.cxx:1114
Bool_t isConstant() const
Check if the "Constant" attribute is set.
Definition: RooAbsArg.h:363
void branchNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with all branch nodes of the arg tree starting with ourself as top node.
Definition: RooAbsArg.cxx:520
void graphVizTree(const char *fileName, const char *delimiter="\n", bool useTitle=false, bool useLatex=false)
Create a GraphViz .dot file visualizing the expression tree headed by this RooAbsArg object.
Definition: RooAbsArg.cxx:2049
virtual void operModeHook()
Definition: RooAbsArg.h:561
Bool_t getTransientAttribute(const Text_t *name) const
Check if a named attribute is set.
Definition: RooAbsArg.cxx:372
std::map< std::string, std::string > _stringAttrib
Definition: RooAbsArg.h:636
Int_t numCaches() const
Return number of registered caches.
Definition: RooAbsArg.cxx:1992
RooAbsArg()
Default constructor.
Definition: RooAbsArg.cxx:125
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Recursively replace all servers with the new servers in newSet.
Definition: RooAbsArg.cxx:1159
void attachDataSet(const RooAbsData &set)
Replace server nodes with names matching the dataset variable names with those data set variables,...
Definition: RooAbsArg.cxx:1573
virtual void printArgs(std::ostream &os) const
Print object arguments, ie its proxies.
Definition: RooAbsArg.cxx:1412
TIteratorToSTLInterface< RefCountList_t::Container_t > RefCountListLegacyIterator_t
Definition: RooAbsArg.h:75
RooAbsArg * findServer(const char *name) const
Return server of this with name name. Returns nullptr if not found.
Definition: RooAbsArg.h:203
void SetName(const char *name)
Set the name of the TNamed.
Definition: RooAbsArg.cxx:2324
void addParameters(RooArgSet &params, const RooArgSet *nset=0, Bool_t stripDisconnected=kTRUE) const
INTERNAL helper function for getParameters()
Definition: RooAbsArg.cxx:587
Bool_t _deleteWatch
Definition: RooAbsArg.h:667
std::vector< RooAbsCache * > _cacheList
Definition: RooAbsArg.h:615
OperMode operMode() const
Query the operation mode of this node.
Definition: RooAbsArg.h:485
Bool_t observableOverlaps(const RooAbsData *dset, const RooAbsArg &testArg) const
Test if any of the dependents of the arg tree (as determined by getObservables) overlaps with those o...
Definition: RooAbsArg.cxx:868
RooAbsCache is the abstract base class for data members of RooAbsArgs that cache other (composite) Ro...
Definition: RooAbsCache.h:27
virtual void operModeHook()
Interface for operation mode changes.
Definition: RooAbsCache.h:46
virtual void findConstantNodes(const RooArgSet &, RooArgSet &, RooLinkedList &)
Interface for constant term node finding calls.
Definition: RooAbsCache.h:52
virtual Bool_t redirectServersHook(const RooAbsCollection &, bool, bool, bool)
Interface for server redirect calls.
Definition: RooAbsCache.h:40
virtual void printCompactTreeHook(std::ostream &, const char *)
Interface for printing of cache guts in tree mode printing.
Definition: RooAbsCache.h:55
virtual void optimizeCacheMode(const RooArgSet &, RooArgSet &, RooLinkedList &)
Interface for processing of cache mode optimization calls.
Definition: RooAbsCache.h:49
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
Int_t getSize() const
void sort(Bool_t reverse=false)
Sort collection using std::sort and name comparison.
RooFIter fwdIterator() const
One-time forward iterator.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Storage_t::size_type size() const
RooAbsArg * first() const
void reserve(Storage_t::size_type count)
void clear()
Clear contents. If the collection is owning, it will also delete the contents.
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
RooAbsCollection * selectByAttrib(const char *name, Bool_t value) const
Create a subset of the current collection, consisting only of those elements with the specified attri...
const char * GetName() const
Returns name of object.
void setName(const char *name)
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsDataStore is the abstract base class for data collection that use a TTree as internal storage m...
virtual const RooArgSet * get(Int_t index) const =0
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:79
virtual const RooArgSet * get() const
Definition: RooAbsData.h:125
RooAbsProxy is the abstact interface for proxy classes.
Definition: RooAbsProxy.h:30
virtual void changeNormSet(const RooArgSet *newNormSet)
Destructor.
Definition: RooAbsProxy.cxx:64
virtual void print(std::ostream &os, Bool_t addContents=kFALSE) const
Print proxy name.
Definition: RooAbsProxy.cxx:74
virtual const char * name() const
Definition: RooAbsProxy.h:40
virtual Bool_t changePointer(const RooAbsCollection &newServerSet, Bool_t nameChange=kFALSE, Bool_t factoryInitMode=kFALSE)=0
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooArgProxy is the abstract interface for RooAbsArg proxy classes.
Definition: RooArgProxy.h:24
Bool_t isShapeServer() const
Definition: RooArgProxy.h:65
Bool_t isValueServer() const
Definition: RooArgProxy.h:61
RooAbsArg * absArg() const
Definition: RooArgProxy.h:37
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:35
RooConstVar represent a constant real-valued object.
Definition: RooConstVar.h:26
RooExpensiveObjectCache is a singleton class that serves as repository for objects that are expensive...
static RooExpensiveObjectCache & instance()
Return reference to singleton instance.
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
RooAbsArg * next()
Return next element or nullptr if at end.
Switches the message service to a different level while the instance is alive.
Definition: RooHelpers.h:42
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:38
TObject * FindObject(const char *name) const
Return pointer to obejct with given name.
RooAbsArg * findArg(const RooAbsArg *) const
Return pointer to object with given name in collection.
virtual void Add(TObject *arg)
Definition: RooLinkedList.h:67
RooListProxy is the concrete proxy for RooArgList objects.
Definition: RooListProxy.h:24
const TNamed * constPtr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:61
static RooNameReg & instance()
Return reference to singleton instance.
Definition: RooNameReg.cxx:51
RooPlotable is a 'mix-in' base class that define the standard RooFit plotting and printing methods.
Definition: RooPrintable.h:25
virtual StyleOption defaultPrintStyle(Option_t *opt) const
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,...
virtual void printValue(std::ostream &os) const
Interface to print value of object.
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
Bool_t isConvolved() const
std::size_t refCount(typename Container_t::const_iterator item) const
Return ref count of item that iterator points to.
const Container_t & containedObjects() const
Direct reference to container of objects held by this list.
void Remove(const T *obj, bool force=false)
Decrease ref count of given object.
Container_t::const_iterator begin() const
Iterator over contained objects.
void Add(T *obj, std::size_t initialCount=1)
Add an object or increase refCount if it is already present.
Container_t::const_iterator end() const
End of contained objects.
bool empty() const
Check if empty.
std::size_t size() const
Number of contained objects (neglecting the ref count).
void RemoveAll(const T *obj)
Remove from list irrespective of ref count.
void reserve(std::size_t amount)
bool containsByNamePtr(const T *obj) const
Check if list contains an item using findByNamePointer().
RooSetProxy is the concrete proxy for RooArgSet objects.
Definition: RooSetProxy.h:23
RooTreeDataStore is a TTree-backed data storage.
RooVectorDataStore uses std::vectors to store data columns.
Bool_t defineSetInternal(const char *name, const RooArgSet &aset)
const RooArgSet * set(const char *name)
Return pointer to previously defined named set with given nmame If no such set is found a null pointe...
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
virtual void SetByteCount(UInt_t cntpos, Bool_t packInVersion=kFALSE)=0
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0
virtual Int_t CheckByteCount(UInt_t startpos, UInt_t bcnt, const TClass *clss)=0
Bool_t IsReading() const
Definition: TBuffer.h:86
virtual UInt_t WriteVersion(const TClass *cl, Bool_t useBcnt=kFALSE)=0
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
TIterator and GenericRooFIter front end with STL back end.
Iterator abstract base class.
Definition: TIterator.h:30
virtual TObject * Next()=0
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
TString fName
Definition: TNamed.h:32
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual void SetNameTitle(const char *name, const char *title)
Set all the TNamed parameters (name and title).
Definition: TNamed.cxx:154
Int_t GetEntriesFast() const
Definition: TObjArray.h:64
virtual void Expand(Int_t newSize)
Expand or shrink the array to newSize elements.
Definition: TObjArray.cxx:387
void Add(TObject *obj)
Definition: TObjArray.h:74
virtual void Compress()
Remove empty slots from array.
Definition: TObjArray.cxx:334
Int_t GetEntries() const
Return the number of objects in array (i.e.
Definition: TObjArray.cxx:523
virtual TObject * FindObject(const char *name) const
Find an object in this collection using its name.
Definition: TObjArray.cxx:415
TIterator * MakeIterator(Bool_t dir=kIterForward) const
Returns an array iterator.
Definition: TObjArray.cxx:649
virtual TObject * Remove(TObject *obj)
Remove object from array.
Definition: TObjArray.cxx:719
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
Mother of all ROOT objects.
Definition: TObject.h:37
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:359
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:130
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
An array of references to TObjects.
Definition: TRefArray.h:39
Basic string class.
Definition: TString.h:136
const char * Data() const
Definition: TString.h:369
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
Bool_t BeginsWith(const char *s, ECaseCompare cmp=kExact) const
Definition: TString.h:615
TString & Append(const char *cs)
Definition: TString.h:564
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
A TTree represents a columnar dataset.
Definition: TTree.h:79
@ Optimization
Definition: RooGlobalFunc.h:61
@ InputArguments
Definition: RooGlobalFunc.h:61
@ Integration
Definition: RooGlobalFunc.h:60
@ LinkStateMgmt
Definition: RooGlobalFunc.h:60
std::string getColonSeparatedNameString(RooArgSet const &argSet)
Create a string with all sorted names of RooArgSet elements separated by colons.
Definition: RooHelpers.cxx:248
static constexpr double s
Definition: first.py:1
Definition: tree.py:1