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