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