Logo ROOT  
Reference Guide
RooDataSet.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitCore *
4  * @(#)root/roofitcore:$Id$
5  * Authors: *
6  * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7  * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8  * *
9  * Copyright (c) 2000-2005, Regents of the University of California *
10  * and Stanford University. All rights reserved. *
11  * *
12  * Redistribution and use in source and binary forms, *
13  * with or without modification, are permitted according to the terms *
14  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15  *****************************************************************************/
16 
17 /**
18 \file RooDataSet.cxx
19 \class RooDataSet
20 \ingroup Roofitcore
21 
22 RooDataSet is a container class to hold unbinned data. The binned equivalent is
23 RooDataHist. In RooDataSet, each data point in N-dimensional space is represented
24 by a RooArgSet of RooRealVar, RooCategory or RooStringVar objects, which can be
25 retrieved using get().
26 
27 Since RooDataSet saves every event, it allows for fits with highest precision. With a large
28 amount of data, however, it could be beneficial to represent them in binned form,
29 i.e., RooDataHist. Binning the data will incur a loss of information, though.
30 RooDataHist on the other hand may suffer from the curse of dimensionality if a high-dimensional
31 problem with a lot of bins on each axis is tackled.
32 
33 ### Inspecting a dataset
34 Inspect a dataset using Print() with the "verbose" option:
35 ```
36 dataset->Print("V");
37 dataset->get(0)->Print("V");
38 dataset->get(1)->Print("V");
39 ...
40 ```
41 
42 ### Plotting data.
43 See RooAbsData::plotOn().
44 
45 
46 ### Storage strategy
47 There are two storage backends:
48 - RooVectorDataStore (default): std::vectors in memory. They are fast, but they
49 cannot be serialised if the dataset exceeds a size of 1 Gb
50 - RooTreeDataStore: Uses a TTree, which can be file backed if a file is opened
51 before creating the dataset. This significantly reduces the memory pressure, as the
52 baskets of the tree can be written to a file, and only the basket that's currently
53 being read stays in RAM.
54  - Enable tree-backed storage similar to this:
55  ```
56  TFile outputFile("filename.root", "RECREATE");
57  RooAbsData::setDefaultStorageType(RooAbsData::Tree);
58  RooDataSet mydata(...);
59  ```
60  - Or convert an existing memory-backed data storage:
61  ```
62  RooDataSet mydata(...);
63 
64  TFile outputFile("filename.root", "RECREATE");
65  mydata.convertToTreeStore();
66  ```
67 
68 For the inverse conversion, see `RooAbsData::convertToVectorStore()`.
69 
70 **/
71 
72 #include "RooDataSet.h"
73 
74 #include "RooPlot.h"
75 #include "RooAbsReal.h"
76 #include "Roo1DTable.h"
77 #include "RooCategory.h"
78 #include "RooFormulaVar.h"
79 #include "RooArgList.h"
80 #include "RooAbsRealLValue.h"
81 #include "RooRealVar.h"
82 #include "RooDataHist.h"
83 #include "RooMsgService.h"
84 #include "RooCmdConfig.h"
85 #include "RooHist.h"
86 #include "RooTreeDataStore.h"
87 #include "RooVectorDataStore.h"
88 #include "RooCompositeDataStore.h"
89 #include "RooSentinel.h"
90 #include "RooTrace.h"
91 #include "RooHelpers.h"
92 
93 #include "TTree.h"
94 #include "TH2.h"
95 #include "TFile.h"
96 #include "TBuffer.h"
97 #include "ROOT/RMakeUnique.hxx"
98 #include "strlcpy.h"
99 #include "snprintf.h"
100 
101 #include <iostream>
102 #include <fstream>
103 
104 
105 using namespace std;
106 
108 
109 #ifndef USEMEMPOOLFORDATASET
111 #else
112 
113 #include "MemPoolForRooSets.h"
114 
117  static auto * memPool = new RooDataSet::MemPool();
118  return memPool;
119 }
120 
121 void RooDataSet::cleanup() {
122  auto pool = memPool();
123  pool->teardown();
124 
125  //The pool will have to leak if it's not empty at this point.
126  if (pool->empty())
127  delete pool;
128 }
129 
130 
131 ////////////////////////////////////////////////////////////////////////////////
132 /// Overloaded new operator guarantees that all RooDataSets allocated with new
133 /// have a unique address, a property that is exploited in several places
134 /// in roofit to quickly index contents on normalization set pointers.
135 /// The memory pool only allocates space for the class itself. The elements
136 /// stored in the set are stored outside the pool.
137 
138 void* RooDataSet::operator new (size_t bytes)
139 {
140  //This will fail if a derived class uses this operator
141  assert(sizeof(RooDataSet) == bytes);
142 
143  return memPool()->allocate(bytes);
144 }
145 
146 
147 
148 ////////////////////////////////////////////////////////////////////////////////
149 /// Memory is owned by pool, we need to do nothing to release it
150 
151 void RooDataSet::operator delete (void* ptr)
152 {
153  // Decrease use count in pool that ptr is on
154  if (memPool()->deallocate(ptr))
155  return;
156 
157  std::cerr << __func__ << " " << ptr << " is not in any of the pools." << std::endl;
158 
159  // Not part of any pool; use global op delete:
160  ::operator delete(ptr);
161 }
162 
163 #endif
164 
165 
166 ////////////////////////////////////////////////////////////////////////////////
167 /// Default constructor for persistence
168 
169 RooDataSet::RooDataSet() : _wgtVar(0)
170 {
172 }
173 
174 
175 
176 
177 
178 ////////////////////////////////////////////////////////////////////////////////
179 /// Construct an unbinned dataset from a RooArgSet defining the dimensions of the data space. Optionally, data
180 /// can be imported at the time of construction.
181 ///
182 /// <table>
183 /// <tr><th> %RooCmdArg <th> Effect
184 /// <tr><td> Import(TTree*) <td> Import contents of given TTree. Only braches of the TTree that have names
185 /// corresponding to those of the RooAbsArgs that define the RooDataSet are
186 /// imported.
187 /// <tr><td> ImportFromFile(const char* fileName, const char* treeName) <td> Import tree with given name from file with given name.
188 /// <tr><td> Import(RooDataSet&)
189 /// <td> Import contents of given RooDataSet. Only observables that are common with the definition of this dataset will be imported
190 /// <tr><td> Index(RooCategory&) <td> Prepare import of datasets into a N+1 dimensional RooDataSet
191 /// where the extra discrete dimension labels the source of the imported histogram.
192 /// <tr><td> Import(const char*, RooDataSet&)
193 /// <td> Import a dataset to be associated with the given state name of the index category
194 /// specified in Index(). If the given state name is not yet defined in the index
195 /// category it will be added on the fly. The import command can be specified multiple times.
196 /// <tr><td> Link(const char*, RooDataSet&) <td> Link contents of supplied RooDataSet to this dataset for given index category state name.
197 /// In this mode, no data is copied and the linked dataset must be remain live for the duration
198 /// of this dataset. Note that link is active for both reading and writing, so modifications
199 /// to the aggregate dataset will also modify its components. Link() and Import() are mutually exclusive.
200 /// <tr><td> OwnLinked() <td> Take ownership of all linked datasets
201 /// <tr><td> Import(map<string,RooDataSet*>&) <td> As above, but allows specification of many imports in a single operation
202 /// <tr><td> Link(map<string,RooDataSet*>&) <td> As above, but allows specification of many links in a single operation
203 /// <tr><td> Cut(const char*) <br>
204 /// Cut(RooFormulaVar&)
205 /// <td> Apply the given cut specification when importing data
206 /// <tr><td> CutRange(const char*) <td> Only accept events in the observable range with the given name
207 /// <tr><td> WeightVar(const char*) <br>
208 /// WeightVar(const RooAbsArg&)
209 /// <td> Interpret the given variable as event weight rather than as observable
210 /// <tr><td> StoreError(const RooArgSet&) <td> Store symmetric error along with value for given subset of observables
211 /// <tr><td> StoreAsymError(const RooArgSet&) <td> Store asymmetric error along with value for given subset of observables
212 /// </table>
213 ///
214 
215 RooDataSet::RooDataSet(const char* name, const char* title, const RooArgSet& vars, const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3,
216  const RooCmdArg& arg4,const RooCmdArg& arg5,const RooCmdArg& arg6,const RooCmdArg& arg7,const RooCmdArg& arg8) :
217  RooAbsData(name,title,RooArgSet(vars,(RooAbsArg*)RooCmdConfig::decodeObjOnTheFly("RooDataSet::RooDataSet", "IndexCat",0,0,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8)))
218 {
219  // Define configuration for this method
220  RooCmdConfig pc(Form("RooDataSet::ctor(%s)",GetName())) ;
221  pc.defineInt("ownLinked","OwnLinked",0) ;
222  pc.defineObject("impTree","ImportTree",0) ;
223  pc.defineObject("impData","ImportData",0) ;
224  pc.defineObject("indexCat","IndexCat",0) ;
225  pc.defineObject("impSliceData","ImportDataSlice",0,0,kTRUE) ; // array
226  pc.defineString("impSliceState","ImportDataSlice",0,"",kTRUE) ; // array
227  pc.defineObject("lnkSliceData","LinkDataSlice",0,0,kTRUE) ; // array
228  pc.defineString("lnkSliceState","LinkDataSlice",0,"",kTRUE) ; // array
229  pc.defineString("cutSpec","CutSpec",0,"") ;
230  pc.defineObject("cutVar","CutVar",0) ;
231  pc.defineString("cutRange","CutRange",0,"") ;
232  pc.defineString("wgtVarName","WeightVarName",0,"") ;
233  pc.defineInt("newWeight1","WeightVarName",0,0) ;
234  pc.defineString("fname","ImportFromFile",0,"") ;
235  pc.defineString("tname","ImportFromFile",1,"") ;
236  pc.defineObject("wgtVar","WeightVar",0) ;
237  pc.defineInt("newWeight2","WeightVar",0,0) ;
238  pc.defineObject("dummy1","ImportDataSliceMany",0) ;
239  pc.defineObject("dummy2","LinkDataSliceMany",0) ;
240  pc.defineSet("errorSet","StoreError",0) ;
241  pc.defineSet("asymErrSet","StoreAsymError",0) ;
242  pc.defineMutex("ImportTree","ImportData","ImportDataSlice","LinkDataSlice","ImportFromFile") ;
243  pc.defineMutex("CutSpec","CutVar") ;
244  pc.defineMutex("WeightVarName","WeightVar") ;
245  pc.defineDependency("ImportDataSlice","IndexCat") ;
246  pc.defineDependency("LinkDataSlice","IndexCat") ;
247  pc.defineDependency("OwnLinked","LinkDataSlice") ;
248 
249 
250  RooLinkedList l ;
251  l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
252  l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
253  l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
254  l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
255 
256  // Process & check varargs
257  pc.process(l) ;
258  if (!pc.ok(kTRUE)) {
259  assert(0) ;
260  return ;
261  }
262 
263  // Extract relevant objects
264  TTree* impTree = static_cast<TTree*>(pc.getObject("impTree")) ;
265  RooDataSet* impData = static_cast<RooDataSet*>(pc.getObject("impData")) ;
266  RooFormulaVar* cutVar = static_cast<RooFormulaVar*>(pc.getObject("cutVar")) ;
267  const char* cutSpec = pc.getString("cutSpec","",kTRUE) ;
268  const char* cutRange = pc.getString("cutRange","",kTRUE) ;
269  const char* wgtVarName = pc.getString("wgtVarName","",kTRUE) ;
270  RooRealVar* wgtVar = static_cast<RooRealVar*>(pc.getObject("wgtVar")) ;
271  const char* impSliceNames = pc.getString("impSliceState","",kTRUE) ;
272  const RooLinkedList& impSliceData = pc.getObjectList("impSliceData") ;
273  const char* lnkSliceNames = pc.getString("lnkSliceState","",kTRUE) ;
274  const RooLinkedList& lnkSliceData = pc.getObjectList("lnkSliceData") ;
275  RooCategory* indexCat = static_cast<RooCategory*>(pc.getObject("indexCat")) ;
276  RooArgSet* errorSet = pc.getSet("errorSet") ;
277  RooArgSet* asymErrorSet = pc.getSet("asymErrSet") ;
278  const char* fname = pc.getString("fname") ;
279  const char* tname = pc.getString("tname") ;
280  Int_t ownLinked = pc.getInt("ownLinked") ;
281  Int_t newWeight = pc.getInt("newWeight1") + pc.getInt("newWeight2") ;
282 
283  // Case 1 --- Link multiple dataset as slices
284  if (lnkSliceNames) {
285 
286  // Make import mapping if index category is specified
287  map<string,RooAbsData*> hmap ;
288  if (indexCat) {
289  char tmp[64000];
290  strlcpy(tmp, lnkSliceNames, 64000);
291  char *token = strtok(tmp, ",");
292  TIterator *hiter = lnkSliceData.MakeIterator();
293  while (token) {
294  hmap[token] = (RooAbsData *)hiter->Next();
295  token = strtok(0, ",");
296  }
297  delete hiter ;
298  }
299 
300  // Lookup name of weight variable if it was specified by object reference
301  if (wgtVar) {
302  // coverity[UNUSED_VALUE]
303  wgtVarName = wgtVar->GetName() ;
304  }
305 
306  appendToDir(this,kTRUE) ;
307 
308  // Initialize RooDataSet with optional weight variable
309  initialize(0) ;
310 
311  map<string,RooAbsDataStore*> storeMap ;
312  RooCategory* icat = (RooCategory*) (indexCat ? _vars.find(indexCat->GetName()) : 0 ) ;
313  if (!icat) {
314  throw std::string("RooDataSet::RooDataSet() ERROR in constructor, cannot find index category") ;
315  }
316  for (map<string,RooAbsData*>::iterator hiter = hmap.begin() ; hiter!=hmap.end() ; ++hiter) {
317  // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
318  if (indexCat && !indexCat->hasLabel(hiter->first)) {
319  indexCat->defineType(hiter->first) ;
320  coutI(InputArguments) << "RooDataSet::ctor(" << GetName() << ") defining state \"" << hiter->first << "\" in index category " << indexCat->GetName() << endl ;
321  }
322  if (icat && !icat->hasLabel(hiter->first)) {
323  icat->defineType(hiter->first) ;
324  }
325  icat->setLabel(hiter->first.c_str()) ;
326  storeMap[icat->getCurrentLabel()]=hiter->second->store() ;
327 
328  // Take ownership of slice if requested
329  if (ownLinked) {
330  addOwnedComponent(hiter->first.c_str(),*hiter->second) ;
331  }
332  }
333 
334  // Create composite datastore
335  _dstore = new RooCompositeDataStore(name,title,_vars,*icat,storeMap) ;
336 
337  } else {
338 
339  if (wgtVar) {
340  wgtVarName = wgtVar->GetName() ;
341  }
342 
343  // Clone weight variable of imported dataset if we are not weighted
344  if (!wgtVar && !wgtVarName && impData && impData->_wgtVar) {
345  _wgtVar = (RooRealVar*) impData->_wgtVar->createFundamental() ;
347  wgtVarName = _wgtVar->GetName() ;
348  }
349 
350  // Create empty datastore
351  RooTreeDataStore* tstore(0) ;
352  RooVectorDataStore* vstore(0) ;
353 
354  if (defaultStorageType==Tree) {
355  tstore = new RooTreeDataStore(name,title,_vars,wgtVarName) ;
356  _dstore = tstore ;
357  } else if (defaultStorageType==Vector) {
358  if (wgtVarName && newWeight) {
359  RooAbsArg* wgttmp = _vars.find(wgtVarName) ;
360  if (wgttmp) {
361  wgttmp->setAttribute("NewWeight") ;
362  }
363  }
364  vstore = new RooVectorDataStore(name,title,_vars,wgtVarName) ;
365  _dstore = vstore ;
366  } else {
367  _dstore = 0 ;
368  }
369 
370 
371  // Make import mapping if index category is specified
372  map<string,RooDataSet*> hmap ;
373  if (indexCat) {
374  TIterator* hiter = impSliceData.MakeIterator() ;
375  for (const auto& token : RooHelpers::tokenise(impSliceNames, ",")) {
376  hmap[token] = (RooDataSet*) hiter->Next() ;
377  }
378  delete hiter ;
379  }
380 
381  // process StoreError requests
382  if (errorSet) {
383  RooArgSet* intErrorSet = (RooArgSet*) _vars.selectCommon(*errorSet) ;
384  intErrorSet->setAttribAll("StoreError") ;
385  TIterator* iter = intErrorSet->createIterator() ;
386  RooAbsArg* arg ;
387  while((arg=(RooAbsArg*)iter->Next())) {
388  arg->attachToStore(*_dstore) ;
389  }
390  delete iter ;
391  delete intErrorSet ;
392  }
393  if (asymErrorSet) {
394  RooArgSet* intAsymErrorSet = (RooArgSet*) _vars.selectCommon(*asymErrorSet) ;
395  intAsymErrorSet->setAttribAll("StoreAsymError") ;
396  TIterator* iter = intAsymErrorSet->createIterator() ;
397  RooAbsArg* arg ;
398  while((arg=(RooAbsArg*)iter->Next())) {
399  arg->attachToStore(*_dstore) ;
400  }
401  delete iter ;
402  delete intAsymErrorSet ;
403  }
404 
405  // Lookup name of weight variable if it was specified by object reference
406  if (wgtVar) {
407  wgtVarName = wgtVar->GetName() ;
408  }
409 
410 
411  appendToDir(this,kTRUE) ;
412 
413  // Initialize RooDataSet with optional weight variable
414  if (wgtVarName && *wgtVarName) {
415  // Use the supplied weight column
416  initialize(wgtVarName) ;
417 
418  } else {
419  if (impData && impData->_wgtVar && vars.find(impData->_wgtVar->GetName())) {
420 
421  // Use the weight column of the source data set
422  initialize(impData->_wgtVar->GetName()) ;
423 
424  } else if (indexCat) {
425 
426  RooDataSet* firstDS = hmap.begin()->second ;
427  if (firstDS->_wgtVar && vars.find(firstDS->_wgtVar->GetName())) {
428  initialize(firstDS->_wgtVar->GetName()) ;
429  } else {
430  initialize(0) ;
431  }
432  } else {
433  initialize(0) ;
434  }
435  }
436 
437  // Import one or more datasets with a cut specification
438  if (cutSpec && *cutSpec) {
439 
440  // Create a RooFormulaVar cut from given cut expression
441  if (indexCat) {
442 
443  // Case 2a --- Import multiple RooDataSets as slices with cutspec
444  RooCategory* icat = (RooCategory*) _vars.find(indexCat->GetName()) ;
445  for (map<string,RooDataSet*>::iterator hiter = hmap.begin() ; hiter!=hmap.end() ; ++hiter) {
446  // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
447  if (!indexCat->hasLabel(hiter->first)) {
448  indexCat->defineType(hiter->first) ;
449  coutI(InputArguments) << "RooDataSet::ctor(" << GetName() << ") defining state \"" << hiter->first << "\" in index category " << indexCat->GetName() << endl ;
450  }
451  if (!icat->hasLabel(hiter->first)) {
452  icat->defineType(hiter->first) ;
453  }
454  icat->setLabel(hiter->first.c_str()) ;
455 
456  RooFormulaVar cutVarTmp(cutSpec,cutSpec,hiter->second->_vars) ;
457  _dstore->loadValues(hiter->second->store(),&cutVarTmp,cutRange) ;
458  }
459 
460  } else if (impData) {
461 
462  // Case 3a --- Import RooDataSet with cutspec
463  RooFormulaVar cutVarTmp(cutSpec,cutSpec,impData->_vars) ;
464  _dstore->loadValues(impData->store(),&cutVarTmp,cutRange);
465  } else if (impTree) {
466 
467  // Case 4a --- Import TTree from memory with cutspec
468  RooFormulaVar cutVarTmp(cutSpec,cutSpec,_vars) ;
469  if (tstore) {
470  tstore->loadValues(impTree,&cutVarTmp,cutRange);
471  } else {
472  RooTreeDataStore tmpstore(name,title,_vars,wgtVarName) ;
473  tmpstore.loadValues(impTree,&cutVarTmp,cutRange) ;
474  _dstore->append(tmpstore) ;
475  }
476  } else if (fname && strlen(fname)) {
477 
478  // Case 5a --- Import TTree from file with cutspec
479  TFile *f = TFile::Open(fname) ;
480  if (!f) {
481  coutE(InputArguments) << "RooDataSet::ctor(" << GetName() << ") ERROR file '" << fname << "' cannot be opened or does not exist" << endl ;
482  throw string(Form("RooDataSet::ctor(%s) ERROR file %s cannot be opened or does not exist",GetName(),fname)) ;
483  }
484  TTree* t = dynamic_cast<TTree*>(f->Get(tname)) ;
485  if (!t) {
486  coutE(InputArguments) << "RooDataSet::ctor(" << GetName() << ") ERROR file '" << fname << "' does not contain a TTree named '" << tname << "'" << endl ;
487  throw string(Form("RooDataSet::ctor(%s) ERROR file %s does not contain a TTree named %s",GetName(),fname,tname)) ;
488  }
489  RooFormulaVar cutVarTmp(cutSpec,cutSpec,_vars) ;
490  if (tstore) {
491  tstore->loadValues(t,&cutVarTmp,cutRange);
492  } else {
493  RooTreeDataStore tmpstore(name,title,_vars,wgtVarName) ;
494  tmpstore.loadValues(t,&cutVarTmp,cutRange) ;
495  _dstore->append(tmpstore) ;
496  }
497  f->Close() ;
498 
499  }
500 
501  // Import one or more datasets with a cut formula
502  } else if (cutVar) {
503 
504  if (indexCat) {
505 
506  // Case 2b --- Import multiple RooDataSets as slices with cutvar
507 
508  RooCategory* icat = (RooCategory*) _vars.find(indexCat->GetName()) ;
509  for (map<string,RooDataSet*>::iterator hiter = hmap.begin() ; hiter!=hmap.end() ; ++hiter) {
510  // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
511  if (!indexCat->hasLabel(hiter->first)) {
512  indexCat->defineType(hiter->first) ;
513  coutI(InputArguments) << "RooDataSet::ctor(" << GetName() << ") defining state \"" << hiter->first << "\" in index category " << indexCat->GetName() << endl ;
514  }
515  if (!icat->hasLabel(hiter->first)) {
516  icat->defineType(hiter->first) ;
517  }
518  icat->setLabel(hiter->first.c_str()) ;
519  _dstore->loadValues(hiter->second->store(),cutVar,cutRange) ;
520  }
521 
522 
523  } else if (impData) {
524  // Case 3b --- Import RooDataSet with cutvar
525  _dstore->loadValues(impData->store(),cutVar,cutRange);
526  } else if (impTree) {
527  // Case 4b --- Import TTree from memory with cutvar
528  if (tstore) {
529  tstore->loadValues(impTree,cutVar,cutRange);
530  } else {
531  RooTreeDataStore tmpstore(name,title,_vars,wgtVarName) ;
532  tmpstore.loadValues(impTree,cutVar,cutRange) ;
533  _dstore->append(tmpstore) ;
534  }
535  } else if (fname && strlen(fname)) {
536  // Case 5b --- Import TTree from file with cutvar
537  TFile *f = TFile::Open(fname) ;
538  if (!f) {
539  coutE(InputArguments) << "RooDataSet::ctor(" << GetName() << ") ERROR file '" << fname << "' cannot be opened or does not exist" << endl ;
540  throw string(Form("RooDataSet::ctor(%s) ERROR file %s cannot be opened or does not exist",GetName(),fname)) ;
541  }
542  TTree* t = dynamic_cast<TTree*>(f->Get(tname)) ;
543  if (!t) {
544  coutE(InputArguments) << "RooDataSet::ctor(" << GetName() << ") ERROR file '" << fname << "' does not contain a TTree named '" << tname << "'" << endl ;
545  throw string(Form("RooDataSet::ctor(%s) ERROR file %s does not contain a TTree named %s",GetName(),fname,tname)) ;
546  }
547  if (tstore) {
548  tstore->loadValues(t,cutVar,cutRange);
549  } else {
550  RooTreeDataStore tmpstore(name,title,_vars,wgtVarName) ;
551  tmpstore.loadValues(t,cutVar,cutRange) ;
552  _dstore->append(tmpstore) ;
553  }
554 
555  f->Close() ;
556  }
557 
558  // Import one or more datasets without cuts
559  } else {
560 
561  if (indexCat) {
562 
563  RooCategory* icat = (RooCategory*) _vars.find(indexCat->GetName()) ;
564  for (map<string,RooDataSet*>::iterator hiter = hmap.begin() ; hiter!=hmap.end() ; ++hiter) {
565  // Define state labels in index category (both in provided indexCat and in internal copy in dataset)
566  if (!indexCat->hasLabel(hiter->first)) {
567  indexCat->defineType(hiter->first) ;
568  coutI(InputArguments) << "RooDataSet::ctor(" << GetName() << ") defining state \"" << hiter->first << "\" in index category " << indexCat->GetName() << endl ;
569  }
570  if (!icat->hasLabel(hiter->first)) {
571  icat->defineType(hiter->first) ;
572  }
573  icat->setLabel(hiter->first.c_str()) ;
574  // Case 2c --- Import multiple RooDataSets as slices
575  _dstore->loadValues(hiter->second->store(),0,cutRange) ;
576  }
577 
578  } else if (impData) {
579  // Case 3c --- Import RooDataSet
580  _dstore->loadValues(impData->store(),0,cutRange);
581 
582  } else if (impTree || (fname && strlen(fname))) {
583  // Case 4c --- Import TTree from memory / file
584  std::unique_ptr<TFile> file;
585 
586  if (impTree == nullptr) {
587  file.reset(TFile::Open(fname));
588  if (!file) {
589  coutE(InputArguments) << "RooDataSet::ctor(" << GetName() << ") ERROR file '" << fname << "' cannot be opened or does not exist" << endl ;
590  throw std::invalid_argument(Form("RooDataSet::ctor(%s) ERROR file %s cannot be opened or does not exist",GetName(),fname)) ;
591  }
592 
593  file->GetObject(tname, impTree);
594  if (!impTree) {
595  coutE(InputArguments) << "RooDataSet::ctor(" << GetName() << ") ERROR file '" << fname << "' does not contain a TTree named '" << tname << "'" << endl ;
596  throw std::invalid_argument(Form("RooDataSet::ctor(%s) ERROR file %s does not contain a TTree named %s",GetName(),fname,tname)) ;
597  }
598  }
599 
600  if (tstore) {
601  tstore->loadValues(impTree,0,cutRange);
602  } else {
603  RooTreeDataStore tmpstore(name,title,_vars,wgtVarName) ;
604  tmpstore.loadValues(impTree,0,cutRange) ;
605  _dstore->append(tmpstore) ;
606  }
607  }
608  }
609 
610  }
612 }
613 
614 
615 
616 ////////////////////////////////////////////////////////////////////////////////
617 /// Constructor of an empty data set from a RooArgSet defining the dimensions
618 /// of the data space.
619 
620 RooDataSet::RooDataSet(const char *name, const char *title, const RooArgSet& vars, const char* wgtVarName) :
621  RooAbsData(name,title,vars)
622 {
623 // cout << "RooDataSet::ctor(" << this << ") storageType = " << ((defaultStorageType==Tree)?"Tree":"Vector") << endl ;
624  _dstore = (defaultStorageType==Tree) ? ((RooAbsDataStore*) new RooTreeDataStore(name,title,_vars,wgtVarName)) :
625  ((RooAbsDataStore*) new RooVectorDataStore(name,title,_vars,wgtVarName)) ;
626 
627  appendToDir(this,kTRUE) ;
628  initialize(wgtVarName) ;
630 }
631 
632 
633 ////////////////////////////////////////////////////////////////////////////////
634 /// Constructor of a data set from (part of) an existing data
635 /// set. The dimensions of the data set are defined by the 'vars'
636 /// RooArgSet, which can be identical to 'dset' dimensions, or a
637 /// subset thereof. The 'cuts' string is an optional RooFormula
638 /// expression and can be used to select the subset of the data
639 /// points in 'dset' to be copied. The cut expression can refer to
640 /// any variable in the source dataset. For cuts involving variables
641 /// other than those contained in the source data set, such as
642 /// intermediate formula objects, use the equivalent constructor
643 /// accepting RooFormulaVar reference as cut specification.
644 ///
645 /// This constructor will internally store the data in a TTree.
646 ///
647 /// For most uses the RooAbsData::reduce() wrapper function, which
648 /// uses this constructor, is the most convenient way to create a
649 /// subset of an existing data
650 ///
651 
652 RooDataSet::RooDataSet(const char *name, const char *title, RooDataSet *dset,
653  const RooArgSet& vars, const char *cuts, const char* wgtVarName) :
654  RooAbsData(name,title,vars)
655 {
656  // Initialize datastore
657  _dstore = new RooTreeDataStore(name,title,_vars,*dset->_dstore,cuts,wgtVarName) ;
658 
659  appendToDir(this,kTRUE) ;
660 
661  if (wgtVarName) {
662  // Use the supplied weight column
663  initialize(wgtVarName) ;
664  } else {
665  if (dset->_wgtVar && vars.find(dset->_wgtVar->GetName())) {
666  // Use the weight column of the source data set
667  initialize(dset->_wgtVar->GetName()) ;
668  } else {
669  initialize(0) ;
670  }
671  }
673 }
674 
675 
676 ////////////////////////////////////////////////////////////////////////////////
677 /// Constructor of a data set from (part of) an existing data
678 /// set. The dimensions of the data set are defined by the 'vars'
679 /// RooArgSet, which can be identical to 'dset' dimensions, or a
680 /// subset thereof. The 'cutVar' formula variable is used to select
681 /// the subset of data points to be copied. For subsets without
682 /// selection on the data points, or involving cuts operating
683 /// exclusively and directly on the data set dimensions, the
684 /// equivalent constructor with a string based cut expression is
685 /// recommended.
686 ///
687 /// This constructor will internally store the data in a TTree.
688 ///
689 /// For most uses the RooAbsData::reduce() wrapper function, which
690 /// uses this constructor, is the most convenient way to create a
691 /// subset of an existing data
692 
693 RooDataSet::RooDataSet(const char *name, const char *title, RooDataSet *dset,
694  const RooArgSet& vars, const RooFormulaVar& cutVar, const char* wgtVarName) :
695  RooAbsData(name,title,vars)
696 {
697  // Initialize datastore
698  _dstore = new RooTreeDataStore(name,title,_vars,*dset->_dstore,cutVar,wgtVarName) ;
699 
700  appendToDir(this,kTRUE) ;
701 
702  if (wgtVarName) {
703  // Use the supplied weight column
704  initialize(wgtVarName) ;
705  } else {
706  if (dset->_wgtVar && vars.find(dset->_wgtVar->GetName())) {
707  // Use the weight column of the source data set
708  initialize(dset->_wgtVar->GetName()) ;
709  } else {
710  initialize(0) ;
711  }
712  }
714 }
715 
716 
717 
718 
719 ////////////////////////////////////////////////////////////////////////////////
720 /// Constructor of a data set from (part of) an ROOT TTRee. The dimensions
721 /// of the data set are defined by the 'vars' RooArgSet. For each dimension
722 /// specified, the TTree must have a branch with the same name. For category
723 /// branches, this branch should contain the numeric index value. Real dimensions
724 /// can be constructed from either 'Double_t' or 'Float_t' tree branches. In the
725 /// latter case, an automatic conversion is applied.
726 ///
727 /// The 'cutVar' formula variable
728 /// is used to select the subset of data points to be copied.
729 /// For subsets without selection on the data points, or involving cuts
730 /// operating exclusively and directly on the data set dimensions, the equivalent
731 /// constructor with a string based cut expression is recommended.
732 
733 RooDataSet::RooDataSet(const char *name, const char *title, TTree *theTree,
734  const RooArgSet& vars, const RooFormulaVar& cutVar, const char* wgtVarName) :
735  RooAbsData(name,title,vars)
736 {
737  // Create tree version of datastore
738  RooTreeDataStore* tstore = new RooTreeDataStore(name,title,_vars,*theTree,cutVar,wgtVarName) ;
739 
740  // Convert to vector datastore if needed
741  if (defaultStorageType==Tree) {
742  _dstore = tstore ;
743  } else if (defaultStorageType==Vector) {
744  RooVectorDataStore* vstore = new RooVectorDataStore(name,title,_vars,wgtVarName) ;
745  _dstore = vstore ;
746  _dstore->append(*tstore) ;
747  delete tstore ;
748  } else {
749  _dstore = 0 ;
750  }
751 
752  appendToDir(this,kTRUE) ;
753  initialize(wgtVarName) ;
755 }
756 
757 
758 
759 ////////////////////////////////////////////////////////////////////////////////
760 /// Constructor of a data set from (part of) a ROOT TTree.
761 ///
762 /// \param[in] name Name of this dataset.
763 /// \param[in] title Title for e.g. plotting.
764 /// \param[in] theTree Tree to be imported.
765 /// \param[in] vars Defines the columns of the data set. For each dimension
766 /// specified, the TTree must have a branch with the same name. For category
767 /// branches, this branch should contain the numeric index value. Real dimensions
768 /// can be constructed from either 'Double_t' or 'Float_t' tree branches. In the
769 /// latter case, an automatic conversion is applied.
770 /// \param[in] cuts Optional RooFormula expression to select the subset of the data points
771 /// to be imported. The cut expression can refer to any variable in `vars`.
772 /// \warning The expression only evaluates variables that are also in `vars`.
773 /// Passing e.g.
774 /// ```
775 /// RooDataSet("data", "data", tree, RooArgSet(x), "x>y")
776 /// ```
777 /// Will load `x` from the tree, but leave `y` at an undefined value.
778 /// If other expressions are needed, such as intermediate formula objects, use
779 /// RooDataSet::RooDataSet(const char*,const char*,TTree*,const RooArgSet&,const RooFormulaVar&,const char*)
780 /// \param[in] wgtVarName Name of the variable in `vars` that represents an event weight.
781 RooDataSet::RooDataSet(const char* name, const char* title, TTree* theTree,
782  const RooArgSet& vars, const char* cuts, const char* wgtVarName) :
783  RooAbsData(name,title,vars)
784 {
785  // Create tree version of datastore
786  RooTreeDataStore* tstore = new RooTreeDataStore(name,title,_vars,*theTree,cuts,wgtVarName);
787 
788  // Convert to vector datastore if needed
789  if (defaultStorageType==Tree) {
790  _dstore = tstore ;
791  } else if (defaultStorageType==Vector) {
792  RooVectorDataStore* vstore = new RooVectorDataStore(name,title,_vars,wgtVarName) ;
793  _dstore = vstore ;
794  _dstore->append(*tstore) ;
795  delete tstore ;
796  } else {
797  _dstore = 0 ;
798  }
799 
800  appendToDir(this,kTRUE) ;
801 
802  initialize(wgtVarName) ;
804 }
805 
806 
807 
808 ////////////////////////////////////////////////////////////////////////////////
809 /// Copy constructor
810 
811 RooDataSet::RooDataSet(RooDataSet const & other, const char* newname) :
812  RooAbsData(other,newname), RooDirItem()
813 {
814  appendToDir(this,kTRUE) ;
815  initialize(other._wgtVar?other._wgtVar->GetName():0) ;
817 }
818 
819 ////////////////////////////////////////////////////////////////////////////////
820 /// Protected constructor for internal use only
821 
822 RooDataSet::RooDataSet(const char *name, const char *title, RooDataSet *dset,
823  const RooArgSet& vars, const RooFormulaVar* cutVar, const char* cutRange,
824  std::size_t nStart, std::size_t nStop, Bool_t copyCache, const char* wgtVarName) :
825  RooAbsData(name,title,vars)
826 {
827  if (defaultStorageType == Tree) {
828  _dstore = new RooTreeDataStore(name, title, *dset->_dstore, _vars, cutVar, cutRange, nStart, nStop,
829  copyCache, wgtVarName);
830  } else {
831  _dstore = new RooVectorDataStore(name, title, *dset->_dstore, _vars, cutVar, cutRange, nStart,
832  nStop, copyCache, wgtVarName);
833  }
834 
836 
837  appendToDir(this, kTRUE);
838  initialize(dset->_wgtVar ? dset->_wgtVar->GetName() : 0);
840 }
841 
842 
843 ////////////////////////////////////////////////////////////////////////////////
844 /// Helper function for constructor that adds optional weight variable to construct
845 /// total set of observables
846 
847 RooArgSet RooDataSet::addWgtVar(const RooArgSet& origVars, const RooAbsArg* wgtVar)
848 {
849  RooArgSet tmp(origVars) ;
850  if (wgtVar) tmp.add(*wgtVar) ;
851  return tmp ;
852 }
853 
854 
855 
856 ////////////////////////////////////////////////////////////////////////////////
857 /// Return a clone of this dataset containing only the cached variables
858 
859 RooAbsData* RooDataSet::cacheClone(const RooAbsArg* newCacheOwner, const RooArgSet* newCacheVars, const char* newName)
860 {
861  RooDataSet* dset = new RooDataSet(newName?newName:GetName(),GetTitle(),this,_vars,(RooFormulaVar*)0,0,0,2000000000,kTRUE,_wgtVar?_wgtVar->GetName():0) ;
862  //if (_wgtVar) dset->setWeightVar(_wgtVar->GetName()) ;
863 
864  RooArgSet* selCacheVars = (RooArgSet*) newCacheVars->selectCommon(dset->_cachedVars) ;
865  dset->attachCache(newCacheOwner, *selCacheVars) ;
866  delete selCacheVars ;
867 
868  return dset ;
869 }
870 
871 
872 
873 ////////////////////////////////////////////////////////////////////////////////
874 /// Return an empty clone of this dataset. If vars is not null, only the variables in vars
875 /// are added to the definition of the empty clone
876 
877 RooAbsData* RooDataSet::emptyClone(const char* newName, const char* newTitle, const RooArgSet* vars, const char* wgtVarName) const
878 {
879  // If variables are given, be sure to include weight variable if it exists and is not included
880  RooArgSet vars2 ;
881  RooRealVar* tmpWgtVar = _wgtVar ;
882  if (wgtVarName && vars && !_wgtVar) {
883  tmpWgtVar = (RooRealVar*) vars->find(wgtVarName) ;
884  }
885 
886  if (vars) {
887  vars2.add(*vars) ;
888  if (_wgtVar && !vars2.find(_wgtVar->GetName())) {
889  vars2.add(*_wgtVar) ;
890  }
891  } else {
892  vars2.add(_vars) ;
893  }
894 
895  RooDataSet* dset = new RooDataSet(newName?newName:GetName(),newTitle?newTitle:GetTitle(),vars2,tmpWgtVar?tmpWgtVar->GetName():0) ;
896  //if (_wgtVar) dset->setWeightVar(_wgtVar->GetName()) ;
897  return dset ;
898 }
899 
900 
901 
902 ////////////////////////////////////////////////////////////////////////////////
903 /// Initialize the dataset. If wgtVarName is not null, interpret the observable
904 /// with that name as event weight
905 
906 void RooDataSet::initialize(const char* wgtVarName)
907 {
909  _varsNoWgt.add(_vars) ;
910  _wgtVar = 0 ;
911  if (wgtVarName) {
912  RooAbsArg* wgt = _varsNoWgt.find(wgtVarName) ;
913  if (!wgt) {
914  coutE(DataHandling) << "RooDataSet::RooDataSet(" << GetName() << "): designated weight variable "
915  << wgtVarName << " not found in set of variables, no weighting will be assigned" << endl ;
916  throw std::invalid_argument("RooDataSet::initialize() weight variable could not be initialised.");
917  } else if (!dynamic_cast<RooRealVar*>(wgt)) {
918  coutE(DataHandling) << "RooDataSet::RooDataSet(" << GetName() << "): designated weight variable "
919  << wgtVarName << " is not of type RooRealVar, no weighting will be assigned" << endl ;
920  throw std::invalid_argument("RooDataSet::initialize() weight variable could not be initialised.");
921  } else {
922  _varsNoWgt.remove(*wgt) ;
923  _wgtVar = (RooRealVar*) wgt ;
924  }
925  }
926 }
927 
928 
929 
930 ////////////////////////////////////////////////////////////////////////////////
931 /// Implementation of RooAbsData virtual method that drives the RooAbsData::reduce() methods
932 
933 RooAbsData* RooDataSet::reduceEng(const RooArgSet& varSubset, const RooFormulaVar* cutVar, const char* cutRange,
934  std::size_t nStart, std::size_t nStop, Bool_t copyCache)
935 {
936  checkInit() ;
937 
938  //cout << "reduceEng varSubset = " << varSubset << " _wgtVar = " << (_wgtVar ? _wgtVar->GetName() : "") << endl;
939 
940  RooArgSet tmp(varSubset) ;
941  if (_wgtVar) {
942  tmp.add(*_wgtVar) ;
943  }
944  RooDataSet* ret = new RooDataSet(GetName(), GetTitle(), this, tmp, cutVar, cutRange, nStart, nStop, copyCache,_wgtVar?_wgtVar->GetName():0) ;
945 
946  // WVE - propagate optional weight variable
947  // check behaviour in plotting.
948  // if (_wgtVar) {
949  // ret->setWeightVar(_wgtVar->GetName()) ;
950  // }
951  return ret ;
952 }
953 
954 
955 
956 ////////////////////////////////////////////////////////////////////////////////
957 /// Destructor
958 
960 {
961  removeFromDir(this) ;
963 }
964 
965 
966 
967 ////////////////////////////////////////////////////////////////////////////////
968 /// Return binned clone of this dataset
969 
970 RooDataHist* RooDataSet::binnedClone(const char* newName, const char* newTitle) const
971 {
972  TString title, name ;
973  if (newName) {
974  name = newName ;
975  } else {
976  name = Form("%s_binned",GetName()) ;
977  }
978  if (newTitle) {
979  title = newTitle ;
980  } else {
981  title = Form("%s_binned",GetTitle()) ;
982  }
983 
984  return new RooDataHist(name,title,*get(),*this) ;
985 }
986 
987 
988 
989 ////////////////////////////////////////////////////////////////////////////////
990 /// Return event weight of current event
991 
993 {
994  return store()->weight() ;
995 }
996 
997 
998 
999 ////////////////////////////////////////////////////////////////////////////////
1000 /// Return squared event weight of current event
1001 
1003 {
1004  return store()->weight()*store()->weight() ;
1005 }
1006 
1007 
1008 ////////////////////////////////////////////////////////////////////////////////
1009 /// \see RooAbsData::getWeightBatch().
1010 RooSpan<const double> RooDataSet::getWeightBatch(std::size_t first, std::size_t len) const {
1011  return _dstore->getWeightBatch(first, len);
1012 }
1013 
1014 
1015 ////////////////////////////////////////////////////////////////////////////////
1016 /// Write information to retrieve data columns into `evalData.spans`.
1017 /// All spans belonging to variables of this dataset are overwritten. Spans to other
1018 /// variables remain intact.
1019 /// \param[out] evalData Store references to all data batches in this struct's `spans`.
1020 /// The key to retrieve an item is the pointer of the variable that owns the data.
1021 /// \param first Index of first event that ends up in the batch.
1022 /// \param len Number of events in each batch.
1023 void RooDataSet::getBatches(RooBatchCompute::RunContext& evalData, std::size_t begin, std::size_t len) const {
1024  for (auto&& batch : store()->getBatches(begin, len).spans) {
1025  evalData.spans[batch.first] = std::move(batch.second);
1026  }
1027 }
1028 
1029 
1030 ////////////////////////////////////////////////////////////////////////////////
1031 
1033 {
1034  store()->weightError(lo,hi,etype) ;
1035 }
1036 
1037 
1038 
1039 ////////////////////////////////////////////////////////////////////////////////
1040 
1042 {
1043  return store()->weightError(etype) ;
1044 }
1045 
1046 
1047 
1048 ////////////////////////////////////////////////////////////////////////////////
1049 /// Return RooArgSet with coordinates of event 'index'
1050 
1051 const RooArgSet* RooDataSet::get(Int_t index) const
1052 {
1053  const RooArgSet* ret = RooAbsData::get(index) ;
1054  return ret ? &_varsNoWgt : 0 ;
1055 }
1056 
1057 
1058 ////////////////////////////////////////////////////////////////////////////////
1059 
1061 {
1062  return store()->sumEntries() ;
1063 
1064  //---------
1065 
1066  // Shortcut for unweighted unselected datasets
1067  if (!isWeighted()) {
1068  return numEntries() ;
1069  }
1070 
1071  // Otherwise sum the weights in the event
1072  Double_t sumw(0), carry(0);
1073  Int_t i ;
1074  for (i=0 ; i<numEntries() ; i++) {
1075  get(i) ;
1076  Double_t y = weight() - carry;
1077  Double_t t = sumw + y;
1078  carry = (t - sumw) - y;
1079  sumw = t;
1080  }
1081 
1082  return sumw ;
1083 }
1084 
1085 
1086 ////////////////////////////////////////////////////////////////////////////////
1087 /// Return the sum of weights in all entries matching cutSpec (if specified)
1088 /// and in named range cutRange (if specified)
1089 
1090 Double_t RooDataSet::sumEntries(const char* cutSpec, const char* cutRange) const
1091 {
1092  // Setup RooFormulaVar for cutSpec if it is present
1093  RooFormula* select = 0 ;
1094  if (cutSpec && strlen(cutSpec) > 0) {
1095  select = new RooFormula("select",cutSpec,*get()) ;
1096  }
1097 
1098  // Shortcut for unweighted unselected datasets
1099  if (!select && !cutRange && !isWeighted()) {
1100  return numEntries() ;
1101  }
1102 
1103  // Otherwise sum the weights in the event
1104  Double_t sumw(0), carry(0);
1105  Int_t i ;
1106  for (i=0 ; i<numEntries() ; i++) {
1107  get(i) ;
1108  if (select && select->eval()==0.) continue ;
1109  if (cutRange && !_vars.allInRange(cutRange)) continue ;
1110  Double_t y = weight() - carry;
1111  Double_t t = sumw + y;
1112  carry = (t - sumw) - y;
1113  sumw = t;
1114  }
1115 
1116  if (select) delete select ;
1117 
1118  return sumw ;
1119 }
1120 
1121 
1122 
1123 
1124 ////////////////////////////////////////////////////////////////////////////////
1125 /// Return true if dataset contains weighted events
1126 
1128 {
1129  return store() ? store()->isWeighted() : false;
1130 }
1131 
1132 
1133 
1134 ////////////////////////////////////////////////////////////////////////////////
1135 /// Returns true if histogram contains bins with entries with a non-integer weight
1136 
1138 {
1139  // Return false if we have no weights
1140  if (!_wgtVar) return kFALSE ;
1141 
1142  // Now examine individual weights
1143  for (int i=0 ; i<numEntries() ; i++) {
1144  get(i) ;
1145  if (fabs(weight()-Int_t(weight()))>1e-10) return kTRUE ;
1146  }
1147  // If sum of weights is less than number of events there are negative (integer) weights
1148  if (sumEntries()<numEntries()) return kTRUE ;
1149 
1150  return kFALSE ;
1151 }
1152 
1153 
1154 
1155 
1156 ////////////////////////////////////////////////////////////////////////////////
1157 /// Return a RooArgSet with the coordinates of the current event
1158 
1160 {
1161  return &_varsNoWgt ;
1162 }
1163 
1164 
1165 
1166 ////////////////////////////////////////////////////////////////////////////////
1167 /// Add a data point, with its coordinates specified in the 'data' argset, to the data set.
1168 /// Any variables present in 'data' but not in the dataset will be silently ignored.
1169 /// \param[in] data Data point.
1170 /// \param[in] wgt Event weight. Defaults to 1. The current value of the weight variable is
1171 /// ignored.
1172 /// \note To obtain weighted events, a variable must be designated `WeightVar` in the constructor.
1173 /// \param[in] wgtError Optional weight error.
1174 /// \note This requires including the weight variable in the set of `StoreError` variables when constructing
1175 /// the dataset.
1176 
1177 void RooDataSet::add(const RooArgSet& data, Double_t wgt, Double_t wgtError)
1178 {
1179  checkInit() ;
1180 
1181  const double oldW = _wgtVar ? _wgtVar->getVal() : 0.;
1182 
1183  _varsNoWgt = data;
1184 
1185  if (_wgtVar) {
1186  _wgtVar->setVal(wgt) ;
1187  if (wgtError!=0.) {
1188  _wgtVar->setError(wgtError) ;
1189  }
1190  } else if ((wgt != 1. || wgtError != 0.) && _errorMsgCount < 5) {
1191  ccoutE(DataHandling) << "An event weight/error was passed but no weight variable was defined"
1192  << " in the dataset '" << GetName() << "'. The weight will be ignored." << std::endl;
1193  ++_errorMsgCount;
1194  }
1195 
1197  && wgtError != 0.
1198  && fabs(wgt*wgt - wgtError)/wgtError > 1.E-15 //Exception for standard wgt^2 errors, which need not be stored.
1199  && _errorMsgCount < 5 && !_wgtVar->getAttribute("StoreError")) {
1200  coutE(DataHandling) << "An event weight error was passed to the RooDataSet '" << GetName()
1201  << "', but the weight variable '" << _wgtVar->GetName()
1202  << "' does not store errors. Check `StoreError` in the RooDataSet constructor." << std::endl;
1203  ++_errorMsgCount;
1204  }
1205 
1206  fill();
1207 
1208  // Restore weight state
1209  if (_wgtVar) {
1210  _wgtVar->setVal(oldW);
1211  _wgtVar->removeError();
1212  }
1213 }
1214 
1215 
1216 
1217 
1218 ////////////////////////////////////////////////////////////////////////////////
1219 /// Add a data point, with its coordinates specified in the 'data' argset, to the data set.
1220 /// Any variables present in 'data' but not in the dataset will be silently ignored.
1221 /// \param[in] indata Data point.
1222 /// \param[in] inweight Event weight. The current value of the weight variable is ignored.
1223 /// \note To obtain weighted events, a variable must be designated `WeightVar` in the constructor.
1224 /// \param[in] weightErrorLo Asymmetric weight error.
1225 /// \param[in] weightErrorHi Asymmetric weight error.
1226 /// \note This requires including the weight variable in the set of `StoreAsymError` variables when constructing
1227 /// the dataset.
1228 
1229 void RooDataSet::add(const RooArgSet& indata, Double_t inweight, Double_t weightErrorLo, Double_t weightErrorHi)
1230 {
1231  checkInit() ;
1232 
1233  const double oldW = _wgtVar ? _wgtVar->getVal() : 0.;
1234 
1235  _varsNoWgt = indata;
1236  if (_wgtVar) {
1237  _wgtVar->setVal(inweight) ;
1238  _wgtVar->setAsymError(weightErrorLo,weightErrorHi) ;
1239  } else if (inweight != 1. && _errorMsgCount < 5) {
1240  ccoutE(DataHandling) << "An event weight was given but no weight variable was defined"
1241  << " in the dataset '" << GetName() << "'. The weight will be ignored." << std::endl;
1242  ++_errorMsgCount;
1243  }
1244 
1246  && _errorMsgCount < 5 && !_wgtVar->getAttribute("StoreAsymError")) {
1247  coutE(DataHandling) << "An event weight error was passed to the RooDataSet '" << GetName()
1248  << "', but the weight variable '" << _wgtVar->GetName()
1249  << "' does not store errors. Check `StoreAsymError` in the RooDataSet constructor." << std::endl;
1250  ++_errorMsgCount;
1251  }
1252 
1253  fill();
1254 
1255  // Restore weight state
1256  if (_wgtVar) {
1257  _wgtVar->setVal(oldW);
1259  }
1260 }
1261 
1262 
1263 
1264 
1265 
1266 ////////////////////////////////////////////////////////////////////////////////
1267 /// Add a data point, with its coordinates specified in the 'data' argset, to the data set.
1268 /// \attention The order and type of the input variables are **assumed** to be the same as
1269 /// for the RooArgSet returned by RooDataSet::get(). Input values will just be written
1270 /// into the internal data columns by ordinal position.
1271 /// \param[in] data Data point.
1272 /// \param[in] wgt Event weight. Defaults to 1. The current value of the weight variable is
1273 /// ignored.
1274 /// \note To obtain weighted events, a variable must be designated `WeightVar` in the constructor.
1275 /// \param[in] wgtError Optional weight error.
1276 /// \note This requires including the weight variable in the set of `StoreError` variables when constructing
1277 /// the dataset.
1278 
1279 void RooDataSet::addFast(const RooArgSet& data, Double_t wgt, Double_t wgtError)
1280 {
1281  checkInit() ;
1282 
1283  const double oldW = _wgtVar ? _wgtVar->getVal() : 0.;
1284 
1286  if (_wgtVar) {
1287  _wgtVar->setVal(wgt) ;
1288  if (wgtError!=0.) {
1289  _wgtVar->setError(wgtError) ;
1290  }
1291  } else if (wgt != 1. && _errorMsgCount < 5) {
1292  ccoutE(DataHandling) << "An event weight was given but no weight variable was defined"
1293  << " in the dataset '" << GetName() << "'. The weight will be ignored." << std::endl;
1294  ++_errorMsgCount;
1295  }
1296 
1297  fill();
1298 
1300  && wgtError != 0. && wgtError != wgt*wgt //Exception for standard weight error, which need not be stored
1301  && _errorMsgCount < 5 && !_wgtVar->getAttribute("StoreError")) {
1302  coutE(DataHandling) << "An event weight error was passed to the RooDataSet '" << GetName()
1303  << "', but the weight variable '" << _wgtVar->GetName()
1304  << "' does not store errors. Check `StoreError` in the RooDataSet constructor." << std::endl;
1305  ++_errorMsgCount;
1306  }
1307  if (_wgtVar && _doWeightErrorCheck) {
1308  _doWeightErrorCheck = false;
1309  }
1310 
1311  if (_wgtVar) {
1312  _wgtVar->setVal(oldW);
1313  _wgtVar->removeError();
1314  }
1315 }
1316 
1317 
1318 
1319 ////////////////////////////////////////////////////////////////////////////////
1320 
1322  RooDataSet* data4, RooDataSet* data5, RooDataSet* data6)
1323 {
1324  checkInit() ;
1325  list<RooDataSet*> dsetList ;
1326  if (data1) dsetList.push_back(data1) ;
1327  if (data2) dsetList.push_back(data2) ;
1328  if (data3) dsetList.push_back(data3) ;
1329  if (data4) dsetList.push_back(data4) ;
1330  if (data5) dsetList.push_back(data5) ;
1331  if (data6) dsetList.push_back(data6) ;
1332  return merge(dsetList) ;
1333 }
1334 
1335 
1336 
1337 ////////////////////////////////////////////////////////////////////////////////
1338 /// Merge columns of supplied data set(s) with this data set. All
1339 /// data sets must have equal number of entries. In case of
1340 /// duplicate columns the column of the last dataset in the list
1341 /// prevails
1342 
1343 Bool_t RooDataSet::merge(list<RooDataSet*>dsetList)
1344 {
1345 
1346  checkInit() ;
1347  // Sanity checks: data sets must have the same size
1348  for (list<RooDataSet*>::iterator iter = dsetList.begin() ; iter != dsetList.end() ; ++iter) {
1349  if (numEntries()!=(*iter)->numEntries()) {
1350  coutE(InputArguments) << "RooDataSet::merge(" << GetName() << ") ERROR: datasets have different size" << endl ;
1351  return kTRUE ;
1352  }
1353  }
1354 
1355  // Extend vars with elements of other dataset
1356  list<RooAbsDataStore*> dstoreList ;
1357  for (list<RooDataSet*>::iterator iter = dsetList.begin() ; iter != dsetList.end() ; ++iter) {
1358  _vars.addClone((*iter)->_vars,kTRUE) ;
1359  dstoreList.push_back((*iter)->store()) ;
1360  }
1361 
1362  // Merge data stores
1363  RooAbsDataStore* mergedStore = _dstore->merge(_vars,dstoreList) ;
1364  mergedStore->SetName(_dstore->GetName()) ;
1365  mergedStore->SetTitle(_dstore->GetTitle()) ;
1366 
1367  // Replace current data store with merged store
1368  delete _dstore ;
1369  _dstore = mergedStore ;
1370 
1372  return kFALSE ;
1373 }
1374 
1375 
1376 ////////////////////////////////////////////////////////////////////////////////
1377 /// Add all data points of given data set to this data set.
1378 /// Observable in 'data' that are not in this dataset
1379 /// with not be transferred
1380 
1382 {
1383  checkInit() ;
1384  _dstore->append(*data._dstore) ;
1385 }
1386 
1387 
1388 
1389 ////////////////////////////////////////////////////////////////////////////////
1390 /// Add a column with the values of the given (function) argument
1391 /// to this dataset. The function value is calculated for each
1392 /// event using the observable values of each event in case the
1393 /// function depends on variables with names that are identical
1394 /// to the observable names in the dataset
1395 
1397 {
1398  checkInit() ;
1399  RooAbsArg* ret = _dstore->addColumn(var,adjustRange) ;
1400  _vars.addOwned(*ret) ;
1402  return ret ;
1403 }
1404 
1405 
1406 ////////////////////////////////////////////////////////////////////////////////
1407 /// Add a column with the values of the given list of (function)
1408 /// argument to this dataset. Each function value is calculated for
1409 /// each event using the observable values of the event in case the
1410 /// function depends on variables with names that are identical to
1411 /// the observable names in the dataset
1412 
1414 {
1415  checkInit() ;
1416  RooArgSet* ret = _dstore->addColumns(varList) ;
1417  _vars.addOwned(*ret) ;
1419  return ret ;
1420 }
1421 
1422 
1423 
1424 ////////////////////////////////////////////////////////////////////////////////
1425 /// Create a TH2F histogram of the distribution of the specified variable
1426 /// using this dataset. Apply any cuts to select which events are used.
1427 /// The variable being plotted can either be contained directly in this
1428 /// dataset, or else be a function of the variables in this dataset.
1429 /// The histogram will be created using RooAbsReal::createHistogram() with
1430 /// the name provided (with our dataset name prepended).
1431 
1432 TH2F* RooDataSet::createHistogram(const RooAbsRealLValue& var1, const RooAbsRealLValue& var2, const char* cuts, const char *name) const
1433 {
1434  checkInit() ;
1435  return createHistogram(var1, var2, var1.getBins(), var2.getBins(), cuts, name);
1436 }
1437 
1438 
1439 
1440 ////////////////////////////////////////////////////////////////////////////////
1441 /// Create a TH2F histogram of the distribution of the specified variable
1442 /// using this dataset. Apply any cuts to select which events are used.
1443 /// The variable being plotted can either be contained directly in this
1444 /// dataset, or else be a function of the variables in this dataset.
1445 /// The histogram will be created using RooAbsReal::createHistogram() with
1446 /// the name provided (with our dataset name prepended).
1447 
1449  Int_t nx, Int_t ny, const char* cuts, const char *name) const
1450 {
1451  checkInit() ;
1452  static Int_t counter(0) ;
1453 
1454  Bool_t ownPlotVarX(kFALSE) ;
1455  // Is this variable in our dataset?
1456  RooAbsReal* plotVarX= (RooAbsReal*)_vars.find(var1.GetName());
1457  if(0 == plotVarX) {
1458  // Is this variable a client of our dataset?
1459  if (!var1.dependsOn(_vars)) {
1460  coutE(InputArguments) << GetName() << "::createHistogram: Argument " << var1.GetName()
1461  << " is not in dataset and is also not dependent on data set" << endl ;
1462  return 0 ;
1463  }
1464 
1465  // Clone derived variable
1466  plotVarX = (RooAbsReal*) var1.Clone() ;
1467  ownPlotVarX = kTRUE ;
1468 
1469  //Redirect servers of derived clone to internal ArgSet representing the data in this set
1470  plotVarX->redirectServers(const_cast<RooArgSet&>(_vars)) ;
1471  }
1472 
1473  Bool_t ownPlotVarY(kFALSE) ;
1474  // Is this variable in our dataset?
1475  RooAbsReal* plotVarY= (RooAbsReal*)_vars.find(var2.GetName());
1476  if(0 == plotVarY) {
1477  // Is this variable a client of our dataset?
1478  if (!var2.dependsOn(_vars)) {
1479  coutE(InputArguments) << GetName() << "::createHistogram: Argument " << var2.GetName()
1480  << " is not in dataset and is also not dependent on data set" << endl ;
1481  return 0 ;
1482  }
1483 
1484  // Clone derived variable
1485  plotVarY = (RooAbsReal*) var2.Clone() ;
1486  ownPlotVarY = kTRUE ;
1487 
1488  //Redirect servers of derived clone to internal ArgSet representing the data in this set
1489  plotVarY->redirectServers(const_cast<RooArgSet&>(_vars)) ;
1490  }
1491 
1492  // Create selection formula if selection cuts are specified
1493  RooFormula* select = 0;
1494  if(0 != cuts && strlen(cuts)) {
1495  select=new RooFormula(cuts,cuts,_vars);
1496  if (!select || !select->ok()) {
1497  delete select;
1498  return 0 ;
1499  }
1500  }
1501 
1502  TString histName(name);
1503  histName.Prepend("_");
1504  histName.Prepend(fName);
1505  histName.Append("_") ;
1506  histName.Append(Form("%08x",counter++)) ;
1507 
1508  // create the histogram
1509  TH2F* histogram=new TH2F(histName.Data(), "Events", nx, var1.getMin(), var1.getMax(),
1510  ny, var2.getMin(), var2.getMax());
1511  if(!histogram) {
1512  coutE(DataHandling) << fName << "::createHistogram: unable to create a new histogram" << endl;
1513  return 0;
1514  }
1515 
1516  // Dump contents
1517  Int_t nevent= numEntries() ;
1518  for(Int_t i=0; i < nevent; ++i)
1519  {
1520  get(i);
1521 
1522  if (select && select->eval()==0) continue ;
1523  histogram->Fill(plotVarX->getVal(), plotVarY->getVal(),weight()) ;
1524  }
1525 
1526  if (ownPlotVarX) delete plotVarX ;
1527  if (ownPlotVarY) delete plotVarY ;
1528  if (select) delete select ;
1529 
1530  return histogram ;
1531 }
1532 
1533 
1534 
1535 
1536 
1537 ////////////////////////////////////////////////////////////////////////////////
1538 /// Special plot method for 'X-Y' datasets used in \f$ \chi^2 \f$ fitting.
1539 /// For general plotting, see RooAbsData::plotOn().
1540 ///
1541 /// These datasets
1542 /// have one observable (X) and have weights (Y) and associated errors.
1543 /// <table>
1544 /// <tr><th> Contents options <th> Effect
1545 /// <tr><td> YVar(RooRealVar& var) <td> Designate specified observable as 'y' variable
1546 /// If not specified, the event weight will be the y variable
1547 /// <tr><th> Histogram drawing options <th> Effect
1548 /// <tr><td> DrawOption(const char* opt) <td> Select ROOT draw option for resulting TGraph object
1549 /// <tr><td> LineStyle(Int_t style) <td> Select line style by ROOT line style code, default is solid
1550 /// <tr><td> LineColor(Int_t color) <td> Select line color by ROOT color code, default is black
1551 /// <tr><td> LineWidth(Int_t width) <td> Select line with in pixels, default is 3
1552 /// <tr><td> MarkerStyle(Int_t style) <td> Select the ROOT marker style, default is 21
1553 /// <tr><td> MarkerColor(Int_t color) <td> Select the ROOT marker color, default is black
1554 /// <tr><td> MarkerSize(Double_t size) <td> Select the ROOT marker size
1555 /// <tr><td> Rescale(Double_t factor) <td> Apply global rescaling factor to histogram
1556 /// <tr><th> Misc. other options <th> Effect
1557 /// <tr><td> Name(const chat* name) <td> Give curve specified name in frame. Useful if curve is to be referenced later
1558 /// <tr><td> Invisible(Bool_t flag) <td> Add curve to frame, but do not display. Useful in combination AddTo()
1559 /// </table>
1560 
1561 RooPlot* RooDataSet::plotOnXY(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
1562  const RooCmdArg& arg3, const RooCmdArg& arg4,
1563  const RooCmdArg& arg5, const RooCmdArg& arg6,
1564  const RooCmdArg& arg7, const RooCmdArg& arg8) const
1565 {
1566  checkInit() ;
1567 
1568  RooLinkedList argList ;
1569  argList.Add((TObject*)&arg1) ; argList.Add((TObject*)&arg2) ;
1570  argList.Add((TObject*)&arg3) ; argList.Add((TObject*)&arg4) ;
1571  argList.Add((TObject*)&arg5) ; argList.Add((TObject*)&arg6) ;
1572  argList.Add((TObject*)&arg7) ; argList.Add((TObject*)&arg8) ;
1573 
1574  // Process named arguments
1575  RooCmdConfig pc(Form("RooDataSet::plotOnXY(%s)",GetName())) ;
1576  pc.defineString("drawOption","DrawOption",0,"P") ;
1577  pc.defineString("histName","Name",0,"") ;
1578  pc.defineInt("lineColor","LineColor",0,-999) ;
1579  pc.defineInt("lineStyle","LineStyle",0,-999) ;
1580  pc.defineInt("lineWidth","LineWidth",0,-999) ;
1581  pc.defineInt("markerColor","MarkerColor",0,-999) ;
1582  pc.defineInt("markerStyle","MarkerStyle",0,8) ;
1583  pc.defineDouble("markerSize","MarkerSize",0,-999) ;
1584  pc.defineInt("fillColor","FillColor",0,-999) ;
1585  pc.defineInt("fillStyle","FillStyle",0,-999) ;
1586  pc.defineInt("histInvisible","Invisible",0,0) ;
1587  pc.defineDouble("scaleFactor","Rescale",0,1.) ;
1588  pc.defineObject("xvar","XVar",0,0) ;
1589  pc.defineObject("yvar","YVar",0,0) ;
1590 
1591 
1592  // Process & check varargs
1593  pc.process(argList) ;
1594  if (!pc.ok(kTRUE)) {
1595  return frame ;
1596  }
1597 
1598  // Extract values from named arguments
1599  const char* drawOptions = pc.getString("drawOption") ;
1600  Int_t histInvisible = pc.getInt("histInvisible") ;
1601  const char* histName = pc.getString("histName",0,kTRUE) ;
1602  Double_t scaleFactor = pc.getDouble("scaleFactor") ;
1603 
1604  RooRealVar* xvar = (RooRealVar*) _vars.find(frame->getPlotVar()->GetName()) ;
1605 
1606  // Determine Y variable (default is weight, if present)
1607  RooRealVar* yvar = (RooRealVar*)(pc.getObject("yvar")) ;
1608 
1609  // Sanity check. XY plotting only applies to weighted datasets if no YVar is specified
1610  if (!_wgtVar && !yvar) {
1611  coutE(InputArguments) << "RooDataSet::plotOnXY(" << GetName() << ") ERROR: no YVar() argument specified and dataset is not weighted" << endl ;
1612  return 0 ;
1613  }
1614 
1615  RooRealVar* dataY = yvar ? (RooRealVar*) _vars.find(yvar->GetName()) : 0 ;
1616  if (yvar && !dataY) {
1617  coutE(InputArguments) << "RooDataSet::plotOnXY(" << GetName() << ") ERROR on YVar() argument, dataset does not contain a variable named " << yvar->GetName() << endl ;
1618  return 0 ;
1619  }
1620 
1621 
1622  // Make RooHist representing XY contents of data
1623  RooHist* graph = new RooHist ;
1624  if (histName) {
1625  graph->SetName(histName) ;
1626  } else {
1627  graph->SetName(Form("hxy_%s",GetName())) ;
1628  }
1629 
1630  for (int i=0 ; i<numEntries() ; i++) {
1631  get(i) ;
1632  Double_t x = xvar->getVal() ;
1633  Double_t exlo = xvar->getErrorLo() ;
1634  Double_t exhi = xvar->getErrorHi() ;
1635  Double_t y,eylo,eyhi ;
1636  if (!dataY) {
1637  y = weight() ;
1638  weightError(eylo,eyhi) ;
1639  } else {
1640  y = dataY->getVal() ;
1641  eylo = dataY->getErrorLo() ;
1642  eyhi = dataY->getErrorHi() ;
1643  }
1644  graph->addBinWithXYError(x,y,-1*exlo,exhi,-1*eylo,eyhi,scaleFactor) ;
1645  }
1646 
1647  // Adjust style options according to named arguments
1648  Int_t lineColor = pc.getInt("lineColor") ;
1649  Int_t lineStyle = pc.getInt("lineStyle") ;
1650  Int_t lineWidth = pc.getInt("lineWidth") ;
1651  Int_t markerColor = pc.getInt("markerColor") ;
1652  Int_t markerStyle = pc.getInt("markerStyle") ;
1653  Size_t markerSize = pc.getDouble("markerSize") ;
1654  Int_t fillColor = pc.getInt("fillColor") ;
1655  Int_t fillStyle = pc.getInt("fillStyle") ;
1656 
1657  if (lineColor!=-999) graph->SetLineColor(lineColor) ;
1658  if (lineStyle!=-999) graph->SetLineStyle(lineStyle) ;
1659  if (lineWidth!=-999) graph->SetLineWidth(lineWidth) ;
1660  if (markerColor!=-999) graph->SetMarkerColor(markerColor) ;
1661  if (markerStyle!=-999) graph->SetMarkerStyle(markerStyle) ;
1662  if (markerSize!=-999) graph->SetMarkerSize(markerSize) ;
1663  if (fillColor!=-999) graph->SetFillColor(fillColor) ;
1664  if (fillStyle!=-999) graph->SetFillStyle(fillStyle) ;
1665 
1666  // Add graph to frame
1667  frame->addPlotable(graph,drawOptions,histInvisible) ;
1668 
1669  return frame ;
1670 }
1671 
1672 
1673 
1674 
1675 ////////////////////////////////////////////////////////////////////////////////
1676 /// Read given list of ascii files, and construct a data set, using the given
1677 /// ArgList as structure definition.
1678 /// \param fileList Multiple file names, comma separated. Each
1679 /// file is optionally prefixed with 'commonPath' if such a path is
1680 /// provided
1681 ///
1682 /// \param varList Specify the dimensions of the dataset to be built.
1683 /// This list describes the order in which these dimensions appear in the
1684 /// ascii files to be read.
1685 /// Each line in the ascii file should contain N white-space separated
1686 /// tokens, with N the number of args in `varList`. Any text beyond
1687 /// N tokens will be ignored with a warning message.
1688 /// (NB: This is the default output of RooArgList::writeToStream())
1689 ///
1690 /// \param verbOpt `Q` be quiet, `D` debug mode (verbose)
1691 ///
1692 /// \param commonPath All filenames in `fileList` will be prefixed with this optional path.
1693 ///
1694 /// \param indexCatName Interpret the data as belonging to category `indexCatName`.
1695 /// When multiple files are read, a RooCategory arg in `varList` can
1696 /// optionally be designated to hold information about the source file
1697 /// of each data point. This feature is enabled by giving the name
1698 /// of the (already existing) category variable in `indexCatName`.
1699 ///
1700 /// \attention If the value of any of the variables on a given line exceeds the
1701 /// fit range associated with that dimension, the entire line will be
1702 /// ignored. A warning message is printed in each case, unless the
1703 /// `Q` verbose option is given. The number of events read and skipped
1704 /// is always summarized at the end.
1705 ///
1706 /// If no further information is given a label name 'fileNNN' will
1707 /// be assigned to each event, where NNN is the sequential number of
1708 /// the source file in `fileList`.
1709 ///
1710 /// Alternatively, it is possible to override the default label names
1711 /// of the index category by specifying them in the fileList string:
1712 /// When instead of `file1.txt,file2.txt` the string
1713 /// `file1.txt:FOO,file2.txt:BAR` is specified, a state named "FOO"
1714 /// is assigned to the index category for each event originating from
1715 /// file1.txt. The labels FOO,BAR may be predefined in the index
1716 /// category via defineType(), but don't have to be.
1717 ///
1718 /// Finally, one can also assign the same label to multiple files,
1719 /// either by specifying `file1.txt:FOO,file2,txt:FOO,file3.txt:BAR`
1720 /// or `file1.txt,file2.txt:FOO,file3.txt:BAR`.
1721 ///
1722 
1723 RooDataSet *RooDataSet::read(const char *fileList, const RooArgList &varList,
1724  const char *verbOpt, const char* commonPath,
1725  const char* indexCatName) {
1726  // Make working copy of variables list
1727  RooArgList variables(varList) ;
1728 
1729  // Append blinding state category to variable list if not already there
1730  Bool_t ownIsBlind(kTRUE) ;
1731  RooAbsArg* blindState = variables.find("blindState") ;
1732  if (!blindState) {
1733  blindState = new RooCategory("blindState","Blinding State") ;
1734  variables.add(*blindState) ;
1735  } else {
1736  ownIsBlind = kFALSE ;
1737  if (blindState->IsA()!=RooCategory::Class()) {
1738  oocoutE((TObject*)0,DataHandling) << "RooDataSet::read: ERROR: variable list already contains"
1739  << "a non-RooCategory blindState member" << endl ;
1740  return 0 ;
1741  }
1742  oocoutW((TObject*)0,DataHandling) << "RooDataSet::read: WARNING: recycling existing "
1743  << "blindState category in variable list" << endl ;
1744  }
1745  RooCategory* blindCat = (RooCategory*) blindState ;
1746 
1747  // Configure blinding state category
1748  blindCat->setAttribute("Dynamic") ;
1749  blindCat->defineType("Normal",0) ;
1750  blindCat->defineType("Blind",1) ;
1751 
1752  // parse the option string
1753  TString opts= verbOpt;
1754  opts.ToLower();
1755  Bool_t verbose= !opts.Contains("q");
1756  Bool_t debug= opts.Contains("d");
1757 
1758  auto data = std::make_unique<RooDataSet>("dataset", fileList, variables);
1759  if (ownIsBlind) { variables.remove(*blindState) ; delete blindState ; }
1760  if(!data) {
1761  oocoutE((TObject*)0,DataHandling) << "RooDataSet::read: unable to create a new dataset"
1762  << endl;
1763  return nullptr;
1764  }
1765 
1766  // Redirect blindCat to point to the copy stored in the data set
1767  blindCat = (RooCategory*) data->_vars.find("blindState") ;
1768 
1769  // Find index category, if requested
1770  RooCategory *indexCat = 0;
1771  //RooCategory *indexCatOrig = 0;
1772  if (indexCatName) {
1773  RooAbsArg* tmp = 0;
1774  tmp = data->_vars.find(indexCatName) ;
1775  if (!tmp) {
1776  oocoutE(data.get(),DataHandling) << "RooDataSet::read: no index category named "
1777  << indexCatName << " in supplied variable list" << endl ;
1778  return nullptr;
1779  }
1780  if (tmp->IsA()!=RooCategory::Class()) {
1781  oocoutE(data.get(),DataHandling) << "RooDataSet::read: variable " << indexCatName
1782  << " is not a RooCategory" << endl ;
1783  return nullptr;
1784  }
1785  indexCat = static_cast<RooCategory*>(tmp);
1786 
1787  // Prevent RooArgSet from attempting to read in indexCat
1788  indexCat->setAttribute("Dynamic") ;
1789  }
1790 
1791 
1792  Int_t outOfRange(0) ;
1793 
1794  // Loop over all names in comma separated list
1795  Int_t fileSeqNum(0);
1796  for (const auto& filename : RooHelpers::tokenise(std::string(fileList), ", ")) {
1797  // Determine index category number, if this option is active
1798  if (indexCat) {
1799 
1800  // Find and detach optional file category name
1801  const char *catname = strchr(filename.c_str(),':');
1802 
1803  if (catname) {
1804  // Use user category name if provided
1805  catname++ ;
1806 
1807  if (indexCat->hasLabel(catname)) {
1808  // Use existing category index
1809  indexCat->setLabel(catname);
1810  } else {
1811  // Register cat name
1812  indexCat->defineType(catname,fileSeqNum) ;
1813  indexCat->setIndex(fileSeqNum) ;
1814  }
1815  } else {
1816  // Assign autogenerated name
1817  char newLabel[128] ;
1818  snprintf(newLabel,128,"file%03d",fileSeqNum) ;
1819  if (indexCat->defineType(newLabel,fileSeqNum)) {
1820  oocoutE(data.get(), DataHandling) << "RooDataSet::read: Error, cannot register automatic type name " << newLabel
1821  << " in index category " << indexCat->GetName() << endl ;
1822  return 0 ;
1823  }
1824  // Assign new category number
1825  indexCat->setIndex(fileSeqNum) ;
1826  }
1827  }
1828 
1829  oocoutI(data.get(), DataHandling) << "RooDataSet::read: reading file " << filename << endl ;
1830 
1831  // Prefix common path
1832  TString fullName(commonPath) ;
1833  fullName.Append(filename) ;
1834  ifstream file(fullName) ;
1835 
1836  if (!file.good()) {
1837  oocoutE(data.get(), DataHandling) << "RooDataSet::read: unable to open '"
1838  << filename << "'. Returning nullptr now." << endl;
1839  return nullptr;
1840  }
1841 
1842  // Double_t value;
1843  Int_t line(0) ;
1844  Bool_t haveBlindString(false) ;
1845 
1846  while(file.good() && !file.eof()) {
1847  line++;
1848  if(debug) oocxcoutD(data.get(),DataHandling) << "reading line " << line << endl;
1849 
1850  // process comment lines
1851  if (file.peek() == '#') {
1852  if(debug) oocxcoutD(data.get(),DataHandling) << "skipping comment on line " << line << endl;
1853  } else {
1854  // Read single line
1855  Bool_t readError = variables.readFromStream(file,kTRUE,verbose) ;
1856  data->_vars = variables ;
1857 
1858  // Stop on read error
1859  if(!file.good()) {
1860  oocoutE(data.get(), DataHandling) << "RooDataSet::read(static): read error at line " << line << endl ;
1861  break;
1862  }
1863 
1864  if (readError) {
1865  outOfRange++ ;
1866  } else {
1867  blindCat->setIndex(haveBlindString) ;
1868  data->fill(); // store this event
1869  }
1870  }
1871 
1872  // Skip all white space (including empty lines).
1873  while (isspace(file.peek())) {
1874  char dummy;
1875  file >> std::noskipws >> dummy >> std::skipws;
1876  }
1877  }
1878 
1879  file.close();
1880 
1881  // get next file name
1882  fileSeqNum++ ;
1883  }
1884 
1885  if (indexCat) {
1886  // Copy dynamically defined types from new data set to indexCat in original list
1887  assert(dynamic_cast<RooCategory*>(variables.find(indexCatName)));
1888  const auto origIndexCat = static_cast<RooCategory*>(variables.find(indexCatName));
1889  for (const auto& nameIdx : *indexCat) {
1890  origIndexCat->defineType(nameIdx.first, nameIdx.second);
1891  }
1892  }
1893  oocoutI(data.get(),DataHandling) << "RooDataSet::read: read " << data->numEntries()
1894  << " events (ignored " << outOfRange << " out of range events)" << endl;
1895 
1896  return data.release();
1897 }
1898 
1899 
1900 
1901 
1902 ////////////////////////////////////////////////////////////////////////////////
1903 /// Write the contents of this dataset to an ASCII file with the specified name.
1904 /// Each event will be written as a single line containing the written values
1905 /// of each observable in the order they were declared in the dataset and
1906 /// separated by whitespaces
1907 
1908 Bool_t RooDataSet::write(const char* filename) const
1909 {
1910  // Open file for writing
1911  ofstream ofs(filename) ;
1912  if (ofs.fail()) {
1913  coutE(DataHandling) << "RooDataSet::write(" << GetName() << ") cannot create file " << filename << endl ;
1914  return kTRUE ;
1915  }
1916 
1917  // Write all lines as arglist in compact mode
1918  coutI(DataHandling) << "RooDataSet::write(" << GetName() << ") writing ASCII file " << filename << endl ;
1919  return write(ofs);
1920 }
1921 
1922 ////////////////////////////////////////////////////////////////////////////////
1923 /// Write the contents of this dataset to the stream.
1924 /// Each event will be written as a single line containing the written values
1925 /// of each observable in the order they were declared in the dataset and
1926 /// separated by whitespaces
1927 
1928 Bool_t RooDataSet::write(ostream & ofs) const {
1929  checkInit();
1930 
1931  for (Int_t i=0; i<numEntries(); ++i) {
1932  get(i)->writeToStream(ofs,kTRUE);
1933  }
1934 
1935  if (ofs.fail()) {
1936  coutW(DataHandling) << "RooDataSet::write(" << GetName() << "): WARNING error(s) have occured in writing" << endl ;
1937  }
1938 
1939  return ofs.fail() ;
1940 }
1941 
1942 
1943 ////////////////////////////////////////////////////////////////////////////////
1944 /// Print info about this dataset to the specified output stream.
1945 ///
1946 /// Standard: number of entries
1947 /// Shape: list of variables we define & were generated with
1948 
1949 void RooDataSet::printMultiline(ostream& os, Int_t contents, Bool_t verbose, TString indent) const
1950 {
1951  checkInit() ;
1953  if (_wgtVar) {
1954  os << indent << " Dataset variable \"" << _wgtVar->GetName() << "\" is interpreted as the event weight" << endl ;
1955  }
1956 }
1957 
1958 
1959 ////////////////////////////////////////////////////////////////////////////////
1960 /// Print value of the dataset, i.e. the sum of weights contained in the dataset
1961 
1962 void RooDataSet::printValue(ostream& os) const
1963 {
1964  os << numEntries() << " entries" ;
1965  if (isWeighted()) {
1966  os << " (" << sumEntries() << " weighted)" ;
1967  }
1968 }
1969 
1970 
1971 
1972 ////////////////////////////////////////////////////////////////////////////////
1973 /// Print argument of dataset, i.e. the observable names
1974 
1975 void RooDataSet::printArgs(ostream& os) const
1976 {
1977  os << "[" ;
1978  TIterator* iter = _varsNoWgt.createIterator() ;
1979  RooAbsArg* arg ;
1980  Bool_t first(kTRUE) ;
1981  while((arg=(RooAbsArg*)iter->Next())) {
1982  if (first) {
1983  first=kFALSE ;
1984  } else {
1985  os << "," ;
1986  }
1987  os << arg->GetName() ;
1988  }
1989  if (_wgtVar) {
1990  os << ",weight:" << _wgtVar->GetName() ;
1991  }
1992  os << "]" ;
1993  delete iter ;
1994 }
1995 
1996 
1997 
1998 ////////////////////////////////////////////////////////////////////////////////
1999 /// Change the name of this dataset into the given name
2000 
2001 void RooDataSet::SetName(const char *name)
2002 {
2003  if (_dir) _dir->GetList()->Remove(this);
2005  if (_dir) _dir->GetList()->Add(this);
2006 }
2007 
2008 
2009 ////////////////////////////////////////////////////////////////////////////////
2010 /// Change the title of this dataset into the given name
2011 
2012 void RooDataSet::SetNameTitle(const char *name, const char* title)
2013 {
2014  if (_dir) _dir->GetList()->Remove(this);
2015  TNamed::SetNameTitle(name,title) ;
2016  if (_dir) _dir->GetList()->Add(this);
2017 }
2018 
2019 
2020 ////////////////////////////////////////////////////////////////////////////////
2021 /// Stream an object of class RooDataSet.
2022 
2023 void RooDataSet::Streamer(TBuffer &R__b)
2024 {
2025  if (R__b.IsReading()) {
2026 
2027  UInt_t R__s, R__c;
2028  Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
2029 
2030  if (R__v>1) {
2031 
2032  // Use new-style streaming for version >1
2033  R__b.ReadClassBuffer(RooDataSet::Class(),this,R__v,R__s,R__c);
2034 
2035  } else {
2036 
2037  // Legacy dataset conversion happens here. Legacy RooDataSet inherits from RooTreeData
2038  // which in turn inherits from RooAbsData. Manually stream RooTreeData contents on
2039  // file here and convert it into a RooTreeDataStore which is installed in the
2040  // new-style RooAbsData base class
2041 
2042  // --- This is the contents of the streamer code of RooTreeData version 1 ---
2043  UInt_t R__s1, R__c1;
2044  Version_t R__v1 = R__b.ReadVersion(&R__s1, &R__c1); if (R__v1) { }
2045 
2046  RooAbsData::Streamer(R__b);
2047  TTree* X_tree(0) ; R__b >> X_tree;
2048  RooArgSet X_truth ; X_truth.Streamer(R__b);
2049  TString X_blindString ; X_blindString.Streamer(R__b);
2050  R__b.CheckByteCount(R__s1, R__c1, TClass::GetClass("RooTreeData"));
2051  // --- End of RooTreeData-v1 streamer
2052 
2053  // Construct RooTreeDataStore from X_tree and complete initialization of new-style RooAbsData
2054  _dstore = new RooTreeDataStore(X_tree,_vars) ;
2055  _dstore->SetName(GetName()) ;
2056  _dstore->SetTitle(GetTitle()) ;
2057  _dstore->checkInit() ;
2058 
2059  // This is the contents of the streamer code of RooDataSet version 1
2060  RooDirItem::Streamer(R__b);
2061  _varsNoWgt.Streamer(R__b);
2062  R__b >> _wgtVar;
2063  R__b.CheckByteCount(R__s, R__c, RooDataSet::IsA());
2064 
2065 
2066  }
2067  } else {
2068  R__b.WriteClassBuffer(RooDataSet::Class(),this);
2069  }
2070 }
2071 
2072 
2073 
2074 ////////////////////////////////////////////////////////////////////////////////
2075 /// Convert vector-based storage to tree-based storage. This implementation overrides the base class
2076 /// implementation because the latter doesn't transfer weights.
2078 {
2079  if (storageType != RooAbsData::Tree) {
2080  RooTreeDataStore *newStore = new RooTreeDataStore(GetName(), GetTitle(), _vars, *_dstore, nullptr, _wgtVar ? _wgtVar->GetName() : nullptr);
2081  delete _dstore;
2082  _dstore = newStore;
2084  }
2085 }
2086 
RooFormulaVar.h
RooAbsArg::Clone
virtual TObject * Clone(const char *newname=0) const
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:84
RooDirItem::removeFromDir
void removeFromDir(TObject *obj)
Remove object from directory it was added to.
Definition: RooDirItem.cxx:43
l
auto * l
Definition: textangle.C:4
RooAbsData::Tree
@ Tree
Definition: RooAbsData.h:248
RooFormula
RooFormula internally uses ROOT's TFormula to compute user-defined expressions of RooAbsArgs.
Definition: RooFormula.h:34
RooDataSet::getBatches
void getBatches(RooBatchCompute::RunContext &evalData, std::size_t first=0, std::size_t len=std::numeric_limits< std::size_t >::max()) const override
Write information to retrieve data columns into evalData.spans.
Definition: RooDataSet.cxx:1023
RooAbsDataStore::weightError
virtual Double_t weightError(RooAbsData::ErrorType etype=RooAbsData::Poisson) const =0
RooLinkedList::MakeIterator
TIterator * MakeIterator(Bool_t forward=kTRUE) const
Create a TIterator for this list.
Definition: RooLinkedList.cxx:747
RooCmdArg
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:27
RooArgSet::addOwned
Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE) override
Add element to an owning set.
Definition: RooArgSet.cxx:274
first
Definition: first.py:1
RooHelpers.h
RooAbsDataStore::cachedVars
const RooArgSet & cachedVars() const
Definition: RooAbsDataStore.h:109
RooRealVar::setVal
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:257
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooCmdConfig.h
e
#define e(i)
Definition: RSha256.hxx:103
RooDataSet::cacheClone
virtual RooAbsData * cacheClone(const RooAbsArg *newCacheOwner, const RooArgSet *newCacheVars, const char *newName=0) override
Return a clone of this dataset containing only the cached variables.
Definition: RooDataSet.cxx:859
TDirectory::GetList
virtual TList * GetList() const
Definition: TDirectory.h:176
RooAbsReal.h
RooDataSet::convertToTreeStore
void convertToTreeStore() override
Convert vector-based storage to tree-based storage.
Definition: RooDataSet.cxx:2077
RooAbsDataStore::isWeighted
virtual Bool_t isWeighted() const =0
Version_t
short Version_t
Definition: RtypesCore.h:65
snprintf
#define snprintf
Definition: civetweb.c:1540
RooArgSet::addClone
RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add clone of specified element to an owning set.
Definition: RooArgSet.cxx:288
RooMsgService.h
f
#define f(i)
Definition: RSha256.hxx:104
TNamed::SetNameTitle
virtual void SetNameTitle(const char *name, const char *title)
Set all the TNamed parameters (name and title).
Definition: TNamed.cxx:154
TNamed::SetName
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
RooDirItem::_dir
TDirectory * _dir
Definition: RooDirItem.h:33
RooAbsData
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:49
RooAbsRealLValue::getMax
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Definition: RooAbsRealLValue.h:89
RooDataSet::MemPool
MemPoolForRooSets< RooDataSet, 5 *150 > MemPool
Definition: RooDataSet.h:168
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:68
oocoutI
#define oocoutI(o, a)
Definition: RooMsgService.h:45
TString::Prepend
TString & Prepend(const char *cs)
Definition: TString.h:661
RooAbsData::printMultiline
void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Interface for detailed printing of object.
Definition: RooAbsData.cxx:802
TH2F
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:251
TString::Data
const char * Data() const
Definition: TString.h:369
RooDataSet::initialize
void initialize(const char *wgtVarName)
Initialize the dataset.
Definition: RooDataSet.cxx:906
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
TNamed::GetTitle
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
TBuffer::ReadClassBuffer
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
RooAbsData::_dstore
RooAbsDataStore * _dstore
External variables cached with this data set.
Definition: RooAbsData.h:294
RooAbsData::fill
virtual void fill()
Definition: RooAbsData.cxx:300
coutE
#define coutE(a)
Definition: RooMsgService.h:33
coutW
#define coutW(a)
Definition: RooMsgService.h:32
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
TTree
A TTree represents a columnar dataset.
Definition: TTree.h:79
RooDirItem
RooDirItem is a utility base class for RooFit objects that are to be attached to ROOT directories.
Definition: RooDirItem.h:22
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooDataSet::printValue
virtual void printValue(std::ostream &os) const override
Print value of the dataset, i.e. the sum of weights contained in the dataset.
Definition: RooDataSet.cxx:1962
RooDataSet::printArgs
virtual void printArgs(std::ostream &os) const override
Print argument of dataset, i.e. the observable names.
Definition: RooDataSet.cxx:1975
TFile::Open
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=ROOT::RCompressionSetting::EDefaults::kUseCompiledDefault, Int_t netopt=0)
Create / open a file.
Definition: TFile.cxx:3997
Int_t
int Int_t
Definition: RtypesCore.h:45
RooRealVar::removeAsymError
void removeAsymError()
Definition: RooRealVar.h:67
RooAbsData::defaultStorageType
static StorageType defaultStorageType
Definition: RooAbsData.h:256
RooAbsCollection::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Definition: RooAbsCollection.cxx:584
TNamed::fName
TString fName
Definition: TNamed.h:32
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
TGeant4Unit::pc
static constexpr double pc
Definition: TGeant4SystemOfUnits.h:130
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
RooCategory::setIndex
virtual Bool_t setIndex(Int_t index, bool printError=true) override
Set value by specifying the index code of the desired state.
Definition: RooCategory.cxx:163
x
Double_t x[n]
Definition: legend1.C:17
RooAbsData::store
RooAbsDataStore * store()
Definition: RooAbsData.h:68
coutI
#define coutI(a)
Definition: RooMsgService.h:30
indent
static void indent(ostringstream &buf, int indent_level)
Definition: TClingCallFunc.cxx:87
RooDataSet::merge
Bool_t merge(RooDataSet *data1, RooDataSet *data2=0, RooDataSet *data3=0, RooDataSet *data4=0, RooDataSet *data5=0, RooDataSet *data6=0)
Definition: RooDataSet.cxx:1321
RooDataSet::append
void append(RooDataSet &data)
Add all data points of given data set to this data set.
Definition: RooDataSet.cxx:1381
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooDataSet::SetName
void SetName(const char *name) override
Change the name of this dataset into the given name.
Definition: RooDataSet.cxx:2001
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
RooAbsData::checkInit
void checkInit() const
Definition: RooAbsData.cxx:2331
RooDataSet::addFast
virtual void addFast(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0)
Add a data point, with its coordinates specified in the 'data' argset, to the data set.
Definition: RooDataSet.cxx:1279
oocoutE
#define oocoutE(o, a)
Definition: RooMsgService.h:48
TTree.h
TBuffer::CheckByteCount
virtual Int_t CheckByteCount(UInt_t startpos, UInt_t bcnt, const TClass *clss)=0
RooDataSet::_wgtVar
RooRealVar * _wgtVar
Definition: RooDataSet.h:164
TString
Basic string class.
Definition: TString.h:136
RooDataSet::cleanup
static void cleanup()
Definition: RooDataSet.cxx:110
RooDataSet::addWgtVar
RooArgSet addWgtVar(const RooArgSet &origVars, const RooAbsArg *wgtVar)
Helper function for constructor that adds optional weight variable to construct total set of observab...
Definition: RooDataSet.cxx:847
RooAbsData::ErrorType
ErrorType
Definition: RooAbsData.h:99
RooAbsCategory::getCurrentLabel
virtual const char * getCurrentLabel() const
Return label string of current state.
Definition: RooAbsCategory.cxx:130
RooDataSet.h
RooTreeDataStore.h
TFile.h
RooCategory::defineType
bool defineType(const std::string &label)
Define a state with given name.
Definition: RooCategory.cxx:208
bool
TIterator
Iterator abstract base class.
Definition: TIterator.h:30
RooAbsRealLValue::getMin
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
Definition: RooAbsRealLValue.h:86
RooCompositeDataStore
RooCompositeDataStore combines several disjunct datasets into one.
Definition: RooCompositeDataStore.h:34
RooCmdConfig
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
Definition: RooCmdConfig.h:27
RooTrace.h
hi
float type_of_call hi(const int &, const int &)
RooDataSet::_errorMsgCount
unsigned short _errorMsgCount
Definition: RooDataSet.h:171
RooDataSet::add
virtual void add(const RooArgSet &row, Double_t weight=1.0, Double_t weightError=0) override
Add a data point, with its coordinates specified in the 'data' argset, to the data set.
Definition: RooDataSet.cxx:1177
RooDataSet::addColumns
virtual RooArgSet * addColumns(const RooArgList &varList)
Add a column with the values of the given list of (function) argument to this dataset.
Definition: RooDataSet.cxx:1413
RooFit::DataHandling
@ DataHandling
Definition: RooGlobalFunc.h:69
RooAbsData::_vars
RooArgSet _vars
Definition: RooAbsData.h:291
ROOT::Math::fabs
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
Definition: UnaryOperators.h:131
RooDataHist
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:39
oocxcoutD
#define oocxcoutD(o, a)
Definition: RooMsgService.h:83
RooPlot::addPlotable
void addPlotable(RooPlotable *plotable, Option_t *drawOptions="", Bool_t invisible=kFALSE, Bool_t refreshNorm=kFALSE)
Add the specified plotable object to our plot.
Definition: RooPlot.cxx:570
RooCategory::setLabel
virtual Bool_t setLabel(const char *label, bool printError=true) override
Set value by specifying the name of the desired state.
Definition: RooCategory.cxx:185
RooFormulaVar
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooAbsDataStore::dirtyProp
Bool_t dirtyProp() const
Definition: RooAbsDataStore.h:114
TBuffer.h
RooAbsDataStore::sumEntries
virtual Double_t sumEntries() const
Definition: RooAbsDataStore.h:79
TRACE_DESTROY
#define TRACE_DESTROY
Definition: RooTrace.h:24
RooFormula::ok
Bool_t ok() const
Definition: RooFormula.h:59
RooRealVar::getErrorHi
Double_t getErrorHi() const
Definition: RooRealVar.h:70
RooDataSet::get
virtual const RooArgSet * get() const override
Return a RooArgSet with the coordinates of the current event.
Definition: RooDataSet.cxx:1159
RooRealVar::removeError
void removeError()
Definition: RooRealVar.h:63
RooDataSet::~RooDataSet
virtual ~RooDataSet()
Destructor.
Definition: RooDataSet.cxx:959
RooDataSet::isNonPoissonWeighted
virtual Bool_t isNonPoissonWeighted() const override
Returns true if histogram contains bins with entries with a non-integer weight.
Definition: RooDataSet.cxx:1137
RooAbsData::Vector
@ Vector
Definition: RooAbsData.h:248
TMVA::variables
void variables(TString dataset, TString fin="TMVA.root", TString dirName="InputVariables_Id", TString title="TMVA Input Variables", Bool_t isRegression=kFALSE, Bool_t useTMVAStyle=kTRUE)
TBuffer::WriteClassBuffer
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
RooRealVar::setError
void setError(Double_t value)
Definition: RooRealVar.h:62
RooDataHist.h
RooDataSet::SetNameTitle
void SetNameTitle(const char *name, const char *title) override
Change the title of this dataset into the given name.
Definition: RooDataSet.cxx:2012
RooAbsData::get
virtual const RooArgSet * get() const
Definition: RooAbsData.h:92
RooArgSet::add
Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add element to non-owning set.
Definition: RooArgSet.cxx:261
RooAbsData::storageType
StorageType storageType
Definition: RooAbsData.h:258
RooAbsCollection::selectCommon
RooAbsCollection * selectCommon(const RooAbsCollection &refColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
Definition: RooAbsCollection.cxx:700
RooLinkedList
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:35
RooPlot.h
RooAbsCollection::createIterator
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:119
RooAbsDataStore::addColumns
virtual RooArgSet * addColumns(const RooArgList &varList)=0
RooPlot
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:44
RooDataSet::printMultiline
void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const override
Print info about this dataset to the specified output stream.
Definition: RooDataSet.cxx:1949
TRACE_CREATE
#define TRACE_CREATE
Definition: RooTrace.h:23
RooAbsDataStore::loadValues
virtual void loadValues(const RooAbsDataStore *tds, const RooFormulaVar *select=0, const char *rangeName=0, std::size_t nStart=0, std::size_t nStop=std::numeric_limits< std::size_t >::max())=0
RooAbsData::numEntries
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Definition: RooAbsData.cxx:307
TClass::GetClass
static TClass * GetClass(const char *name, Bool_t load=kTRUE, Bool_t silent=kFALSE)
Static method returning pointer to TClass of the specified class name.
Definition: TClass.cxx:2946
RooCategory.h
RooLinkedList::Add
virtual void Add(TObject *arg)
Definition: RooLinkedList.h:62
y
Double_t y[n]
Definition: legend1.C:17
RooDataSet::_varsNoWgt
RooArgSet _varsNoWgt
Definition: RooDataSet.h:163
RooAbsDataStore::addColumn
virtual RooAbsArg * addColumn(RooAbsArg &var, Bool_t adjustRange=kTRUE)=0
RooDataSet::_doWeightErrorCheck
bool _doWeightErrorCheck
Counter to silence error messages when filling dataset.
Definition: RooDataSet.h:172
oocoutW
#define oocoutW(o, a)
Definition: RooMsgService.h:47
RooRealVar.h
TNamed::SetTitle
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
line
TLine * line
Definition: entrylistblock_figure1.C:235
RooAbsDataStore::merge
virtual RooAbsDataStore * merge(const RooArgSet &allvars, std::list< RooAbsDataStore * > dstoreList)=0
RooDirItem::appendToDir
void appendToDir(TObject *obj, Bool_t forceMemoryResident=kFALSE)
Append object to directory.
Definition: RooDirItem.cxx:55
RooHist
A RooHist is a graphical representation of binned data based on the TGraphAsymmErrors class.
Definition: RooHist.h:27
RooSentinel::activate
static void activate()
Install atexit handler that calls CleanupRooFitAtExit() on program termination.
Definition: RooSentinel.cxx:57
TH2.h
TBuffer::ReadVersion
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0
RooDataSet::sumEntries
virtual Double_t sumEntries() const override
Return effective number of entries in dataset, i.e., sum all weights.
Definition: RooDataSet.cxx:1060
RooDataSet::RooDataSet
RooDataSet()
Default constructor for persistence.
Definition: RooDataSet.cxx:169
TFile
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition: TFile.h:54
RooArgSet::writeToStream
virtual void writeToStream(std::ostream &os, Bool_t compact, const char *section=0) const
Write the contents of the argset in ASCII form to given stream.
Definition: RooArgSet.cxx:559
RooTreeDataStore
RooTreeDataStore is a TTree-backed data storage.
Definition: RooTreeDataStore.h:32
TIterator::Next
virtual TObject * Next()=0
TH2::Fill
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:294
RooAbsDataStore::weight
virtual Double_t weight() const =0
MemPoolForRooSets
Memory pool for RooArgSet and RooDataSet.
Definition: MemPoolForRooSets.h:34
unsigned int
RooHist.h
ccoutE
#define ccoutE(a)
Definition: RooMsgService.h:41
RooDataSet::isWeighted
virtual Bool_t isWeighted() const override
Return true if dataset contains weighted events.
Definition: RooDataSet.cxx:1127
RooAbsData::addOwnedComponent
void addOwnedComponent(const char *idxlabel, RooAbsData &data)
Definition: RooAbsData.cxx:2306
RooHelpers::tokenise
std::vector< std::string > tokenise(const std::string &str, const std::string &delims, bool returnEmptyToken=true)
Tokenise the string by splitting at the characters in delims.
Definition: RooHelpers.cxx:62
RooDataSet::read
static RooDataSet * read(const char *filename, const RooArgList &variables, const char *opts="", const char *commonPath="", const char *indexCatName=0)
Read given list of ascii files, and construct a data set, using the given ArgList as structure defini...
Definition: RooDataSet.cxx:1723
RooDataSet::weightSquared
virtual Double_t weightSquared() const override
Return squared event weight of current event.
Definition: RooDataSet.cxx:1002
TBuffer::IsReading
Bool_t IsReading() const
Definition: TBuffer.h:86
RooAbsDataStore::checkInit
virtual void checkInit() const
Definition: RooAbsDataStore.h:116
Double_t
double Double_t
Definition: RtypesCore.h:59
MemPoolForRooSets.h
RooAbsDataStore::append
virtual void append(RooAbsDataStore &other)=0
TList::Remove
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:822
Roo1DTable.h
RooCategory
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:27
RooDataSet::plotOnXY
virtual RooPlot * plotOnXY(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Special plot method for 'X-Y' datasets used in fitting.
Definition: RooDataSet.cxx:1561
file
Definition: file.py:1
RooDataSet::memPool
static MemPool * memPool()
graph
Definition: graph.py:1
RooVectorDataStore
RooVectorDataStore uses std::vectors to store data columns.
Definition: RooVectorDataStore.h:36
TList::Add
virtual void Add(TObject *obj)
Definition: TList.h:87
TObject
Mother of all ROOT objects.
Definition: TObject.h:37
RooDataSet::binnedClone
RooDataHist * binnedClone(const char *newName=0, const char *newTitle=0) const
Return binned clone of this dataset.
Definition: RooDataSet.cxx:970
RooAbsArg::setAttribute
void setAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:291
RooAbsCategory::hasLabel
bool hasLabel(const std::string &label) const
Check if a state with name label exists.
Definition: RooAbsCategory.h:64
RooAbsReal::createFundamental
RooAbsArg * createFundamental(const char *newname=0) const
Create a RooRealVar fundamental object with our properties.
Definition: RooAbsReal.cxx:3390
name
char name[80]
Definition: TGX11.cxx:110
RooAbsArg::redirectServers
Bool_t redirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t isRecursionStep=kFALSE)
Replace all direct servers of this object with the new servers in newServerList.
Definition: RooAbsArg.cxx:959
RooAbsArg::dependsOn
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:764
RooDataSet::emptyClone
virtual RooAbsData * emptyClone(const char *newName=0, const char *newTitle=0, const RooArgSet *vars=0, const char *wgtVarName=0) const override
Return an empty clone of this dataset.
Definition: RooDataSet.cxx:877
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
RooDataSet
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooTreeDataStore::loadValues
void loadValues(const TTree *t, const RooFormulaVar *select=0, const char *rangeName=0, Int_t nStart=0, Int_t nStop=2000000000)
Load values from tree 't' into this data collection, optionally selecting events using the RooFormula...
Definition: RooTreeDataStore.cxx:455
RooAbsRealLValue::getBins
virtual Int_t getBins(const char *name=0) const
Get number of bins of currently defined range.
Definition: RooAbsRealLValue.h:83
RMakeUnique.hxx
RooRealVar::setAsymError
void setAsymError(Double_t lo, Double_t hi)
Definition: RooRealVar.h:68
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooRealVar::getErrorLo
Double_t getErrorLo() const
Definition: RooRealVar.h:69
RooVectorDataStore.h
RooAbsDataStore
RooAbsDataStore is the abstract base class for data collection that use a TTree as internal storage m...
Definition: RooAbsDataStore.h:34
RooAbsData::attachCache
virtual void attachCache(const RooAbsArg *newOwner, const RooArgSet &cachedVars)
Internal method – Attach dataset copied with cache contents to copied instances of functions.
Definition: RooAbsData.cxx:347
RooDataSet::reduceEng
RooAbsData * reduceEng(const RooArgSet &varSubset, const RooFormulaVar *cutVar, const char *cutRange=0, std::size_t nStart=0, std::size_t nStop=std::numeric_limits< std::size_t >::max(), Bool_t copyCache=kTRUE) override
Implementation of RooAbsData virtual method that drives the RooAbsData::reduce() methods.
Definition: RooDataSet.cxx:933
RooAbsCollection::allInRange
Bool_t allInRange(const char *rangeSpec) const
Return true if all contained object report to have their value inside the specified range.
Definition: RooAbsCollection.cxx:1197
Class
void Class()
Definition: Class.C:29
TString::ToLower
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1125
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooAbsCollection::removeAll
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
Definition: RooAbsCollection.cxx:643
RooFormula::eval
Double_t eval(const RooArgSet *nset=0) const
Evalute all parameters/observables, and then evaluate formula.
Definition: RooFormula.cxx:343
RooAbsRealLValue
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Definition: RooAbsRealLValue.h:31
RooArgList.h
RooSentinel.h
RooAbsRealLValue.h
RooAbsArg::attachToStore
void attachToStore(RooAbsDataStore &store)
Attach this argument to the data store such that it reads data from there.
Definition: RooAbsArg.cxx:2204
RooPlot::getPlotVar
RooAbsRealLValue * getPlotVar() const
Definition: RooPlot.h:139
RooBatchCompute::RunContext::spans
std::unordered_map< const RooAbsReal *, RooSpan< const double > > spans
Once an object has computed its value(s), the span pointing to the results is registered here.
Definition: RunContext.h:51
RooBatchCompute::RunContext
Data that has to be passed around when evaluating functions / PDFs.
Definition: RunContext.h:32
RooCompositeDataStore.h
RooDataSet::createHistogram
TH1 * createHistogram(const char *name, const RooAbsRealLValue &xvar, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Calls createHistogram(const char *name, const RooAbsRealLValue& xvar, const RooLinkedList& argList) c...
Definition: RooAbsData.cxx:629
RooSpan
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
RooDataSet::weightError
virtual void weightError(Double_t &lo, Double_t &hi, ErrorType etype=SumW2) const override
Return asymmetric error on weight. (Dummy implementation returning zero)
Definition: RooDataSet.cxx:1032
RooAbsData::_cachedVars
RooArgSet _cachedVars
Definition: RooAbsData.h:292
RooDataSet::getWeightBatch
virtual RooSpan< const double > getWeightBatch(std::size_t first, std::size_t len) const override
Definition: RooDataSet.cxx:1010
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
RooAbsCollection::assignFast
void assignFast(const RooAbsCollection &other, Bool_t setValDirty=kTRUE)
Functional equivalent of operator=() but assumes this and other collection have same layout.
Definition: RooAbsCollection.cxx:356
RooDataSet::weight
virtual Double_t weight() const override
Return event weight of current event.
Definition: RooDataSet.cxx:992
RooDataSet::addColumn
virtual RooAbsArg * addColumn(RooAbsArg &var, Bool_t adjustRange=kTRUE)
Add a column with the values of the given (function) argument to this dataset.
Definition: RooDataSet.cxx:1396
RooAbsCollection::setAttribAll
void setAttribAll(const Text_t *name, Bool_t value=kTRUE)
Set given attribute in each element of the collection by calling each elements setAttribute() functio...
Definition: RooAbsCollection.cxx:662
int
Size_t
float Size_t
Definition: RtypesCore.h:87
RooAbsDataStore::getWeightBatch
virtual RooSpan< const double > getWeightBatch(std::size_t first, std::size_t len) const =0
RooDataSet::write
Bool_t write(const char *filename) const
Write the contents of this dataset to an ASCII file with the specified name.
Definition: RooDataSet.cxx:1908