Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooWorkspace.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 * *
8 * Copyright (c) 2000-2005, Regents of the University of California *
9 * and Stanford University. All rights reserved. *
10 * *
11 * Redistribution and use in source and binary forms, *
12 * with or without modification, are permitted according to the terms *
13 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
14 *****************************************************************************/
15
16/**
17\file RooWorkspace.cxx
18\class RooWorkspace
19\ingroup Roofitcore
20
21Persistable container for RooFit projects. A workspace
22can contain and own variables, p.d.f.s, functions and datasets. All objects
23that live in the workspace are owned by the workspace. The `import()` method
24enforces consistency of objects upon insertion into the workspace (e.g. no
25duplicate object with the same name are allowed) and makes sure all objects
26in the workspace are connected to each other. Easy accessor methods like
27`pdf()`, `var()` and `data()` allow to refer to the contents of the workspace by
28object name. The entire RooWorkspace can be saved into a ROOT TFile and organises
29the consistent streaming of its contents without duplication.
30If a RooWorkspace contains custom classes, i.e. classes not in the
31ROOT distribution, portability of workspaces can be enhanced by
32storing the source code of those classes in the workspace as well.
33This process is also organized by the workspace through the
34`importClassCode()` method.
35
36### Seemingly random crashes when reading large workspaces
37When reading or loading workspaces with deeply nested PDFs, one can encounter
38ouf-of-memory errors if the stack size is too small. This manifests in crashes
39at seemingly random locations, or in the process silently ending.
40Unfortunately, ROOT neither recover from this situation, nor warn or give useful
41instructions. When suspecting to have run out of stack memory, check
42```
43ulimit -s
44```
45and try reading again.
46**/
47
48#include <RooWorkspace.h>
49
50#include <RooAbsData.h>
51#include <RooAbsPdf.h>
52#include <RooAbsStudy.h>
53#include <RooCategory.h>
54#include <RooCmdConfig.h>
55#include <RooConstVar.h>
56#include <RooFactoryWSTool.h>
57#include <RooLinkedListIter.h>
58#include <RooMsgService.h>
59#include <RooPlot.h>
60#include <RooRandom.h>
61#include <RooRealVar.h>
62#include <RooResolutionModel.h>
63#include <RooTObjWrap.h>
64#include <RooWorkspaceHandle.h>
65
66#include "TBuffer.h"
67#include "TInterpreter.h"
68#include "TClassTable.h"
69#include "TBaseClass.h"
70#include "TSystem.h"
71#include "TRegexp.h"
72#include "TROOT.h"
73#include "TFile.h"
74#include "TH1.h"
75#include "TClass.h"
76#include "strlcpy.h"
77
78#ifdef ROOFIT_LEGACY_EVAL_BACKEND
80#endif
81
82#include "ROOT/StringUtils.hxx"
83
84#include <map>
85#include <sstream>
86#include <string>
87#include <iostream>
88#include <fstream>
89#include <cstring>
90
91namespace {
92
93// Infer from a RooArgSet name whether this set is used internally by
94// RooWorkspace to cache things.
95bool isCacheSet(std::string const& setName) {
96 // Check if the setName starts with CACHE_.
97 return setName.rfind("CACHE_", 0) == 0;
98}
99
100} // namespace
101
102using std::string, std::list, std::map, std::vector, std::ifstream, std::ofstream, std::fstream, std::make_unique;
103
104
105////////////////////////////////////////////////////////////////////////////////
106
107
108////////////////////////////////////////////////////////////////////////////////
109
110
113string RooWorkspace::_classFileExportDir = ".wscode.%s.%s" ;
115
116
117////////////////////////////////////////////////////////////////////////////////
118/// Add `dir` to search path for class declaration (header) files. This is needed
119/// to find class headers custom classes are imported into the workspace.
121{
122 _classDeclDirList.push_back(dir) ;
123}
124
125
126////////////////////////////////////////////////////////////////////////////////
127/// Add `dir` to search path for class implementation (.cxx) files. This is needed
128/// to find class headers custom classes are imported into the workspace.
130{
131 _classImplDirList.push_back(dir) ;
132}
133
134
135////////////////////////////////////////////////////////////////////////////////
136/// Specify the name of the directory in which embedded source
137/// code is unpacked and compiled. The specified string may contain
138/// one '%s' token which will be substituted by the workspace name
139
141{
142 if (dir) {
143 _classFileExportDir = dir ;
144 } else {
145 _classFileExportDir = ".wscode.%s.%s" ;
146 }
147}
148
149
150////////////////////////////////////////////////////////////////////////////////
151/// If flag is true, source code of classes not the ROOT distribution
152/// is automatically imported if on object of such a class is imported
153/// in the workspace
154
159
160
161
162////////////////////////////////////////////////////////////////////////////////
163/// Default constructor
164
166{
167}
168
169
170
171////////////////////////////////////////////////////////////////////////////////
172/// Construct empty workspace with given name and title
173
174RooWorkspace::RooWorkspace(const char* name, const char* title) :
175 TNamed(name,title?title:name), _classes(this)
176{
177}
178
179////////////////////////////////////////////////////////////////////////////////
180/// Construct empty workspace with given name and option to export reference to
181/// all workspace contents to a CINT namespace with the same name.
182
183RooWorkspace::RooWorkspace(const char* name, bool /*doCINTExport*/) :
184 TNamed(name,name), _classes(this)
185{
186}
187
188
189////////////////////////////////////////////////////////////////////////////////
190/// Workspace copy constructor
191
193 TNamed(other), _uuid(other._uuid), _classes(other._classes,this)
194{
195 // Copy owned nodes
196 other._allOwnedNodes.snapshot(_allOwnedNodes,true) ;
197
198 // Copy datasets
199 for(TObject *data2 : other._dataList) _dataList.Add(data2->Clone());
200
201 // Copy snapshots
202 for(auto * snap : static_range_cast<RooArgSet*>(other._snapshots)) {
203 auto snapClone = new RooArgSet;
204 snap->snapshot(*snapClone);
205 snapClone->setName(snap->GetName()) ;
207 }
208
209 // Copy named sets
210 for (map<string,RooArgSet>::const_iterator iter3 = other._namedSets.begin() ; iter3 != other._namedSets.end() ; ++iter3) {
211 // Make RooArgSet with equivalent content of this workspace
212 _namedSets[iter3->first].add(*std::unique_ptr<RooArgSet>{_allOwnedNodes.selectCommon(iter3->second)});
213 }
214
215 // Copy generic objects
216 for(TObject * gobj : other._genObjects) {
217 _genObjects.Add(gobj->Clone());
218 }
219
220 for(TObject * gobj : allGenericObjects()) {
221 if (auto handle = dynamic_cast<RooWorkspaceHandle*>(gobj)) {
222 handle->ReplaceWS(this);
223 }
224 }
225}
226
227
228/// TObject::Clone() needs to be overridden.
230{
231 auto out = new RooWorkspace{*this};
232 if(newname && std::string(newname) != GetName()) {
233 out->SetName(newname);
234 }
235 return out;
236}
237
238
239////////////////////////////////////////////////////////////////////////////////
240/// Workspace destructor
241
243{
244 // Delete contents
245 _dataList.Delete() ;
246 if (_dir) {
247 delete _dir ;
248 }
250
251 // WVE named sets too?
252
254
256 _views.Delete();
258
259}
260
261
262////////////////////////////////////////////////////////////////////////////////
263/// Import a RooAbsArg or RooAbsData set from a workspace in a file. Filespec should be constructed as "filename:wspacename:objectname"
264/// The arguments will be passed to the relevant import() or import(RooAbsData&, ...) import calls
265/// \note From python, use `Import()`, since `import` is a reserved keyword.
267 const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
268 const RooCmdArg& arg4, const RooCmdArg& arg5, const RooCmdArg& arg6,
269 const RooCmdArg& arg7, const RooCmdArg& arg8, const RooCmdArg& arg9)
270{
271 // Parse file/workspace/objectname specification
272 std::vector<std::string> tokens = ROOT::Split(fileSpec, ":");
273
274 // Check that parsing was successful
275 if (tokens.size() != 3) {
276 std::ostringstream stream;
277 for (const auto& token : tokens) {
278 stream << "\n\t" << token;
279 }
280 coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR in file specification, expecting 'filename:wsname:objname', but '" << fileSpec << "' given."
281 << "\nTokens read are:" << stream.str() << std::endl;
282 return true ;
283 }
284
285 const std::string& filename = tokens[0];
286 const std::string& wsname = tokens[1];
287 const std::string& objname = tokens[2];
288
289 // Check that file can be opened
290 std::unique_ptr<TFile> f{TFile::Open(filename.c_str())};
291 if (f==nullptr) {
292 coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR opening file " << filename << std::endl ;
293 return false;
294 }
295
296 // That that file contains workspace
297 RooWorkspace* w = dynamic_cast<RooWorkspace*>(f->Get(wsname.c_str())) ;
298 if (w==nullptr) {
299 coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR: No object named " << wsname << " in file " << filename
300 << " or object is not a RooWorkspace" << std::endl ;
301 return false;
302 }
303
304 // Check that workspace contains object and forward to appropriate import method
305 RooAbsArg* warg = w->arg(objname.c_str()) ;
306 if (warg) {
307 bool ret = import(*warg,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9) ;
308 return ret ;
309 }
310 RooAbsData* wdata = w->data(objname.c_str()) ;
311 if (wdata) {
312 bool ret = import(*wdata,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9) ;
313 return ret ;
314 }
315
316 coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR: No RooAbsArg or RooAbsData object named " << objname
317 << " in workspace " << wsname << " in file " << filename << std::endl ;
318 return true ;
319}
320
321
322////////////////////////////////////////////////////////////////////////////////
323/// Import multiple RooAbsArg objects into workspace. For details on arguments see documentation
324/// of import() method for single RooAbsArg
325/// \note From python, use `Import()`, since `import` is a reserved keyword.
327 const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
328 const RooCmdArg& arg4, const RooCmdArg& arg5, const RooCmdArg& arg6,
329 const RooCmdArg& arg7, const RooCmdArg& arg8, const RooCmdArg& arg9)
330{
331 bool ret(false) ;
332 for(RooAbsArg * oneArg : args) {
334 }
335 return ret ;
336}
337
338
339
340////////////////////////////////////////////////////////////////////////////////
341/// Import a RooAbsArg object, e.g. function, p.d.f or variable into the workspace. This import function clones the input argument and will
342/// own the clone. If a composite object is offered for import, e.g. a p.d.f with parameters and observables, the
343/// complete tree of objects is imported. If any of the _variables_ of a composite object (parameters/observables) are already
344/// in the workspace the imported p.d.f. is connected to the already existing variables. If any of the _function_ objects (p.d.f, formulas)
345/// to be imported already exists in the workspace an error message is printed and the import of the entire tree of objects is cancelled.
346/// Several optional arguments can be provided to modify the import procedure.
347///
348/// <table>
349/// <tr><th> Accepted arguments
350/// <tr><td> `RenameConflictNodes(const char* suffix)` <td> Add suffix to branch node name if name conflicts with existing node in workspace
351/// <tr><td> `RenameAllNodes(const char* suffix)` <td> Add suffix to all branch node names including top level node.
352/// <tr><td> `RenameAllVariables(const char* suffix)` <td> Add suffix to all variables of objects being imported.
353/// <tr><td> `RenameAllVariablesExcept(const char* suffix, const char* exceptionList)` <td> Add suffix to all variables names, except ones listed
354/// <tr><td> `RenameVariable(const char* inputName, const char* outputName)` <td> Rename a single variable as specified upon import.
355/// <tr><td> `RecycleConflictNodes()` <td> If any of the function objects to be imported already exist in the name space, connect the
356/// imported expression to the already existing nodes.
357/// \attention Use with care! If function definitions do not match, this alters the definition of your function upon import
358///
359/// <tr><td> `Silence()` <td> Do not issue any info message
360/// </table>
361///
362/// The RenameConflictNodes, RenameNodes and RecycleConflictNodes arguments are mutually exclusive. The RenameVariable argument can be repeated
363/// as often as necessary to rename multiple variables. Alternatively, a single RenameVariable argument can be given with
364/// two comma separated lists.
365/// \note From python, use `Import()`, since `import` is a reserved keyword.
367 const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
368 const RooCmdArg& arg4, const RooCmdArg& arg5, const RooCmdArg& arg6,
369 const RooCmdArg& arg7, const RooCmdArg& arg8, const RooCmdArg& arg9)
370{
371 RooLinkedList args ;
372 args.Add((TObject*)&arg1) ;
373 args.Add((TObject*)&arg2) ;
374 args.Add((TObject*)&arg3) ;
375 args.Add((TObject*)&arg4) ;
376 args.Add((TObject*)&arg5) ;
377 args.Add((TObject*)&arg6) ;
378 args.Add((TObject*)&arg7) ;
379 args.Add((TObject*)&arg8) ;
380 args.Add((TObject*)&arg9) ;
381
382 // Select the pdf-specific commands
383 RooCmdConfig pc("RooWorkspace::import(" + std::string(GetName()) + ")");
384
385 pc.defineString("conflictSuffix","RenameConflictNodes",0) ;
386 pc.defineInt("renameConflictOrig","RenameConflictNodes",0,0) ;
387 pc.defineString("allSuffix","RenameAllNodes",0) ;
388 pc.defineString("allVarsSuffix","RenameAllVariables",0) ;
389 pc.defineString("allVarsExcept","RenameAllVariables",1) ;
390 pc.defineString("varChangeIn","RenameVar",0,"",true) ;
391 pc.defineString("varChangeOut","RenameVar",1,"",true) ;
392 pc.defineString("factoryTag","FactoryTag",0) ;
393 pc.defineInt("useExistingNodes","RecycleConflictNodes",0,0) ;
394 pc.defineInt("silence","Silence",0,0) ;
395 pc.defineInt("noRecursion","NoRecursion",0,0) ;
396 pc.defineMutex("RenameConflictNodes","RenameAllNodes") ;
397 pc.defineMutex("RenameConflictNodes","RecycleConflictNodes") ;
398 pc.defineMutex("RenameAllNodes","RecycleConflictNodes") ;
399 pc.defineMutex("RenameVariable","RenameAllVariables") ;
400
401 // Process and check varargs
402 pc.process(args) ;
403 if (!pc.ok(true)) {
404 return true ;
405 }
406
407 // Decode renaming logic into suffix string and boolean for conflictOnly mode
408 const char* suffixC = pc.getString("conflictSuffix") ;
409 const char* suffixA = pc.getString("allSuffix") ;
410 const char* suffixV = pc.getString("allVarsSuffix") ;
411 const char* exceptVars = pc.getString("allVarsExcept") ;
412 const char* varChangeIn = pc.getString("varChangeIn") ;
413 const char* varChangeOut = pc.getString("varChangeOut") ;
414 bool renameConflictOrig = pc.getInt("renameConflictOrig") ;
415 Int_t useExistingNodes = pc.getInt("useExistingNodes") ;
416 Int_t silence = pc.getInt("silence") ;
417 Int_t noRecursion = pc.getInt("noRecursion") ;
418
419
420 // Turn zero length strings into null pointers
421 if (suffixC && strlen(suffixC)==0) suffixC = nullptr ;
422 if (suffixA && strlen(suffixA)==0) suffixA = nullptr ;
423
424 bool conflictOnly = suffixA ? false : true ;
425 const char* suffix = suffixA ? suffixA : suffixC ;
426
427 // Process any change in variable names
428 std::map<string,string> varMap ;
429 if (strlen(varChangeIn)>0) {
430
431 // Parse comma separated lists into map<string,string>
432 const std::vector<std::string> tokIn = ROOT::Split(varChangeIn, ", ", /*skipEmpty= */ true);
433 const std::vector<std::string> tokOut = ROOT::Split(varChangeOut, ", ", /*skipEmpty= */ true);
434 for (unsigned int i=0; i < tokIn.size(); ++i) {
435 varMap.insert(std::make_pair(tokIn[i], tokOut[i]));
436 }
437
438 assert(tokIn.size() == tokOut.size());
439 }
440
441 // Process RenameAllVariables argument if specified
442 // First convert exception list if provided
443 std::set<string> exceptVarNames ;
444 if (exceptVars && strlen(exceptVars)) {
445 const std::vector<std::string> toks = ROOT::Split(exceptVars, ", ", /*skipEmpty= */ true);
446 exceptVarNames.insert(toks.begin(), toks.end());
447 }
448
449 if (suffixV != nullptr && strlen(suffixV)>0) {
450 std::unique_ptr<RooArgSet> vars{inArg.getVariables()};
451 for (const auto v : *vars) {
452 if (exceptVarNames.find(v->GetName())==exceptVarNames.end()) {
453 varMap[v->GetName()] = Form("%s_%s",v->GetName(),suffixV) ;
454 }
455 }
456 }
457
458 // Scan for overlaps with current contents
459 RooAbsArg* wsarg = _allOwnedNodes.find(inArg.GetName()) ;
460
461 // Check for factory specification match
462 const char* tagIn = inArg.getStringAttribute("factory_tag") ;
463 const char* tagWs = wsarg ? wsarg->getStringAttribute("factory_tag") : nullptr ;
464 bool factoryMatch = (tagIn && tagWs && !strcmp(tagIn,tagWs)) ;
465 if (factoryMatch) {
466 ((RooAbsArg&)inArg).setAttribute("RooWorkspace::Recycle") ;
467 }
468
469 if (!suffix && wsarg && !useExistingNodes && !(inArg.isFundamental() && !varMap[inArg.GetName()].empty())) {
470 if (!factoryMatch) {
471 if (wsarg!=&inArg) {
472 coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR importing object named " << inArg.GetName()
473 << ": another instance with same name already in the workspace and no conflict resolution protocol specified" << std::endl ;
474 return true ;
475 } else {
476 if (!silence) {
477 coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") Object " << inArg.GetName() << " is already in workspace!" << std::endl ;
478 }
479 return true ;
480 }
481 } else {
482 if(!silence) {
483 coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") Recycling existing object " << inArg.GetName() << " created with identical factory specification" << std::endl ;
484 }
485 }
486 }
487
488 // Make list of conflicting nodes
491 if (noRecursion) {
492 branchSet.add(inArg) ;
493 } else {
494 inArg.branchNodeServerList(&branchSet) ;
495 }
496
497 for (const auto branch : branchSet) {
499 if (wsbranch && wsbranch!=branch && !branch->getAttribute("RooWorkspace::Recycle") && !useExistingNodes) {
500 conflictNodes.add(*branch) ;
501 }
502 }
503
504 // Terminate here if there are conflicts and no resolution protocol
505 if (!conflictNodes.empty() && !suffix && !useExistingNodes) {
506 coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR object named " << inArg.GetName() << ": component(s) "
507 << conflictNodes << " already in the workspace and no conflict resolution protocol specified" << std::endl ;
508 return true ;
509 }
510
511 // Now create a working copy of the incoming object tree
513 cloneSet.useHashMapForFind(true); // Accelerate finding
515 RooAbsArg* cloneTop = cloneSet.find(inArg.GetName()) ;
516
517 // Mark all nodes for renaming if we are not in conflictOnly mode
518 if (!conflictOnly) {
519 conflictNodes.removeAll() ;
521 }
522
523 // Mark nodes that are to be renamed with special attribute
524 string topName2 = cloneTop->GetName() ;
525 if (!renameConflictOrig) {
526 // Mark all nodes to be imported for renaming following conflict resolution protocol
527 for (const auto cnode : conflictNodes) {
528 RooAbsArg* cnode2 = cloneSet.find(cnode->GetName()) ;
529 string origName = cnode2->GetName() ;
530 cnode2->SetName(Form("%s_%s",cnode2->GetName(),suffix)) ;
531 cnode2->SetTitle(Form("%s (%s)",cnode2->GetTitle(),suffix)) ;
532 string tag = Form("ORIGNAME:%s",origName.c_str()) ;
533 cnode2->setAttribute(tag.c_str()) ;
534 if (!cnode2->getStringAttribute("origName")) {
535 cnode2->setStringAttribute("origName",origName.c_str());
536 }
537
538 // Save name of new top level node for later use
539 if (cnode2==cloneTop) {
540 topName2 = cnode2->GetName() ;
541 }
542
543 if (!silence) {
544 coutI(ObjectHandling) << "RooWorkspace::import(" << GetName()
545 << ") Resolving name conflict in workspace by changing name of imported node "
546 << origName << " to " << cnode2->GetName() << std::endl ;
547 }
548 }
549 } else {
550
551 // Rename all nodes already in the workspace to 'clear the way' for the imported nodes
552 for (const auto cnode : conflictNodes) {
553
554 string origName = cnode->GetName() ;
556 if (wsnode) {
557
558 if (!wsnode->getStringAttribute("origName")) {
559 wsnode->setStringAttribute("origName",wsnode->GetName()) ;
560 }
561
562 if (!_allOwnedNodes.find(Form("%s_%s",cnode->GetName(),suffix))) {
563 wsnode->SetName(Form("%s_%s",cnode->GetName(),suffix)) ;
564 wsnode->SetTitle(Form("%s (%s)",cnode->GetTitle(),suffix)) ;
565 } else {
566 // Name with suffix already taken, add additional suffix
567 for (unsigned int n=1; true; ++n) {
568 string newname = Form("%s_%s_%d",cnode->GetName(),suffix,n) ;
569 if (!_allOwnedNodes.find(newname.c_str())) {
570 wsnode->SetName(newname.c_str()) ;
571 wsnode->SetTitle(Form("%s (%s %d)",cnode->GetTitle(),suffix,n)) ;
572 break ;
573 }
574 }
575 }
576 if (!silence) {
577 coutI(ObjectHandling) << "RooWorkspace::import(" << GetName()
578 << ") Resolving name conflict in workspace by changing name of original node "
579 << origName << " to " << wsnode->GetName() << std::endl ;
580 }
581 } else {
582 coutW(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") Internal error: expected to find existing node "
583 << origName << " to be renamed, but didn't find it..." << std::endl ;
584 }
585
586 }
587 }
588
589 // Process any change in variable names
590 if (strlen(varChangeIn)>0 || (suffixV && strlen(suffixV)>0)) {
591
592 // Process all changes in variable names
593 for (const auto cnode : cloneSet) {
594
595 if (varMap.find(cnode->GetName())!=varMap.end()) {
596 string origName = cnode->GetName() ;
597 cnode->SetName(varMap[cnode->GetName()].c_str()) ;
598 string tag = Form("ORIGNAME:%s",origName.c_str()) ;
599 cnode->setAttribute(tag.c_str()) ;
600 if (!cnode->getStringAttribute("origName")) {
601 cnode->setStringAttribute("origName",origName.c_str()) ;
602 }
603
604 if (!silence) {
605 coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") Changing name of variable "
606 << origName << " to " << cnode->GetName() << " on request" << std::endl ;
607 }
608
609 if (cnode==cloneTop) {
610 topName2 = cnode->GetName() ;
611 }
612
613 }
614 }
615 }
616
617 // Now clone again with renaming effective
619 cloneSet2.useHashMapForFind(true); // Faster finding
621 RooAbsArg* cloneTop2 = cloneSet2.find(topName2.c_str()) ;
622
623 // Make final check list of conflicting nodes
626 for (const auto branch2 : branchSet2) {
627 if (_allOwnedNodes.find(branch2->GetName())) {
628 conflictNodes2.add(*branch2) ;
629 }
630 }
631
632 // Terminate here if there are conflicts and no resolution protocol
633 if (!conflictNodes2.empty()) {
634 coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR object named " << inArg.GetName() << ": component(s) "
635 << conflictNodes2 << " cause naming conflict after conflict resolution protocol was executed" << std::endl ;
636 return true ;
637 }
638
639 // Perform any auxiliary imports at this point
640 for (const auto node : cloneSet2) {
641 if (node->importWorkspaceHook(*this)) {
642 coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR object named " << node->GetName()
643 << " has an error in importing in one or more of its auxiliary objects, aborting" << std::endl ;
644 return true ;
645 }
646 }
647
650 for (const auto node : cloneSet2) {
651 if (_autoClass) {
652 if (!_classes.autoImportClass(node->IsA())) {
653 coutW(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") WARNING: problems import class code of object "
654 << node->ClassName() << "::" << node->GetName() << ", reading of workspace will require external definition of class" << std::endl ;
655 }
656 }
657
658 // Point expensiveObjectCache to copy in this workspace
659 RooExpensiveObjectCache& oldCache = node->expensiveObjectCache() ;
660 node->setExpensiveObjectCache(_eocache) ;
661 _eocache.importCacheObjects(oldCache,node->GetName(),true) ;
662
663 // Check if node is already in workspace (can only happen for variables or identical instances, unless RecycleConflictNodes is specified)
664 RooAbsArg* wsnode = _allOwnedNodes.find(node->GetName()) ;
665
666 if (wsnode) {
667 // Do not import node, add not to list of nodes that require reconnection
668 if (!silence && useExistingNodes) {
669 coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") using existing copy of " << node->ClassName()
670 << "::" << node->GetName() << " for import of " << cloneTop2->ClassName() << "::"
671 << cloneTop2->GetName() << std::endl ;
672 }
673 recycledNodes.add(*_allOwnedNodes.find(node->GetName())) ;
674
675 // Delete clone of incoming node
676 nodesToBeDeleted.addOwned(std::unique_ptr<RooAbsArg>{node});
677
678 //cout << "WV: recycling existing node " << existingNode << " = " << existingNode->GetName() << " for imported node " << node << std::endl ;
679
680 } else {
681 // Import node
682 if (!silence) {
683 coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") importing " << node->ClassName() << "::"
684 << node->GetName() << std::endl ;
685 }
686 _allOwnedNodes.addOwned(std::unique_ptr<RooAbsArg>{node});
687 node->setWorkspace(*this);
688 if (_openTrans) {
689 _sandboxNodes.add(*node) ;
690 } else {
691 if (_dir && node->IsA() != RooConstVar::Class()) {
692 _dir->InternalAppend(node) ;
693 }
694 }
695 }
696 }
697
698 // Reconnect any nodes that need to be
699 if (!recycledNodes.empty()) {
700 for (const auto node : cloneSet2) {
701 node->redirectServers(recycledNodes) ;
702 }
703 }
704
705 cloneSet2.releaseOwnership() ;
706
707 return false ;
708}
709
710
711
712////////////////////////////////////////////////////////////////////////////////
713/// Import a dataset (RooDataSet or RooDataHist) into the workspace. The workspace will contain a copy of the data.
714/// The dataset and its variables can be renamed upon insertion with the options below
715///
716/// <table>
717/// <tr><th> Accepted arguments
718/// <tr><td> `Rename(const char* suffix)` <td> Rename dataset upon insertion
719/// <tr><td> `RenameVariable(const char* inputName, const char* outputName)` <td> Change names of observables in dataset upon insertion
720/// <tr><td> `Silence` <td> Be quiet, except in case of errors
721/// \note From python, use `Import()`, since `import` is a reserved keyword.
723 const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
724 const RooCmdArg& arg4, const RooCmdArg& arg5, const RooCmdArg& arg6,
725 const RooCmdArg& arg7, const RooCmdArg& arg8, const RooCmdArg& arg9)
726
727{
728
729 RooLinkedList args ;
730 args.Add((TObject*)&arg1) ;
731 args.Add((TObject*)&arg2) ;
732 args.Add((TObject*)&arg3) ;
733 args.Add((TObject*)&arg4) ;
734 args.Add((TObject*)&arg5) ;
735 args.Add((TObject*)&arg6) ;
736 args.Add((TObject*)&arg7) ;
737 args.Add((TObject*)&arg8) ;
738 args.Add((TObject*)&arg9) ;
739
740 // Select the pdf-specific commands
741 RooCmdConfig pc(Form("RooWorkspace::import(%s)",GetName())) ;
742
743 pc.defineString("dsetName","Rename",0,"") ;
744 pc.defineString("varChangeIn","RenameVar",0,"",true) ;
745 pc.defineString("varChangeOut","RenameVar",1,"",true) ;
746 pc.defineInt("embedded","Embedded",0,0) ;
747 pc.defineInt("silence","Silence",0,0) ;
748
749 // Process and check varargs
750 pc.process(args) ;
751 if (!pc.ok(true)) {
752 return true ;
753 }
754
755 // Decode renaming logic into suffix string and boolean for conflictOnly mode
756 const char* dsetName = pc.getString("dsetName") ;
757 const char* varChangeIn = pc.getString("varChangeIn") ;
758 const char* varChangeOut = pc.getString("varChangeOut") ;
759 bool embedded = pc.getInt("embedded") ;
760 Int_t silence = pc.getInt("silence") ;
761
762 if (!silence)
763 coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") importing dataset " << inData.GetName() << std::endl ;
764
765 // Transform empty string into null pointer
766 if (dsetName && strlen(dsetName)==0) {
767 dsetName=nullptr ;
768 }
769
771 if (dataList.size() > 50 && dataList.getHashTableSize() == 0) {
772 // When the workspaces get larger, traversing the linked list becomes a bottleneck:
773 dataList.setHashTableSize(200);
774 }
775
776 // Check that no dataset with target name already exists
777 if (dsetName && dataList.FindObject(dsetName)) {
778 coutE(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") ERROR dataset with name " << dsetName << " already exists in workspace, import aborted" << std::endl ;
779 return true ;
780 }
781 if (!dsetName && dataList.FindObject(inData.GetName())) {
782 coutE(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") ERROR dataset with name " << inData.GetName() << " already exists in workspace, import aborted" << std::endl ;
783 return true ;
784 }
785
786 // Rename dataset if required
787 RooAbsData* clone ;
788 if (dsetName) {
789 if (!silence)
790 coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") changing name of dataset from " << inData.GetName() << " to " << dsetName << std::endl ;
791 clone = static_cast<RooAbsData*>(inData.Clone(dsetName)) ;
792 } else {
793 clone = static_cast<RooAbsData*>(inData.Clone(inData.GetName())) ;
794 }
795
796
797 // Process any change in variable names
798 if (strlen(varChangeIn)>0) {
799 // Parse comma separated lists of variable name changes
800 const std::vector<std::string> tokIn = ROOT::Split(varChangeIn, ",");
801 const std::vector<std::string> tokOut = ROOT::Split(varChangeOut, ",");
802 for (unsigned int i=0; i < tokIn.size(); ++i) {
803 if (!silence)
804 coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") changing name of dataset observable " << tokIn[i] << " to " << tokOut[i] << std::endl ;
805 clone->changeObservableName(tokIn[i].c_str(), tokOut[i].c_str());
806 }
807 }
808
809 // Now import the dataset observables, unless dataset is embedded
810 if (!embedded) {
811 for(RooAbsArg* carg : *clone->get()) {
812 if (!arg(carg->GetName())) {
813 import(*carg) ;
814 }
815 }
816 }
817
818 dataList.Add(clone) ;
819 if (_dir) {
820 _dir->InternalAppend(clone) ;
821 }
822
823 // Set expensive object cache of dataset internal buffers to that of workspace
824 for(RooAbsArg* carg : *clone->get()) {
825 carg->setExpensiveObjectCache(expensiveObjectCache()) ;
826 }
827
828
829 return false ;
830}
831
832
833
834
835////////////////////////////////////////////////////////////////////////////////
836/// Define a named RooArgSet with given constituents. If importMissing is true, any constituents
837/// of aset that are not in the workspace will be imported, otherwise an error is returned
838/// for missing components
839
841{
842 // Check if set was previously defined, if so print warning
843 map<string,RooArgSet>::iterator i = _namedSets.find(name) ;
844 if (i!=_namedSets.end()) {
845 coutW(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") WARNING redefining previously defined named set " << name << std::endl ;
846 }
847
849
850 // Check all constituents of provided set
851 for (RooAbsArg* sarg : aset) {
852 // If missing, either import or report error
853 if (!arg(sarg->GetName())) {
854 if (importMissing) {
855 import(*sarg) ;
856 } else {
857 coutE(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") ERROR set constituent \"" << sarg->GetName()
858 << "\" is not in workspace and importMissing option is disabled" << std::endl ;
859 return true ;
860 }
861 }
862 wsargs.add(*arg(sarg->GetName())) ;
863 }
864
865
866 // Install named set
867 _namedSets[name].removeAll() ;
868 _namedSets[name].add(wsargs) ;
869
870 return false ;
871}
872
873//_____________________________________________________________________________
875{
876 // Define a named RooArgSet with given constituents. If importMissing is true, any constituents
877 // of aset that are not in the workspace will be imported, otherwise an error is returned
878 // for missing components
879
880 // Check if set was previously defined, if so print warning
881 map<string, RooArgSet>::iterator i = _namedSets.find(name);
882 if (i != _namedSets.end()) {
883 coutW(InputArguments) << "RooWorkspace::defineSet(" << GetName()
884 << ") WARNING redefining previously defined named set " << name << std::endl;
885 }
886
887 // Install named set
888 _namedSets[name].removeAll();
889 _namedSets[name].add(aset);
890
891 return false;
892}
893
894////////////////////////////////////////////////////////////////////////////////
895/// Define a named set in the workspace through a comma separated list of
896/// names of objects already in the workspace
897
898bool RooWorkspace::defineSet(const char* name, const char* contentList)
899{
900 // Check if set was previously defined, if so print warning
901 map<string,RooArgSet>::iterator i = _namedSets.find(name) ;
902 if (i!=_namedSets.end()) {
903 coutW(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") WARNING redefining previously defined named set " << name << std::endl ;
904 }
905
907
908 // Check all constituents of provided set
909 for (const std::string& token : ROOT::Split(contentList, ",")) {
910 // If missing, either import or report error
911 if (!arg(token.c_str())) {
912 coutE(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") ERROR proposed set constituent \"" << token
913 << "\" is not in workspace" << std::endl ;
914 return true ;
915 }
916 wsargs.add(*arg(token.c_str())) ;
917 }
918
919 // Install named set
920 _namedSets[name].removeAll() ;
921 _namedSets[name].add(wsargs) ;
922
923 return false ;
924}
925
926
927
928
929////////////////////////////////////////////////////////////////////////////////
930/// Define a named set in the workspace through a comma separated list of
931/// names of objects already in the workspace
932
933bool RooWorkspace::extendSet(const char* name, const char* newContents)
934{
936
937 // Check all constituents of provided set
938 for (const std::string& token : ROOT::Split(newContents, ",")) {
939 // If missing, either import or report error
940 if (!arg(token.c_str())) {
941 coutE(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") ERROR proposed set constituent \"" << token
942 << "\" is not in workspace" << std::endl ;
943 return true ;
944 }
945 wsargs.add(*arg(token.c_str())) ;
946 }
947
948 // Extend named set
949 _namedSets[name].add(wsargs,true) ;
950
951 return false ;
952}
953
954
955
956////////////////////////////////////////////////////////////////////////////////
957/// Return pointer to previously defined named set with given nmame
958/// If no such set is found a null pointer is returned
959
961{
962 std::map<string,RooArgSet>::iterator i = _namedSets.find(name.c_str());
963 return (i!=_namedSets.end()) ? &(i->second) : nullptr;
964}
965
966
967
968
969////////////////////////////////////////////////////////////////////////////////
970/// Rename set to a new name
971
972bool RooWorkspace::renameSet(const char* name, const char* newName)
973{
974 // First check if set exists
975 if (!set(name)) {
976 coutE(InputArguments) << "RooWorkspace::renameSet(" << GetName() << ") ERROR a set with name " << name
977 << " does not exist" << std::endl ;
978 return true ;
979 }
980
981 // Check if no set exists with new name
982 if (set(newName)) {
983 coutE(InputArguments) << "RooWorkspace::renameSet(" << GetName() << ") ERROR a set with name " << newName
984 << " already exists" << std::endl ;
985 return true ;
986 }
987
988 // Copy entry under 'name' to 'newName'
990
991 // Remove entry under old name
992 _namedSets.erase(name) ;
993
994 return false ;
995}
996
997
998
999
1000////////////////////////////////////////////////////////////////////////////////
1001/// Remove a named set from the workspace
1002
1004{
1005 // First check if set exists
1006 if (!set(name)) {
1007 coutE(InputArguments) << "RooWorkspace::removeSet(" << GetName() << ") ERROR a set with name " << name
1008 << " does not exist" << std::endl ;
1009 return true ;
1010 }
1011
1012 // Remove set with given name
1013 _namedSets.erase(name) ;
1014
1015 return false ;
1016}
1017
1018
1019
1020
1021////////////////////////////////////////////////////////////////////////////////
1022/// Open an import transaction operations. Returns true if successful, false
1023/// if there is already an ongoing transaction
1024
1026{
1027 // Check that there was no ongoing transaction
1028 if (_openTrans) {
1029 return false ;
1030 }
1031
1032 // Open transaction
1033 _openTrans = true ;
1034 return true ;
1035}
1036
1037
1038
1039
1040////////////////////////////////////////////////////////////////////////////////
1041/// Cancel an ongoing import transaction. All objects imported since startTransaction()
1042/// will be removed and the transaction will be terminated. Return true if cancel operation
1043/// succeeds, return false if there was no open transaction
1044
1046{
1047 // Check that there is an ongoing transaction
1048 if (!_openTrans) {
1049 return false ;
1050 }
1051
1052 // Delete all objects in the sandbox
1053 for(RooAbsArg * tmpArg : _sandboxNodes) {
1055 }
1057
1058 // Mark transaction as finished
1059 _openTrans = false ;
1060
1061 return true ;
1062}
1063
1065{
1066 // Commit an ongoing import transaction. Returns true if commit succeeded,
1067 // return false if there was no ongoing transaction
1068
1069 // Check that there is an ongoing transaction
1070 if (!_openTrans) {
1071 return false ;
1072 }
1073
1074 // Publish sandbox nodes in directory and/or CINT if requested
1075 for(RooAbsArg* sarg : _sandboxNodes) {
1076 if (_dir && sarg->IsA() != RooConstVar::Class()) {
1078 }
1079 }
1080
1081 // Remove all committed objects from the sandbox
1083
1084 // Mark transaction as finished
1085 _openTrans = false ;
1086
1087 return true ;
1088}
1089
1090
1091
1092
1093////////////////////////////////////////////////////////////////////////////////
1094
1099
1100
1101
1102////////////////////////////////////////////////////////////////////////////////
1103/// Import code of all classes in the workspace that have a class name
1104/// that matches pattern 'pat' and which are not found to be part of
1105/// the standard ROOT distribution. If doReplace is true any existing
1106/// class code saved in the workspace is replaced
1107
1109{
1110 bool ret(true) ;
1111
1112 TRegexp re(pat,true) ;
1113 for (RooAbsArg * carg : _allOwnedNodes) {
1114 TString className = carg->ClassName() ;
1115 if (className.Index(re)>=0 && !_classes.autoImportClass(carg->IsA(),doReplace)) {
1116 coutW(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") WARNING: problems import class code of object "
1117 << carg->ClassName() << "::" << carg->GetName() << ", reading of workspace will require external definition of class" << std::endl ;
1118 ret = false ;
1119 }
1120 }
1121
1122 return ret ;
1123}
1124
1125
1126
1127
1128
1129////////////////////////////////////////////////////////////////////////////////
1130/// Save snapshot of values and attributes (including "Constant") of given parameters.
1131/// \param[in] name Name of the snapshot.
1132/// \param[in] paramNames Comma-separated list of parameter names to be snapshot.
1134{
1135 return saveSnapshot(name,argSet(paramNames),false) ;
1136}
1137
1138
1139
1140
1141
1142////////////////////////////////////////////////////////////////////////////////
1143/// Save snapshot of values and attributes (including "Constant") of parameters 'params'.
1144/// If importValues is FALSE, the present values from the object in the workspace are
1145/// saved. If importValues is TRUE, the values of the objects passed in the 'params'
1146/// argument are saved
1147
1149{
1152 auto snapshot = new RooArgSet;
1153 actualParams.snapshot(*snapshot);
1154
1155 snapshot->setName(name.c_str()) ;
1156
1157 if (importValues) {
1158 snapshot->assign(params) ;
1159 }
1160
1161 if (std::unique_ptr<RooArgSet> oldSnap{static_cast<RooArgSet*>(_snapshots.FindObject(name.c_str()))}) {
1162 coutI(ObjectHandling) << "RooWorkspace::saveSnapshot(" << GetName() << ") replacing previous snapshot with name " << name << std::endl ;
1163 _snapshots.Remove(oldSnap.get()) ;
1164 }
1165
1166 _snapshots.Add(snapshot) ;
1167
1168 return true ;
1169}
1170
1171
1172
1173
1174////////////////////////////////////////////////////////////////////////////////
1175/// Load the values and attributes of the parameters in the snapshot saved with
1176/// the given name
1177
1179{
1180 RooArgSet* snap = static_cast<RooArgSet*>(_snapshots.find(name)) ;
1181 if (!snap) {
1182 coutE(ObjectHandling) << "RooWorkspace::loadSnapshot(" << GetName() << ") no snapshot with name " << name << " is available" << std::endl ;
1183 return false ;
1184 }
1185
1188 actualParams.assign(*snap) ;
1189
1190 return true ;
1191}
1192
1193
1194////////////////////////////////////////////////////////////////////////////////
1195/// Return the RooArgSet containing a snapshot of variables contained in the workspace
1196///
1197/// Note that the variables of the objects in the snapshots are **copies** of the
1198/// variables in the workspace. To load the values of a snapshot in the workspace
1199/// variables, use loadSnapshot() instead.
1200
1202{
1203 return static_cast<RooArgSet*>(_snapshots.find(name));
1204}
1205
1206
1207////////////////////////////////////////////////////////////////////////////////
1208/// Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found
1209
1211{
1212 return dynamic_cast<RooAbsPdf*>(_allOwnedNodes.find(name.c_str())) ;
1213}
1214
1215
1216////////////////////////////////////////////////////////////////////////////////
1217/// Retrieve function (RooAbsReal) with given name. Note that all RooAbsPdfs are also RooAbsReals. A null pointer is returned if not found.
1218
1220{
1221 return dynamic_cast<RooAbsReal*>(_allOwnedNodes.find(name.c_str())) ;
1222}
1223
1224
1225////////////////////////////////////////////////////////////////////////////////
1226/// Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found
1227
1229{
1230 return dynamic_cast<RooRealVar*>(_allOwnedNodes.find(name.c_str())) ;
1231}
1232
1233
1234////////////////////////////////////////////////////////////////////////////////
1235/// Retrieve discrete variable (RooCategory) with given name. A null pointer is returned if not found
1236
1238{
1239 return dynamic_cast<RooCategory*>(_allOwnedNodes.find(name.c_str())) ;
1240}
1241
1242
1243////////////////////////////////////////////////////////////////////////////////
1244/// Retrieve discrete function (RooAbsCategory) with given name. A null pointer is returned if not found
1245
1247{
1248 return dynamic_cast<RooAbsCategory*>(_allOwnedNodes.find(name.c_str())) ;
1249}
1250
1251
1252
1253////////////////////////////////////////////////////////////////////////////////
1254/// Return RooAbsArg with given name. A null pointer is returned if none is found.
1255
1257{
1258 return _allOwnedNodes.find(name.c_str()) ;
1259}
1260
1261
1262
1263////////////////////////////////////////////////////////////////////////////////
1264/// Return set of RooAbsArgs matching to given list of names
1265
1267{
1268 RooArgSet ret ;
1269
1270 for (const std::string& token : ROOT::Split(nameList, ",")) {
1271 RooAbsArg* oneArg = arg(token.c_str()) ;
1272 if (oneArg) {
1273 ret.add(*oneArg) ;
1274 } else {
1275 std::stringstream ss;
1276 ss << " RooWorkspace::argSet(" << GetName() << ") no RooAbsArg named \"" << token << "\" in workspace" ;
1277 const std::string errorMsg = ss.str();
1278 coutE(InputArguments) << errorMsg << std::endl;
1279 throw std::runtime_error(errorMsg);
1280 }
1281 }
1282 return ret ;
1283}
1284
1285
1286
1287////////////////////////////////////////////////////////////////////////////////
1288/// Return fundamental (i.e. non-derived) RooAbsArg with given name. Fundamental types
1289/// are e.g. RooRealVar, RooCategory. A null pointer is returned if none is found.
1290
1292{
1293 RooAbsArg* tmp = arg(name) ;
1294 if (!tmp) {
1295 return nullptr;
1296 }
1297 return tmp->isFundamental() ? tmp : nullptr;
1298}
1299
1300
1301
1302////////////////////////////////////////////////////////////////////////////////
1303/// Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found
1304
1306{
1307 return static_cast<RooAbsData*>(_dataList.FindObject(name.c_str())) ;
1308}
1309
1310
1311////////////////////////////////////////////////////////////////////////////////
1312/// Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found
1313
1315{
1316 return static_cast<RooAbsData*>(_embeddedDataList.FindObject(name.c_str())) ;
1317}
1318
1319
1320
1321
1322////////////////////////////////////////////////////////////////////////////////
1323/// Return set with all variable objects
1324
1326{
1327 RooArgSet ret ;
1328
1329 // Split list of components in pdfs, functions and variables
1330 for(RooAbsArg* parg : _allOwnedNodes) {
1331 if (parg->IsA()->InheritsFrom(RooRealVar::Class())) {
1332 ret.add(*parg) ;
1333 }
1334 }
1335
1336 return ret ;
1337}
1338
1339
1340////////////////////////////////////////////////////////////////////////////////
1341/// Return set with all category objects
1342
1344{
1345 RooArgSet ret ;
1346
1347 // Split list of components in pdfs, functions and variables
1348 for(RooAbsArg* parg : _allOwnedNodes) {
1349 if (parg->IsA()->InheritsFrom(RooCategory::Class())) {
1350 ret.add(*parg) ;
1351 }
1352 }
1353
1354 return ret ;
1355}
1356
1357
1358
1359////////////////////////////////////////////////////////////////////////////////
1360/// Return set with all function objects
1361
1363{
1364 RooArgSet ret ;
1365
1366 // Split list of components in pdfs, functions and variables
1367 for(RooAbsArg* parg : _allOwnedNodes) {
1368 if (parg->IsA()->InheritsFrom(RooAbsReal::Class()) &&
1369 !parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
1370 !parg->IsA()->InheritsFrom(RooConstVar::Class()) &&
1371 !parg->IsA()->InheritsFrom(RooRealVar::Class())) {
1372 ret.add(*parg) ;
1373 }
1374 }
1375
1376 return ret ;
1377}
1378
1379
1380////////////////////////////////////////////////////////////////////////////////
1381/// Return set with all category function objects
1382
1384{
1385 RooArgSet ret ;
1386
1387 // Split list of components in pdfs, functions and variables
1388 for(RooAbsArg* parg : _allOwnedNodes) {
1389 if (parg->IsA()->InheritsFrom(RooAbsCategory::Class()) &&
1390 !parg->IsA()->InheritsFrom(RooCategory::Class())) {
1391 ret.add(*parg) ;
1392 }
1393 }
1394 return ret ;
1395}
1396
1397
1398
1399////////////////////////////////////////////////////////////////////////////////
1400/// Return set with all resolution model objects
1401
1403{
1404 RooArgSet ret ;
1405
1406 // Split list of components in pdfs, functions and variables
1407 for(RooAbsArg* parg : _allOwnedNodes) {
1408 if (parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
1409 if (!(static_cast<RooResolutionModel*>(parg))->isConvolved()) {
1410 ret.add(*parg) ;
1411 }
1412 }
1413 }
1414 return ret ;
1415}
1416
1417
1418////////////////////////////////////////////////////////////////////////////////
1419/// Return set with all probability density function objects
1420
1422{
1423 RooArgSet ret ;
1424
1425 // Split list of components in pdfs, functions and variables
1426 for(RooAbsArg* parg : _allOwnedNodes) {
1427 if (parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
1428 !parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
1429 ret.add(*parg) ;
1430 }
1431 }
1432 return ret ;
1433}
1434
1435
1436
1437////////////////////////////////////////////////////////////////////////////////
1438/// Return list of all dataset in the workspace
1439
1440std::list<RooAbsData*> RooWorkspace::allData() const
1441{
1442 std::list<RooAbsData*> ret ;
1444 ret.push_back(dat) ;
1445 }
1446 return ret ;
1447}
1448
1449
1450////////////////////////////////////////////////////////////////////////////////
1451/// Return list of all dataset in the workspace
1452
1453std::list<RooAbsData*> RooWorkspace::allEmbeddedData() const
1454{
1455 std::list<RooAbsData*> ret ;
1457 ret.push_back(dat) ;
1458 }
1459 return ret ;
1460}
1461
1462
1463
1464////////////////////////////////////////////////////////////////////////////////
1465/// Return list of all generic objects in the workspace
1466
1467std::list<TObject*> RooWorkspace::allGenericObjects() const
1468{
1469 std::list<TObject*> ret ;
1470 for(TObject * gobj : _genObjects) {
1471
1472 // If found object is wrapper, return payload
1473 if (gobj->IsA()==RooTObjWrap::Class()) {
1474 ret.push_back((static_cast<RooTObjWrap*>(gobj))->obj()) ;
1475 } else {
1476 ret.push_back(gobj) ;
1477 }
1478 }
1479 return ret ;
1480}
1481
1482
1483namespace {
1484
1485std::string findFileInPath(std::string const &file, std::list<std::string> const &dirList)
1486{
1487 // Check list of additional paths
1488 for (std::string const &diter : dirList) {
1489
1490 char *cpath = gSystem->ConcatFileName(diter.c_str(), file.c_str());
1491 std::string path = cpath;
1492 delete[] cpath;
1493 if (!gSystem->AccessPathName(path.c_str())) {
1494 // found file
1495 return path;
1496 }
1497 }
1498 return "";
1499}
1500
1501} // namespace
1502
1503
1504////////////////////////////////////////////////////////////////////////////////
1505/// Import code of class 'tc' into the repository. If code is already in repository it is only imported
1506/// again if doReplace is false. The names and location of the source files is determined from the information
1507/// in TClass. If no location is found in the TClass information, the files are searched in the workspace
1508/// search path, defined by addClassDeclImportDir() and addClassImplImportDir() for declaration and implementation
1509/// files respectively. If files cannot be found, abort with error status, otherwise update the internal
1510/// class-to-file map and import the contents of the files, if they are not imported yet.
1511
1513{
1514
1515 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") request to import code of class " << tc->GetName() << std::endl ;
1516
1517 // *** PHASE 1 *** Check if file needs to be imported, or is in ROOT distribution, and check if it can be persisted
1518
1519 // Check if we already have the class (i.e. it is in the classToFile map)
1520 if (!doReplace && _c2fmap.find(tc->GetName())!=_c2fmap.end()) {
1521 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") code of class " << tc->GetName() << " already imported, skipping" << std::endl ;
1522 return true ;
1523 }
1524
1525 // Check if class is listed in a ROOTMAP file - if so we can skip it because it is in the root distribution
1526 const char* mapEntry = gInterpreter->GetClassSharedLibs(tc->GetName()) ;
1527 if (mapEntry && strlen(mapEntry)>0) {
1528 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") code of class " << tc->GetName() << " is in ROOT distribution, skipping " << std::endl ;
1529 return true ;
1530 }
1531
1532 // Retrieve file names through ROOT TClass interface
1533 string implfile = tc->GetImplFileName() ;
1534 string declfile = tc->GetDeclFileName() ;
1535
1536 // Check that file names are not empty
1537 if (implfile.empty() || declfile.empty()) {
1538 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") ERROR: cannot retrieve code file names for class "
1539 << tc->GetName() << " through ROOT TClass interface, unable to import code" << std::endl ;
1540 return false ;
1541 }
1542
1543 // Check if header filename is found in ROOT distribution, if so, do not import class
1544 TString rootsys = gSystem->Getenv("ROOTSYS") ;
1545 if (TString(implfile.c_str()).Index(rootsys)>=0) {
1546 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") code of class " << tc->GetName() << " is in ROOT distribution, skipping " << std::endl ;
1547 return true ;
1548 }
1549
1550 // Require that class meets technical criteria to be persistable (i.e it has a default constructor)
1551 // (We also need a default constructor of abstract classes, but cannot check that through is interface
1552 // as TClass::HasDefaultCtor only returns true for callable default constructors)
1553 if (!(tc->Property() & kIsAbstract) && !tc->HasDefaultConstructor()) {
1554 oocoutW(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() << ") WARNING cannot import class "
1555 << tc->GetName() << " : it cannot be persisted because it doesn't have a default constructor. Please fix " << std::endl ;
1556 return false ;
1557 }
1558
1559
1560 // *** PHASE 2 *** Check if declaration and implementation files can be located
1561
1562 std::string declpath;
1563 std::string implpath;
1564
1565 // Check if header file can be found in specified location
1566 // If not, scan through list of 'class declaration' paths in RooWorkspace
1567 if (gSystem->AccessPathName(declfile.c_str())) {
1568
1570
1571 // Header file cannot be found anywhere, warn user and abort operation
1572 if (declpath.empty()) {
1573 oocoutW(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() << ") WARNING Cannot access code of class "
1574 << tc->GetName() << " because header file " << declfile << " is not found in current directory nor in $ROOTSYS" ;
1575 if (!_classDeclDirList.empty()) {
1576 ooccoutW(_wspace,ObjectHandling) << ", nor in the search path " ;
1578
1579 while(diter!= RooWorkspace::_classDeclDirList.end()) {
1580
1582 ooccoutW(_wspace,ObjectHandling) << "," ;
1583 }
1584 ooccoutW(_wspace,ObjectHandling) << diter->c_str() ;
1585 ++diter ;
1586 }
1587 }
1588 ooccoutW(_wspace,ObjectHandling) << ". To fix this problem, add the required directory to the search "
1589 << "path using RooWorkspace::addClassDeclImportDir(const char* dir)" << std::endl ;
1590
1591 return false ;
1592 }
1593 }
1594
1595
1596 // Check if implementation file can be found in specified location
1597 // If not, scan through list of 'class implementation' paths in RooWorkspace
1598 if (gSystem->AccessPathName(implfile.c_str())) {
1599
1601
1602 // Implementation file cannot be found anywhere, warn user and abort operation
1603 if (implpath.empty()) {
1604 oocoutW(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() << ") WARNING Cannot access code of class "
1605 << tc->GetName() << " because implementation file " << implfile << " is not found in current directory nor in $ROOTSYS" ;
1606 if (!_classDeclDirList.empty()) {
1607 ooccoutW(_wspace,ObjectHandling) << ", nor in the search path " ;
1609
1610 while(iiter!= RooWorkspace::_classImplDirList.end()) {
1611
1613 ooccoutW(_wspace,ObjectHandling) << "," ;
1614 }
1615 ooccoutW(_wspace,ObjectHandling) << iiter->c_str() ;
1616 ++iiter ;
1617 }
1618 }
1619 ooccoutW(_wspace,ObjectHandling) << ". To fix this problem add the required directory to the search "
1620 << "path using RooWorkspace::addClassImplImportDir(const char* dir)" << std::endl;
1621 return false;
1622 }
1623 }
1624
1625 char buf[64000];
1626
1627 // *** Phase 3 *** Prepare to import code from files into STL string buffer
1628 //
1629 // Code storage is organized in two linked maps
1630 //
1631 // _fmap contains stl strings with code, indexed on declaration file name
1632 //
1633 // _c2fmap contains list of declaration file names and list of base classes
1634 // and is indexed on class name
1635 //
1636 // Phase 3 is skipped if fmap already contains an entry with given filebasename
1637
1638 const std::string declfilename = !declpath.empty() ? gSystem->BaseName(declpath.c_str())
1639 : gSystem->BaseName(declfile.c_str());
1640
1641 // Split in base and extension
1642 int dotpos2 = strrchr(declfilename.c_str(),'.') - declfilename.c_str() ;
1643 string declfilebase = declfilename.substr(0,dotpos2) ;
1644 string declfileext = declfilename.substr(dotpos2+1) ;
1645
1647
1648 // If file has not been stored yet, enter stl strings with implementation and declaration in file map
1649 if (_fmap.find(declfilebase) == _fmap.end()) {
1650
1651 // Open declaration file
1652 std::fstream fdecl(!declpath.empty() ? declpath.c_str() : declfile.c_str());
1653
1654 // Abort import if declaration file cannot be opened
1655 if (!fdecl) {
1656 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1657 << ") ERROR opening declaration file " << declfile << std::endl ;
1658 return false ;
1659 }
1660
1661 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1662 << ") importing code of class " << tc->GetName()
1663 << " from " << (!implpath.empty() ? implpath.c_str() : implfile.c_str())
1664 << " and " << (!declpath.empty() ? declpath.c_str() : declfile.c_str()) << std::endl ;
1665
1666
1667 // Read entire file into an stl string
1668 string decl ;
1669 while(fdecl.getline(buf,1023)) {
1670
1671 // Look for include state of self
1672 bool processedInclude = false ;
1673 char* extincfile = nullptr ;
1674
1675 // Look for include of declaration file corresponding to this implementation file
1676 if (strstr(buf,"#include")) {
1677 // Process #include statements here
1678 char tmp[64000];
1679 strlcpy(tmp, buf, 64000);
1680 bool stdinclude = strchr(buf, '<');
1681 strtok(tmp, " <\"");
1682 char *incfile = strtok(nullptr, " <>\"");
1683
1684 if (!stdinclude) {
1685 // check if it lives in $ROOTSYS/include
1686 TString hpath = gSystem->Getenv("ROOTSYS");
1687 hpath += "/include/";
1688 hpath += incfile;
1689 if (gSystem->AccessPathName(hpath.Data())) {
1690 oocoutI(_wspace, ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1691 << ") scheduling include file " << incfile << " for import" << std::endl;
1692 extraHeaders.push_back(incfile);
1694 processedInclude = true;
1695 }
1696 }
1697 }
1698
1699 if (processedInclude) {
1700 decl += "// external include file below retrieved from workspace code storage\n" ;
1701 decl += Form("#include \"%s\"\n",extincfile) ;
1702 } else {
1703 decl += buf ;
1704 decl += '\n' ;
1705 }
1706 }
1707
1708 // Open implementation file
1709 fstream fimpl(!implpath.empty() ? implpath.c_str() : implfile.c_str()) ;
1710
1711 // Abort import if implementation file cannot be opened
1712 if (!fimpl) {
1713 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1714 << ") ERROR opening implementation file " << implfile << std::endl ;
1715 return false ;
1716 }
1717
1718
1719 // Import entire implementation file into stl string
1720 string impl ;
1721 while(fimpl.getline(buf,1023)) {
1722 // Process #include statements here
1723
1724 // Look for include state of self
1725 bool foundSelfInclude=false ;
1726 bool processedInclude = false ;
1727 char* extincfile = nullptr ;
1728
1729 // Look for include of declaration file corresponding to this implementation file
1730 if (strstr(buf,"#include")) {
1731 // Process #include statements here
1732 char tmp[64000];
1733 strlcpy(tmp, buf, 64000);
1734 bool stdinclude = strchr(buf, '<');
1735 strtok(tmp, " <\"");
1736 char *incfile = strtok(nullptr, " <>\"");
1737
1738 if (strstr(incfile, declfilename.c_str())) {
1739 foundSelfInclude = true;
1740 }
1741
1742 if (!stdinclude && !foundSelfInclude) {
1743 // check if it lives in $ROOTSYS/include
1744 TString hpath = gSystem->Getenv("ROOTSYS");
1745 hpath += "/include/";
1746 hpath += incfile;
1747
1748 if (gSystem->AccessPathName(hpath.Data())) {
1749 oocoutI(_wspace, ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1750 << ") scheduling include file " << incfile << " for import" << std::endl;
1751 extraHeaders.push_back(incfile);
1753 processedInclude = true;
1754 }
1755 }
1756 }
1757
1758 // Explicitly rewrite include of own declaration file to string
1759 // any directory prefixes, copy all other lines verbatim in stl string
1760 if (foundSelfInclude) {
1761 // If include of self is found, substitute original include
1762 // which may have directory structure with a plain include
1763 impl += "// class declaration include file below retrieved from workspace code storage\n" ;
1764 impl += Form("#include \"%s.%s\"\n",declfilebase.c_str(),declfileext.c_str()) ;
1765 } else if (processedInclude) {
1766 impl += "// external include file below retrieved from workspace code storage\n" ;
1767 impl += Form("#include \"%s\"\n",extincfile) ;
1768 } else {
1769 impl += buf ;
1770 impl += '\n' ;
1771 }
1772 }
1773
1774 // Create entry in file map
1775 _fmap[declfilebase]._hfile = decl ;
1776 _fmap[declfilebase]._cxxfile = impl ;
1777 _fmap[declfilebase]._hext = declfileext ;
1778
1779 // Process extra includes now
1780 for (list<string>::iterator ehiter = extraHeaders.begin() ; ehiter != extraHeaders.end() ; ++ehiter ) {
1781 if (_ehmap.find(*ehiter) == _ehmap.end()) {
1782
1783 ExtraHeader eh ;
1784 eh._hname = ehiter->c_str() ;
1785 fstream fehdr(ehiter->c_str()) ;
1786 string ehimpl ;
1787 char buf2[1024] ;
1788 while(fehdr.getline(buf2,1023)) {
1789
1790 // Look for include of declaration file corresponding to this implementation file
1791 if (strstr(buf2,"#include")) {
1792 // Process #include statements here
1793 char tmp[64000];
1794 strlcpy(tmp, buf2, 64000);
1795 bool stdinclude = strchr(buf, '<');
1796 strtok(tmp, " <\"");
1797 char *incfile = strtok(nullptr, " <>\"");
1798
1799 if (!stdinclude) {
1800 // check if it lives in $ROOTSYS/include
1801 TString hpath = gSystem->Getenv("ROOTSYS");
1802 hpath += "/include/";
1803 hpath += incfile;
1804 if (gSystem->AccessPathName(hpath.Data())) {
1805 oocoutI(_wspace, ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1806 << ") scheduling recursive include file " << incfile << " for import"
1807 << std::endl;
1808 extraHeaders.push_back(incfile);
1809 }
1810 }
1811 }
1812
1813 ehimpl += buf2;
1814 ehimpl += '\n';
1815 }
1816 eh._hfile = ehimpl.c_str();
1817
1818 _ehmap[ehiter->c_str()] = eh;
1819 }
1820 }
1821
1822 } else {
1823
1824 // Inform that existing file entry is being recycled because it already contained class code
1825 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName()
1826 << ") code of class " << tc->GetName()
1827 << " was already imported from " << (!implpath.empty() ? implpath : implfile)
1828 << " and " << (!declpath.empty() ? declpath.c_str() : declfile.c_str()) << std::endl;
1829
1830 }
1831
1832
1833 // *** PHASE 4 *** Import stl strings with code into workspace
1834 //
1835 // If multiple classes are declared in a single code unit, there will be
1836 // multiple _c2fmap entries all pointing to the same _fmap entry.
1837
1838 // Make list of all immediate base classes of this class
1840 TList* bl = tc->GetListOfBases() ;
1841 std::list<TClass*> bases ;
1842 for(auto * base : static_range_cast<TBaseClass*>(*bl)) {
1843 if (baseNameList.Length()>0) {
1844 baseNameList += "," ;
1845 }
1846 baseNameList += base->GetClassPointer()->GetName() ;
1847 bases.push_back(base->GetClassPointer()) ;
1848 }
1849
1850 // Map class name to above _fmap entries, along with list of base classes
1851 // in _c2fmap
1852 _c2fmap[tc->GetName()]._baseName = baseNameList ;
1853 _c2fmap[tc->GetName()]._fileBase = declfilebase ;
1854
1855 // Recursive store all base classes.
1856 for(TClass* bclass : bases) {
1858 }
1859
1860 return true ;
1861}
1862
1863
1864////////////////////////////////////////////////////////////////////////////////
1865/// Create transient TDirectory representation of this workspace. This directory
1866/// will appear as a subdirectory of the directory that contains the workspace
1867/// and will have the name of the workspace suffixed with "Dir". The TDirectory
1868/// interface is read-only. Any attempt to insert objects into the workspace
1869/// directory representation will result in an error message. Note that some
1870/// ROOT object like TH1 automatically insert themselves into the current directory
1871/// when constructed. This will give error messages when done in a workspace
1872/// directory.
1873
1875{
1876 if (_dir) return true ;
1877
1878 std::string title= "TDirectory representation of RooWorkspace " + std::string(GetName());
1879 _dir = new WSDir(GetName(),title.c_str(),this) ;
1880
1881 for (RooAbsArg * darg : _allOwnedNodes) {
1882 if (darg->IsA() != RooConstVar::Class()) {
1884 }
1885 }
1886
1887 return true ;
1888}
1889
1890
1891
1892////////////////////////////////////////////////////////////////////////////////
1893/// Import a clone of a generic TObject into workspace generic object container. Imported
1894/// object can be retrieved by name through the obj() method. The object is cloned upon
1895/// importation and the input argument does not need to live beyond the import call
1896///
1897/// Returns true if an error has occurred.
1898
1900{
1901 // First check if object with given name already exists
1902 std::unique_ptr<TObject> oldObj{_genObjects.FindObject(object.GetName())};
1903 if (oldObj && !replaceExisting) {
1904 coutE(InputArguments) << "RooWorkspace::import(" << GetName() << ") generic object with name "
1905 << object.GetName() << " is already in workspace and replaceExisting flag is set to false" << std::endl ;
1906 return true ;
1907 }
1908
1909 // Grab the current state of the directory Auto-Add
1910 ROOT::DirAutoAdd_t func = object.IsA()->GetDirectoryAutoAdd();
1911 object.IsA()->SetDirectoryAutoAdd(nullptr);
1912 bool tmp = RooPlot::setAddDirectoryStatus(false) ;
1913
1914 if (oldObj) {
1915 _genObjects.Replace(oldObj.get(),object.Clone()) ;
1916 } else {
1917 _genObjects.Add(object.Clone()) ;
1918 }
1919
1920 // Reset the state of the directory Auto-Add
1921 object.IsA()->SetDirectoryAutoAdd(func);
1923
1924 return false ;
1925}
1926
1927
1928
1929
1930////////////////////////////////////////////////////////////////////////////////
1931/// Import a clone of a generic TObject into workspace generic object container.
1932/// The imported object will be stored under the given alias name rather than its
1933/// own name. Imported object can be retrieved its alias name through the obj() method.
1934/// The object is cloned upon importation and the input argument does not need to live beyond the import call
1935/// This method is mostly useful for importing objects that do not have a settable name such as TMatrix
1936///
1937/// Returns true if an error has occurred.
1938
1939bool RooWorkspace::import(TObject const& object, const char* aliasName, bool replaceExisting)
1940{
1941 // First check if object with given name already exists
1942 std::unique_ptr<TObject> oldObj{_genObjects.FindObject(aliasName)};
1943 if (oldObj && !replaceExisting) {
1944 coutE(InputArguments) << "RooWorkspace::import(" << GetName() << ") generic object with name "
1945 << aliasName << " is already in workspace and replaceExisting flag is set to false" << std::endl ;
1946 return true ;
1947 }
1948
1949 TH1::AddDirectory(false) ;
1950 auto wrapper = new RooTObjWrap(object.Clone()) ;
1951 TH1::AddDirectory(true) ;
1952 wrapper->setOwning(true) ;
1953 wrapper->SetName(aliasName) ;
1954 wrapper->SetTitle(aliasName) ;
1955
1956 if (oldObj) {
1958 } else {
1960 }
1961 return false ;
1962}
1963
1964
1965
1966
1967////////////////////////////////////////////////////////////////////////////////
1968/// Insert RooStudyManager module
1969
1971{
1972 RooAbsStudy* clone = static_cast<RooAbsStudy*>(study.Clone()) ;
1973 _studyMods.Add(clone) ;
1974 return false ;
1975}
1976
1977
1978
1979
1980////////////////////////////////////////////////////////////////////////////////
1981/// Remove all RooStudyManager modules
1982
1987
1988
1989
1990
1991////////////////////////////////////////////////////////////////////////////////
1992/// Return any type of object (RooAbsArg, RooAbsData or generic object) with given name)
1993
1995{
1996 // Try RooAbsArg first
1997 TObject* ret = arg(name) ;
1998 if (ret) return ret ;
1999
2000 // Then try RooAbsData
2001 ret = data(name) ;
2002 if (ret) return ret ;
2003
2004 // Finally try generic object store
2005 return genobj(name) ;
2006}
2007
2008
2009
2010////////////////////////////////////////////////////////////////////////////////
2011/// Return generic object with given name
2012
2014{
2015 // Find object by name
2016 TObject* gobj = _genObjects.FindObject(name.c_str()) ;
2017
2018 // Exit here if not found
2019 if (!gobj) return nullptr;
2020
2021 // If found object is wrapper, return payload
2022 if (gobj->IsA()==RooTObjWrap::Class()) return (static_cast<RooTObjWrap*>(gobj))->obj() ;
2023
2024 return gobj ;
2025}
2026
2027
2028
2029////////////////////////////////////////////////////////////////////////////////
2030
2031bool RooWorkspace::cd(const char* path)
2032{
2033 makeDir() ;
2034 return _dir->cd(path) ;
2035}
2036
2037
2038
2039////////////////////////////////////////////////////////////////////////////////
2040/// Save this current workspace into given file
2041
2042bool RooWorkspace::writeToFile(const char* fileName, bool recreate)
2043{
2044 TFile f(fileName,recreate?"RECREATE":"UPDATE") ;
2045 Write() ;
2046 return false ;
2047}
2048
2049
2050
2051////////////////////////////////////////////////////////////////////////////////
2052/// Return instance to factory tool
2053
2055{
2056 if (_factory) {
2057 return *_factory;
2058 }
2059 cxcoutD(ObjectHandling) << "INFO: Creating RooFactoryWSTool associated with this workspace" << std::endl ;
2061 return *_factory;
2062}
2063
2064
2065
2066
2067////////////////////////////////////////////////////////////////////////////////
2068/// Short-hand function for `factory()->process(expr);`
2069///
2070/// \copydoc RooFactoryWSTool::process(const char*)
2075
2076
2077
2078
2079////////////////////////////////////////////////////////////////////////////////
2080/// Print contents of the workspace
2081
2083{
2084 bool treeMode(false) ;
2085 bool verbose(false);
2086 if (TString(opts).Contains("t")) {
2087 treeMode=true ;
2088 }
2089 if (TString(opts).Contains("v")) {
2090 verbose = true;
2091 }
2092
2093 std::cout << std::endl << "RooWorkspace(" << GetName() << ") " << GetTitle() << " contents" << std::endl << std::endl ;
2094
2097 RooArgSet varSet ;
2101
2102
2103 // Split list of components in pdfs, functions and variables
2104 for(RooAbsArg* parg : _allOwnedNodes) {
2105
2106 //---------------
2107
2108 if (treeMode) {
2109
2110 // In tree mode, only add nodes with no clients to the print lists
2111
2112 if (parg->IsA()->InheritsFrom(RooAbsPdf::Class())) {
2113 if (!parg->hasClients()) {
2114 pdfSet.add(*parg) ;
2115 }
2116 }
2117
2118 if (parg->IsA()->InheritsFrom(RooAbsReal::Class()) &&
2119 !parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
2120 !parg->IsA()->InheritsFrom(RooConstVar::Class()) &&
2121 !parg->IsA()->InheritsFrom(RooRealVar::Class())) {
2122 if (!parg->hasClients()) {
2123 funcSet.add(*parg) ;
2124 }
2125 }
2126
2127
2128 if (parg->IsA()->InheritsFrom(RooAbsCategory::Class()) &&
2129 !parg->IsA()->InheritsFrom(RooCategory::Class())) {
2130 if (!parg->hasClients()) {
2131 catfuncSet.add(*parg) ;
2132 }
2133 }
2134
2135 } else {
2136
2137 if (parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
2138 if ((static_cast<RooResolutionModel*>(parg))->isConvolved()) {
2139 convResoSet.add(*parg) ;
2140 } else {
2141 resoSet.add(*parg) ;
2142 }
2143 }
2144
2145 if (parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
2146 !parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
2147 pdfSet.add(*parg) ;
2148 }
2149
2150 if (parg->IsA()->InheritsFrom(RooAbsReal::Class()) &&
2151 !parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
2152 !parg->IsA()->InheritsFrom(RooConstVar::Class()) &&
2153 !parg->IsA()->InheritsFrom(RooRealVar::Class())) {
2154 funcSet.add(*parg) ;
2155 }
2156
2157 if (parg->IsA()->InheritsFrom(RooAbsCategory::Class()) &&
2158 !parg->IsA()->InheritsFrom(RooCategory::Class())) {
2159 catfuncSet.add(*parg) ;
2160 }
2161 }
2162
2163 if (parg->IsA()->InheritsFrom(RooRealVar::Class())) {
2164 varSet.add(*parg) ;
2165 }
2166
2167 if (parg->IsA()->InheritsFrom(RooCategory::Class())) {
2168 varSet.add(*parg) ;
2169 }
2170
2171 }
2172
2173
2174 RooFit::MsgLevel oldLevel = RooMsgService::instance().globalKillBelow() ;
2175 RooMsgService::instance().setGlobalKillBelow(RooFit::WARNING) ;
2176
2177 if (!varSet.empty()) {
2178 varSet.sort() ;
2179 std::cout << "variables" << std::endl ;
2180 std::cout << "---------" << std::endl ;
2181 std::cout << varSet << std::endl ;
2182 std::cout << std::endl ;
2183 }
2184
2185 if (!pdfSet.empty()) {
2186 std::cout << "p.d.f.s" << std::endl ;
2187 std::cout << "-------" << std::endl ;
2188 pdfSet.sort() ;
2189 for(RooAbsArg* parg : pdfSet) {
2190 if (treeMode) {
2191 parg->printComponentTree() ;
2192 } else {
2193 parg->Print() ;
2194 }
2195 }
2196 std::cout << std::endl ;
2197 }
2198
2199 if (!treeMode) {
2200 if (!resoSet.empty()) {
2201 std::cout << "analytical resolution models" << std::endl ;
2202 std::cout << "----------------------------" << std::endl ;
2203 resoSet.sort() ;
2204 for(RooAbsArg* parg : resoSet) {
2205 parg->Print() ;
2206 }
2207 std::cout << std::endl ;
2208 }
2209 }
2210
2211 if (!funcSet.empty()) {
2212 std::cout << "functions" << std::endl ;
2213 std::cout << "--------" << std::endl ;
2214 funcSet.sort() ;
2215 for(RooAbsArg * parg : funcSet) {
2216 if (treeMode) {
2217 parg->printComponentTree() ;
2218 } else {
2219 parg->Print() ;
2220 }
2221 }
2222 std::cout << std::endl ;
2223 }
2224
2225 if (!catfuncSet.empty()) {
2226 std::cout << "category functions" << std::endl ;
2227 std::cout << "------------------" << std::endl ;
2228 catfuncSet.sort() ;
2229 for(RooAbsArg* parg : catfuncSet) {
2230 if (treeMode) {
2231 parg->printComponentTree() ;
2232 } else {
2233 parg->Print() ;
2234 }
2235 }
2236 std::cout << std::endl ;
2237 }
2238
2239 if (!_dataList.empty()) {
2240 std::cout << "datasets" << std::endl ;
2241 std::cout << "--------" << std::endl ;
2243 std::cout << data2->ClassName() << "::" << data2->GetName() << *data2->get() << std::endl;
2244 }
2245 std::cout << std::endl ;
2246 }
2247
2248 if (!_embeddedDataList.empty()) {
2249 std::cout << "embedded datasets (in pdfs and functions)" << std::endl ;
2250 std::cout << "-----------------------------------------" << std::endl ;
2252 std::cout << data2->ClassName() << "::" << data2->GetName() << *data2->get() << std::endl ;
2253 }
2254 std::cout << std::endl ;
2255 }
2256
2257 if (!_snapshots.empty()) {
2258 std::cout << "parameter snapshots" << std::endl ;
2259 std::cout << "-------------------" << std::endl ;
2261 std::cout << snap->GetName() << " = (" ;
2262 bool first(true) ;
2263 for(RooAbsArg* a : *snap) {
2264 if (first) { first=false ; } else { std::cout << "," ; }
2265 std::cout << a->GetName() << "=" ;
2266 a->printValue(std::cout) ;
2267 if (a->isConstant()) {
2268 std::cout << "[C]" ;
2269 }
2270 }
2271 std::cout << ")" << std::endl ;
2272 }
2273 std::cout << std::endl ;
2274 }
2275
2276
2277 if (!_namedSets.empty()) {
2278 std::cout << "named sets" << std::endl ;
2279 std::cout << "----------" << std::endl ;
2280 for (map<string,RooArgSet>::const_iterator it = _namedSets.begin() ; it != _namedSets.end() ; ++it) {
2281 if (verbose || !isCacheSet(it->first)) {
2282 std::cout << it->first << ":" << it->second << std::endl;
2283 }
2284 }
2285
2286 std::cout << std::endl ;
2287 }
2288
2289
2290 if (!_genObjects.empty()) {
2291 std::cout << "generic objects" << std::endl ;
2292 std::cout << "---------------" << std::endl ;
2293 for(TObject* gobj : _genObjects) {
2294 if (gobj->IsA()==RooTObjWrap::Class()) {
2295 std::cout << (static_cast<RooTObjWrap*>(gobj))->obj()->ClassName() << "::" << gobj->GetName() << std::endl ;
2296 } else {
2297 std::cout << gobj->ClassName() << "::" << gobj->GetName() << std::endl ;
2298 }
2299 }
2300 std::cout << std::endl ;
2301
2302 }
2303
2304 if (!_studyMods.empty()) {
2305 std::cout << "study modules" << std::endl ;
2306 std::cout << "-------------" << std::endl ;
2307 for(TObject* smobj : _studyMods) {
2308 std::cout << smobj->ClassName() << "::" << smobj->GetName() << std::endl ;
2309 }
2310 std::cout << std::endl ;
2311
2312 }
2313
2314 if (!_classes.listOfClassNames().empty()) {
2315 std::cout << "embedded class code" << std::endl ;
2316 std::cout << "-------------------" << std::endl ;
2317 std::cout << _classes.listOfClassNames() << std::endl ;
2318 std::cout << std::endl ;
2319 }
2320
2321 if (!_eocache.empty()) {
2322 std::cout << "embedded precalculated expensive components" << std::endl ;
2323 std::cout << "-------------------------------------------" << std::endl ;
2324 _eocache.print() ;
2325 }
2326
2327 RooMsgService::instance().setGlobalKillBelow(oldLevel) ;
2328
2329 return ;
2330}
2331
2332
2333////////////////////////////////////////////////////////////////////////////////
2334/// Custom streamer for the workspace. Stream contents of workspace
2335/// and code repository. When reading, read code repository first
2336/// and compile missing classes before proceeding with streaming
2337/// of workspace contents to avoid errors.
2338
2340{
2342
2343 // Stream an object of class RooWorkspace::CodeRepo.
2344 if (R__b.IsReading()) {
2345
2346 UInt_t R__s;
2347 UInt_t R__c;
2348 Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
2349
2350 // Stream contents of ClassFiles map
2351 Int_t count(0);
2352 R__b >> count;
2353 while (count--) {
2354 TString name;
2355 name.Streamer(R__b);
2356 _fmap[name]._hext.Streamer(R__b);
2357 _fmap[name]._hfile.Streamer(R__b);
2358 _fmap[name]._cxxfile.Streamer(R__b);
2359 }
2360
2361 // Stream contents of ClassRelInfo map
2362 count = 0;
2363 R__b >> count;
2364 while (count--) {
2365 TString name;
2366 name.Streamer(R__b);
2367 _c2fmap[name]._baseName.Streamer(R__b);
2368 _c2fmap[name]._fileBase.Streamer(R__b);
2369 }
2370
2371 if (R__v == 2) {
2372
2373 count = 0;
2374 R__b >> count;
2375 while (count--) {
2376 TString name;
2377 name.Streamer(R__b);
2378 _ehmap[name]._hname.Streamer(R__b);
2379 _ehmap[name]._hfile.Streamer(R__b);
2380 }
2381 }
2382
2383 R__b.CheckByteCount(R__s, R__c, thisClass::IsA());
2384
2385 // Instantiate any classes that are not defined in current session
2386 _compiledOK = !compileClasses();
2387
2388 } else {
2389
2390 UInt_t R__c;
2391 R__c = R__b.WriteVersion(thisClass::IsA(), true);
2392
2393 // Stream contents of ClassFiles map
2394 UInt_t count = _fmap.size();
2395 R__b << count;
2396 map<TString, ClassFiles>::iterator iter = _fmap.begin();
2397 while (iter != _fmap.end()) {
2398 TString key_copy(iter->first);
2399 key_copy.Streamer(R__b);
2400 iter->second._hext.Streamer(R__b);
2401 iter->second._hfile.Streamer(R__b);
2402 iter->second._cxxfile.Streamer(R__b);
2403
2404 ++iter;
2405 }
2406
2407 // Stream contents of ClassRelInfo map
2408 count = _c2fmap.size();
2409 R__b << count;
2410 map<TString, ClassRelInfo>::iterator iter2 = _c2fmap.begin();
2411 while (iter2 != _c2fmap.end()) {
2412 TString key_copy(iter2->first);
2413 key_copy.Streamer(R__b);
2414 iter2->second._baseName.Streamer(R__b);
2415 iter2->second._fileBase.Streamer(R__b);
2416 ++iter2;
2417 }
2418
2419 // Stream contents of ExtraHeader map
2420 count = _ehmap.size();
2421 R__b << count;
2422 map<TString, ExtraHeader>::iterator iter3 = _ehmap.begin();
2423 while (iter3 != _ehmap.end()) {
2424 TString key_copy(iter3->first);
2425 key_copy.Streamer(R__b);
2426 iter3->second._hname.Streamer(R__b);
2427 iter3->second._hfile.Streamer(R__b);
2428 ++iter3;
2429 }
2430
2431 R__b.SetByteCount(R__c, true);
2432 }
2433}
2434
2435
2436////////////////////////////////////////////////////////////////////////////////
2437/// Stream an object of class RooWorkspace. This is a standard ROOT streamer for the
2438/// I/O part. This custom function exists to detach all external client links
2439/// from the payload prior to writing the payload so that these client links
2440/// are not persisted. (Client links occur if external function objects use
2441/// objects contained in the workspace as input)
2442/// After the actual writing, these client links are restored.
2443
2445{
2446 if (R__b.IsReading()) {
2447
2448 R__b.ReadClassBuffer(RooWorkspace::Class(), this);
2449
2450 // Perform any pass-2 schema evolution here
2451 for (RooAbsArg *node : _allOwnedNodes) {
2452 node->ioStreamerPass2();
2453 }
2455
2456 // Make expensive object cache of all objects point to intermal copy.
2457 // Somehow this doesn't work OK automatically
2458 for (RooAbsArg *node : _allOwnedNodes) {
2459 node->setExpensiveObjectCache(_eocache);
2460 node->setWorkspace(*this);
2461#ifdef ROOFIT_LEGACY_EVAL_BACKEND
2462 if (dynamic_cast<RooAbsOptTestStatistic *>(node)) {
2463 RooAbsOptTestStatistic *tmp = static_cast<RooAbsOptTestStatistic *>(node);
2464 if (tmp->isSealed() && tmp->sealNotice() && strlen(tmp->sealNotice()) > 0) {
2465 std::cout << "RooWorkspace::Streamer(" << GetName() << ") " << node->ClassName() << "::" << node->GetName()
2466 << " : " << tmp->sealNotice() << std::endl;
2467 }
2468 }
2469#endif
2470 }
2471
2472 for(TObject * gobj : allGenericObjects()) {
2473 if (auto handle = dynamic_cast<RooWorkspaceHandle*>(gobj)) {
2474 handle->ReplaceWS(this);
2475 }
2476 }
2477
2478 } else {
2479
2480 // Make lists of external clients of WS objects, and remove those links temporarily
2481
2485
2487
2488 // Loop over client list of this arg
2489 std::vector<RooAbsArg *> clientsTmp{tmparg->_clientList.begin(), tmparg->_clientList.end()};
2490 for (auto client : clientsTmp) {
2491 if (!_allOwnedNodes.containsInstance(*client)) {
2492
2493 const auto refCount = tmparg->_clientList.refCount(client);
2494 auto &bufferVec = extClients[tmparg];
2495
2496 bufferVec.insert(bufferVec.end(), refCount, client);
2497 tmparg->_clientList.Remove(client, true);
2498 }
2499 }
2500
2501 // Loop over value client list of this arg
2502 clientsTmp.assign(tmparg->_clientListValue.begin(), tmparg->_clientListValue.end());
2503 for (auto vclient : clientsTmp) {
2505 cxcoutD(ObjectHandling) << "RooWorkspace::Streamer(" << GetName() << ") element " << tmparg->GetName()
2506 << " has external value client link to " << vclient << " (" << vclient->GetName()
2507 << ") with ref count " << tmparg->_clientListValue.refCount(vclient) << std::endl;
2508
2509 const auto refCount = tmparg->_clientListValue.refCount(vclient);
2511
2512 bufferVec.insert(bufferVec.end(), refCount, vclient);
2513 tmparg->_clientListValue.Remove(vclient, true);
2514 }
2515 }
2516
2517 // Loop over shape client list of this arg
2518 clientsTmp.assign(tmparg->_clientListShape.begin(), tmparg->_clientListShape.end());
2519 for (auto sclient : clientsTmp) {
2521 cxcoutD(ObjectHandling) << "RooWorkspace::Streamer(" << GetName() << ") element " << tmparg->GetName()
2522 << " has external shape client link to " << sclient << " (" << sclient->GetName()
2523 << ") with ref count " << tmparg->_clientListShape.refCount(sclient) << std::endl;
2524
2525 const auto refCount = tmparg->_clientListShape.refCount(sclient);
2527
2528 bufferVec.insert(bufferVec.end(), refCount, sclient);
2529 tmparg->_clientListShape.Remove(sclient, true);
2530 }
2531 }
2532 }
2533
2534 R__b.WriteClassBuffer(RooWorkspace::Class(), this);
2535
2536 // Reinstate clients here
2537
2538 for (auto &iterx : extClients) {
2539 for (auto client : iterx.second) {
2540 iterx.first->_clientList.Add(client);
2541 }
2542 }
2543
2544 for (auto &iterx : extValueClients) {
2545 for (auto client : iterx.second) {
2546 iterx.first->_clientListValue.Add(client);
2547 }
2548 }
2549
2550 for (auto &iterx : extShapeClients) {
2551 for (auto client : iterx.second) {
2552 iterx.first->_clientListShape.Add(client);
2553 }
2554 }
2555 }
2556}
2557
2558
2559
2560
2561////////////////////////////////////////////////////////////////////////////////
2562/// Return STL string with last of class names contained in the code repository
2563
2565{
2566 string ret ;
2567 map<TString,ClassRelInfo>::const_iterator iter = _c2fmap.begin() ;
2568 while(iter!=_c2fmap.end()) {
2569 if (!ret.empty()) {
2570 ret += ", " ;
2571 }
2572 ret += iter->first ;
2573 ++iter ;
2574 }
2575
2576 return ret ;
2577}
2578
2579namespace {
2580UInt_t crc32(const char* data, ULong_t sz, UInt_t crc)
2581{
2582 // update CRC32 with new data
2583
2584 // use precomputed table, rather than computing it on the fly
2585 static const UInt_t crctab[256] = { 0x00000000,
2586 0x04c11db7, 0x09823b6e, 0x0d4326d9, 0x130476dc, 0x17c56b6b,
2587 0x1a864db2, 0x1e475005, 0x2608edb8, 0x22c9f00f, 0x2f8ad6d6,
2588 0x2b4bcb61, 0x350c9b64, 0x31cd86d3, 0x3c8ea00a, 0x384fbdbd,
2589 0x4c11db70, 0x48d0c6c7, 0x4593e01e, 0x4152fda9, 0x5f15adac,
2590 0x5bd4b01b, 0x569796c2, 0x52568b75, 0x6a1936c8, 0x6ed82b7f,
2591 0x639b0da6, 0x675a1011, 0x791d4014, 0x7ddc5da3, 0x709f7b7a,
2592 0x745e66cd, 0x9823b6e0, 0x9ce2ab57, 0x91a18d8e, 0x95609039,
2593 0x8b27c03c, 0x8fe6dd8b, 0x82a5fb52, 0x8664e6e5, 0xbe2b5b58,
2594 0xbaea46ef, 0xb7a96036, 0xb3687d81, 0xad2f2d84, 0xa9ee3033,
2595 0xa4ad16ea, 0xa06c0b5d, 0xd4326d90, 0xd0f37027, 0xddb056fe,
2596 0xd9714b49, 0xc7361b4c, 0xc3f706fb, 0xceb42022, 0xca753d95,
2597 0xf23a8028, 0xf6fb9d9f, 0xfbb8bb46, 0xff79a6f1, 0xe13ef6f4,
2598 0xe5ffeb43, 0xe8bccd9a, 0xec7dd02d, 0x34867077, 0x30476dc0,
2599 0x3d044b19, 0x39c556ae, 0x278206ab, 0x23431b1c, 0x2e003dc5,
2600 0x2ac12072, 0x128e9dcf, 0x164f8078, 0x1b0ca6a1, 0x1fcdbb16,
2601 0x018aeb13, 0x054bf6a4, 0x0808d07d, 0x0cc9cdca, 0x7897ab07,
2602 0x7c56b6b0, 0x71159069, 0x75d48dde, 0x6b93dddb, 0x6f52c06c,
2603 0x6211e6b5, 0x66d0fb02, 0x5e9f46bf, 0x5a5e5b08, 0x571d7dd1,
2604 0x53dc6066, 0x4d9b3063, 0x495a2dd4, 0x44190b0d, 0x40d816ba,
2605 0xaca5c697, 0xa864db20, 0xa527fdf9, 0xa1e6e04e, 0xbfa1b04b,
2606 0xbb60adfc, 0xb6238b25, 0xb2e29692, 0x8aad2b2f, 0x8e6c3698,
2607 0x832f1041, 0x87ee0df6, 0x99a95df3, 0x9d684044, 0x902b669d,
2608 0x94ea7b2a, 0xe0b41de7, 0xe4750050, 0xe9362689, 0xedf73b3e,
2609 0xf3b06b3b, 0xf771768c, 0xfa325055, 0xfef34de2, 0xc6bcf05f,
2610 0xc27dede8, 0xcf3ecb31, 0xcbffd686, 0xd5b88683, 0xd1799b34,
2611 0xdc3abded, 0xd8fba05a, 0x690ce0ee, 0x6dcdfd59, 0x608edb80,
2612 0x644fc637, 0x7a089632, 0x7ec98b85, 0x738aad5c, 0x774bb0eb,
2613 0x4f040d56, 0x4bc510e1, 0x46863638, 0x42472b8f, 0x5c007b8a,
2614 0x58c1663d, 0x558240e4, 0x51435d53, 0x251d3b9e, 0x21dc2629,
2615 0x2c9f00f0, 0x285e1d47, 0x36194d42, 0x32d850f5, 0x3f9b762c,
2616 0x3b5a6b9b, 0x0315d626, 0x07d4cb91, 0x0a97ed48, 0x0e56f0ff,
2617 0x1011a0fa, 0x14d0bd4d, 0x19939b94, 0x1d528623, 0xf12f560e,
2618 0xf5ee4bb9, 0xf8ad6d60, 0xfc6c70d7, 0xe22b20d2, 0xe6ea3d65,
2619 0xeba91bbc, 0xef68060b, 0xd727bbb6, 0xd3e6a601, 0xdea580d8,
2620 0xda649d6f, 0xc423cd6a, 0xc0e2d0dd, 0xcda1f604, 0xc960ebb3,
2621 0xbd3e8d7e, 0xb9ff90c9, 0xb4bcb610, 0xb07daba7, 0xae3afba2,
2622 0xaafbe615, 0xa7b8c0cc, 0xa379dd7b, 0x9b3660c6, 0x9ff77d71,
2623 0x92b45ba8, 0x9675461f, 0x8832161a, 0x8cf30bad, 0x81b02d74,
2624 0x857130c3, 0x5d8a9099, 0x594b8d2e, 0x5408abf7, 0x50c9b640,
2625 0x4e8ee645, 0x4a4ffbf2, 0x470cdd2b, 0x43cdc09c, 0x7b827d21,
2626 0x7f436096, 0x7200464f, 0x76c15bf8, 0x68860bfd, 0x6c47164a,
2627 0x61043093, 0x65c52d24, 0x119b4be9, 0x155a565e, 0x18197087,
2628 0x1cd86d30, 0x029f3d35, 0x065e2082, 0x0b1d065b, 0x0fdc1bec,
2629 0x3793a651, 0x3352bbe6, 0x3e119d3f, 0x3ad08088, 0x2497d08d,
2630 0x2056cd3a, 0x2d15ebe3, 0x29d4f654, 0xc5a92679, 0xc1683bce,
2631 0xcc2b1d17, 0xc8ea00a0, 0xd6ad50a5, 0xd26c4d12, 0xdf2f6bcb,
2632 0xdbee767c, 0xe3a1cbc1, 0xe760d676, 0xea23f0af, 0xeee2ed18,
2633 0xf0a5bd1d, 0xf464a0aa, 0xf9278673, 0xfde69bc4, 0x89b8fd09,
2634 0x8d79e0be, 0x803ac667, 0x84fbdbd0, 0x9abc8bd5, 0x9e7d9662,
2635 0x933eb0bb, 0x97ffad0c, 0xafb010b1, 0xab710d06, 0xa6322bdf,
2636 0xa2f33668, 0xbcb4666d, 0xb8757bda, 0xb5365d03, 0xb1f740b4
2637 };
2638
2639 crc = ~crc;
2640 while (sz--) crc = (crc << 8) ^ UInt_t(*data++) ^ crctab[crc >> 24];
2641
2642 return ~crc;
2643}
2644
2645UInt_t crc32(const char* data)
2646{
2647 // Calculate crc32 checksum on given string
2648 unsigned long sz = strlen(data);
2649 switch (strlen(data)) {
2650 case 0:
2651 return 0;
2652 case 1:
2653 return data[0];
2654 case 2:
2655 return (data[0] << 8) | data[1];
2656 case 3:
2657 return (data[0] << 16) | (data[1] << 8) | data[2];
2658 case 4:
2659 return (data[0] << 24) | (data[1] << 16) | (data[2] << 8) | data[3];
2660 default:
2661 return crc32(data + 4, sz - 4, (data[0] << 24) | (data[1] << 16) |
2662 (data[2] << 8) | data[3]);
2663 }
2664}
2665
2666}
2667
2668////////////////////////////////////////////////////////////////////////////////
2669/// For all classes in the workspace for which no class definition is
2670/// found in the ROOT class table extract source code stored in code
2671/// repository into temporary directory set by
2672/// setClassFileExportDir(), compile classes and link them with
2673/// current ROOT session. If a compilation error occurs print
2674/// instructions for user how to fix errors and recover workspace and
2675/// abort import procedure.
2676
2678{
2679 bool haveDir=false ;
2680
2681 // Retrieve name of directory in which to export code files
2682 string dirName = Form(_classFileExportDir.c_str(),_wspace->uuid().AsString(),_wspace->GetName()) ;
2683
2684 bool writeExtraHeaders(false) ;
2685
2686 // Process all class entries in repository
2687 map<TString,ClassRelInfo>::iterator iter = _c2fmap.begin() ;
2688 while(iter!=_c2fmap.end()) {
2689
2690 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() now processing class " << iter->first.Data() << std::endl ;
2691
2692 // If class is already known, don't load
2693 if (gClassTable->GetDict(iter->first.Data())) {
2694 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Embedded class "
2695 << iter->first << " already in ROOT class table, skipping" << std::endl ;
2696 ++iter ;
2697 continue ;
2698 }
2699
2700 // Check that export directory exists
2701 if (!haveDir) {
2702
2703 // If not, make local directory to extract files
2704 if (!gSystem->AccessPathName(dirName.c_str())) {
2705 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() reusing code export directory " << dirName.c_str()
2706 << " to extract coded embedded in workspace" << std::endl ;
2707 } else {
2708 if (gSystem->MakeDirectory(dirName.c_str())==0) {
2709 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() creating code export directory " << dirName.c_str()
2710 << " to extract coded embedded in workspace" << std::endl ;
2711 } else {
2712 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR creating code export directory " << dirName.c_str()
2713 << " to extract coded embedded in workspace" << std::endl ;
2714 return false ;
2715 }
2716 }
2717 haveDir=true ;
2718
2719 }
2720
2721 // First write any extra header files
2722 if (!writeExtraHeaders) {
2724
2725 map<TString,ExtraHeader>::iterator extraIter = _ehmap.begin() ;
2726 while(extraIter!=_ehmap.end()) {
2727
2728 // Check if identical declaration file (header) is already written
2729 bool needEHWrite=true ;
2730 string fdname = Form("%s/%s",dirName.c_str(),extraIter->second._hname.Data()) ;
2731 ifstream ifdecl(fdname.c_str()) ;
2732 if (ifdecl) {
2733 TString contents ;
2734 char buf[64000];
2735 while (ifdecl.getline(buf, 64000)) {
2736 contents += buf;
2737 contents += "\n";
2738 }
2739 UInt_t crcFile = crc32(contents.Data());
2740 UInt_t crcWS = crc32(extraIter->second._hfile.Data());
2741 needEHWrite = (crcFile != crcWS);
2742 }
2743
2744 // Write declaration file if required
2745 if (needEHWrite) {
2746 oocoutI(_wspace, ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Extracting extra header file "
2747 << fdname << std::endl;
2748
2749 // Extra headers may contain non-existing path - create first to be sure
2751
2752 ofstream fdecl(fdname.c_str());
2753 if (!fdecl) {
2754 oocoutE(_wspace, ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR opening file " << fdname
2755 << " for writing" << std::endl;
2756 return false;
2757 }
2758 fdecl << extraIter->second._hfile.Data();
2759 fdecl.close();
2760 }
2761 ++extraIter;
2762 }
2763 }
2764
2765
2766 // Navigate from class to file
2767 ClassFiles& cfinfo = _fmap[iter->second._fileBase] ;
2768
2769 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() now processing file with base " << iter->second._fileBase << std::endl ;
2770
2771 // If file is already processed, skip to next class
2772 if (cfinfo._extracted) {
2773 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() file with base name " << iter->second._fileBase
2774 << " has already been extracted, skipping to next class" << std::endl ;
2775 continue ;
2776 }
2777
2778 // Check if identical declaration file (header) is already written
2779 bool needDeclWrite=true ;
2780 string fdname = Form("%s/%s.%s",dirName.c_str(),iter->second._fileBase.Data(),cfinfo._hext.Data()) ;
2781 ifstream ifdecl(fdname.c_str()) ;
2782 if (ifdecl) {
2783 TString contents ;
2784 char buf[64000];
2785 while (ifdecl.getline(buf, 64000)) {
2786 contents += buf;
2787 contents += "\n";
2788 }
2789 UInt_t crcFile = crc32(contents.Data()) ;
2790 UInt_t crcWS = crc32(cfinfo._hfile.Data()) ;
2792 }
2793
2794 // Write declaration file if required
2795 if (needDeclWrite) {
2796 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Extracting declaration code of class " << iter->first << ", file " << fdname << std::endl ;
2797 ofstream fdecl(fdname.c_str()) ;
2798 if (!fdecl) {
2799 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR opening file "
2800 << fdname << " for writing" << std::endl ;
2801 return false ;
2802 }
2803 fdecl << cfinfo._hfile ;
2804 fdecl.close() ;
2805 }
2806
2807 // Check if identical implementation file is already written
2808 bool needImplWrite=true ;
2809 string finame = Form("%s/%s.cxx",dirName.c_str(),iter->second._fileBase.Data()) ;
2810 ifstream ifimpl(finame.c_str()) ;
2811 if (ifimpl) {
2812 TString contents ;
2813 char buf[64000];
2814 while (ifimpl.getline(buf, 64000)) {
2815 contents += buf;
2816 contents += "\n";
2817 }
2818 UInt_t crcFile = crc32(contents.Data()) ;
2819 UInt_t crcWS = crc32(cfinfo._cxxfile.Data()) ;
2821 }
2822
2823 // Write implementation file if required
2824 if (needImplWrite) {
2825 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Extracting implementation code of class " << iter->first << ", file " << finame << std::endl ;
2826 ofstream fimpl(finame.c_str()) ;
2827 if (!fimpl) {
2828 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR opening file"
2829 << finame << " for writing" << std::endl ;
2830 return false ;
2831 }
2832 fimpl << cfinfo._cxxfile ;
2833 fimpl.close() ;
2834 }
2835
2836 // Mark this file as extracted
2837 cfinfo._extracted = true ;
2838 oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() marking code unit " << iter->second._fileBase << " as extracted" << std::endl ;
2839
2840 // Compile class
2841 oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Compiling code unit " << iter->second._fileBase.Data() << " to define class " << iter->first << std::endl ;
2842 bool ok = gSystem->CompileMacro(finame.c_str(),"k") ;
2843
2844 if (!ok) {
2845 oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR compiling class " << iter->first.Data() << ", to fix this you can do the following: " << std::endl
2846 << " 1) Fix extracted source code files in directory " << dirName.c_str() << "/" << std::endl
2847 << " 2) In clean ROOT session compiled fixed classes by hand using '.x " << dirName.c_str() << "/ClassName.cxx+'" << std::endl
2848 << " 3) Reopen file with RooWorkspace with broken source code in UPDATE mode. Access RooWorkspace to force loading of class" << std::endl
2849 << " Broken instances in workspace will _not_ be compiled, instead precompiled fixed instances will be used." << std::endl
2850 << " 4) Reimport fixed code in workspace using 'RooWorkspace::importClassCode(\"*\",true)' method, Write() updated workspace to file and close file" << std::endl
2851 << " 5) Reopen file in clean ROOT session to confirm that problems are fixed" << std::endl ;
2852 return false ;
2853 }
2854
2855 ++iter ;
2856 }
2857
2858 return true ;
2859}
2860
2861
2862
2863////////////////////////////////////////////////////////////////////////////////
2864/// Internal access to TDirectory append method
2865
2870
2871
2872////////////////////////////////////////////////////////////////////////////////
2873/// Overload TDirectory interface method to prohibit insertion of objects in read-only directory workspace representation
2874
2876{
2877 if (dynamic_cast<RooAbsArg*>(obj) || dynamic_cast<RooAbsData*>(obj)) {
2878 coutE(ObjectHandling) << "RooWorkspace::WSDir::Add(" << GetName() << ") ERROR: Directory is read-only representation of a RooWorkspace, use RooWorkspace::import() to add objects" << std::endl ;
2879 } else {
2880 InternalAppend(obj) ;
2881 }
2882}
2883
2884
2885////////////////////////////////////////////////////////////////////////////////
2886/// Overload TDirectory interface method to prohibit insertion of objects in read-only directory workspace representation
2887
2889{
2890 if (dynamic_cast<RooAbsArg*>(obj) || dynamic_cast<RooAbsData*>(obj)) {
2891 coutE(ObjectHandling) << "RooWorkspace::WSDir::Add(" << GetName() << ") ERROR: Directory is read-only representation of a RooWorkspace, use RooWorkspace::import() to add objects" << std::endl ;
2892 } else {
2893 InternalAppend(obj) ;
2894 }
2895}
2896
2897
2898////////////////////////////////////////////////////////////////////////////////
2899/// If one of the TObject we have a referenced to is deleted, remove the
2900/// reference.
2901
2903{
2905 if (removedObj == _dir) _dir = nullptr;
2906
2908
2915
2916 std::vector<std::string> invalidSets;
2917
2918 for(auto &c : _namedSets) {
2919 auto const& setName = c.first;
2920 auto& set = c.second;
2921 std::size_t oldSize = set.size();
2922 set.RecursiveRemove(removedObj);
2923 // If the set is used internally by RooFit to cache parameters or
2924 // constraints, it is invalidated by object removal. We will keep track
2925 // of its name to remove the cache set later.
2926 if(set.size() < oldSize && isCacheSet(setName)) {
2927 invalidSets.emplace_back(setName);
2928 }
2929 }
2930
2931 // Remove the sets that got invalidated by the object removal
2932 for(std::string const& setName : invalidSets) {
2933 removeSet(setName.c_str());
2934 }
2935
2936 _eocache.RecursiveRemove(removedObj); // RooExpensiveObjectCache
2937}
#define f(i)
Definition RSha256.hxx:104
#define c(i)
Definition RSha256.hxx:101
#define a(i)
Definition RSha256.hxx:99
#define coutI(a)
#define cxcoutD(a)
#define oocoutW(o, a)
#define oocxcoutD(o, a)
#define coutW(a)
#define oocoutE(o, a)
#define oocoutI(o, a)
#define coutE(a)
#define ooccoutW(o, a)
short Version_t
Class version identifier (short)
Definition RtypesCore.h:79
unsigned long ULong_t
Unsigned long integer 4 bytes (unsigned long). Size depends on architecture.
Definition RtypesCore.h:69
unsigned int UInt_t
Unsigned integer 4 bytes (unsigned int)
Definition RtypesCore.h:60
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
R__EXTERN TClassTable * gClassTable
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
@ kIsAbstract
Definition TDictionary.h:71
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char filename
char name[80]
Definition TGX11.cxx:110
#define gInterpreter
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2496
R__EXTERN TSystem * gSystem
Definition TSystem.h:572
const_iterator begin() const
const_iterator end() const
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:76
static void ioStreamerPass2Finalize()
Method called by workspace container to finalize schema evolution issues that cannot be handled in a ...
A space to attach TBranches.
static TClass * Class()
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
void RecursiveRemove(TObject *obj) override
If one of the TObject we have a referenced to is deleted, remove the reference.
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
void sort(bool reverse=false)
Sort collection using std::sort and name comparison.
RooAbsArg * find(const char *name) const
Find object with given name in list.
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
virtual const RooArgSet * get() const
Definition RooAbsData.h:101
virtual bool changeObservableName(const char *from, const char *to)
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:32
static TClass * Class()
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:63
static TClass * Class()
Abstract base class for RooStudyManager modules.
Definition RooAbsStudy.h:33
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition RooArgSet.h:159
bool containsInstance(const RooAbsArg &var) const override
Check if this exact instance is in this collection.
Definition RooArgSet.h:132
RooArgSet * selectCommon(const RooAbsCollection &refColl) const
Use RooAbsCollection::selecCommon(), but return as RooArgSet.
Definition RooArgSet.h:154
Object to represent discrete states.
Definition RooCategory.h:28
static TClass * Class()
Named container for two doubles, two integers two object points and three string pointers that can be...
Definition RooCmdArg.h:26
Configurable parser for RooCmdArg named arguments.
void defineMutex(const char *head, Args_t &&... tail)
Define arguments where any pair is mutually exclusive.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
bool ok(bool verbose) const
Return true of parsing was successful.
const char * getString(const char *name, const char *defaultValue="", bool convEmptyToNull=false) const
Return string property registered with name 'name'.
bool defineString(const char *name, const char *argName, int stringNum, const char *defValue="", bool appendMode=false)
Define double property name 'name' mapped to double in slot 'stringNum' in RooCmdArg with name argNam...
bool defineInt(const char *name, const char *argName, int intNum, int defValue=0)
Define integer property name 'name' mapped to integer in slot 'intNum' in RooCmdArg with name argName...
int getInt(const char *name, int defaultValue=0) const
Return integer property registered with name 'name'.
static TClass * Class()
Singleton class that serves as repository for objects that are expensive to calculate.
void importCacheObjects(RooExpensiveObjectCache &other, const char *ownerName, bool verbose=false)
Implementation detail of the RooWorkspace.
RooAbsArg * process(const char *expr)
Create a RooFit object from the given expression.
Collection class for internal use, storing a collection of RooAbsArg pointers in a doubly linked list...
bool empty() const
void RecursiveRemove(TObject *obj) override
If one of the TObject we have a referenced to is deleted, remove the reference.
bool Replace(const TObject *oldArg, const TObject *newArg)
Replace object 'oldArg' in collection with new object 'newArg'.
void Delete(Option_t *o=nullptr) override
Remove all elements in collection and delete all elements NB: Collection does not own elements,...
TObject * find(const char *name) const
Return pointer to object with given name in collection.
virtual void Add(TObject *arg)
TObject * FindObject(const char *name) const override
Return pointer to object with given name.
virtual bool Remove(TObject *arg)
Remove object from collection.
static RooMsgService & instance()
Return reference to singleton instance.
static bool setAddDirectoryStatus(bool flag)
Configure whether new instances of RooPlot will add themselves to gDirectory.
Definition RooPlot.cxx:77
Variable that can be changed from the outside.
Definition RooRealVar.h:37
static TClass * Class()
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
static TClass * Class()
The RooStringView is a wrapper around a C-style string that can also be constructed from a std::strin...
An interface to set and retrieve a workspace.
std::map< TString, ExtraHeader > _ehmap
RooWorkspace * _wspace
void Streamer(TBuffer &) override
Custom streamer for the workspace.
std::string listOfClassNames() const
Return STL string with last of class names contained in the code repository.
bool autoImportClass(TClass *tc, bool doReplace=false)
Import code of class 'tc' into the repository.
bool compileClasses()
For all classes in the workspace for which no class definition is found in the ROOT class table extra...
std::map< TString, ClassRelInfo > _c2fmap
std::map< TString, ClassFiles > _fmap
void InternalAppend(TObject *obj)
Internal access to TDirectory append method.
void Add(TObject *, bool) override
Overload TDirectory interface method to prohibit insertion of objects in read-only directory workspac...
TClass * IsA() const override
void Append(TObject *, bool) override
Overload TDirectory interface method to prohibit insertion of objects in read-only directory workspac...
Persistable container for RooFit projects.
RooExpensiveObjectCache _eocache
Cache for expensive objects.
TObject * obj(RooStringView name) const
Return any type of object (RooAbsArg, RooAbsData or generic object) with given name)
RooLinkedList _genObjects
List of generic objects.
static bool _autoClass
static std::list< std::string > _classDeclDirList
const RooArgSet * getSnapshot(const char *name) const
Return the RooArgSet containing a snapshot of variables contained in the workspace.
static void addClassDeclImportDir(const char *dir)
Add dir to search path for class declaration (header) files.
void Print(Option_t *opts=nullptr) const override
Print contents of the workspace.
RooLinkedList _dataList
List of owned datasets.
RooAbsCategory * catfunc(RooStringView name) const
Retrieve discrete function (RooAbsCategory) with given name. A null pointer is returned if not found.
WSDir * _dir
! Transient ROOT directory representation of workspace
static void addClassImplImportDir(const char *dir)
Add dir to search path for class implementation (.cxx) files.
RooAbsPdf * pdf(RooStringView name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
std::map< std::string, RooArgSet > _namedSets
Map of named RooArgSets.
RooAbsData * embeddedData(RooStringView name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found.
RooCategory * cat(RooStringView name) const
Retrieve discrete variable (RooCategory) with given name. A null pointer is returned if not found.
void clearStudies()
Remove all RooStudyManager modules.
bool renameSet(const char *name, const char *newName)
Rename set to a new name.
std::unique_ptr< RooFactoryWSTool > _factory
! Factory tool associated with workspace
RooArgSet allVars() const
Return set with all variable objects.
RooArgSet argSet(RooStringView nameList) const
Return set of RooAbsArgs matching to given list of names.
bool writeToFile(const char *fileName, bool recreate=true)
Save this current workspace into given file.
const RooArgSet * set(RooStringView name)
Return pointer to previously defined named set with given nmame If no such set is found a null pointe...
bool cd(const char *path=nullptr)
RooArgSet allCats() const
Return set with all category objects.
void RecursiveRemove(TObject *obj) override
If one of the TObject we have a referenced to is deleted, remove the reference.
RooAbsArg * fundArg(RooStringView name) const
Return fundamental (i.e.
RooLinkedList _views
List of model views.
bool commitTransaction()
~RooWorkspace() override
Workspace destructor.
bool cancelTransaction()
Cancel an ongoing import transaction.
bool startTransaction()
Open an import transaction operations.
TObject * Clone(const char *newname="") const override
TObject::Clone() needs to be overridden.
RooArgSet allResolutionModels() const
Return set with all resolution model objects.
RooLinkedList _snapshots
List of parameter snapshots.
bool saveSnapshot(RooStringView, const char *paramNames)
Save snapshot of values and attributes (including "Constant") of given parameters.
RooArgSet allPdfs() const
Return set with all probability density function objects.
void Streamer(TBuffer &) override
Stream an object of class RooWorkspace.
TObject * genobj(RooStringView name) const
Return generic object with given name.
std::list< RooAbsData * > allData() const
Return list of all dataset in the workspace.
RooLinkedList _studyMods
List if StudyManager modules.
std::list< TObject * > allGenericObjects() const
Return list of all generic objects in the workspace.
static void setClassFileExportDir(const char *dir=nullptr)
Specify the name of the directory in which embedded source code is unpacked and compiled.
bool importClassCode(const char *pat="*", bool doReplace=false)
Import code of all classes in the workspace that have a class name that matches pattern 'pat' and whi...
bool makeDir()
Create transient TDirectory representation of this workspace.
RooArgSet allCatFunctions() const
Return set with all category function objects.
static std::string _classFileExportDir
static std::list< std::string > _classImplDirList
RooAbsReal * function(RooStringView name) const
Retrieve function (RooAbsReal) with given name. Note that all RooAbsPdfs are also RooAbsReals....
RooAbsArg * arg(RooStringView name) const
Return RooAbsArg with given name. A null pointer is returned if none is found.
RooWorkspace()
Default constructor.
bool removeSet(const char *name)
Remove a named set from the workspace.
static TClass * Class()
CodeRepo _classes
RooArgSet allFunctions() const
Return set with all function objects.
RooFactoryWSTool & factory()
Return instance to factory tool.
bool extendSet(const char *name, const char *newContents)
Define a named set in the workspace through a comma separated list of names of objects already in the...
RooExpensiveObjectCache & expensiveObjectCache()
RooArgSet _sandboxNodes
! Sandbox for incoming objects in a transaction
bool defineSetInternal(const char *name, const RooArgSet &aset)
bool _openTrans
! Is there a transaction open?
RooRealVar * var(RooStringView name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found.
bool addStudy(RooAbsStudy &study)
Insert RooStudyManager module.
static void autoImportClassCode(bool flag)
If flag is true, source code of classes not the ROOT distribution is automatically imported if on obj...
RooLinkedList _embeddedDataList
List of owned datasets that are embedded in pdfs.
RooArgSet _allOwnedNodes
List of owned pdfs and components.
RooAbsData * data(RooStringView name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found.
std::list< RooAbsData * > allEmbeddedData() const
Return list of all dataset in the workspace.
bool loadSnapshot(const char *name)
Load the values and attributes of the parameters in the snapshot saved with the given name.
bool defineSet(const char *name, const RooArgSet &aset, bool importMissing=false)
Define a named RooArgSet with given constituents.
bool import(const RooAbsArg &arg, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={})
Import a RooAbsArg object, e.g.
Buffer base class used for serializing objects.
Definition TBuffer.h:43
static DictFuncPtr_t GetDict(const char *cname)
Given the class name returns the Dictionary() function of a class (uses hash of name).
TClass instances represent classes, structs and namespaces in the ROOT type system.
Definition TClass.h:84
Bool_t cd() override
Change current directory to "this" directory.
virtual void Append(TObject *obj, Bool_t replace=kFALSE)
Append object to this directory.
A ROOT file is an on-disk file, usually with extension .root, that stores objects in a file-system-li...
Definition TFile.h:131
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition TFile.cxx:3765
static void AddDirectory(Bool_t add=kTRUE)
Sets the flag controlling the automatic add of histograms in memory.
Definition TH1.cxx:1264
A doubly linked list.
Definition TList.h:38
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:50
Mother of all ROOT objects.
Definition TObject.h:41
virtual void RecursiveRemove(TObject *obj)
Recursively remove this object from a list.
Definition TObject.cxx:679
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:227
virtual Int_t Write(const char *name=nullptr, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
Definition TObject.cxx:965
Regular expression class.
Definition TRegexp.h:31
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
Ssiz_t Index(const char *pat, Ssiz_t i=0, ECaseCompare cmp=kExact) const
Definition TString.h:659
virtual const char * Getenv(const char *env)
Get environment variable.
Definition TSystem.cxx:1678
virtual char * ConcatFileName(const char *dir, const char *name)
Concatenate a directory and a file name.
Definition TSystem.cxx:1084
virtual int MakeDirectory(const char *name)
Make a directory.
Definition TSystem.cxx:838
virtual Bool_t AccessPathName(const char *path, EAccessMode mode=kFileExists)
Returns FALSE if one can access a file using the specified access mode.
Definition TSystem.cxx:1309
virtual const char * BaseName(const char *pathname)
Base name of a file name. Base name of /user/root is root.
Definition TSystem.cxx:946
virtual int CompileMacro(const char *filename, Option_t *opt="", const char *library_name="", const char *build_dir="", UInt_t dirmode=0)
This method compiles and loads a shared library containing the code from the file "filename".
Definition TSystem.cxx:2849
virtual TString GetDirName(const char *pathname)
Return the directory name in pathname.
Definition TSystem.cxx:1044
const Int_t n
Definition legend1.C:16
void(* DirAutoAdd_t)(void *, TDirectory *)
Definition Rtypes.h:120
std::vector< std::string > Split(std::string_view str, std::string_view delims, bool skipEmpty=false)
Splits a string at each character in delims.
MsgLevel
Verbosity level for RooMsgService::StreamConfig in RooMsgService.