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