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