RooWorkspace.cxx

Go to the documentation of this file.
00001 /*****************************************************************************
00002  * Project: RooFit                                                           *
00003  * Package: RooFitCore                                                       *
00004  * @(#)root/roofitcore:$Id: RooWorkspace.cxx 37380 2010-12-07 21:26:55Z wouter $
00005  * Authors:                                                                  *
00006  *   WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu       *
00007  *                                                                           *
00008  * Copyright (c) 2000-2005, Regents of the University of California          *
00009  *                          and Stanford University. All rights reserved.    *
00010  *                                                                           *
00011  * Redistribution and use in source and binary forms,                        *
00012  * with or without modification, are permitted according to the terms        *
00013  * listed in LICENSE (http://roofit.sourceforge.net/license.txt)             *
00014  *****************************************************************************/
00015 
00016 //////////////////////////////////////////////////////////////////////////////
00017 //
00018 // BEGIN_HTML
00019 // The RooWorkspace is a persistable container for RooFit projects. A workspace
00020 // can contain and own variables, p.d.f.s, functions and datasets. All objects
00021 // that live in the workspace are owned by the workspace. The import() method
00022 // enforces consistency of objects upon insertion into the workspace (e.g. no
00023 // duplicate object with the same name are allowed) and makes sure all objects
00024 // in the workspace are connected to each other. Easy accessor methods like
00025 // pdf(), var() and data() allow to refer to the contents of the workspace by
00026 // object name. The entire RooWorkspace can be saved into a ROOT TFile and organises
00027 // the consistent streaming of its contents without duplication.
00028 // <p>
00029 // If a RooWorkspace contains custom classes, i.e. classes not in the 
00030 // ROOT distribution, portability of workspaces can be enhanced by
00031 // storing the source code of those classes in the workspace as well.
00032 // This process is also organized by the workspace through the
00033 // importClassCode() method.
00034 // END_HTML
00035 //
00036 
00037 #include "RooFit.h"
00038 #include "RooWorkspace.h"
00039 #include "RooAbsPdf.h"
00040 #include "RooRealVar.h"
00041 #include "RooCategory.h"
00042 #include "RooAbsData.h"
00043 #include "RooCmdConfig.h"
00044 #include "RooMsgService.h"
00045 #include "RooConstVar.h"
00046 #include "RooResolutionModel.h"
00047 #include "TInterpreter.h"
00048 #include "TClassTable.h"
00049 #include "TBaseClass.h"
00050 #include "TSystem.h"
00051 #include "TRegexp.h"
00052 #include "RooFactoryWSTool.h"
00053 #include "RooAbsStudy.h"
00054 #include "RooTObjWrap.h"
00055 #include "RooAbsOptTestStatistic.h"
00056 #include "TROOT.h"
00057 #include "TFile.h"
00058 #include "TH1.h"
00059 #include "Api.h"
00060 #include <map>
00061 #include <string>
00062 #include <list>
00063 #include <set>
00064 
00065 using namespace std ;
00066 
00067 
00068 #if ROOT_VERSION_CODE <= ROOT_VERSION(5,19,02)
00069 #include "Api.h"
00070 #endif
00071 
00072 
00073 #include "TClass.h"
00074 #include "Riostream.h"
00075 #include <string.h>
00076 #include <assert.h>
00077 
00078 ClassImp(RooWorkspace)
00079 ;
00080 
00081 //_____________________________________________________________________________
00082 ClassImp(RooWorkspace::CodeRepo)
00083 ;
00084 
00085 //_____________________________________________________________________________
00086 ClassImp(RooWorkspace::WSDir)
00087 ;
00088 
00089 list<string> RooWorkspace::_classDeclDirList ;
00090 list<string> RooWorkspace::_classImplDirList ;
00091 string RooWorkspace::_classFileExportDir = ".wscode.%s.%s" ;
00092 Bool_t RooWorkspace::_autoClass = kFALSE ;
00093 
00094 
00095 //_____________________________________________________________________________
00096 void RooWorkspace::addClassDeclImportDir(const char* dir) 
00097 {
00098   // Add 'dir' to search path for class declaration (header) files, when
00099   // attempting to import class code with importClassClode()
00100 
00101   _classDeclDirList.push_back(dir) ;
00102 }
00103 
00104 
00105 //_____________________________________________________________________________
00106 void RooWorkspace::addClassImplImportDir(const char* dir) 
00107 {
00108   // Add 'dir' to search path for class implementation (.cxx) files, when
00109   // attempting to import class code with importClassClode()
00110 
00111   _classImplDirList.push_back(dir) ;
00112 }
00113 
00114 
00115 //_____________________________________________________________________________
00116 void RooWorkspace::setClassFileExportDir(const char* dir) 
00117 {
00118   // Specify the name of the directory in which embedded source
00119   // code is unpacked and compiled. The specified string may contain
00120   // one '%s' token which will be substituted by the workspace name
00121 
00122   if (dir) {
00123     _classFileExportDir = dir ;
00124   } else {
00125     _classFileExportDir = ".wscode.%s.%s" ;
00126   }
00127 }
00128 
00129 
00130 //_____________________________________________________________________________
00131 void RooWorkspace::autoImportClassCode(Bool_t flag) 
00132 {
00133   // If flag is true, source code of classes not the the ROOT distribution
00134   // is automatically imported if on object of such a class is imported
00135   // in the workspace
00136   _autoClass = flag ; 
00137 }
00138 
00139 
00140 
00141 //_____________________________________________________________________________
00142 RooWorkspace::RooWorkspace() : _classes(this), _dir(0), _factory(0), _doExport(kFALSE), _openTrans(kFALSE)
00143 {
00144   // Default constructor
00145 }
00146 
00147 
00148 
00149 //_____________________________________________________________________________
00150 RooWorkspace::RooWorkspace(const char* name, const char* title) : 
00151   TNamed(name,title?title:name), _classes(this), _dir(0), _factory(0), _doExport(kFALSE), _openTrans(kFALSE)
00152 {
00153   // Construct empty workspace with given name and title
00154 }
00155 
00156 
00157 RooWorkspace::RooWorkspace(const char* name, Bool_t doCINTExport)  : 
00158   TNamed(name,name), _classes(this), _dir(0), _factory(0), _doExport(kFALSE), _openTrans(kFALSE)
00159 {
00160   // Construct empty workspace with given name and option to export reference to all workspace contents to a CINT namespace with the same name
00161   if (doCINTExport) {
00162     exportToCint(name) ;
00163   }
00164 }
00165 
00166 
00167 //_____________________________________________________________________________
00168 RooWorkspace::RooWorkspace(const RooWorkspace& other) : 
00169   TNamed(other), _uuid(other._uuid), _classes(this), _dir(0), _factory(0), _doExport(kFALSE), _openTrans(kFALSE)
00170 {
00171   // Workspace copy constructor
00172 
00173   // Copy owned nodes
00174   other._allOwnedNodes.snapshot(_allOwnedNodes,kTRUE) ;
00175 
00176   // Copy datasets
00177   TIterator* iter = other._dataList.MakeIterator() ;
00178   TObject* data2 ;
00179   while((data2=iter->Next())) {
00180     _dataList.Add(data2->Clone()) ;
00181   }
00182   delete iter ;
00183 
00184   // Copy snapshots
00185   TIterator* iter2 = other._snapshots.MakeIterator() ;
00186   RooArgSet* snap ;
00187   while((snap=(RooArgSet*)iter2->Next())) {
00188     RooArgSet* snapClone = (RooArgSet*) snap->snapshot() ;
00189     snapClone->setName(snap->GetName()) ;
00190     _snapshots.Add(snapClone) ;
00191   }
00192   delete iter2 ;
00193 
00194   // Copy named sets
00195   for (map<string,RooArgSet>::const_iterator iter3 = other._namedSets.begin() ; iter3 != other._namedSets.end() ; ++iter3) {
00196     // Make RooArgSet with equivalent content of this workspace
00197     RooArgSet* tmp = (RooArgSet*) _allOwnedNodes.selectCommon(iter3->second) ;
00198     _namedSets[iter3->first].add(*tmp) ;
00199     delete tmp ;
00200   }
00201 
00202   // Copy generic objects
00203   TIterator* iter4 = other._genObjects.MakeIterator() ;
00204   TObject* gobj ;
00205   while((gobj=iter4->Next())) {
00206     _genObjects.Add(gobj->Clone()) ;
00207   }
00208   delete iter4 ;
00209   
00210 }
00211 
00212 
00213 
00214 //_____________________________________________________________________________
00215 RooWorkspace::~RooWorkspace() 
00216 {
00217   // Workspace destructor
00218 
00219   // Delete references to variables that were declared in CINT 
00220   if (_doExport) {
00221     unExport() ;
00222   }
00223 
00224   // Delete contents
00225   _dataList.Delete() ;
00226   if (_dir) {
00227     delete _dir ;
00228   }
00229   _snapshots.Delete() ;
00230 
00231   // WVE named sets too?
00232 
00233   _genObjects.Delete() ;
00234 }
00235 
00236 
00237 //_____________________________________________________________________________
00238 Bool_t RooWorkspace::import(const char* fileSpec, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3) 
00239 {
00240   // Import a RooAbsArg or RooAbsData set from a workspace in a file. Filespec should be constructed as "filename:wspacename:objectname"
00241   // The arguments will be passed on to the relevant RooAbsArg& or RooAbsData& import call
00242 
00243   // Parse file/workspace/objectname specification
00244   char buf[1024] ;
00245   strlcpy(buf,fileSpec,1024) ;
00246   char* filename = strtok(buf,":") ;
00247   char* wsname = strtok(0,":") ;
00248   char* objname = strtok(0,":") ;
00249 
00250   // Check that parsing was successful
00251   if (!filename||!wsname||!objname) {
00252     coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR in file specification, expecting for 'filename:wsname:objname'" << endl ;
00253     return kTRUE ;
00254   }
00255 
00256   // Check that file can be opened
00257   TFile* f = TFile::Open(filename) ;
00258   if (f==0) {
00259     coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR opening file " << filename << endl ;
00260     return 0 ;
00261   }
00262 
00263   // That that file contains workspace
00264   RooWorkspace* w = dynamic_cast<RooWorkspace*>(f->Get(wsname)) ;
00265   if (w==0) {
00266     coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR: No object named " << wsname << " in file " << filename 
00267                           << " or object is not a RooWorkspace" << endl ;
00268     return 0 ;
00269   }
00270 
00271   // Check that workspace contains object and forward to appropriate import method
00272   RooAbsArg* warg = w->arg(objname) ;
00273   if (warg) {
00274     Bool_t ret = import(*warg,arg1,arg2,arg3) ;
00275     delete f ;
00276     return ret ;    
00277   }
00278   RooAbsData* wdata = w->data(objname) ;
00279   if (wdata) {
00280     Bool_t ret = import(*wdata,arg1,arg2,arg3) ;
00281     delete f ;
00282     return ret ;    
00283   }
00284 
00285   coutE(InputArguments) << "RooWorkspace(" << GetName() << ") ERROR: No RooAbsArg or RooAbsData object named " << objname 
00286                         << " in workspace " << wsname << " in file " << filename << endl ;
00287   return kTRUE ;  
00288 }
00289 
00290 
00291 //_____________________________________________________________________________
00292 Bool_t RooWorkspace::import(const RooArgSet& args, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3) 
00293 {
00294   // Import multiple RooAbsArg objects into workspace. For details on arguments see documentation
00295   // of import() method for single RooAbsArg
00296 
00297   TIterator* iter = args.createIterator() ;
00298   RooAbsArg* oneArg ;
00299   Bool_t ret(kFALSE) ;
00300   while((oneArg=(RooAbsArg*)iter->Next())) {
00301     ret |= import(*oneArg,arg1,arg2,arg3) ;
00302   }
00303   return ret ;
00304 }
00305 
00306 
00307 
00308 //_____________________________________________________________________________
00309 Bool_t RooWorkspace::import(const RooAbsArg& inArg, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3) 
00310 {
00311   //  Import a RooAbsArg object, e.g. function, p.d.f or variable into the workspace. This import function clones the input argument and will
00312   //  own the clone. If a composite object is offered for import, e.g. a p.d.f with parameters and observables, the
00313   //  complete tree of objects is imported. If any of the _variables_ of a composite object (parameters/observables) are already 
00314   //  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) 
00315   //  to be imported already exists in the workspace an error message is printed and the import of the entire tree of objects is cancelled. 
00316   //  Several optional arguments can be provided to modify the import procedure.
00317   //
00318   //  Accepted arguments
00319   //  -------------------------------
00320   //  RenameConflictNodes(const char* suffix) -- Add suffix to branch node name if name conflicts with existing node in workspace
00321   //  RenameAllNodes(const char* suffix) -- Add suffix to all branch node names including top level node
00322   //  RenameAllVariables(const char* suffix) -- Add suffix to all variables names
00323   //  RenameAllVariablesExcept(const char* suffix, const char* exceptionList) -- Add suffix to all variables names, except ones listed
00324   //  RenameVariable(const char* inputName, const char* outputName) -- Rename variable as specified upon import.
00325   //  RecycleConflictNodes() -- If any of the function objects to be imported already exist in the name space, connect the
00326   //                            imported expression to the already existing nodes. WARNING: use with care! If function definitions
00327   //                            do not match, this alters the definition of your function upon import
00328   //  Silence() -- Do not issue any info message
00329   //
00330   //  The RenameConflictNodes, RenameNodes and RecycleConflictNodes arguments are mutually exclusive. The RenameVariable argument can be repeated
00331   //  as often as necessary to rename multiple variables. Alternatively, a single RenameVariable argument can be given with
00332   //  two comma separated lists.
00333 
00334   RooLinkedList args ;
00335   args.Add((TObject*)&arg1) ;
00336   args.Add((TObject*)&arg2) ;
00337   args.Add((TObject*)&arg3) ;
00338 
00339   // Select the pdf-specific commands 
00340   RooCmdConfig pc(Form("RooWorkspace::import(%s)",GetName())) ;
00341 
00342   pc.defineString("conflictSuffix","RenameConflictNodes",0) ;
00343   pc.defineString("allSuffix","RenameAllNodes",0) ;
00344   pc.defineString("allVarsSuffix","RenameAllVariables",0) ;
00345   pc.defineString("allVarsExcept","RenameAllVariables",1) ;
00346   pc.defineString("varChangeIn","RenameVar",0,"",kTRUE) ;
00347   pc.defineString("varChangeOut","RenameVar",1,"",kTRUE) ;
00348   pc.defineString("factoryTag","FactoryTag",0) ;
00349   pc.defineInt("useExistingNodes","RecycleConflictNodes",0,0) ;
00350   pc.defineInt("silence","Silence",0,0) ;
00351   pc.defineMutex("RenameConflictNodes","RenameAllNodes") ;
00352   pc.defineMutex("RenameConflictNodes","RecycleConflictNodes") ;
00353   pc.defineMutex("RenameAllNodes","RecycleConflictNodes") ;
00354   pc.defineMutex("RenameVariable","RenameAllVariables") ;
00355 
00356   // Process and check varargs 
00357   pc.process(args) ;
00358   if (!pc.ok(kTRUE)) {
00359     return kTRUE ;
00360   }
00361 
00362   // Decode renaming logic into suffix string and boolean for conflictOnly mode
00363   const char* suffixC = pc.getString("conflictSuffix") ;
00364   const char* suffixA = pc.getString("allSuffix") ;
00365   const char* suffixV = pc.getString("allVarsSuffix") ;
00366   const char* exceptVars = pc.getString("allVarsExcept") ;
00367   const char* varChangeIn = pc.getString("varChangeIn") ;
00368   const char* varChangeOut = pc.getString("varChangeOut") ;
00369   Int_t useExistingNodes = pc.getInt("useExistingNodes") ;
00370   Int_t silence = pc.getInt("silence") ;
00371 
00372   // Turn zero length strings into null pointers 
00373   if (suffixC && strlen(suffixC)==0) suffixC = 0 ;
00374   if (suffixA && strlen(suffixA)==0) suffixA = 0 ;
00375 
00376   Bool_t conflictOnly = suffixA ? kFALSE : kTRUE ;
00377   const char* suffix = suffixA ? suffixA : suffixC ;
00378 
00379   // Process any change in variable names 
00380   map<string,string> varMap ;
00381   if (strlen(varChangeIn)>0) {
00382     
00383     // Parse comma separated lists into map<string,string>
00384     char tmp[1024] ;
00385     strlcpy(tmp,varChangeIn,1024) ;
00386     list<string> tmpIn,tmpOut ;
00387     char* ptr = strtok(tmp,",") ;
00388     while (ptr) {
00389       tmpIn.push_back(ptr) ;
00390       ptr = strtok(0,",") ;
00391     }
00392     strlcpy(tmp,varChangeOut,1024) ;
00393     ptr = strtok(tmp,",") ;
00394     while (ptr) {
00395       tmpOut.push_back(ptr) ;
00396       ptr = strtok(0,",") ;
00397     }    
00398     list<string>::iterator iin = tmpIn.begin() ;
00399     list<string>::iterator iout = tmpOut.begin() ;
00400     for (;iin!=tmpIn.end() ; ++iin,++iout) {
00401       varMap[*iin]=*iout ;
00402     }       
00403   }
00404 
00405   // Process RenameAllVariables argument if specified  
00406   // First convert exception list if provided
00407   std::set<string> exceptVarNames ;
00408   char tmp[1024] ;
00409   if (exceptVars && strlen(exceptVars)) {
00410     strlcpy(tmp,exceptVars,1024) ;
00411     char* ptr = strtok(tmp,",") ;
00412     while(ptr) {
00413       exceptVarNames.insert(ptr) ;
00414       ptr = strtok(0,",") ;
00415     }
00416   }
00417 
00418   if (suffixV != 0 && strlen(suffixV)>0) {
00419     RooArgSet* vars = inArg.getVariables() ;
00420     TIterator* iter = vars->createIterator() ;
00421     RooAbsArg* v ;
00422     while((v=(RooAbsArg*)iter->Next())) {
00423       if (exceptVarNames.find(v->GetName())==exceptVarNames.end()) {
00424         varMap[v->GetName()] = Form("%s_%s",v->GetName(),suffixV) ;
00425       }
00426     }
00427     delete iter ;
00428     delete vars ;
00429   }
00430   
00431   // Scan for overlaps with current contents
00432   RooAbsArg* wsarg = _allOwnedNodes.find(inArg.GetName()) ;
00433 
00434   // Check for factory specification match
00435   const char* tagIn = inArg.getStringAttribute("factory_tag") ;
00436   const char* tagWs = wsarg ? wsarg->getStringAttribute("factory_tag") : 0 ;
00437   Bool_t factoryMatch = (tagIn && tagWs && !strcmp(tagIn,tagWs)) ;
00438   if (factoryMatch) {
00439     ((RooAbsArg&)inArg).setAttribute("RooWorkspace::Recycle") ;
00440   }
00441 
00442   if (!suffix && wsarg && !useExistingNodes && !(inArg.isFundamental() && varMap[inArg.GetName()]!="")) {
00443     if (!factoryMatch) {
00444       if (wsarg!=&inArg) {
00445         coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR importing object named " << inArg.GetName() 
00446                               << ": another instance with same name already in the workspace and no conflict resolution protocol specified" << endl ;
00447         return kTRUE ;    
00448       } else {
00449         if (!silence) {
00450           coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") Object " << inArg.GetName() << " is already in workspace!" << endl ;
00451         }
00452         return kTRUE ;    
00453       }
00454     } else {
00455       coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") Recycling existing object " << inArg.GetName() << " created with identical factory specification" << endl ;
00456     }
00457   }
00458 
00459   // Make list of conflicting nodes
00460   RooArgSet conflictNodes ;
00461   RooArgSet branchSet ;
00462   inArg.branchNodeServerList(&branchSet) ;
00463   TIterator* iter = branchSet.createIterator() ;
00464   RooAbsArg* branch ;
00465   while ((branch=(RooAbsArg*)iter->Next())) {
00466     RooAbsArg* wsbranch = _allOwnedNodes.find(branch->GetName()) ;
00467     if (wsbranch && wsbranch!=branch && !branch->getAttribute("RooWorkspace::Recycle") && !useExistingNodes) {
00468       conflictNodes.add(*branch) ;
00469     }
00470   }
00471   delete iter ;
00472   
00473   // Terminate here if there are conflicts and no resolution protocol
00474   if (conflictNodes.getSize()>0 && !suffix && !useExistingNodes) {
00475       coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR object named " << inArg.GetName() << ": component(s) " 
00476            << conflictNodes << " already in the workspace and no conflict resolution protocol specified" << endl ;      
00477       return kTRUE ;
00478   }
00479     
00480   // Now create a working copy of the incoming object tree
00481   RooArgSet* cloneSet = (RooArgSet*) RooArgSet(inArg).snapshot(kTRUE) ;
00482   RooAbsArg* cloneTop = cloneSet->find(inArg.GetName()) ;
00483 
00484   // Mark all nodes for renaming if we are not in conflictOnly mode
00485   if (!conflictOnly) {
00486     conflictNodes.removeAll() ;
00487     conflictNodes.add(branchSet) ;
00488   }
00489 
00490   // Mark nodes that are to be renamed with special attribute
00491   TIterator* citer = conflictNodes.createIterator() ;
00492   string topName2 = cloneTop->GetName() ;
00493   RooAbsArg* cnode ;
00494   while ((cnode=(RooAbsArg*)citer->Next())) {
00495     RooAbsArg* cnode2 = cloneSet->find(cnode->GetName()) ;
00496     string origName = cnode2->GetName() ;
00497     cnode2->SetName(Form("%s_%s",cnode2->GetName(),suffix)) ;
00498     cnode2->SetTitle(Form("%s (%s)",cnode2->GetTitle(),suffix)) ;
00499     string tag = Form("ORIGNAME:%s",origName.c_str()) ;
00500     cnode2->setAttribute(tag.c_str()) ;
00501 
00502     // Save name of new top level node for later use
00503     if (cnode2==cloneTop) {
00504       topName2 = cnode2->GetName() ;      
00505     }
00506 
00507     if (!silence) {
00508       coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() 
00509                             << ") Resolving name conflict in workspace by changing name of imported node  " 
00510                             << origName << " to " << cnode2->GetName() << endl ;
00511     }
00512   }  
00513   delete citer ;
00514 
00515   // Process any change in variable names 
00516   if (strlen(varChangeIn)>0 || (suffixV && strlen(suffixV)>0)) {
00517     
00518     // Process all changes in variable names
00519     TIterator* cliter = cloneSet->createIterator() ;
00520     while ((cnode=(RooAbsArg*)cliter->Next())) {
00521       
00522       if (varMap.find(cnode->GetName())!=varMap.end()) {        
00523         string origName = cnode->GetName() ;
00524         cnode->SetName(varMap[cnode->GetName()].c_str()) ;
00525         string tag = Form("ORIGNAME:%s",origName.c_str()) ;
00526         cnode->setAttribute(tag.c_str()) ;
00527         if (!silence) {
00528           coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") Changing name of variable " 
00529                                 << origName << " to " << cnode->GetName() << " on request" << endl ;
00530         }
00531 
00532         if (cnode==cloneTop) {
00533           topName2 = cnode->GetName() ;
00534         }
00535 
00536       }    
00537     }
00538     delete cliter ;
00539   }
00540   
00541   // Now clone again with renaming effective
00542   RooArgSet* cloneSet2 = (RooArgSet*) RooArgSet(*cloneTop).snapshot(kTRUE) ;
00543   RooAbsArg* cloneTop2 = cloneSet2->find(topName2.c_str()) ;
00544 
00545   // Make final check list of conflicting nodes
00546   RooArgSet conflictNodes2 ;
00547   RooArgSet branchSet2 ;
00548   inArg.branchNodeServerList(&branchSet) ;
00549   TIterator* iter2 = branchSet2.createIterator() ;
00550   RooAbsArg* branch2 ;
00551   while ((branch2=(RooAbsArg*)iter2->Next())) {
00552     if (_allOwnedNodes.find(branch2->GetName())) {
00553       conflictNodes2.add(*branch2) ;
00554     }
00555   }
00556   delete iter2 ;
00557 
00558   // Terminate here if there are conflicts and no resolution protocol
00559   if (conflictNodes2.getSize()) {
00560     coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR object named " << inArg.GetName() << ": component(s) " 
00561                           << conflictNodes2 << " cause naming conflict after conflict resolution protocol was executed" << endl ;      
00562     return kTRUE ;
00563   }
00564     
00565   // Print a message for each imported node
00566   iter = cloneSet2->createIterator() ;
00567   
00568   // Perform any auxiliary imports at this point
00569   RooAbsArg* node ;
00570   while((node=(RooAbsArg*)iter->Next())) {
00571     if (node->importWorkspaceHook(*this)) {
00572       coutE(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") ERROR object named " << node->GetName() 
00573                             << " has an error in importing in one or more of its auxiliary objects, aborting" << endl ;
00574       return kTRUE ;
00575     }
00576   }
00577   iter->Reset() ;
00578 
00579   RooArgSet recycledNodes ;
00580   RooArgSet nodesToBeDeleted ;
00581   while((node=(RooAbsArg*)iter->Next())) {
00582 
00583     if (_autoClass) {
00584       if (!_classes.autoImportClass(node->IsA())) {
00585         coutW(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") WARNING: problems import class code of object " 
00586                               << node->IsA()->GetName() << "::" << node->GetName() << ", reading of workspace will require external definition of class" << endl ;
00587       }
00588     }
00589 
00590     // Point expensiveObjectCache to copy in this workspace
00591     RooExpensiveObjectCache& oldCache = node->expensiveObjectCache() ;
00592     node->setExpensiveObjectCache(_eocache) ;    
00593     _eocache.importCacheObjects(oldCache,node->GetName(),kTRUE) ;
00594 
00595     // Check if node is already in workspace (can only happen for variables or identical instances, unless RecycleConflictNodes is specified)
00596     RooAbsArg* wsnode = _allOwnedNodes.find(node->GetName()) ;
00597 
00598     if (wsnode) {
00599       // Do not import node, add not to list of nodes that require reconnection
00600       if (!silence && useExistingNodes) {
00601         coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") using existing copy of " << node->IsA()->GetName() 
00602                               << "::" << node->GetName() << " for import of " << cloneTop2->IsA()->GetName() << "::" 
00603                               << cloneTop2->GetName() << endl ;      
00604       }
00605       recycledNodes.add(*_allOwnedNodes.find(node->GetName())) ;
00606 
00607       // Delete clone of incoming node
00608       nodesToBeDeleted.addOwned(*node) ;
00609       
00610     } else {
00611       // Import node
00612       if (!silence) {
00613         coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") importing " << node->IsA()->GetName() << "::" 
00614                               << node->GetName() << endl ;
00615       }
00616       _allOwnedNodes.addOwned(*node) ;
00617       if (_openTrans) {
00618         _sandboxNodes.add(*node) ;
00619       } else {
00620         if (_dir && node->IsA() != RooConstVar::Class()) {      
00621           _dir->InternalAppend(node) ;
00622         }
00623         if (_doExport && node->IsA() != RooConstVar::Class()) {
00624           exportObj(node) ;
00625         }
00626       }
00627     }
00628   }
00629 
00630   // Release working copy
00631   delete cloneSet ;
00632 
00633   // Reconnect any nodes that need to be
00634   if (recycledNodes.getSize()>0) {
00635     iter->Reset() ;
00636     while((node=(RooAbsArg*)iter->Next())) {
00637       node->redirectServers(recycledNodes) ;
00638     }
00639   }
00640   delete iter ;
00641 
00642   cloneSet2->releaseOwnership() ;
00643   delete cloneSet2 ;  
00644 
00645   return kFALSE ;
00646 }
00647 
00648 
00649 
00650 //_____________________________________________________________________________
00651 Bool_t RooWorkspace::import(RooAbsData& inData, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3) 
00652 {
00653   //  Import a dataset (RooDataSet or RooDataHist) into the work space. The workspace will contain a copy of the data
00654   //  The dataset and its variables can be renamed upon insertion with the options below
00655   //
00656   //  Accepted arguments
00657   //  -------------------------------
00658   //  Rename(const char* suffix) -- Rename dataset upon insertion
00659   //  RenameVariable(const char* inputName, const char* outputName) -- Change names of observables in dataset upon insertion
00660 
00661   coutI(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") importing dataset " << inData.GetName() << endl ;
00662 
00663   RooLinkedList args ;
00664   args.Add((TObject*)&arg1) ;
00665   args.Add((TObject*)&arg2) ;
00666   args.Add((TObject*)&arg3) ;
00667 
00668   // Select the pdf-specific commands 
00669   RooCmdConfig pc(Form("RooWorkspace::import(%s)",GetName())) ;
00670 
00671   pc.defineString("dsetName","Rename",0,"") ;
00672   pc.defineString("varChangeIn","RenameVar",0,"",kTRUE) ;
00673   pc.defineString("varChangeOut","RenameVar",1,"",kTRUE) ;
00674 
00675   // Process and check varargs 
00676   pc.process(args) ;
00677   if (!pc.ok(kTRUE)) {
00678     return kTRUE ;
00679   }
00680 
00681   // Decode renaming logic into suffix string and boolean for conflictOnly mode
00682   const char* dsetName = pc.getString("dsetName") ;
00683   const char* varChangeIn = pc.getString("varChangeIn") ;
00684   const char* varChangeOut = pc.getString("varChangeOut") ;
00685 
00686   // Transform emtpy string into null pointer
00687   if (dsetName && strlen(dsetName)==0) {
00688     dsetName=0 ;
00689   }
00690 
00691   // Check that no dataset with target name already exists
00692   if (dsetName && _dataList.FindObject(dsetName)) {
00693     coutE(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") ERROR dataset with name " << dsetName << " already exists in workspace, import aborted" << endl ;
00694     return kTRUE ;
00695   }
00696   if (!dsetName && _dataList.FindObject(inData.GetName())) {
00697     coutE(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") ERROR dataset with name " << inData.GetName() << " already exists in workspace, import aborted" << endl ;
00698     return kTRUE ;
00699   }
00700 
00701   // Rename dataset if required
00702   RooAbsData* clone ;
00703   if (dsetName) {
00704     coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") changing name of dataset from  " << inData.GetName() << " to " << dsetName << endl ;
00705     clone = (RooAbsData*) inData.Clone(dsetName) ;
00706   } else {
00707     clone = (RooAbsData*) inData.Clone(inData.GetName()) ;
00708   }
00709 
00710 
00711   // Process any change in variable names 
00712   if (strlen(varChangeIn)>0) {
00713     
00714     // Parse comma separated lists of variable name changes
00715     char tmp[1024] ;
00716     strlcpy(tmp,varChangeIn,1024) ;
00717     list<string> tmpIn,tmpOut ;
00718     char* ptr = strtok(tmp,",") ;
00719     while (ptr) {
00720       tmpIn.push_back(ptr) ;
00721       ptr = strtok(0,",") ;
00722     }
00723     strlcpy(tmp,varChangeOut,1024) ;
00724     ptr = strtok(tmp,",") ;
00725     while (ptr) {
00726       tmpOut.push_back(ptr) ;
00727       ptr = strtok(0,",") ;
00728     }    
00729     list<string>::iterator iin = tmpIn.begin() ;
00730     list<string>::iterator iout = tmpOut.begin() ;
00731 
00732     for (; iin!=tmpIn.end() ; ++iin,++iout) {
00733       coutI(ObjectHandling) << "RooWorkSpace::import(" << GetName() << ") changing name of dataset observable " << *iin << " to " << *iout << endl ;
00734       clone->changeObservableName(iin->c_str(),iout->c_str()) ;
00735     }
00736   }
00737 
00738   // Now import the dataset observables
00739   TIterator* iter = clone->get()->createIterator() ;
00740   RooAbsArg* carg ;
00741   while((carg=(RooAbsArg*)iter->Next())) {
00742     if (!arg(carg->GetName())) {
00743       import(*carg) ;
00744     }
00745   }
00746   delete iter ;
00747     
00748   _dataList.Add(clone) ;
00749   if (_dir) {   
00750     _dir->InternalAppend(clone) ;
00751   }
00752   if (_doExport) {
00753     exportObj(clone) ;
00754   }
00755   return kFALSE ;
00756 }
00757 
00758 
00759 
00760 
00761 //_____________________________________________________________________________
00762 Bool_t RooWorkspace::defineSet(const char* name, const RooArgSet& aset, Bool_t importMissing) 
00763 {
00764   // Define a named RooArgSet with given constituents. If importMissing is true, any constituents
00765   // of aset that are not in the workspace will be imported, otherwise an error is returned
00766   // for missing components
00767 
00768   // Check if set was previously defined, if so print warning
00769   map<string,RooArgSet>::iterator i = _namedSets.find(name) ;
00770   if (i!=_namedSets.end()) {
00771     coutW(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") WARNING redefining previously defined named set " << name << endl ;
00772   }
00773 
00774   RooArgSet wsargs ;
00775 
00776   // Check all constituents of provided set
00777   TIterator* iter = aset.createIterator() ;
00778   RooAbsArg* sarg ;
00779   while((sarg=(RooAbsArg*)iter->Next())) {
00780     // If missing, either import or report error
00781     if (!arg(sarg->GetName())) {
00782       if (importMissing) {
00783         import(*sarg) ;
00784       } else {
00785         coutE(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") ERROR set constituent " << sarg->GetName() 
00786                               << " is not in workspace and importMissing option is disabled" << endl ;
00787         return kTRUE ;
00788       }
00789     }
00790     wsargs.add(*arg(sarg->GetName())) ;    
00791   }
00792   delete iter ;
00793 
00794   // Install named set
00795   _namedSets[name].removeAll() ;
00796   _namedSets[name].add(wsargs) ;
00797    
00798   return kFALSE ;
00799 }
00800 
00801 
00802 
00803 
00804 //_____________________________________________________________________________
00805 Bool_t RooWorkspace::defineSet(const char* name, const char* contentList) 
00806 {
00807   // Define a named set in the work space through a comma separated list of
00808   // names of objects already in the workspace
00809 
00810   // Check if set was previously defined, if so print warning
00811   map<string,RooArgSet>::iterator i = _namedSets.find(name) ;
00812   if (i!=_namedSets.end()) {
00813     coutW(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") WARNING redefining previously defined named set " << name << endl ;
00814   }
00815 
00816   RooArgSet wsargs ;
00817 
00818   // Check all constituents of provided set
00819   char buf[1024] ;
00820   strlcpy(buf,contentList,1024) ;
00821   char* token = strtok(buf,",") ;
00822   while(token) {
00823     // If missing, either import or report error
00824     if (!arg(token)) {
00825       coutE(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") ERROR proposed set constituent " << token
00826                             << " is not in workspace" << endl ;
00827       return kTRUE ;
00828     }
00829     wsargs.add(*arg(token)) ;    
00830     token = strtok(0,",") ;
00831   }
00832 
00833   // Install named set
00834   _namedSets[name].removeAll() ;
00835   _namedSets[name].add(wsargs) ;
00836    
00837   return kFALSE ;
00838 }
00839 
00840 
00841 
00842 
00843 //_____________________________________________________________________________
00844 Bool_t RooWorkspace::extendSet(const char* name, const char* newContents) 
00845 {
00846   // Define a named set in the work space through a comma separated list of
00847   // names of objects already in the workspace
00848 
00849   RooArgSet wsargs ;
00850 
00851   // Check all constituents of provided set
00852   char buf[1024] ;
00853   strlcpy(buf,newContents,1024) ;
00854   char* token = strtok(buf,",") ;
00855   while(token) {
00856     // If missing, either import or report error
00857     if (!arg(token)) {
00858       coutE(InputArguments) << "RooWorkspace::defineSet(" << GetName() << ") ERROR proposed set constituent " << token
00859                             << " is not in workspace" << endl ;
00860       return kTRUE ;
00861     }
00862     wsargs.add(*arg(token)) ;    
00863     token = strtok(0,",") ;
00864   }
00865 
00866   // Extend named set
00867   _namedSets[name].add(wsargs,kTRUE) ;
00868    
00869   return kFALSE ;
00870 }
00871 
00872 
00873 
00874 //_____________________________________________________________________________
00875 const RooArgSet* RooWorkspace::set(const char* name) 
00876 {
00877   // Return pointer to previously defined named set with given nmame
00878   // If no such set is found a null pointer is returned
00879 
00880   map<string,RooArgSet>::iterator i = _namedSets.find(name) ;
00881   return (i!=_namedSets.end()) ? &(i->second) : 0 ;
00882 }
00883 
00884 
00885 
00886 
00887 //_____________________________________________________________________________
00888 Bool_t RooWorkspace::startTransaction() 
00889 {
00890   // Open an import transaction operations. Returns kTRUE if successful, kFALSE
00891   // if there is already an ongoing transaction
00892 
00893   // Check that there was no ongoing transaction
00894   if (_openTrans) {
00895     return kFALSE ;
00896   } 
00897 
00898   // Open transaction
00899   _openTrans = kTRUE ;
00900   return kTRUE ;
00901 }
00902 
00903 Bool_t RooWorkspace::cancelTransaction() 
00904 {
00905   // Cancel an ongoing import transaction. All objects imported since startTransaction()
00906   // will be removed and the transaction will be terminated. Return kTRUE if cancel operation
00907   // succeeds, return kFALSE if there was no open transaction
00908   
00909   // Check that there is an ongoing transaction
00910   if (!_openTrans) {
00911     return kFALSE ;
00912   }
00913 
00914   // Delete all objects in the sandbox
00915   TIterator* iter = _sandboxNodes.createIterator() ;
00916   RooAbsArg* tmpArg ;
00917   while((tmpArg=(RooAbsArg*)iter->Next())) {
00918     _allOwnedNodes.remove(*tmpArg) ;
00919   }
00920   delete iter ;
00921   _sandboxNodes.removeAll() ;
00922   
00923   // Mark transaction as finished
00924   _openTrans = kFALSE ;
00925 
00926   return kTRUE ;
00927 }
00928 
00929 Bool_t RooWorkspace::commitTransaction() 
00930 {
00931   // Commit an ongoing import transaction. Returns kTRUE if commit succeeded,
00932   // return kFALSE if there was no ongoing transaction
00933 
00934   // Check that there is an ongoing transaction
00935   if (!_openTrans) {
00936     return kFALSE ;
00937   }
00938 
00939   // Publish sandbox nodes in directory and/or CINT if requested 
00940   TIterator* iter = _sandboxNodes.createIterator() ;
00941   RooAbsArg* sarg ;
00942   while((sarg=(RooAbsArg*)iter->Next())) {
00943     if (_dir && sarg->IsA() != RooConstVar::Class()) {  
00944       _dir->InternalAppend(sarg) ;
00945     }
00946     if (_doExport && sarg->IsA() != RooConstVar::Class()) {
00947       exportObj(sarg) ;
00948     }    
00949   }
00950   delete iter ;
00951 
00952   // Remove all commited objects from the sandbox
00953   _sandboxNodes.removeAll() ;
00954 
00955   // Mark transaction as finished
00956   _openTrans = kFALSE ;
00957 
00958   return kTRUE ;
00959 }
00960 
00961 
00962 
00963 
00964 //_____________________________________________________________________________
00965 Bool_t RooWorkspace::importClassCode(TClass* theClass, Bool_t doReplace) 
00966 {
00967   return _classes.autoImportClass(theClass,doReplace) ;
00968 }
00969 
00970 
00971 
00972 //_____________________________________________________________________________
00973 Bool_t RooWorkspace::importClassCode(const char* pat, Bool_t doReplace)  
00974 {
00975   // Inport code of all classes in the workspace that have a class name
00976   // that matches pattern 'pat' and which are not found to be part of
00977   // the standard ROOT distribution. If doReplace is true any existing
00978   // class code saved in the workspace is replaced
00979 
00980   Bool_t ret(kTRUE) ;
00981 
00982   TRegexp re(pat,kTRUE) ;
00983   TIterator* iter = componentIterator() ;
00984   RooAbsArg* carg ;
00985   while((carg=(RooAbsArg*)iter->Next())) {
00986     TString className = carg->IsA()->GetName() ;
00987     if (className.Index(re)>=0 && !_classes.autoImportClass(carg->IsA(),doReplace)) {
00988       coutW(ObjectHandling) << "RooWorkspace::import(" << GetName() << ") WARNING: problems import class code of object " 
00989                             << carg->IsA()->GetName() << "::" << carg->GetName() << ", reading of workspace will require external definition of class" << endl ;
00990       ret = kFALSE ;
00991     }
00992   }  
00993   delete iter ;
00994   return ret ;
00995 }
00996 
00997 
00998 
00999 
01000 
01001 //_____________________________________________________________________________
01002 Bool_t RooWorkspace::saveSnapshot(const char* name, const char* paramNames) 
01003 {
01004   // Save snapshot of values and attributes (including "Constant") of parameters 'params'
01005   // If importValues is FALSE, the present values from the object in the workspace are
01006   // saved. If importValues is TRUE, the values of the objects passed in the 'params'
01007   // argument are saved
01008 
01009   return saveSnapshot(name,argSet(paramNames),kFALSE) ;
01010 }
01011 
01012 
01013 
01014 
01015 
01016 //_____________________________________________________________________________
01017 Bool_t RooWorkspace::saveSnapshot(const char* name, const RooArgSet& params, Bool_t importValues) 
01018 {
01019   // Save snapshot of values and attributes (including "Constant") of parameters 'params'
01020   // If importValues is FALSE, the present values from the object in the workspace are
01021   // saved. If importValues is TRUE, the values of the objects passed in the 'params'
01022   // argument are saved
01023 
01024   RooArgSet* actualParams = (RooArgSet*) _allOwnedNodes.selectCommon(params) ;
01025   RooArgSet* snapshot = (RooArgSet*) actualParams->snapshot() ;
01026   delete actualParams ;
01027 
01028   snapshot->setName(name) ;
01029 
01030   if (importValues) {
01031     *snapshot = params ;
01032   }
01033 
01034   RooArgSet* oldSnap = (RooArgSet*) _snapshots.FindObject(name) ;
01035   if (oldSnap) {
01036     coutI(ObjectHandling) << "RooWorkspace::saveSnaphot(" << GetName() << ") replacing previous snapshot with name " << name << endl ;
01037     _snapshots.Remove(oldSnap) ;
01038     delete oldSnap ;
01039   }
01040 
01041   _snapshots.Add(snapshot) ;
01042 
01043   return kTRUE ;
01044 }
01045 
01046 
01047 
01048 
01049 //_____________________________________________________________________________
01050 Bool_t RooWorkspace::loadSnapshot(const char* name) 
01051 {
01052   // Load the values and attributes of the parameters in the snapshot saved with
01053   // the given name
01054 
01055   RooArgSet* snap = (RooArgSet*) _snapshots.find(name) ;
01056   if (!snap) {
01057     coutE(ObjectHandling) << "RooWorkspace::loadSnapshot(" << GetName() << ") no snapshot with name " << name << " is available" << endl ;
01058     return kFALSE ;
01059   }
01060 
01061   RooArgSet* actualParams = (RooArgSet*) _allOwnedNodes.selectCommon(*snap) ;
01062   *actualParams = *snap ;
01063   delete actualParams ;
01064 
01065   return kTRUE ;
01066 }
01067 
01068 
01069 
01070 // //_____________________________________________________________________________
01071 // RooAbsPdf* RooWorkspace::joinPdf(const char* jointPdfName, const char* indexName, const char* inputMapping) 
01072 // {
01073 //   // Join given list of p.d.f.s into a simultaneous p.d.f with given name. If the named index category
01074 //   // does not exist, it is created. 
01075 //   //
01076 //   //  Example : joinPdf("simPdf","expIndex","A=pdfA,B=pdfB") ;
01077 //   //
01078 //   //            will return a RooSimultaneous named 'simPdf' with index category 'expIndex' with
01079 //   //            state names A and B. Pdf 'pdfA' will be associated with state A and pdf 'pdfB'
01080 //   //            will be associated with state B
01081 //   //
01082 //   return 0 ;
01083 // }
01084 
01085 // //_____________________________________________________________________________
01086 // RooAbsData* RooWorkspace::joinData(const char* jointDataName, const char* indexName, const char* inputMapping) 
01087 // {
01088 //   // Join given list of dataset into a joint dataset for use with a simultaneous pdf
01089 //   // (as e.g. created by joingPdf"
01090 //   //
01091 //   //  Example : joinData("simData","expIndex","A=dataA,B=dataB") ;
01092 //   //
01093 //   //            will return a RooDataSet named 'simData' that consist of the entries of both
01094 //   //            dataA and dataB. An extra category column 'expIndex' is added that labels
01095 //   //            each entry with state 'A' and 'B' to indicate the originating dataset
01096 //   return 0 ;
01097 // }
01098 
01099 
01100 //_____________________________________________________________________________
01101 RooAbsPdf* RooWorkspace::pdf(const char* name) const
01102 { 
01103   // Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found
01104 
01105   return dynamic_cast<RooAbsPdf*>(_allOwnedNodes.find(name)) ; 
01106 }
01107 
01108 
01109 //_____________________________________________________________________________
01110 RooAbsReal* RooWorkspace::function(const char* name) const 
01111 { 
01112   // Retrieve function (RooAbsReal) with given name. Note that all RooAbsPdfs are also RooAbsReals. A null pointer is returned if not found.
01113 
01114   return dynamic_cast<RooAbsReal*>(_allOwnedNodes.find(name)) ; 
01115 }
01116 
01117 
01118 //_____________________________________________________________________________
01119 RooRealVar* RooWorkspace::var(const char* name) const
01120 { 
01121   // Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found
01122 
01123   return dynamic_cast<RooRealVar*>(_allOwnedNodes.find(name)) ; 
01124 }
01125 
01126 
01127 //_____________________________________________________________________________
01128 RooCategory* RooWorkspace::cat(const char* name) const
01129 { 
01130   // Retrieve discrete variable (RooCategory) with given name. A null pointer is returned if not found
01131 
01132   return dynamic_cast<RooCategory*>(_allOwnedNodes.find(name)) ; 
01133 }
01134 
01135 
01136 //_____________________________________________________________________________
01137 RooAbsCategory* RooWorkspace::catfunc(const char* name) const
01138 {
01139   // Retrieve discrete function (RooAbsCategory) with given name. A null pointer is returned if not found
01140   return dynamic_cast<RooAbsCategory*>(_allOwnedNodes.find(name)) ; 
01141 }
01142 
01143 
01144 
01145 //_____________________________________________________________________________
01146 RooAbsArg* RooWorkspace::arg(const char* name) const
01147 {
01148   // Return RooAbsArg with given name. A null pointer is returned if none is found.
01149   return _allOwnedNodes.find(name) ;
01150 }
01151 
01152 
01153 
01154 //_____________________________________________________________________________
01155 RooArgSet RooWorkspace::argSet(const char* nameList) const
01156 {
01157   // Return set of RooAbsArgs matching to given list of names
01158   RooArgSet ret ;
01159 
01160   char tmp[1024] ;
01161   strlcpy(tmp,nameList,1024) ;
01162   char* token = strtok(tmp,",") ;
01163   while(token) {
01164     RooAbsArg* oneArg = arg(token) ;
01165     if (oneArg) {
01166       ret.add(*oneArg) ;
01167     } else {
01168       coutE(InputArguments) << " RooWorkspace::argSet(" << GetName() << ") no RooAbsArg named " << token << " in workspace" << endl ;
01169     }
01170     token = strtok(0,",") ; 
01171   }
01172   return ret ;
01173 }
01174 
01175 
01176 
01177 //_____________________________________________________________________________
01178 RooAbsArg* RooWorkspace::fundArg(const char* name) const
01179 {
01180   // Return fundamental (i.e. non-derived) RooAbsArg with given name. Fundamental types
01181   // are e.g. RooRealVar, RooCategory. A null pointer is returned if none is found.
01182   RooAbsArg* tmp = arg(name) ;
01183   if (!tmp) {
01184     return 0 ;
01185   }
01186   return tmp->isFundamental() ? tmp : 0 ;
01187 }
01188 
01189 
01190 
01191 //_____________________________________________________________________________
01192 RooAbsData* RooWorkspace::data(const char* name) const
01193 {
01194   // Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found
01195 
01196   return (RooAbsData*)_dataList.FindObject(name) ;
01197 }
01198 
01199 
01200 
01201 
01202 //_____________________________________________________________________________
01203 RooArgSet RooWorkspace::allVars() const
01204 {
01205   // Return set with all variable objects
01206   RooArgSet ret ;
01207 
01208   // Split list of components in pdfs, functions and variables
01209   TIterator* iter = _allOwnedNodes.createIterator() ;
01210   RooAbsArg* parg ;
01211   while((parg=(RooAbsArg*)iter->Next())) {
01212     if (parg->IsA()->InheritsFrom(RooRealVar::Class())) {
01213       ret.add(*parg) ;
01214     }
01215   }
01216   delete iter ;
01217 
01218   return ret ;
01219 }
01220 
01221 
01222 //_____________________________________________________________________________
01223 RooArgSet RooWorkspace::allCats() const
01224 {
01225   // Return set with all category objects
01226   RooArgSet ret ;
01227 
01228   // Split list of components in pdfs, functions and variables
01229   TIterator* iter = _allOwnedNodes.createIterator() ;
01230   RooAbsArg* parg ;
01231   while((parg=(RooAbsArg*)iter->Next())) {
01232     if (parg->IsA()->InheritsFrom(RooCategory::Class())) {
01233       ret.add(*parg) ;
01234     }
01235   }
01236   delete iter ;
01237 
01238   return ret ;
01239 }
01240 
01241 
01242 
01243 //_____________________________________________________________________________
01244 RooArgSet RooWorkspace::allFunctions() const
01245 {
01246   // Return set with all function objects
01247   RooArgSet ret ;
01248 
01249   // Split list of components in pdfs, functions and variables
01250   TIterator* iter = _allOwnedNodes.createIterator() ;
01251   RooAbsArg* parg ;
01252   while((parg=(RooAbsArg*)iter->Next())) {
01253     if (parg->IsA()->InheritsFrom(RooAbsReal::Class()) && 
01254         !parg->IsA()->InheritsFrom(RooAbsPdf::Class()) && 
01255         !parg->IsA()->InheritsFrom(RooConstVar::Class()) && 
01256         !parg->IsA()->InheritsFrom(RooRealVar::Class())) {
01257       ret.add(*parg) ;
01258     }
01259   }
01260 
01261   return ret ;
01262 }
01263 
01264 
01265 //_____________________________________________________________________________
01266 RooArgSet RooWorkspace::allCatFunctions() const
01267 {
01268   // Return set with all category function objects
01269   RooArgSet ret ;
01270 
01271   // Split list of components in pdfs, functions and variables
01272   TIterator* iter = _allOwnedNodes.createIterator() ;
01273   RooAbsArg* parg ;
01274   while((parg=(RooAbsArg*)iter->Next())) {  
01275     if (parg->IsA()->InheritsFrom(RooAbsCategory::Class()) && 
01276         !parg->IsA()->InheritsFrom(RooCategory::Class())) {
01277       ret.add(*parg) ;
01278     }
01279   }
01280   return ret ;
01281 }
01282 
01283 
01284 
01285 //_____________________________________________________________________________
01286 RooArgSet RooWorkspace::allResolutionModels() const
01287 {
01288   // Return set with all resolution model objects
01289   RooArgSet ret ;
01290 
01291   // Split list of components in pdfs, functions and variables
01292   TIterator* iter = _allOwnedNodes.createIterator() ;
01293   RooAbsArg* parg ;
01294   while((parg=(RooAbsArg*)iter->Next())) {  
01295     if (parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
01296       if (!((RooResolutionModel*)parg)->isConvolved()) {
01297         ret.add(*parg) ;
01298       }
01299     }
01300   }
01301   return ret ;
01302 }
01303 
01304 
01305 //_____________________________________________________________________________
01306 RooArgSet RooWorkspace::allPdfs() const
01307 {
01308   // Return set with all probability density function objects
01309   RooArgSet ret ;
01310 
01311   // Split list of components in pdfs, functions and variables
01312   TIterator* iter = _allOwnedNodes.createIterator() ;
01313   RooAbsArg* parg ;
01314   while((parg=(RooAbsArg*)iter->Next())) {  
01315     if (parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
01316         !parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
01317       ret.add(*parg) ;
01318     }
01319   }
01320   return ret ;
01321 }
01322 
01323 
01324 
01325 //_____________________________________________________________________________
01326 list<RooAbsData*> RooWorkspace::allData() const 
01327 {
01328   // Return list of all dataset in the workspace
01329 
01330   list<RooAbsData*> ret ;
01331   TIterator* iter = _dataList.MakeIterator() ;
01332   RooAbsData* dat ;
01333   while((dat=(RooAbsData*)iter->Next())) {
01334     ret.push_back(dat) ;
01335   }
01336   delete iter ;
01337   return ret ;
01338 }
01339 
01340 
01341 
01342 //_____________________________________________________________________________
01343 list<TObject*> RooWorkspace::allGenericObjects() const 
01344 {
01345   // Return list of all generic objects in the workspace
01346 
01347   list<TObject*> ret ;
01348   TIterator* iter = _genObjects.MakeIterator() ;
01349   TObject* gobj ;
01350   while((gobj=(RooAbsData*)iter->Next())) {    
01351 
01352     // If found object is wrapper, return payload
01353     if (gobj->IsA()==RooTObjWrap::Class()) {
01354       ret.push_back(((RooTObjWrap*)gobj)->obj()) ;
01355     } else {
01356       ret.push_back(gobj) ;
01357     }
01358   }
01359   delete iter ;
01360   return ret ;
01361 }
01362 
01363 
01364 
01365 
01366 //_____________________________________________________________________________
01367 Bool_t RooWorkspace::CodeRepo::autoImportClass(TClass* tc, Bool_t doReplace) 
01368 {
01369   // Import code of class 'tc' into the repository. If code is already in repository it is only imported
01370   // again if doReplace is false. The names and location of the source files is determined from the information
01371   // in TClass. If no location is found in the TClass information, the files are searched in the workspace
01372   // search path, defined by addClassDeclImportDir() and addClassImplImportDir() for declaration and implementation
01373   // files respectively. If files cannot be found, abort with error status, otherwise update the internal
01374   // class-to-file map and import the contents of the files, if they are not imported yet.
01375 
01376 
01377   oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") request to import code of class " << tc->GetName() << endl ;
01378 
01379   // *** PHASE 1 *** Check if file needs to be imported, or is in ROOT distribution, and check if it can be persisted
01380 
01381   // Check if we already have the class (i.e. it is in the classToFile map)
01382   if (!doReplace && _c2fmap.find(tc->GetName())!=_c2fmap.end()) {
01383     oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") code of class " << tc->GetName() << " already imported, skipping" << endl ;
01384     return kTRUE ;
01385   }
01386 
01387   // Retrieve file names through ROOT TClass interface
01388   string implfile = tc->GetImplFileName() ;
01389   string declfile = tc->GetDeclFileName() ;
01390 
01391   // Check that file names are not empty
01392   if (implfile.empty() || declfile.empty()) {
01393     oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") ERROR: cannot retrieve code file names for class " 
01394                                    << tc->GetName() << " through ROOT TClass interface, unable to import code" << endl ;
01395     return kFALSE ;
01396   }
01397 
01398   // Check if header filename is found in ROOT distribution, if so, do not import class
01399   TString rootsys = gSystem->Getenv("ROOTSYS") ;
01400   if (TString(implfile.c_str()).Index(rootsys)>=0) {
01401     oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo(" << _wspace->GetName() << ") code of class " << tc->GetName() << " is in ROOT distribution, skipping " << endl ;
01402     return kTRUE ;
01403   }
01404   const char* implpath=0 ;
01405 
01406   // Require that class meets technical criteria to be persistable (i.e it has a default ctor)
01407   // (We also need a default ctor of abstract classes, but cannot check that through is interface
01408   //  as TClass::HasDefaultCtor only returns true for callable default ctors)
01409 #if ROOT_VERSION_CODE <= ROOT_VERSION(5,19,02)
01410   if (!(tc->GetClassInfo()->Property()&G__BIT_ISABSTRACT) && !tc->HasDefaultConstructor()) {
01411 #else
01412   if (!(gCint->ClassInfo_Property(tc->GetClassInfo())&G__BIT_ISABSTRACT) && !tc->HasDefaultConstructor()) {
01413 #endif
01414     oocoutW(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() << ") WARNING cannot import class " 
01415                                     << tc->GetName() << " : it cannot be persisted because it doesn't have a default constructor. Please fix " << endl ;
01416     return kFALSE ;      
01417   }
01418 
01419 
01420   // *** PHASE 2 *** Check if declaration and implementation files can be located 
01421 
01422   char* declpath = 0 ;
01423 
01424   // Check if header file can be found in specified location
01425   // If not, scan through list of 'class declaration' paths in RooWorkspace
01426   if (gSystem->AccessPathName(declfile.c_str())) {
01427 
01428     // Check list of additional declaration paths
01429     list<string>::iterator diter = RooWorkspace::_classDeclDirList.begin() ;
01430 
01431     while(diter!= RooWorkspace::_classDeclDirList.end()) {
01432       
01433       declpath = gSystem->ConcatFileName(diter->c_str(),declfile.c_str()) ;      
01434       if (!gSystem->AccessPathName(declpath)) {
01435         // found declaration file
01436         break ;
01437       }
01438       // cleanup and continue ;
01439       delete[] declpath ;
01440       declpath=0 ;
01441 
01442       ++diter ;
01443     }
01444     
01445     // Header file cannot be found anywhere, warn user and abort operation
01446     if (!declpath) {
01447       oocoutW(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() << ") WARNING Cannot access code of class " 
01448                                       << tc->GetName() << " because header file " << declfile << " is not found in current directory nor in $ROOTSYS" ;
01449       if (_classDeclDirList.size()>0) {
01450         ooccoutW(_wspace,ObjectHandling) << ", nor in the search path " ;
01451         diter = RooWorkspace::_classDeclDirList.begin() ;
01452 
01453         while(diter!= RooWorkspace::_classDeclDirList.end()) {
01454 
01455           if (diter!=RooWorkspace::_classDeclDirList.begin()) {
01456             ooccoutW(_wspace,ObjectHandling) << "," ;
01457           }
01458           ooccoutW(_wspace,ObjectHandling) << diter->c_str() ;
01459           ++diter ;
01460         }
01461       }
01462       ooccoutW(_wspace,ObjectHandling) << ". To fix this problem add the required directory to the search "
01463                                        << "path using RooWorkspace::addClassDeclDir(const char* dir)" << endl ;
01464       
01465       return kFALSE ;
01466     }
01467   }
01468 
01469   
01470   // Check if implementation file can be found in specified location
01471   // If not, scan through list of 'class implementation' paths in RooWorkspace
01472   if (gSystem->AccessPathName(implfile.c_str())) {
01473 
01474     // Check list of additional declaration paths
01475     list<string>::iterator iiter = RooWorkspace::_classImplDirList.begin() ;
01476 
01477     while(iiter!= RooWorkspace::_classImplDirList.end()) {
01478       
01479       implpath = gSystem->ConcatFileName(iiter->c_str(),implfile.c_str()) ;      
01480       if (!gSystem->AccessPathName(implpath)) {
01481         // found implementation file
01482         break ;
01483       }
01484       // cleanup and continue ;
01485       delete[] implpath ;
01486       implpath=0 ;
01487 
01488       ++iiter ;
01489     }
01490      
01491     // Implementation file cannot be found anywhere, warn user and abort operation
01492     if (!implpath) {
01493       oocoutW(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() << ") WARNING Cannot access code of class " 
01494                                       << tc->GetName() << " because implementation file " << implfile << " is not found in current directory nor in $ROOTSYS" ;
01495       if (_classDeclDirList.size()>0) {
01496         ooccoutW(_wspace,ObjectHandling) << ", nor in the search path " ;
01497         iiter = RooWorkspace::_classImplDirList.begin() ;
01498 
01499         while(iiter!= RooWorkspace::_classImplDirList.end()) {
01500 
01501           if (iiter!=RooWorkspace::_classImplDirList.begin()) {
01502             ooccoutW(_wspace,ObjectHandling) << "," ;
01503           }
01504           ooccoutW(_wspace,ObjectHandling) << iiter->c_str() ;
01505           ++iiter ;
01506         }
01507       }
01508       ooccoutW(_wspace,ObjectHandling) << ". To fix this problem add the required directory to the search "
01509                                        << "path using RooWorkspace::addClassImplDir(const char* dir)" << endl ;    
01510       return kFALSE ;
01511     }
01512   }
01513   
01514   char buf[1024] ;
01515 
01516   // *** Phase 3 *** Prepare to import code from files into STL string buffer
01517   //
01518   // Code storage is organized in two linked maps
01519   //
01520   // _fmap contains stl strings with code, indexed on declaration file name
01521   //
01522   // _c2fmap contains list of declaration file names and list of base classes
01523   //                  and is indexed on class name
01524   //
01525   // Phase 3 is skipped if fmap already contains an entry with given filebasename
01526 
01527   string declfilename = declpath?gSystem->BaseName(declpath):gSystem->BaseName(declfile.c_str()) ;
01528 
01529   // Split in base and extension
01530   int dotpos2 = strrchr(declfilename.c_str(),'.') - declfilename.c_str() ;
01531   string declfilebase = declfilename.substr(0,dotpos2) ;
01532   string declfileext = declfilename.substr(dotpos2+1) ;
01533 
01534   // If file has not beed stored yet, enter stl strings with implementation and declaration in file map
01535   if (_fmap.find(declfilebase) == _fmap.end()) {
01536 
01537     // Open declaration file
01538     fstream fdecl(declpath?declpath:declfile.c_str()) ;
01539     
01540     // Abort import if declaration file cannot be opened
01541     if (!fdecl) {
01542       oocoutE(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() 
01543                                       << ") ERROR opening declaration file " <<  declfile << endl ;
01544       return kFALSE ;      
01545     }
01546     
01547     oocoutI(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() 
01548                                     << ") importing code of class " << tc->GetName() 
01549                                     << " from " << (implpath?implpath:implfile.c_str()) 
01550                                     << " and " << (declpath?declpath:declfile.c_str()) << endl ;
01551     
01552     
01553     // Read entire file into an stl string
01554     string decl ;
01555     while(fdecl.getline(buf,1023)) {
01556       decl += buf ;
01557       decl += '\n' ;
01558     }
01559     
01560     // Open implementation file
01561     fstream fimpl(implpath?implpath:implfile.c_str()) ;
01562     
01563     // Abort import if implementation file cannot be opened
01564     if (!fimpl) {
01565       oocoutE(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() 
01566                                       << ") ERROR opening implementation file " <<  implfile << endl ;
01567       return kFALSE ;      
01568     }
01569     
01570     
01571     // Import entire implentation file into stl string
01572     string impl ;
01573     while(fimpl.getline(buf,1023)) {
01574       // Process #include statements here
01575       
01576       // Look for include state of self
01577       Bool_t foundSelfInclude=kFALSE ;
01578       
01579       // Look for include of declaration file corresponding to this implementation file
01580       if (strstr(buf,"#include")) {
01581         // Process #include statements here
01582         char tmp[1024] ;
01583         strlcpy(tmp,buf,1024) ;
01584         strtok(tmp," <\"") ;
01585         char* incfile = strtok(0," <\"") ;
01586         
01587         if (strstr(incfile,declfilename.c_str())) {
01588           foundSelfInclude=kTRUE ;
01589         }
01590       } 
01591       
01592       // Explicitly rewrite include of own declaration file to string
01593       // any directory prefixes, copy all other lines verbatim in stl string
01594       if (foundSelfInclude) {
01595         // If include of self is found, substitute original include 
01596         // which may have directory structure with a plain include
01597         impl += Form("#include \"%s.%s\"\n",declfilebase.c_str(),declfileext.c_str()) ;
01598       } else {
01599         impl += buf ;
01600         impl += '\n' ;
01601       }
01602     }
01603     
01604         
01605     // Create entry in file map
01606     _fmap[declfilebase]._hfile = decl ;
01607     _fmap[declfilebase]._cxxfile = impl ;   
01608     _fmap[declfilebase]._hext = declfileext ;
01609 
01610   } else {
01611 
01612     // Inform that existing file entry is being recycled because it already contained class code
01613     oocoutI(_wspace,ObjectHandling) << "RooWorkspace::autoImportClass(" << _wspace->GetName() 
01614                                     << ") code of class " << tc->GetName() 
01615                                     << " was already imported from " << (implpath?implpath:implfile.c_str()) 
01616                                     << " and " << (declpath?declpath:declfile.c_str()) << endl ;
01617     
01618   }
01619 
01620 
01621   // *** PHASE 4 *** Import stl strings with code into workspace 
01622   //
01623   // If multiple classes are declared in a single code unit, there will be
01624   // multiple _c2fmap entries all pointing to the same _fmap entry.  
01625   
01626   // Make list of all immediate base classes of this class
01627   TString baseNameList ;
01628   TList* bl = tc->GetListOfBases() ;
01629   TIterator* iter = bl->MakeIterator() ;
01630   TBaseClass* base ;
01631   list<TClass*> bases ;
01632   while((base=(TBaseClass*)iter->Next())) {
01633     if (baseNameList.Length()>0) {
01634       baseNameList += "," ;
01635     }
01636     baseNameList += base->GetClassPointer()->GetName() ;
01637     bases.push_back(base->GetClassPointer()) ;
01638   }
01639   
01640   // Map class name to above _fmap entries, along with list of base classes
01641   // in _c2fmap
01642   _c2fmap[tc->GetName()]._baseName = baseNameList ;
01643   _c2fmap[tc->GetName()]._fileBase = declfilebase ;
01644   
01645   // Recursive store all base classes.
01646   list<TClass*>::iterator biter = bases.begin() ;
01647   while(biter!=bases.end()) {
01648     autoImportClass(*biter,doReplace) ;
01649     ++biter ;
01650   }
01651 
01652   // Cleanup 
01653   if (implpath) {
01654     delete[] implpath ;
01655   }
01656   if (declpath) {
01657     delete[] declpath ;
01658   }
01659 
01660   return kTRUE ;
01661 }
01662 
01663 
01664 //_____________________________________________________________________________
01665 Bool_t RooWorkspace::makeDir() 
01666 {
01667   // Create transient TDirectory representation of this workspace. This directory
01668   // will appear as a subdirectory of the directory that contains the workspace
01669   // and will have the name of the workspace suffixed with "Dir". The TDirectory
01670   // interface is read-only. Any attempt to insert objects into the workspace
01671   // directory representation will result in an error message. Note that some
01672   // ROOT object like TH1 automatically insert themselves into the current directory
01673   // when constructed. This will give error messages when done in a workspace
01674   // directory.
01675 
01676   if (_dir) return kTRUE ;
01677 
01678   TString title= Form("TDirectory representation of RooWorkspace %s",GetName()) ;
01679   _dir = new WSDir(GetName(),title.Data(),this) ;
01680   
01681   TIterator* iter = componentIterator() ;
01682   RooAbsArg* darg ;
01683   while((darg=(RooAbsArg*)iter->Next())) {
01684     if (darg->IsA() != RooConstVar::Class()) {
01685       _dir->InternalAppend(darg) ;
01686     }
01687   }
01688   
01689   return kTRUE ;
01690 }
01691  
01692  
01693 
01694 //_____________________________________________________________________________
01695 Bool_t RooWorkspace::import(TObject& object, Bool_t replaceExisting) 
01696 {
01697   // Import a clone of a generic TObject into workspace generic object container. Imported
01698   // object can be retrieved by name through the obj() method. The object is cloned upon
01699   // importation and the input argument does not need to live beyond the import call
01700   // 
01701   // Returns kTRUE if an error has occurred.
01702   
01703   // First check if object with given name already exists
01704   TObject* oldObj = _genObjects.FindObject(object.GetName()) ;
01705   if (oldObj && !replaceExisting) {
01706     coutE(InputArguments) << "RooWorkspace::import(" << GetName() << ") generic object with name " 
01707                           << object.GetName() << " is already in workspace and replaceExisting flag is set to false" << endl ;
01708     return kTRUE ;
01709   }  
01710   TH1::AddDirectory(kFALSE) ;
01711   if (oldObj) {
01712     _genObjects.Replace(oldObj,object.Clone()) ;
01713     delete oldObj ;
01714   } else {
01715     _genObjects.Add(object.Clone()) ;
01716   }
01717   TH1::AddDirectory(kTRUE) ;
01718   return kFALSE ;
01719 }
01720 
01721 
01722 
01723 
01724 //_____________________________________________________________________________
01725 Bool_t RooWorkspace::import(TObject& object, const char* aliasName, Bool_t replaceExisting) 
01726 {
01727   // Import a clone of a generic TObject into workspace generic object container. 
01728   // The imported object will be stored under the given alias name rather than its
01729   // own name. Imported object can be retrieved its alias name through the obj() method. 
01730   // The object is cloned upon importation and the input argument does not need to live beyond the import call
01731   // This method is mostly useful for importing objects that do not have a settable name such as TMatrix
01732   // 
01733   // Returns kTRUE if an error has occurred.
01734   
01735   // First check if object with given name already exists
01736   TObject* oldObj = _genObjects.FindObject(object.GetName()) ;
01737   if (oldObj && !replaceExisting) {
01738     coutE(InputArguments) << "RooWorkspace::import(" << GetName() << ") generic object with name " 
01739                           << object.GetName() << " is already in workspace and replaceExisting flag is set to false" << endl ;
01740     return kTRUE ;
01741   }  
01742   
01743   TH1::AddDirectory(kFALSE) ;
01744   RooTObjWrap* wrapper = new RooTObjWrap(object.Clone()) ;
01745   TH1::AddDirectory(kTRUE) ;
01746   wrapper->setOwning(kTRUE) ;
01747   wrapper->SetName(aliasName) ;
01748   wrapper->SetTitle(aliasName) ;
01749     
01750   if (oldObj) {
01751     _genObjects.Replace(oldObj,wrapper) ;
01752     delete oldObj ;
01753   } else {
01754     _genObjects.Add(wrapper) ;
01755   }
01756   return kFALSE ;
01757 }
01758 
01759 
01760 
01761 
01762 //_____________________________________________________________________________
01763 Bool_t RooWorkspace::addStudy(RooAbsStudy& study) 
01764 {
01765   // Insert RooStudyManager module
01766   RooAbsStudy* clone = (RooAbsStudy*) study.Clone() ;
01767   _studyMods.Add(clone) ;
01768   return kFALSE ;
01769 }
01770 
01771 
01772 
01773 
01774 //_____________________________________________________________________________
01775 void RooWorkspace::clearStudies() 
01776 {
01777   // Remove all RooStudyManager modules
01778   _studyMods.Delete() ;
01779 }
01780 
01781 
01782 
01783 
01784 //_____________________________________________________________________________
01785 TObject* RooWorkspace::obj(const char* name) const
01786 {
01787   // Return any type of object (RooAbsArg, RooAbsData or generic object) with given name)
01788 
01789   // Try RooAbsArg first
01790   TObject* ret = arg(name) ;
01791   if (ret) return ret ;
01792 
01793   // Then try RooAbsData
01794   ret = data(name) ;
01795   if (ret) return ret ;
01796 
01797   // Finally try generic object store
01798   return genobj(name) ;
01799 }
01800 
01801 
01802 
01803 //_____________________________________________________________________________
01804 TObject* RooWorkspace::genobj(const char* name)  const
01805 {
01806   // Return generic object with given name
01807 
01808   // Find object by name
01809   TObject* gobj = _genObjects.FindObject(name) ;
01810 
01811   // Exit here if not found
01812   if (!gobj) return 0 ;
01813 
01814   // If found object is wrapper, return payload
01815   if (gobj->IsA()==RooTObjWrap::Class()) return ((RooTObjWrap*)gobj)->obj() ;
01816 
01817   return gobj ;
01818 }
01819 
01820 
01821 
01822 //_____________________________________________________________________________
01823 Bool_t RooWorkspace::cd(const char* path) 
01824 {
01825   makeDir() ;
01826   return _dir->cd(path) ;
01827 }
01828 
01829 
01830 
01831 //_____________________________________________________________________________
01832 Bool_t RooWorkspace::writeToFile(const char* fileName, Bool_t recreate) 
01833 {
01834   // Save this current workspace into given file
01835   TFile f(fileName,recreate?"RECREATE":"UPDATE") ;
01836   Write() ;
01837   return kFALSE ;
01838 }
01839  
01840  
01841 
01842 //_____________________________________________________________________________
01843 RooFactoryWSTool& RooWorkspace::factory()
01844 {
01845   // Return instance to factory tool 
01846 
01847   if (_factory) {
01848     return *_factory ;  
01849   }
01850   cxcoutD(ObjectHandling) << "INFO: Creating RooFactoryWSTool associated with this workspace" << endl ;
01851   _factory = new RooFactoryWSTool(*this) ;
01852   return *_factory  ;
01853 }
01854 
01855 
01856 
01857 
01858 //_____________________________________________________________________________
01859 RooAbsArg* RooWorkspace::factory(const char* expr)
01860 {
01861   // Short-hand function for factory()->process(expr) ;
01862   return factory().process(expr) ;
01863 }
01864 
01865 
01866 
01867 
01868 //_____________________________________________________________________________
01869 void RooWorkspace::Print(Option_t* opts) const 
01870 {
01871   // Print contents of the workspace 
01872 
01873   Bool_t treeMode(kFALSE) ;
01874   if (TString(opts).Contains("t")) {
01875     treeMode=kTRUE ;
01876   }
01877 
01878   cout << endl << "RooWorkspace(" << GetName() << ") " << GetTitle() << " contents" << endl << endl  ;
01879   
01880   RooAbsArg* parg ;
01881 
01882   RooArgSet pdfSet ;
01883   RooArgSet funcSet ;
01884   RooArgSet varSet ;
01885   RooArgSet catfuncSet ;
01886   RooArgSet convResoSet ;
01887   RooArgSet resoSet ;
01888 
01889 
01890   // Split list of components in pdfs, functions and variables
01891   TIterator* iter = _allOwnedNodes.createIterator() ;
01892   while((parg=(RooAbsArg*)iter->Next())) {
01893 
01894     //---------------
01895 
01896     if (treeMode) {
01897       
01898       // In tree mode, only add nodes with no clients to the print lists
01899 
01900       if (parg->IsA()->InheritsFrom(RooAbsPdf::Class())) {
01901         if (!parg->hasClients()) {
01902           pdfSet.add(*parg) ;
01903         }
01904       }
01905       
01906       if (parg->IsA()->InheritsFrom(RooAbsReal::Class()) && 
01907           !parg->IsA()->InheritsFrom(RooAbsPdf::Class()) && 
01908           !parg->IsA()->InheritsFrom(RooConstVar::Class()) && 
01909           !parg->IsA()->InheritsFrom(RooRealVar::Class())) {
01910         if (!parg->hasClients()) {
01911           funcSet.add(*parg) ;
01912         }
01913       }
01914       
01915       
01916       if (parg->IsA()->InheritsFrom(RooAbsCategory::Class()) && 
01917           !parg->IsA()->InheritsFrom(RooCategory::Class())) {
01918         if (!parg->hasClients()) {      
01919           catfuncSet.add(*parg) ;
01920         }
01921       }
01922 
01923     } else {
01924 
01925       if (parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
01926         if (((RooResolutionModel*)parg)->isConvolved()) {
01927           convResoSet.add(*parg) ;
01928         } else {
01929           resoSet.add(*parg) ;
01930         }
01931       }
01932       
01933       if (parg->IsA()->InheritsFrom(RooAbsPdf::Class()) &&
01934           !parg->IsA()->InheritsFrom(RooResolutionModel::Class())) {
01935         pdfSet.add(*parg) ;
01936       }
01937       
01938       if (parg->IsA()->InheritsFrom(RooAbsReal::Class()) && 
01939           !parg->IsA()->InheritsFrom(RooAbsPdf::Class()) && 
01940           !parg->IsA()->InheritsFrom(RooConstVar::Class()) && 
01941           !parg->IsA()->InheritsFrom(RooRealVar::Class())) {
01942         funcSet.add(*parg) ;
01943       }
01944       
01945       if (parg->IsA()->InheritsFrom(RooAbsCategory::Class()) && 
01946           !parg->IsA()->InheritsFrom(RooCategory::Class())) {
01947         catfuncSet.add(*parg) ;
01948       }
01949     }
01950 
01951     if (parg->IsA()->InheritsFrom(RooRealVar::Class())) {
01952       varSet.add(*parg) ;
01953     }
01954 
01955     if (parg->IsA()->InheritsFrom(RooCategory::Class())) {
01956       varSet.add(*parg) ;
01957     }
01958 
01959   }
01960   delete iter ;
01961 
01962 
01963   RooFit::MsgLevel oldLevel = RooMsgService::instance().globalKillBelow() ;
01964   RooMsgService::instance().setGlobalKillBelow(RooFit::WARNING) ;
01965 
01966   if (varSet.getSize()>0) {
01967     varSet.sort() ;
01968     cout << "variables" << endl ;
01969     cout << "---------" << endl ;
01970     cout << varSet << endl ;
01971     cout << endl ;
01972   }
01973 
01974   if (pdfSet.getSize()>0) {
01975     cout << "p.d.f.s" << endl ;
01976     cout << "-------" << endl ;
01977     pdfSet.sort() ;
01978     iter = pdfSet.createIterator() ;
01979     while((parg=(RooAbsArg*)iter->Next())) {
01980       if (treeMode) {
01981         parg->printComponentTree() ;
01982       } else {
01983         parg->Print() ;
01984       }
01985     }
01986     delete iter ;
01987     cout << endl ;
01988   }
01989 
01990   if (!treeMode) {
01991     if (resoSet.getSize()>0) {
01992       cout << "analytical resolution models" << endl ;
01993       cout << "----------------------------" << endl ;
01994       resoSet.sort() ;
01995       iter = resoSet.createIterator() ;
01996       while((parg=(RooAbsArg*)iter->Next())) {
01997         parg->Print() ;
01998       }
01999       delete iter ;
02000       //     iter = convResoSet.createIterator() ;
02001       //     while((parg=(RooAbsArg*)iter->Next())) {
02002       //       parg->Print() ;
02003       //     }
02004       //     delete iter ;
02005       cout << endl ;
02006     }
02007   }
02008 
02009   if (funcSet.getSize()>0) {
02010     cout << "functions" << endl ;
02011     cout << "--------" << endl ;
02012     funcSet.sort() ;
02013     iter = funcSet.createIterator() ;
02014     while((parg=(RooAbsArg*)iter->Next())) {
02015       if (treeMode) {
02016         parg->printComponentTree() ;
02017       } else {
02018         parg->Print() ;
02019       }
02020     }
02021     delete iter ;
02022     cout << endl ;
02023   }
02024 
02025   if (catfuncSet.getSize()>0) {
02026     cout << "category functions" << endl ;
02027     cout << "------------------" << endl ;
02028     catfuncSet.sort() ;
02029     iter = catfuncSet.createIterator() ;
02030     while((parg=(RooAbsArg*)iter->Next())) {
02031       if (treeMode) {
02032         parg->printComponentTree() ;
02033       } else {
02034         parg->Print() ;
02035       }
02036     }
02037     delete iter ;
02038     cout << endl ;
02039   }
02040 
02041   if (_dataList.GetSize()>0) {
02042     cout << "datasets" << endl ;
02043     cout << "--------" << endl ;
02044     iter = _dataList.MakeIterator() ;
02045     RooAbsData* data2 ;
02046     while((data2=(RooAbsData*)iter->Next())) {
02047       cout << data2->IsA()->GetName() << "::" << data2->GetName() << *data2->get() << endl ;
02048     }
02049     delete iter ;
02050     cout << endl ;
02051   }
02052 
02053   if (_snapshots.GetSize()>0) {
02054     cout << "parameter snapshots" << endl ;
02055     cout << "-------------------" << endl ;
02056     iter = _snapshots.MakeIterator() ;
02057     RooArgSet* snap ;
02058     while((snap=(RooArgSet*)iter->Next())) {
02059       cout << snap->GetName() << " = (" ;
02060       TIterator* aiter = snap->createIterator() ;
02061       RooAbsArg* a ;
02062       Bool_t first(kTRUE) ;
02063       while((a=(RooAbsArg*)aiter->Next())) {
02064         if (first) { first=kFALSE ; } else { cout << "," ; }
02065         cout << a->GetName() << "=" ; 
02066         a->printValue(cout) ;
02067         if (a->isConstant()) {
02068           cout << "[C]" ;
02069         }
02070       }
02071       cout << ")" << endl ;
02072       delete aiter ;
02073     }
02074     delete iter ;
02075     cout << endl ;
02076   }
02077 
02078 
02079   if (_namedSets.size()>0) {
02080     cout << "named sets" << endl ;
02081     cout << "----------" << endl ;
02082     for (map<string,RooArgSet>::const_iterator it = _namedSets.begin() ; it != _namedSets.end() ; it++) {
02083       cout << it->first << ":" << it->second << endl ;
02084     }
02085     
02086     cout << endl ;
02087   }
02088 
02089  
02090   if (_genObjects.GetSize()>0) {
02091     cout << "generic objects" << endl ;
02092     cout << "---------------" << endl ;
02093     iter = _genObjects.MakeIterator() ;
02094     TObject* gobj ;
02095     while((gobj=(TObject*)iter->Next())) {
02096       if (gobj->IsA()==RooTObjWrap::Class()) {
02097         cout << ((RooTObjWrap*)gobj)->obj()->IsA()->GetName() << "::" << gobj->GetName() << endl ;
02098       } else {
02099         cout << gobj->IsA()->GetName() << "::" << gobj->GetName() << endl ;
02100       }
02101     }
02102     delete iter ;
02103     cout << endl ;
02104     
02105   }
02106 
02107   if (_studyMods.GetSize()>0) {
02108     cout << "study modules" << endl ;
02109     cout << "-------------" << endl ;
02110     iter = _studyMods.MakeIterator() ;
02111     TObject* smobj ;
02112     while((smobj=(TObject*)iter->Next())) {
02113       cout << smobj->IsA()->GetName() << "::" << smobj->GetName() << endl ;
02114     }
02115     delete iter ;
02116     cout << endl ;
02117     
02118   }
02119 
02120   if (_classes.listOfClassNames().size()>0) {
02121     cout << "embedded class code" << endl ;
02122     cout << "-------------------" << endl ;    
02123     cout << _classes.listOfClassNames() << endl ;
02124     cout << endl ;
02125   }
02126 
02127   if (_eocache.size()>0) {
02128     cout << "embedded precalculated expensive components" << endl ;
02129     cout << "-------------------------------------------" << endl ;
02130     _eocache.print() ;
02131   }
02132 
02133   RooMsgService::instance().setGlobalKillBelow(oldLevel) ;
02134 
02135   return ;
02136 }
02137 
02138 
02139 //_____________________________________________________________________________
02140 void RooWorkspace::CodeRepo::Streamer(TBuffer &R__b)
02141 {
02142   // Custom streamer for the workspace. Stream contents of workspace
02143   // and code repository. When reading, read code repository first
02144   // and compile missing classes before proceeding with streaming
02145   // of workspace contents to avoid errors.
02146 
02147   typedef ::RooWorkspace::CodeRepo thisClass;
02148 
02149    // Stream an object of class RooWorkspace::CodeRepo.
02150    if (R__b.IsReading()) {
02151 
02152      UInt_t R__s, R__c;
02153      R__b.ReadVersion(&R__s, &R__c); 
02154 
02155      // Stream contents of ClassFiles map
02156      Int_t count(0) ;
02157      R__b >> count ;
02158      while(count--) {
02159        TString name ;
02160        name.Streamer(R__b) ;       
02161        _fmap[name]._hext.Streamer(R__b) ;
02162        _fmap[name]._hfile.Streamer(R__b) ;
02163        _fmap[name]._cxxfile.Streamer(R__b) ;    
02164      }     
02165  
02166      // Stream contents of ClassRelInfo map
02167      count=0 ;
02168      R__b >> count ;
02169      while(count--) {
02170        TString name,bname,fbase ;
02171        name.Streamer(R__b) ;       
02172        _c2fmap[name]._baseName.Streamer(R__b) ;
02173        _c2fmap[name]._fileBase.Streamer(R__b) ;
02174      }     
02175      R__b.CheckByteCount(R__s, R__c, thisClass::IsA());
02176 
02177      // Instantiate any classes that are not defined in current session
02178      _compiledOK = !compileClasses() ;
02179 
02180    } else {
02181      
02182      UInt_t R__c;
02183      R__c = R__b.WriteVersion(thisClass::IsA(), kTRUE);
02184      
02185      // Stream contents of ClassFiles map
02186      UInt_t count = _fmap.size() ;
02187      R__b << count ;
02188      map<TString,ClassFiles>::iterator iter = _fmap.begin() ;
02189      while(iter!=_fmap.end()) {       
02190        TString key_copy(iter->first) ;
02191        key_copy.Streamer(R__b) ;
02192        iter->second._hext.Streamer(R__b) ;
02193        iter->second._hfile.Streamer(R__b);
02194        iter->second._cxxfile.Streamer(R__b);
02195 
02196        ++iter ;
02197      }
02198      
02199      // Stream contents of ClassRelInfo map
02200      count = _c2fmap.size() ;
02201      R__b << count ;
02202      map<TString,ClassRelInfo>::iterator iter2 = _c2fmap.begin() ;
02203      while(iter2!=_c2fmap.end()) {
02204        TString key_copy(iter2->first) ;
02205        key_copy.Streamer(R__b) ;
02206        iter2->second._baseName.Streamer(R__b) ;
02207        iter2->second._fileBase.Streamer(R__b);
02208        ++iter2 ;
02209      }
02210 
02211      R__b.SetByteCount(R__c, kTRUE);
02212      
02213    }
02214 }
02215 
02216 
02217 //______________________________________________________________________________
02218 void RooWorkspace::Streamer(TBuffer &R__b)
02219 {
02220   // Stream an object of class RooWorkspace. This is a standard ROOT streamer for the
02221   // I/O part. This custom function exists to detach all external client links
02222   // from the payload prior to writing the payload so that these client links
02223   // are not persisted. (Client links occur if external function objects use 
02224   // objects contained in the workspace as input)
02225   // After the actual writing, these client links are restored.
02226 
02227    if (R__b.IsReading()) {
02228 
02229       R__b.ReadClassBuffer(RooWorkspace::Class(),this);
02230 
02231       // Make expensive object cache of all objects point to intermal copy.
02232       // Somehow this doesn't work OK automatically
02233       TIterator* iter = _allOwnedNodes.createIterator() ;
02234       RooAbsArg* node ;
02235       while((node=(RooAbsArg*)iter->Next())) {
02236         node->setExpensiveObjectCache(_eocache) ;
02237         if (node->IsA()->InheritsFrom(RooAbsOptTestStatistic::Class())) {
02238           RooAbsOptTestStatistic* tmp = (RooAbsOptTestStatistic*) node ;
02239           if (tmp->isSealed() && tmp->sealNotice() && strlen(tmp->sealNotice())>0) {
02240             cout << "RooWorkspace::Streamer(" << GetName() << ") " << node->IsA()->GetName() << "::" << node->GetName() << " : " << tmp->sealNotice() << endl ;
02241           }
02242         }
02243       }
02244       delete iter ;
02245 
02246 
02247    } else {
02248 
02249 
02250      // Make lists of external clients of WS objects, and remove those links temporarily
02251 
02252      map<RooAbsArg*,list<RooAbsArg*> > extClients, extValueClients, extShapeClients ;
02253 
02254      TIterator* iter = _allOwnedNodes.createIterator() ;
02255      RooAbsArg* tmparg ;
02256      while((tmparg=(RooAbsArg*)iter->Next())) {
02257 
02258        // Loop over client list of this arg
02259        TIterator* clientIter = tmparg->_clientList.MakeIterator() ;
02260        RooAbsArg* client ;
02261        while((client=(RooAbsArg*)clientIter->Next())) {
02262          if (!_allOwnedNodes.containsInstance(*client)) {
02263            while(tmparg->_clientList.refCount(client)>0) {
02264              tmparg->_clientList.Remove(client) ;
02265              extClients[tmparg].push_back(client) ;
02266            }
02267          }
02268        }
02269        delete clientIter ;
02270 
02271        // Loop over value client list of this arg
02272        TIterator* vclientIter = tmparg->_clientListValue.MakeIterator() ;
02273        RooAbsArg* vclient ;
02274        while((vclient=(RooAbsArg*)vclientIter->Next())) {
02275          if (!_allOwnedNodes.containsInstance(*vclient)) {
02276            cxcoutD(ObjectHandling) << "RooWorkspace::Streamer(" << GetName() << ") element " << tmparg->GetName() 
02277                                    << " has external value client link to " << vclient << " (" << vclient->GetName() << ") with ref count " << tmparg->_clientListValue.refCount(vclient) << endl ;
02278            while(tmparg->_clientListValue.refCount(vclient)>0) {
02279              tmparg->_clientListValue.Remove(vclient) ;
02280              extValueClients[tmparg].push_back(vclient) ;
02281            }
02282          }
02283        }
02284        delete vclientIter ;
02285 
02286        // Loop over shape client list of this arg
02287        TIterator* sclientIter = tmparg->_clientListShape.MakeIterator() ;
02288        RooAbsArg* sclient ;
02289        while((sclient=(RooAbsArg*)sclientIter->Next())) {
02290          if (!_allOwnedNodes.containsInstance(*sclient)) {
02291            cxcoutD(ObjectHandling) << "RooWorkspace::Streamer(" << GetName() << ") element " << tmparg->GetName() 
02292                                  << " has external shape client link to " << sclient << " (" << sclient->GetName() << ") with ref count " << tmparg->_clientListShape.refCount(sclient) << endl ;
02293            while(tmparg->_clientListShape.refCount(sclient)>0) {
02294              tmparg->_clientListShape.Remove(sclient) ;
02295              extShapeClients[tmparg].push_back(sclient) ;
02296            }
02297          }
02298        }
02299        delete sclientIter ;
02300 
02301      }
02302      delete iter ;
02303 
02304      R__b.WriteClassBuffer(RooWorkspace::Class(),this);
02305 
02306      // Reinstate clients here
02307 
02308      
02309      for (map<RooAbsArg*,list<RooAbsArg*> >::iterator iterx = extClients.begin() ; iterx!=extClients.end() ; iterx++) {
02310        for (list<RooAbsArg*>::iterator citer = iterx->second.begin() ; citer!=iterx->second.end() ; citer++) {
02311          iterx->first->_clientList.Add(*citer) ;
02312        }
02313      }
02314 
02315      for (map<RooAbsArg*,list<RooAbsArg*> >::iterator iterx = extValueClients.begin() ; iterx!=extValueClients.end() ; iterx++) {
02316        for (list<RooAbsArg*>::iterator citer = iterx->second.begin() ; citer!=iterx->second.end() ; citer++) {
02317          iterx->first->_clientListValue.Add(*citer) ;
02318        }
02319      }
02320 
02321      for (map<RooAbsArg*,list<RooAbsArg*> >::iterator iterx = extShapeClients.begin() ; iterx!=extShapeClients.end() ; iterx++) {
02322        for (list<RooAbsArg*>::iterator citer = iterx->second.begin() ; citer!=iterx->second.end() ; citer++) {
02323          iterx->first->_clientListShape.Add(*citer) ;
02324        }
02325      }
02326 
02327    }
02328 }
02329 
02330 
02331 
02332 
02333 //_____________________________________________________________________________
02334 std::string RooWorkspace::CodeRepo::listOfClassNames() const 
02335 {
02336   // Return STL string with last of class names contained in the code repository
02337 
02338   string ret ;
02339   map<TString,ClassRelInfo>::const_iterator iter = _c2fmap.begin() ;
02340   while(iter!=_c2fmap.end()) {
02341     if (ret.size()>0) {
02342       ret += ", " ;
02343     }
02344     ret += iter->first ;    
02345     ++iter ;
02346   }  
02347   
02348   return ret ;
02349 }
02350 
02351 
02352 
02353 //_____________________________________________________________________________
02354 Bool_t RooWorkspace::CodeRepo::compileClasses() 
02355 {
02356   // For all classes in the workspace for which no class definition is
02357   // found in the ROOT class table extract source code stored in code
02358   // repository into temporary directory set by
02359   // setClassFileExportDir(), compile classes and link them with
02360   // current ROOT session. If a compilation error occurs print
02361   // instructions for user how to fix errors and recover workspace and
02362   // abort import procedure.
02363 
02364   Bool_t haveDir=kFALSE ;
02365 
02366   // Retrieve name of directory in which to export code files
02367   string dirName = Form(_classFileExportDir.c_str(),_wspace->uuid().AsString(),_wspace->GetName()) ;
02368 
02369   // Process all class entries in repository
02370   map<TString,ClassRelInfo>::iterator iter = _c2fmap.begin() ;
02371   while(iter!=_c2fmap.end()) {
02372 
02373     oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() now processing class " << iter->first.Data() << endl ;
02374 
02375     // If class is already known, don't load
02376     if (gClassTable->GetDict(iter->first.Data())) {
02377       oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Embedded class " 
02378                                       << iter->first << " already in ROOT class table, skipping" << endl ;
02379       ++iter ;
02380       continue ;
02381     }
02382 
02383     // Check that export directory exists
02384     if (!haveDir) {
02385 
02386       // If not, make local directory to extract files 
02387       if (!gSystem->AccessPathName(dirName.c_str())) {
02388         oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() reusing code export directory " << dirName.c_str() 
02389                                         << " to extract coded embedded in workspace" << endl ;
02390       } else {
02391         if (gSystem->MakeDirectory(dirName.c_str())==0) { 
02392           oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() creating code export directory " << dirName.c_str() 
02393                                           << " to extract coded embedded in workspace" << endl ;
02394         } else {
02395           oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR creating code export directory " << dirName.c_str() 
02396                                           << " to extract coded embedded in workspace" << endl ;
02397           return kFALSE ;
02398         }
02399       }
02400       haveDir=kTRUE ;
02401     }
02402 
02403     // Navigate from class to file
02404     ClassFiles& cfinfo = _fmap[iter->second._fileBase] ;
02405 
02406     oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() now processing file with base " << iter->second._fileBase << endl ;
02407     
02408     // If file is already processed, skip to next class
02409     if (cfinfo._extracted) {
02410       oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() file with base name " << iter->second._fileBase 
02411                                          << " has already been extracted, skipping to next class" << endl ;
02412       continue ;
02413     }
02414 
02415     // Check if identical declaration file (header) is already written
02416     Bool_t needDeclWrite=kTRUE ;
02417     string fdname = Form("%s/%s.%s",dirName.c_str(),iter->second._fileBase.Data(),cfinfo._hext.Data()) ;
02418     ifstream ifdecl(fdname.c_str()) ;
02419     if (ifdecl) {
02420       TString contents ;
02421       char buf[1024] ;
02422       while(ifdecl.getline(buf,1024)) {
02423         contents += buf ;
02424         contents += "\n" ;
02425       }      
02426       UInt_t crcFile = RooAbsArg::crc32(contents.Data()) ;
02427       UInt_t crcWS   = RooAbsArg::crc32(cfinfo._hfile.Data()) ;
02428       needDeclWrite = (crcFile!=crcWS) ;
02429     }
02430 
02431     // Write declaration file if required 
02432     if (needDeclWrite) {
02433       oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Extracting declaration code of class " << iter->first << ", file " << fdname << endl ;
02434       ofstream fdecl(fdname.c_str()) ;
02435       if (!fdecl) {
02436         oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR opening file" 
02437                                         << fdname << " for writing" << endl ;
02438         return kFALSE ;
02439       }
02440       fdecl << cfinfo._hfile ;
02441       fdecl.close() ;
02442     }
02443 
02444     // Check if identical implementation file is already written
02445     Bool_t needImplWrite=kTRUE ;
02446     string finame = Form("%s/%s.cxx",dirName.c_str(),iter->second._fileBase.Data()) ;
02447     ifstream ifimpl(finame.c_str()) ;
02448     if (ifimpl) {
02449       TString contents ;
02450       char buf[1024] ;
02451       while(ifimpl.getline(buf,1024)) {
02452         contents += buf ;
02453         contents += "\n" ;
02454       }      
02455       UInt_t crcFile = RooAbsArg::crc32(contents.Data()) ;
02456       UInt_t crcWS   = RooAbsArg::crc32(cfinfo._cxxfile.Data()) ;
02457       needImplWrite = (crcFile!=crcWS) ;
02458     }
02459 
02460     // Write implementation file if required
02461     if (needImplWrite) {
02462       oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Extracting implementation code of class " << iter->first << ", file " << finame << endl ;
02463       ofstream fimpl(finame.c_str()) ;
02464       if (!fimpl) {
02465         oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR opening file" 
02466                                         << finame << " for writing" << endl ;
02467         return kFALSE ;
02468       }
02469       fimpl << cfinfo._cxxfile ;
02470       fimpl.close() ;
02471     }
02472 
02473     // Mark this file as extracted
02474     cfinfo._extracted = kTRUE ;
02475     oocxcoutD(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() marking code unit  " << iter->second._fileBase << " as extracted" << endl ;
02476 
02477     // Compile class
02478     oocoutI(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() Compiling code unit " << iter->second._fileBase.Data() << " to define class " << iter->first << endl ;
02479     Bool_t ok = gSystem->CompileMacro(finame.c_str(),"k") ;
02480     
02481     if (!ok) {
02482       oocoutE(_wspace,ObjectHandling) << "RooWorkspace::CodeRepo::compileClasses() ERROR compiling class " << iter->first.Data() << ", to fix this you can do the following: " << endl 
02483                                       << "  1) Fix extracted source code files in directory " << dirName.c_str() << "/" << endl 
02484                                       << "  2) In clean ROOT session compiled fixed classes by hand using '.x " << dirName.c_str() << "/ClassName.cxx+'" << endl
02485                                       << "  3) Reopen file with RooWorkspace with broken source code in UPDATE mode. Access RooWorkspace to force loading of class" << endl
02486                                       << "     Broken instances in workspace will _not_ be compiled, instead precompiled fixed instances will be used." << endl
02487                                       << "  4) Reimport fixed code in workspace using 'RooWorkspace::importClassCode(\"*\",kTRUE)' method, Write() updated workspace to file and close file" << endl
02488                                       << "  5) Reopen file in clean ROOT session to confirm that problems are fixed" << endl ;
02489         return kFALSE ;
02490     }
02491     
02492     ++iter ;
02493   }
02494 
02495   return kTRUE ;
02496 }
02497 
02498 
02499 
02500 //_____________________________________________________________________________
02501 void RooWorkspace::WSDir::InternalAppend(TObject* obj) 
02502 {
02503   // Internal access to TDirectory append method
02504 
02505 #if ROOT_VERSION_CODE <= ROOT_VERSION(5,19,02)
02506   TDirectory::Append(obj) ;
02507 #else
02508   TDirectory::Append(obj,kFALSE) ;
02509 #endif
02510 
02511 }
02512 
02513 
02514 //_____________________________________________________________________________
02515 #if ROOT_VERSION_CODE <= ROOT_VERSION(5,19,02)
02516 void RooWorkspace::WSDir::Add(TObject* obj) 
02517 #else
02518 void RooWorkspace::WSDir::Add(TObject* obj,Bool_t) 
02519 #endif
02520 {
02521   // Overload TDirectory interface method to prohibit insertion of objects in read-only directory workspace representation
02522   if (dynamic_cast<RooAbsArg*>(obj) || dynamic_cast<RooAbsData*>(obj)) {
02523     coutE(ObjectHandling) << "RooWorkspace::WSDir::Add(" << GetName() << ") ERROR: Directory is read-only representation of a RooWorkspace, use RooWorkspace::import() to add objects" << endl ;
02524   } else {
02525     InternalAppend(obj) ;
02526   }
02527 } 
02528 
02529 
02530 //_____________________________________________________________________________
02531 #if ROOT_VERSION_CODE <= ROOT_VERSION(5,19,02)
02532 void RooWorkspace::WSDir::Append(TObject* obj) 
02533 #else
02534 void RooWorkspace::WSDir::Append(TObject* obj,Bool_t) 
02535 #endif
02536 {
02537   // Overload TDirectory interface method to prohibit insertion of objects in read-only directory workspace representation
02538   if (dynamic_cast<RooAbsArg*>(obj) || dynamic_cast<RooAbsData*>(obj)) {
02539     coutE(ObjectHandling) << "RooWorkspace::WSDir::Add(" << GetName() << ") ERROR: Directory is read-only representation of a RooWorkspace, use RooWorkspace::import() to add objects" << endl ;
02540   } else {
02541     InternalAppend(obj) ;
02542   }
02543 }
02544 
02545 
02546 
02547 //_____________________________________________________________________________
02548 void RooWorkspace::exportToCint(const char* nsname) 
02549 {
02550   // Activate export of workspace symbols to CINT in a namespace with given name. If no name
02551   // is given the namespace will have the same name as the workspace
02552 
02553   // If export is already active, do nothing
02554   if (_doExport) {
02555     coutE(ObjectHandling) << "RooWorkspace::exportToCint(" << GetName() << ") WARNING: repeated calls to exportToCint() have no effect" << endl ;
02556     return ;
02557   }
02558 
02559 
02560   // Set flag so that future import to workspace are automatically exported to CINT
02561   _doExport = kTRUE ;
02562 
02563   // If no name is provided choose name of workspace
02564   if (!nsname) nsname = GetName() ;
02565   _exportNSName = nsname ;
02566 
02567   coutI(ObjectHandling) << "RooWorkspace::exportToCint(" << GetName() 
02568                         << ") INFO: references to all objects in this workspace will be created in CINT in 'namespace " << _exportNSName << "'" << endl ;
02569 
02570   // Export present contents of workspace to CINT
02571   TIterator* iter = _allOwnedNodes.createIterator() ;
02572   TObject* wobj ;
02573   while((wobj=iter->Next())) {
02574     exportObj(wobj) ;
02575   }  
02576   delete iter ;
02577   iter = _dataList.MakeIterator() ;
02578   while((wobj=iter->Next())) {
02579     exportObj(wobj) ;
02580   }  
02581   delete iter ;
02582 }
02583 
02584 
02585 //_____________________________________________________________________________
02586 void RooWorkspace::exportObj(TObject* wobj)
02587 {
02588   // Export reference to given workspace object to CINT
02589 
02590   // Do nothing if export flag is not set
02591   if (!_doExport) return ;
02592 
02593   // Do not export RooConstVars
02594   if (wobj->IsA() == RooConstVar::Class()) {
02595     return ;
02596   }
02597 
02598 
02599   // Determine if object name is a valid C++ identifier name
02600 
02601   // Do not export objects that have names that are not valid C++ identifiers
02602   if (!isValidCPPID(wobj->GetName())) {
02603     cxcoutD(ObjectHandling) << "RooWorkspace::exportObj(" << GetName() << ") INFO: Workspace object name " << wobj->GetName() << " is not a valid C++ identifier and is not exported to CINT" << endl ;
02604     return ;
02605   }
02606 
02607   // Declare correctly typed reference to object in CINT in the namespace associated with this workspace
02608   string cintExpr = Form("namespace %s { %s& %s = *(%s *)%p ; }",_exportNSName.c_str(),wobj->IsA()->GetName(),wobj->GetName(),wobj->IsA()->GetName(),wobj) ;
02609   gROOT->ProcessLine(cintExpr.c_str()) ;  
02610 }
02611 
02612 
02613 
02614 //_____________________________________________________________________________
02615 Bool_t RooWorkspace::isValidCPPID(const char* name)   
02616 {
02617   // Return true if given name is a valid C++ identifier name
02618 
02619   string oname(name) ;
02620   if (isdigit(oname[0])) {
02621     return kFALSE ;
02622   } else {
02623     for (UInt_t i=0 ; i<oname.size() ; i++) {
02624       char c = oname[i] ;
02625       if (!isalnum(c) && (c!='_')) {
02626         return kFALSE ;
02627       }    
02628     }
02629   }
02630   return kTRUE ;
02631 }
02632 
02633 //_____________________________________________________________________________
02634 void RooWorkspace::unExport()
02635 {
02636   // Delete exported reference in CINT namespace 
02637   char buf[1024] ;
02638   TIterator* iter = _allOwnedNodes.createIterator() ;
02639   TObject* wobj ;
02640   while((wobj=iter->Next())) {
02641     if (isValidCPPID(wobj->GetName())) {
02642       strlcpy(buf,Form("%s::%s",_exportNSName.c_str(),wobj->GetName()),1024) ;
02643       G__deletevariable(buf) ;
02644     }
02645   }
02646   delete iter ;
02647 }
02648 

Generated on Tue Jul 5 15:07:45 2011 for ROOT_528-00b_version by  doxygen 1.5.1