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