Logo ROOT  
Reference Guide
RooAbsData.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 RooAbsData.cxx
19 \class RooAbsData
20 \ingroup Roofitcore
21 
22 RooAbsData is the common abstract base class for binned and unbinned
23 datasets. The abstract interface defines plotting and tabulating entry
24 points for its contents and provides an iterator over its elements
25 (bins for binned data sets, data points for unbinned datasets).
26 **/
27 
28 #include "RooAbsData.h"
29 
30 #include "TBuffer.h"
31 #include "TClass.h"
32 #include "TMath.h"
33 #include "TTree.h"
34 
35 #include "RooFormulaVar.h"
36 #include "RooCmdConfig.h"
37 #include "RooAbsRealLValue.h"
38 #include "RooMsgService.h"
39 #include "RooMultiCategory.h"
40 #include "Roo1DTable.h"
41 #include "RooAbsDataStore.h"
42 #include "RooVectorDataStore.h"
43 #include "RooTreeDataStore.h"
44 #include "RooDataHist.h"
45 #include "RooCompositeDataStore.h"
46 #include "RooCategory.h"
47 #include "RooTrace.h"
48 #include "RooUniformBinning.h"
49 
50 #include "RooRealVar.h"
51 #include "RooGlobalFunc.h"
52 #include "RooPlot.h"
53 #include "RooCurve.h"
54 #include "RooHist.h"
55 #include "RooHelpers.h"
56 
57 #include "TMatrixDSym.h"
58 #include "TPaveText.h"
59 #include "TH1.h"
60 #include "TH2.h"
61 #include "TH3.h"
62 #include "Math/Util.h"
63 
64 #include <iostream>
65 #include <memory>
66 
67 
68 using namespace std;
69 
71 ;
72 
73 static std::map<RooAbsData*,int> _dcc ;
74 
76 
77 ////////////////////////////////////////////////////////////////////////////////
78 
80 {
81  if (RooAbsData::Composite == s) {
82  cout << "Composite storage is not a valid *default* storage type." << endl;
83  } else {
84  defaultStorageType = s;
85  }
86 }
87 
88 ////////////////////////////////////////////////////////////////////////////////
89 
91 {
92  return defaultStorageType;
93 }
94 
95 ////////////////////////////////////////////////////////////////////////////////
96 
98 {
99  _dcc[data]++ ;
100  //cout << "RooAbsData(" << data << ") claim incremented to " << _dcc[data] << endl ;
101 }
102 
103 ////////////////////////////////////////////////////////////////////////////////
104 /// If return value is true variables can be deleted
105 
107 {
108  if (_dcc[data]>0) {
109  _dcc[data]-- ;
110  }
111 
112  //cout << "RooAbsData(" << data << ") claim decremented to " << _dcc[data] << endl ;
113  return (_dcc[data]==0) ;
114 }
115 
116 ////////////////////////////////////////////////////////////////////////////////
117 /// Default constructor
118 
120 {
121  claimVars(this) ;
122  _dstore = 0 ;
123  storageType = defaultStorageType;
124 
125  RooTrace::create(this) ;
126 }
127 
128 ////////////////////////////////////////////////////////////////////////////////
129 /// Constructor from a set of variables. Only fundamental elements of vars
130 /// (RooRealVar,RooCategory etc) are stored as part of the dataset
131 
132 RooAbsData::RooAbsData(const char *name, const char *title, const RooArgSet& vars, RooAbsDataStore* dstore) :
133  TNamed(name,title), _vars("Dataset Variables"), _cachedVars("Cached Variables"), _dstore(dstore)
134 {
135  if (dynamic_cast<RooTreeDataStore *>(dstore)) {
137  } else if (dynamic_cast<RooVectorDataStore *>(dstore)) {
139  } else {
141  }
142  // cout << "created dataset " << this << endl ;
143  claimVars(this);
144 
145  // clone the fundamentals of the given data set into internal buffer
146  for (const auto var : vars) {
147  if (!var->isFundamental()) {
148  coutE(InputArguments) << "RooAbsDataStore::initialize(" << GetName()
149  << "): Data set cannot contain non-fundamental types, ignoring " << var->GetName()
150  << endl;
151  throw std::invalid_argument(std::string("Only fundamental variables can be placed into datasets. This is violated for ") + var->GetName());
152  } else {
153  _vars.addClone(*var);
154  }
155  }
156 
157  // reconnect any parameterized ranges to internal dataset observables
158  for (auto var : _vars) {
159  var->attachDataSet(*this);
160  }
161 
162  RooTrace::create(this);
163 }
164 
165 ////////////////////////////////////////////////////////////////////////////////
166 /// Copy constructor
167 
168 RooAbsData::RooAbsData(const RooAbsData& other, const char* newname) :
169  TNamed(newname?newname:other.GetName(),other.GetTitle()),
170  RooPrintable(other), _vars(),
171  _cachedVars("Cached Variables")
172 {
173  //cout << "created dataset " << this << endl ;
174  claimVars(this) ;
175  _vars.addClone(other._vars) ;
176 
177  // reconnect any parameterized ranges to internal dataset observables
178  for (const auto var : _vars) {
179  var->attachDataSet(*this) ;
180  }
181 
182 
183  if (other._ownedComponents.size()>0) {
184 
185  // copy owned components here
186 
187  map<string,RooAbsDataStore*> smap ;
188  for (auto& itero : other._ownedComponents) {
189  RooAbsData* dclone = (RooAbsData*) itero.second->Clone();
190  _ownedComponents[itero.first] = dclone;
191  smap[itero.first] = dclone->store();
192  }
193 
194  RooCategory* idx = (RooCategory*) _vars.find(*((RooCompositeDataStore*)other.store())->index()) ;
195  _dstore = new RooCompositeDataStore(newname?newname:other.GetName(),other.GetTitle(),_vars,*idx,smap) ;
197 
198  } else {
199 
200  // Convert to vector store if default is vector
201  _dstore = other._dstore->clone(_vars,newname?newname:other.GetName()) ;
202  storageType = other.storageType;
203  }
204 
205  RooTrace::create(this) ;
206 }
207 
209  TNamed::operator=(other);
211 
212  claimVars(this);
213  _vars.Clear();
214  _vars.addClone(other._vars);
215 
216  // reconnect any parameterized ranges to internal dataset observables
217  for (const auto var : _vars) {
218  var->attachDataSet(*this) ;
219  }
220 
221 
222  if (other._ownedComponents.size()>0) {
223 
224  // copy owned components here
225 
226  map<string,RooAbsDataStore*> smap ;
227  for (auto& itero : other._ownedComponents) {
228  RooAbsData* dclone = (RooAbsData*) itero.second->Clone();
229  _ownedComponents[itero.first] = dclone;
230  smap[itero.first] = dclone->store();
231  }
232 
233  RooCategory* idx = (RooCategory*) _vars.find(*((RooCompositeDataStore*)other.store())->index()) ;
234  _dstore = new RooCompositeDataStore(GetName(), GetTitle(), _vars, *idx, smap);
236 
237  } else {
238 
239  // Convert to vector store if default is vector
240  _dstore = other._dstore->clone(_vars);
241  storageType = other.storageType;
242  }
243 
244  return *this;
245 }
246 
247 ////////////////////////////////////////////////////////////////////////////////
248 /// Destructor
249 
251 {
252  if (releaseVars(this)) {
253  // will cause content to be deleted subsequently in dtor
254  } else {
256  }
257 
258  // delete owned contents.
259  delete _dstore ;
260 
261  // Delete owned dataset components
262  for(map<std::string,RooAbsData*>::iterator iter = _ownedComponents.begin() ; iter!= _ownedComponents.end() ; ++iter) {
263  delete iter->second ;
264  }
265 
266  RooTrace::destroy(this) ;
267 }
268 
269 ////////////////////////////////////////////////////////////////////////////////
270 /// Convert tree-based storage to vector-based storage
271 
273 {
274  if (storageType == RooAbsData::Tree) {
276  delete _dstore;
277  _dstore = newStore;
279  }
280 }
281 
282 ////////////////////////////////////////////////////////////////////////////////
283 
284 Bool_t RooAbsData::changeObservableName(const char* from, const char* to)
285 {
286  Bool_t ret = _dstore->changeObservableName(from,to) ;
287 
288  RooAbsArg* tmp = _vars.find(from) ;
289  if (tmp) {
290  tmp->SetName(to) ;
291  }
292  return ret ;
293 }
294 
295 ////////////////////////////////////////////////////////////////////////////////
296 
298 {
299  _dstore->fill() ;
300 }
301 
302 ////////////////////////////////////////////////////////////////////////////////
303 
305 {
306  return nullptr != _dstore ? _dstore->numEntries() : 0;
307 }
308 
309 ////////////////////////////////////////////////////////////////////////////////
310 
312 {
313  _dstore->reset() ;
314 }
315 
316 ////////////////////////////////////////////////////////////////////////////////
317 
318 const RooArgSet* RooAbsData::get(Int_t index) const
319 {
320  checkInit() ;
321  return _dstore->get(index) ;
322 }
323 
324 ////////////////////////////////////////////////////////////////////////////////
325 /// Internal method -- Cache given set of functions with data
326 
327 void RooAbsData::cacheArgs(const RooAbsArg* cacheOwner, RooArgSet& varSet, const RooArgSet* nset, Bool_t skipZeroWeights)
328 {
329  _dstore->cacheArgs(cacheOwner,varSet,nset,skipZeroWeights) ;
330 }
331 
332 ////////////////////////////////////////////////////////////////////////////////
333 /// Internal method -- Remove cached function values
334 
336 {
337  _dstore->resetCache() ;
339 }
340 
341 ////////////////////////////////////////////////////////////////////////////////
342 /// Internal method -- Attach dataset copied with cache contents to copied instances of functions
343 
344 void RooAbsData::attachCache(const RooAbsArg* newOwner, const RooArgSet& cachedVars)
345 {
346  _dstore->attachCache(newOwner, cachedVars) ;
347 }
348 
349 ////////////////////////////////////////////////////////////////////////////////
350 
351 void RooAbsData::setArgStatus(const RooArgSet& set, Bool_t active)
352 {
353  _dstore->setArgStatus(set,active) ;
354 }
355 
356 ////////////////////////////////////////////////////////////////////////////////
357 /// Control propagation of dirty flags from observables in dataset
358 
360 {
361  _dstore->setDirtyProp(flag) ;
362 }
363 
364 ////////////////////////////////////////////////////////////////////////////////
365 /// Create a reduced copy of this dataset. The caller takes ownership of the returned dataset
366 ///
367 /// The following optional named arguments are accepted
368 /// <table>
369 /// <tr><td> `SelectVars(const RooArgSet& vars)` <td> Only retain the listed observables in the output dataset
370 /// <tr><td> `Cut(const char* expression)` <td> Only retain event surviving the given cut expression
371 /// <tr><td> `Cut(const RooFormulaVar& expr)` <td> Only retain event surviving the given cut formula
372 /// <tr><td> `CutRange(const char* name)` <td> Only retain events inside range with given name. Multiple CutRange
373 /// arguments may be given to select multiple ranges
374 /// <tr><td> `EventRange(int lo, int hi)` <td> Only retain events with given sequential event numbers
375 /// <tr><td> `Name(const char* name)` <td> Give specified name to output dataset
376 /// <tr><td> `Title(const char* name)` <td> Give specified title to output dataset
377 /// </table>
378 
379 RooAbsData* RooAbsData::reduce(const RooCmdArg& arg1,const RooCmdArg& arg2,const RooCmdArg& arg3,const RooCmdArg& arg4,
380  const RooCmdArg& arg5,const RooCmdArg& arg6,const RooCmdArg& arg7,const RooCmdArg& arg8)
381 {
382  // Define configuration for this method
383  RooCmdConfig pc(Form("RooAbsData::reduce(%s)",GetName())) ;
384  pc.defineString("name","Name",0,"") ;
385  pc.defineString("title","Title",0,"") ;
386  pc.defineString("cutRange","CutRange",0,"") ;
387  pc.defineString("cutSpec","CutSpec",0,"") ;
388  pc.defineObject("cutVar","CutVar",0,0) ;
389  pc.defineInt("evtStart","EventRange",0,0) ;
390  pc.defineInt("evtStop","EventRange",1,std::numeric_limits<int>::max()) ;
391  pc.defineObject("varSel","SelectVars",0,0) ;
392  pc.defineMutex("CutVar","CutSpec") ;
393 
394  // Process & check varargs
395  pc.process(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8) ;
396  if (!pc.ok(kTRUE)) {
397  return 0 ;
398  }
399 
400  // Extract values from named arguments
401  const char* cutRange = pc.getString("cutRange",0,kTRUE) ;
402  const char* cutSpec = pc.getString("cutSpec",0,kTRUE) ;
403  RooFormulaVar* cutVar = static_cast<RooFormulaVar*>(pc.getObject("cutVar",0)) ;
404  Int_t nStart = pc.getInt("evtStart",0) ;
405  Int_t nStop = pc.getInt("evtStop",std::numeric_limits<int>::max()) ;
406  RooArgSet* varSet = static_cast<RooArgSet*>(pc.getObject("varSel")) ;
407  const char* name = pc.getString("name",0,kTRUE) ;
408  const char* title = pc.getString("title",0,kTRUE) ;
409 
410  // Make sure varSubset doesn't contain any variable not in this dataset
411  RooArgSet varSubset ;
412  if (varSet) {
413  varSubset.add(*varSet) ;
414  for (const auto arg : varSubset) {
415  if (!_vars.find(arg->GetName())) {
416  coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
417  << arg->GetName() << " not in dataset, ignored" << endl ;
418  varSubset.remove(*arg) ;
419  }
420  }
421  } else {
422  varSubset.add(*get()) ;
423  }
424 
425  RooAbsData* ret = 0 ;
426  if (cutSpec) {
427 
428  RooFormulaVar cutVarTmp(cutSpec,cutSpec,*get()) ;
429  ret = reduceEng(varSubset,&cutVarTmp,cutRange,nStart,nStop,kFALSE) ;
430 
431  } else if (cutVar) {
432 
433  ret = reduceEng(varSubset,cutVar,cutRange,nStart,nStop,kFALSE) ;
434 
435  } else {
436 
437  ret = reduceEng(varSubset,0,cutRange,nStart,nStop,kFALSE) ;
438 
439  }
440 
441  if (!ret) return 0 ;
442 
443  if (name) {
444  ret->SetName(name) ;
445  }
446  if (title) {
447  ret->SetTitle(title) ;
448  }
449 
450  return ret ;
451 }
452 
453 ////////////////////////////////////////////////////////////////////////////////
454 /// Create a subset of the data set by applying the given cut on the data points.
455 /// The cut expression can refer to any variable in the data set. For cuts involving
456 /// other variables, such as intermediate formula objects, use the equivalent
457 /// reduce method specifying the as a RooFormulVar reference.
458 
460 {
461  RooFormulaVar cutVar(cut,cut,*get()) ;
462  return reduceEng(*get(),&cutVar,0,0,std::numeric_limits<std::size_t>::max(),kFALSE) ;
463 }
464 
465 ////////////////////////////////////////////////////////////////////////////////
466 /// Create a subset of the data set by applying the given cut on the data points.
467 /// The 'cutVar' formula variable is used to select the subset of data points to be
468 /// retained in the reduced data collection.
469 
471 {
472  return reduceEng(*get(),&cutVar,0,0,std::numeric_limits<std::size_t>::max(),kFALSE) ;
473 }
474 
475 ////////////////////////////////////////////////////////////////////////////////
476 /// Create a subset of the data set by applying the given cut on the data points
477 /// and reducing the dimensions to the specified set.
478 ///
479 /// The cut expression can refer to any variable in the data set. For cuts involving
480 /// other variables, such as intermediate formula objects, use the equivalent
481 /// reduce method specifying the as a RooFormulVar reference.
482 
483 RooAbsData* RooAbsData::reduce(const RooArgSet& varSubset, const char* cut)
484 {
485  // Make sure varSubset doesn't contain any variable not in this dataset
486  RooArgSet varSubset2(varSubset) ;
487  for (const auto arg : varSubset) {
488  if (!_vars.find(arg->GetName())) {
489  coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
490  << arg->GetName() << " not in dataset, ignored" << endl ;
491  varSubset2.remove(*arg) ;
492  }
493  }
494 
495  if (cut && strlen(cut)>0) {
496  RooFormulaVar cutVar(cut, cut, *get(), false);
497  return reduceEng(varSubset2,&cutVar,0,0,std::numeric_limits<std::size_t>::max(),kFALSE) ;
498  }
499  return reduceEng(varSubset2,0,0,0,std::numeric_limits<std::size_t>::max(),kFALSE) ;
500 }
501 
502 ////////////////////////////////////////////////////////////////////////////////
503 /// Create a subset of the data set by applying the given cut on the data points
504 /// and reducing the dimensions to the specified set.
505 ///
506 /// The 'cutVar' formula variable is used to select the subset of data points to be
507 /// retained in the reduced data collection.
508 
509 RooAbsData* RooAbsData::reduce(const RooArgSet& varSubset, const RooFormulaVar& cutVar)
510 {
511  // Make sure varSubset doesn't contain any variable not in this dataset
512  RooArgSet varSubset2(varSubset) ;
513  TIterator* iter = varSubset.createIterator() ;
514  RooAbsArg* arg ;
515  while((arg=(RooAbsArg*)iter->Next())) {
516  if (!_vars.find(arg->GetName())) {
517  coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
518  << arg->GetName() << " not in dataset, ignored" << endl ;
519  varSubset2.remove(*arg) ;
520  }
521  }
522  delete iter ;
523 
524  return reduceEng(varSubset2,&cutVar,0,0,std::numeric_limits<std::size_t>::max(),kFALSE) ;
525 }
526 
527 ////////////////////////////////////////////////////////////////////////////////
528 /// Return error on current weight (dummy implementation returning zero)
529 
531 {
532  return 0 ;
533 }
534 
535 ////////////////////////////////////////////////////////////////////////////////
536 /// Return asymmetric error on weight. (Dummy implementation returning zero)
537 
539 {
540  lo=0 ; hi=0 ;
541 }
542 
543 
544 RooPlot* RooAbsData::plotOn(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
545  const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
546  const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8) const
547 {
548  RooLinkedList l ;
549  l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
550  l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
551  l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
552  l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
553  return plotOn(frame,l) ;
554 }
555 
556 ////////////////////////////////////////////////////////////////////////////////
557 /// Create and fill a ROOT histogram TH1,TH2 or TH3 with the values of this dataset for the variables with given names
558 /// The range of each observable that is histogrammed is always automatically calculated from the distribution in
559 /// the dataset. The number of bins can be controlled using the [xyz]bins parameters. For a greater degree of control
560 /// use the createHistogram() method below with named arguments
561 ///
562 /// The caller takes ownership of the returned histogram
563 
564 TH1 *RooAbsData::createHistogram(const char* varNameList, Int_t xbins, Int_t ybins, Int_t zbins) const
565 {
566  // Parse list of variable names
567  const auto varNames = RooHelpers::tokenise(varNameList, ",:");
568  RooLinkedList argList;
569  RooRealVar* vars[3] = {nullptr, nullptr, nullptr};
570 
571  for (unsigned int i = 0; i < varNames.size(); ++i) {
572  if (i >= 3) {
573  coutW(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): Can only create 3-dimensional histograms. Variable "
574  << i << " " << varNames[i] << " unused." << std::endl;
575  continue;
576  }
577 
578  vars[i] = static_cast<RooRealVar*>( get()->find(varNames[i].data()) );
579  if (!vars[i]) {
580  coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << ") ERROR: dataset does not contain an observable named " << varNames[i] << std::endl;
581  return nullptr;
582  }
583  }
584 
585  if (!vars[0]) {
586  coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): No variable to be histogrammed in list '" << varNameList << "'" << std::endl;
587  return nullptr;
588  }
589 
590  if (xbins<=0 || !vars[0]->hasMax() || !vars[0]->hasMin() ) {
591  argList.Add(RooFit::AutoBinning(xbins==0?vars[0]->numBins():abs(xbins)).Clone()) ;
592  } else {
593  argList.Add(RooFit::Binning(xbins).Clone()) ;
594  }
595 
596  if (vars[1]) {
597  if (ybins<=0 || !vars[1]->hasMax() || !vars[1]->hasMin() ) {
598  argList.Add(RooFit::YVar(*vars[1],RooFit::AutoBinning(ybins==0?vars[1]->numBins():abs(ybins))).Clone()) ;
599  } else {
600  argList.Add(RooFit::YVar(*vars[1],RooFit::Binning(ybins)).Clone()) ;
601  }
602  }
603  if (vars[2]) {
604  if (zbins<=0 || !vars[2]->hasMax() || !vars[2]->hasMin() ) {
605  argList.Add(RooFit::ZVar(*vars[2],RooFit::AutoBinning(zbins==0?vars[2]->numBins():abs(zbins))).Clone()) ;
606  } else {
607  argList.Add(RooFit::ZVar(*vars[2],RooFit::Binning(zbins)).Clone()) ;
608  }
609  }
610 
611 
612  // Call implementation function
613  TH1* result = createHistogram(GetName(), *vars[0], argList);
614 
615  // Delete temporary list of RooCmdArgs
616  argList.Delete() ;
617 
618  return result ;
619 }
620 
621 
623  const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3, const RooCmdArg& arg4,
624  const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8) const
625 {
626  RooLinkedList l ;
627  l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
628  l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
629  l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
630  l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
631 
632  return createHistogram(name,xvar,l) ;
633 }
634 
635 ////////////////////////////////////////////////////////////////////////////////
636 ///
637 /// This function accepts the following arguments
638 ///
639 /// \param[in] name Name of the ROOT histogram
640 /// \param[in] xvar Observable to be mapped on x axis of ROOT histogram
641 /// \return Histogram now owned by user.
642 ///
643 /// <table>
644 /// <tr><td> `AutoBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin, set binning to nbins
645 /// <tr><td> `AutoSymBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin,
646 /// with additional constraint that mean of data is in center of range, set binning to nbins
647 /// <tr><td> `Binning(const char* name)` <td> Apply binning with given name to x axis of histogram
648 /// <tr><td> `Binning(RooAbsBinning& binning)` <td> Apply specified binning to x axis of histogram
649 /// <tr><td> `Binning(int nbins, double lo, double hi)` <td> Apply specified binning to x axis of histogram
650 ///
651 /// <tr><td> `YVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on y axis of ROOT histogram
652 /// <tr><td> `ZVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on z axis of ROOT histogram
653 /// </table>
654 ///
655 /// The YVar() and ZVar() arguments can be supplied with optional Binning() Auto(Sym)Range() arguments to control the binning of the Y and Z axes, e.g.
656 /// ```
657 /// createHistogram("histo",x,Binning(-1,1,20), YVar(y,Binning(-1,1,30)), ZVar(z,Binning("zbinning")))
658 /// ```
659 ///
660 /// The caller takes ownership of the returned histogram
661 
662 TH1 *RooAbsData::createHistogram(const char *name, const RooAbsRealLValue& xvar, const RooLinkedList& argListIn) const
663 {
664  RooLinkedList argList(argListIn) ;
665 
666  // Define configuration for this method
667  RooCmdConfig pc(Form("RooAbsData::createHistogram(%s)",GetName())) ;
668  pc.defineString("cutRange","CutRange",0,"",kTRUE) ;
669  pc.defineString("cutString","CutSpec",0,"") ;
670  pc.defineObject("yvar","YVar",0,0) ;
671  pc.defineObject("zvar","ZVar",0,0) ;
672  pc.allowUndefined() ;
673 
674  // Process & check varargs
675  pc.process(argList) ;
676  if (!pc.ok(kTRUE)) {
677  return 0 ;
678  }
679 
680  const char* cutSpec = pc.getString("cutString",0,kTRUE) ;
681  const char* cutRange = pc.getString("cutRange",0,kTRUE) ;
682 
683  RooArgList vars(xvar) ;
684  RooAbsArg* yvar = static_cast<RooAbsArg*>(pc.getObject("yvar")) ;
685  if (yvar) {
686  vars.add(*yvar) ;
687  }
688  RooAbsArg* zvar = static_cast<RooAbsArg*>(pc.getObject("zvar")) ;
689  if (zvar) {
690  vars.add(*zvar) ;
691  }
692 
693  pc.stripCmdList(argList,"CutRange,CutSpec") ;
694 
695  // Swap Auto(Sym)RangeData with a Binning command
696  RooLinkedList ownedCmds ;
697  RooCmdArg* autoRD = (RooCmdArg*) argList.find("AutoRangeData") ;
698  if (autoRD) {
699  Double_t xmin,xmax ;
700  if (!getRange((RooRealVar&)xvar,xmin,xmax,autoRD->getDouble(0),autoRD->getInt(0))) {
701  RooCmdArg* bincmd = (RooCmdArg*) RooFit::Binning(autoRD->getInt(1),xmin,xmax).Clone() ;
702  ownedCmds.Add(bincmd) ;
703  argList.Replace(autoRD,bincmd) ;
704  }
705  }
706 
707  if (yvar) {
708  RooCmdArg* autoRDY = (RooCmdArg*) ((RooCmdArg*)argList.find("YVar"))->subArgs().find("AutoRangeData") ;
709  if (autoRDY) {
710  Double_t ymin,ymax ;
711  if (!getRange((RooRealVar&)(*yvar),ymin,ymax,autoRDY->getDouble(0),autoRDY->getInt(0))) {
712  RooCmdArg* bincmd = (RooCmdArg*) RooFit::Binning(autoRDY->getInt(1),ymin,ymax).Clone() ;
713  //ownedCmds.Add(bincmd) ;
714  ((RooCmdArg*)argList.find("YVar"))->subArgs().Replace(autoRDY,bincmd) ;
715  }
716  delete autoRDY ;
717  }
718  }
719 
720  if (zvar) {
721  RooCmdArg* autoRDZ = (RooCmdArg*) ((RooCmdArg*)argList.find("ZVar"))->subArgs().find("AutoRangeData") ;
722  if (autoRDZ) {
723  Double_t zmin,zmax ;
724  if (!getRange((RooRealVar&)(*zvar),zmin,zmax,autoRDZ->getDouble(0),autoRDZ->getInt(0))) {
725  RooCmdArg* bincmd = (RooCmdArg*) RooFit::Binning(autoRDZ->getInt(1),zmin,zmax).Clone() ;
726  //ownedCmds.Add(bincmd) ;
727  ((RooCmdArg*)argList.find("ZVar"))->subArgs().Replace(autoRDZ,bincmd) ;
728  }
729  delete autoRDZ ;
730  }
731  }
732 
733 
734  TH1* histo = xvar.createHistogram(name,argList) ;
735  fillHistogram(histo,vars,cutSpec,cutRange) ;
736 
737  ownedCmds.Delete() ;
738 
739  return histo ;
740 }
741 
742 ////////////////////////////////////////////////////////////////////////////////
743 /// Construct table for product of categories in catSet
744 
745 Roo1DTable* RooAbsData::table(const RooArgSet& catSet, const char* cuts, const char* opts) const
746 {
747  RooArgSet catSet2 ;
748 
749  string prodName("(") ;
750  TIterator* iter = catSet.createIterator() ;
751  RooAbsArg* arg ;
752  while((arg=(RooAbsArg*)iter->Next())) {
753  if (dynamic_cast<RooAbsCategory*>(arg)) {
754  RooAbsCategory* varsArg = dynamic_cast<RooAbsCategory*>(_vars.find(arg->GetName())) ;
755  if (varsArg != 0) catSet2.add(*varsArg) ;
756  else catSet2.add(*arg) ;
757  if (prodName.length()>1) {
758  prodName += " x " ;
759  }
760  prodName += arg->GetName() ;
761  } else {
762  coutW(InputArguments) << "RooAbsData::table(" << GetName() << ") non-RooAbsCategory input argument " << arg->GetName() << " ignored" << endl ;
763  }
764  }
765  prodName += ")" ;
766  delete iter ;
767 
768  RooMultiCategory tmp(prodName.c_str(),prodName.c_str(),catSet2) ;
769  return table(tmp,cuts,opts) ;
770 }
771 
772 ////////////////////////////////////////////////////////////////////////////////
773 /// Print name of dataset
774 
775 void RooAbsData::printName(ostream& os) const
776 {
777  os << GetName() ;
778 }
779 
780 ////////////////////////////////////////////////////////////////////////////////
781 /// Print title of dataset
782 
783 void RooAbsData::printTitle(ostream& os) const
784 {
785  os << GetTitle() ;
786 }
787 
788 ////////////////////////////////////////////////////////////////////////////////
789 /// Print class name of dataset
790 
791 void RooAbsData::printClassName(ostream& os) const
792 {
793  os << IsA()->GetName() ;
794 }
795 
796 ////////////////////////////////////////////////////////////////////////////////
797 
798 void RooAbsData::printMultiline(ostream& os, Int_t contents, Bool_t verbose, TString indent) const
799 {
800  _dstore->printMultiline(os,contents,verbose,indent) ;
801 }
802 
803 ////////////////////////////////////////////////////////////////////////////////
804 /// Define default print options, for a given print style
805 
807 {
808  return kName|kClassName|kArgs|kValue ;
809 }
810 
811 ////////////////////////////////////////////////////////////////////////////////
812 /// Calculate standardized moment.
813 ///
814 /// \param[in] var Variable to be used for calculating the moment.
815 /// \param[in] order Order of the moment.
816 /// \param[in] cutSpec If specified, the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
817 /// \param[in] cutRange If specified, calculate inside the range named 'cutRange' (also applies cut spec)
818 /// \return \f$ \frac{\left< \left( X - \left< X \right> \right)^n \right>}{\sigma^n} \f$, where n = order.
819 
820 Double_t RooAbsData::standMoment(const RooRealVar &var, Double_t order, const char* cutSpec, const char* cutRange) const
821 {
822  // Hardwire invariant answer for first and second moment
823  if (order==1) return 0 ;
824  if (order==2) return 1 ;
825 
826  return moment(var,order,cutSpec,cutRange) / TMath::Power(sigma(var,cutSpec,cutRange),order) ;
827 }
828 
829 ////////////////////////////////////////////////////////////////////////////////
830 /// Calculate moment of requested order.
831 ///
832 /// \param[in] var Variable to be used for calculating the moment.
833 /// \param[in] order Order of the moment.
834 /// \param[in] cutSpec If specified, the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
835 /// \param[in] cutRange If specified, calculate inside the range named 'cutRange' (also applies cut spec)
836 /// \return \f$ \left< \left( X - \left< X \right> \right)^n \right> \f$ of order \f$n\f$.
837 ///
838 
839 Double_t RooAbsData::moment(const RooRealVar& var, Double_t order, const char* cutSpec, const char* cutRange) const
840 {
841  Double_t offset = order>1 ? moment(var,1,cutSpec,cutRange) : 0 ;
842  return moment(var,order,offset,cutSpec,cutRange) ;
843 
844 }
845 
846 ////////////////////////////////////////////////////////////////////////////////
847 /// Return the 'order'-ed moment of observable 'var' in this dataset. If offset is non-zero it is subtracted
848 /// from the values of 'var' prior to the moment calculation. If cutSpec and/or cutRange are specified
849 /// the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
850 /// and/or are inside the range named 'cutRange'
851 
852 Double_t RooAbsData::moment(const RooRealVar& var, Double_t order, Double_t offset, const char* cutSpec, const char* cutRange) const
853 {
854  // Lookup variable in dataset
855  auto arg = _vars.find(var.GetName());
856  if (!arg) {
857  coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") ERROR: unknown variable: " << var.GetName() << endl ;
858  return 0;
859  }
860 
861  auto varPtr = dynamic_cast<const RooRealVar*>(arg);
862  // Check if found variable is of type RooRealVar
863  if (!varPtr) {
864  coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") ERROR: variable " << var.GetName() << " is not of type RooRealVar" << endl ;
865  return 0;
866  }
867 
868  // Check if dataset is not empty
869  if(sumEntries(cutSpec, cutRange) == 0.) {
870  coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") WARNING: empty dataset" << endl ;
871  return 0;
872  }
873 
874  // Setup RooFormulaVar for cutSpec if it is present
875  std::unique_ptr<RooFormula> select;
876  if (cutSpec) {
877  select.reset(new RooFormula("select",cutSpec,*get()));
878  }
879 
880 
881  // Calculate requested moment
883  for(Int_t index= 0; index < numEntries(); index++) {
884  const RooArgSet* vars = get(index) ;
885  if (select && select->eval()==0) continue ;
886  if (cutRange && vars->allInRange(cutRange)) continue ;
887 
888  sum += weight() * TMath::Power(varPtr->getVal() - offset,order);
889  }
890 
891  return sum/sumEntries(cutSpec, cutRange);
892 }
893 
894 ////////////////////////////////////////////////////////////////////////////////
895 /// Internal method to check if given RooRealVar maps to a RooRealVar in this dataset
896 
897 RooRealVar* RooAbsData::dataRealVar(const char* methodname, const RooRealVar& extVar) const
898 {
899  // Lookup variable in dataset
900  RooRealVar *xdata = (RooRealVar*) _vars.find(extVar.GetName());
901  if(!xdata) {
902  coutE(InputArguments) << "RooDataSet::" << methodname << "(" << GetName() << ") ERROR: variable : " << extVar.GetName() << " is not in data" << endl ;
903  return 0;
904  }
905  // Check if found variable is of type RooRealVar
906  if (!dynamic_cast<RooRealVar*>(xdata)) {
907  coutE(InputArguments) << "RooDataSet::" << methodname << "(" << GetName() << ") ERROR: variable : " << extVar.GetName() << " is not of type RooRealVar in data" << endl ;
908  return 0;
909  }
910  return xdata;
911 }
912 
913 ////////////////////////////////////////////////////////////////////////////////
914 /// Internal method to calculate single correlation and covariance elements
915 
916 Double_t RooAbsData::corrcov(const RooRealVar &x, const RooRealVar &y, const char* cutSpec, const char* cutRange, Bool_t corr) const
917 {
918  // Lookup variable in dataset
919  RooRealVar *xdata = dataRealVar(corr?"correlation":"covariance",x) ;
920  RooRealVar *ydata = dataRealVar(corr?"correlation":"covariance",y) ;
921  if (!xdata||!ydata) return 0 ;
922 
923  // Check if dataset is not empty
924  if(sumEntries(cutSpec, cutRange) == 0.) {
925  coutW(InputArguments) << "RooDataSet::" << (corr?"correlation":"covariance") << "(" << GetName() << ") WARNING: empty dataset, returning zero" << endl ;
926  return 0;
927  }
928 
929  // Setup RooFormulaVar for cutSpec if it is present
930  RooFormula* select = cutSpec ? new RooFormula("select",cutSpec,*get()) : 0 ;
931 
932  // Calculate requested moment
933  Double_t xysum(0),xsum(0),ysum(0),x2sum(0),y2sum(0);
934  const RooArgSet* vars ;
935  for(Int_t index= 0; index < numEntries(); index++) {
936  vars = get(index) ;
937  if (select && select->eval()==0) continue ;
938  if (cutRange && vars->allInRange(cutRange)) continue ;
939 
940  xysum += weight()*xdata->getVal()*ydata->getVal() ;
941  xsum += weight()*xdata->getVal() ;
942  ysum += weight()*ydata->getVal() ;
943  if (corr) {
944  x2sum += weight()*xdata->getVal()*xdata->getVal() ;
945  y2sum += weight()*ydata->getVal()*ydata->getVal() ;
946  }
947  }
948 
949  // Normalize entries
950  xysum/=sumEntries(cutSpec, cutRange) ;
951  xsum/=sumEntries(cutSpec, cutRange) ;
952  ysum/=sumEntries(cutSpec, cutRange) ;
953  if (corr) {
954  x2sum/=sumEntries(cutSpec, cutRange) ;
955  y2sum/=sumEntries(cutSpec, cutRange) ;
956  }
957 
958  // Cleanup
959  if (select) delete select ;
960 
961  // Return covariance or correlation as requested
962  if (corr) {
963  return (xysum-xsum*ysum)/(sqrt(x2sum-(xsum*xsum))*sqrt(y2sum-(ysum*ysum))) ;
964  } else {
965  return (xysum-xsum*ysum);
966  }
967 }
968 
969 ////////////////////////////////////////////////////////////////////////////////
970 /// Return covariance matrix from data for given list of observables
971 
972 TMatrixDSym* RooAbsData::corrcovMatrix(const RooArgList& vars, const char* cutSpec, const char* cutRange, Bool_t corr) const
973 {
974  RooArgList varList ;
975  TIterator* iter = vars.createIterator() ;
976  RooRealVar* var ;
977  while((var=(RooRealVar*)iter->Next())) {
978  RooRealVar* datavar = dataRealVar("covarianceMatrix",*var) ;
979  if (!datavar) {
980  delete iter ;
981  return 0 ;
982  }
983  varList.add(*datavar) ;
984  }
985  delete iter ;
986 
987 
988  // Check if dataset is not empty
989  if(sumEntries(cutSpec, cutRange) == 0.) {
990  coutW(InputArguments) << "RooDataSet::covariance(" << GetName() << ") WARNING: empty dataset, returning zero" << endl ;
991  return 0;
992  }
993 
994  // Setup RooFormulaVar for cutSpec if it is present
995  RooFormula* select = cutSpec ? new RooFormula("select",cutSpec,*get()) : 0 ;
996 
997  iter = varList.createIterator() ;
998  TIterator* iter2 = varList.createIterator() ;
999 
1000  TMatrixDSym xysum(varList.getSize()) ;
1001  vector<double> xsum(varList.getSize()) ;
1002  vector<double> x2sum(varList.getSize()) ;
1003 
1004  // Calculate <x_i> and <x_i y_j>
1005  for(Int_t index= 0; index < numEntries(); index++) {
1006  const RooArgSet* dvars = get(index) ;
1007  if (select && select->eval()==0) continue ;
1008  if (cutRange && dvars->allInRange(cutRange)) continue ;
1009 
1010  RooRealVar* varx, *vary ;
1011  iter->Reset() ;
1012  Int_t ix=0,iy=0 ;
1013  while((varx=(RooRealVar*)iter->Next())) {
1014  xsum[ix] += weight()*varx->getVal() ;
1015  if (corr) {
1016  x2sum[ix] += weight()*varx->getVal()*varx->getVal() ;
1017  }
1018 
1019  *iter2=*iter ; iy=ix ;
1020  vary=varx ;
1021  while(vary) {
1022  xysum(ix,iy) += weight()*varx->getVal()*vary->getVal() ;
1023  xysum(iy,ix) = xysum(ix,iy) ;
1024  iy++ ;
1025  vary=(RooRealVar*)iter2->Next() ;
1026  }
1027  ix++ ;
1028  }
1029 
1030  }
1031 
1032  // Normalize sums
1033  for (Int_t ix=0 ; ix<varList.getSize() ; ix++) {
1034  xsum[ix] /= sumEntries(cutSpec, cutRange) ;
1035  if (corr) {
1036  x2sum[ix] /= sumEntries(cutSpec, cutRange) ;
1037  }
1038  for (Int_t iy=0 ; iy<varList.getSize() ; iy++) {
1039  xysum(ix,iy) /= sumEntries(cutSpec, cutRange) ;
1040  }
1041  }
1042 
1043  // Calculate covariance matrix
1044  TMatrixDSym* C = new TMatrixDSym(varList.getSize()) ;
1045  for (Int_t ix=0 ; ix<varList.getSize() ; ix++) {
1046  for (Int_t iy=0 ; iy<varList.getSize() ; iy++) {
1047  (*C)(ix,iy) = xysum(ix,iy)-xsum[ix]*xsum[iy] ;
1048  if (corr) {
1049  (*C)(ix,iy) /= sqrt((x2sum[ix]-(xsum[ix]*xsum[ix]))*(x2sum[iy]-(xsum[iy]*xsum[iy]))) ;
1050  }
1051  }
1052  }
1053 
1054  if (select) delete select ;
1055  delete iter ;
1056  delete iter2 ;
1057 
1058  return C ;
1059 }
1060 
1061 ////////////////////////////////////////////////////////////////////////////////
1062 /// Create a RooRealVar containing the mean of observable 'var' in
1063 /// this dataset. If cutSpec and/or cutRange are specified the
1064 /// moment is calculated on the subset of the data which pass the C++
1065 /// cut specification expression 'cutSpec' and/or are inside the
1066 /// range named 'cutRange'
1067 
1068 RooRealVar* RooAbsData::meanVar(const RooRealVar &var, const char* cutSpec, const char* cutRange) const
1069 {
1070  // Create a new variable with appropriate strings. The error is calculated as
1071  // RMS/Sqrt(N) which is generally valid.
1072 
1073  // Create holder variable for mean
1074  TString name(var.GetName()),title("Mean of ") ;
1075  name.Append("Mean");
1076  title.Append(var.GetTitle());
1077  RooRealVar *meanv= new RooRealVar(name,title,0) ;
1078  meanv->setConstant(kFALSE) ;
1079 
1080  // Adjust plot label
1081  TString label("<") ;
1082  label.Append(var.getPlotLabel());
1083  label.Append(">");
1084  meanv->setPlotLabel(label.Data());
1085 
1086  // fill in this variable's value and error
1087  Double_t meanVal=moment(var,1,0,cutSpec,cutRange) ;
1088  Double_t N(sumEntries(cutSpec,cutRange)) ;
1089 
1090  Double_t rmsVal= sqrt(moment(var,2,meanVal,cutSpec,cutRange)*N/(N-1));
1091  meanv->setVal(meanVal) ;
1092  meanv->setError(N > 0 ? rmsVal/sqrt(N) : 0);
1093 
1094  return meanv;
1095 }
1096 
1097 ////////////////////////////////////////////////////////////////////////////////
1098 /// Create a RooRealVar containing the RMS of observable 'var' in
1099 /// this dataset. If cutSpec and/or cutRange are specified the
1100 /// moment is calculated on the subset of the data which pass the C++
1101 /// cut specification expression 'cutSpec' and/or are inside the
1102 /// range named 'cutRange'
1103 
1104 RooRealVar* RooAbsData::rmsVar(const RooRealVar &var, const char* cutSpec, const char* cutRange) const
1105 {
1106  // Create a new variable with appropriate strings. The error is calculated as
1107  // RMS/(2*Sqrt(N)) which is only valid if the variable has a Gaussian distribution.
1108 
1109  // Create RMS value holder
1110  TString name(var.GetName()),title("RMS of ") ;
1111  name.Append("RMS");
1112  title.Append(var.GetTitle());
1113  RooRealVar *rms= new RooRealVar(name,title,0) ;
1114  rms->setConstant(kFALSE) ;
1115 
1116  // Adjust plot label
1117  TString label(var.getPlotLabel());
1118  label.Append("_{RMS}");
1119  rms->setPlotLabel(label);
1120 
1121  // Fill in this variable's value and error
1122  Double_t meanVal(moment(var,1,0,cutSpec,cutRange)) ;
1123  Double_t N(sumEntries(cutSpec, cutRange));
1124  Double_t rmsVal= sqrt(moment(var,2,meanVal,cutSpec,cutRange)*N/(N-1));
1125  rms->setVal(rmsVal) ;
1126  rms->setError(rmsVal/sqrt(2*N));
1127 
1128  return rms;
1129 }
1130 
1131 ////////////////////////////////////////////////////////////////////////////////
1132 /// Add a box with statistics information to the specified frame. By default a box with the
1133 /// event count, mean and rms of the plotted variable is added.
1134 ///
1135 /// The following optional named arguments are accepted
1136 /// <table>
1137 /// <tr><td> `What(const char* whatstr)` <td> Controls what is printed: "N" = count, "M" is mean, "R" is RMS.
1138 /// <tr><td> `Format(const char* optStr)` <td> \deprecated Classing parameter formatting options, provided for backward compatibility
1139 ///
1140 /// <tr><td> `Format(const char* what,...)` <td> Parameter formatting options.
1141 /// <table>
1142 /// <tr><td> const char* what <td> Controls what is shown:
1143 /// - "N" adds name
1144 /// - "E" adds error
1145 /// - "A" shows asymmetric error
1146 /// - "U" shows unit
1147 /// - "H" hides the value
1148 /// <tr><td> `FixedPrecision(int n)` <td> Controls precision, set fixed number of digits
1149 /// <tr><td> `AutoPrecision(int n)` <td> Controls precision. Number of shown digits is calculated from error + n specified additional digits (1 is sensible default)
1150 /// <tr><td> `VerbatimName(Bool_t flag)` <td> Put variable name in a \\verb+ + clause.
1151 /// </table>
1152 /// <tr><td> `Label(const chat* label)` <td> Add header label to parameter box
1153 /// <tr><td> `Layout(Double_t xmin, Double_t xmax, Double_t ymax)` <td> Specify relative position of left,right side of box and top of box. Position of
1154 /// bottom of box is calculated automatically from number lines in box
1155 /// <tr><td> `Cut(const char* expression)` <td> Apply given cut expression to data when calculating statistics
1156 /// <tr><td> `CutRange(const char* rangeName)` <td> Only consider events within given range when calculating statistics. Multiple
1157 /// CutRange() argument may be specified to combine ranges.
1158 ///
1159 /// </table>
1160 
1161 RooPlot* RooAbsData::statOn(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
1162  const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
1163  const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
1164 {
1165  // Stuff all arguments in a list
1166  RooLinkedList cmdList;
1167  cmdList.Add(const_cast<RooCmdArg*>(&arg1)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg2)) ;
1168  cmdList.Add(const_cast<RooCmdArg*>(&arg3)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg4)) ;
1169  cmdList.Add(const_cast<RooCmdArg*>(&arg5)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg6)) ;
1170  cmdList.Add(const_cast<RooCmdArg*>(&arg7)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg8)) ;
1171 
1172  // Select the pdf-specific commands
1173  RooCmdConfig pc(Form("RooTreeData::statOn(%s)",GetName())) ;
1174  pc.defineString("what","What",0,"MNR") ;
1175  pc.defineString("label","Label",0,"") ;
1176  pc.defineDouble("xmin","Layout",0,0.65) ;
1177  pc.defineDouble("xmax","Layout",1,0.99) ;
1178  pc.defineInt("ymaxi","Layout",0,Int_t(0.95*10000)) ;
1179  pc.defineString("formatStr","Format",0,"NELU") ;
1180  pc.defineInt("sigDigit","Format",0,2) ;
1181  pc.defineInt("dummy","FormatArgs",0,0) ;
1182  pc.defineString("cutRange","CutRange",0,"",kTRUE) ;
1183  pc.defineString("cutString","CutSpec",0,"") ;
1184  pc.defineMutex("Format","FormatArgs") ;
1185 
1186  // Process and check varargs
1187  pc.process(cmdList) ;
1188  if (!pc.ok(kTRUE)) {
1189  return frame ;
1190  }
1191 
1192  const char* label = pc.getString("label") ;
1193  Double_t xmin = pc.getDouble("xmin") ;
1194  Double_t xmax = pc.getDouble("xmax") ;
1195  Double_t ymax = pc.getInt("ymaxi") / 10000. ;
1196  const char* formatStr = pc.getString("formatStr") ;
1197  Int_t sigDigit = pc.getInt("sigDigit") ;
1198  const char* what = pc.getString("what") ;
1199 
1200  const char* cutSpec = pc.getString("cutString",0,kTRUE) ;
1201  const char* cutRange = pc.getString("cutRange",0,kTRUE) ;
1202 
1203  if (pc.hasProcessed("FormatArgs")) {
1204  RooCmdArg* formatCmd = static_cast<RooCmdArg*>(cmdList.FindObject("FormatArgs")) ;
1205  return statOn(frame,what,label,0,0,xmin,xmax,ymax,cutSpec,cutRange,formatCmd) ;
1206  } else {
1207  return statOn(frame,what,label,sigDigit,formatStr,xmin,xmax,ymax,cutSpec,cutRange) ;
1208  }
1209 }
1210 
1211 ////////////////////////////////////////////////////////////////////////////////
1212 /// Implementation back-end of statOn() method with named arguments
1213 
1214 RooPlot* RooAbsData::statOn(RooPlot* frame, const char* what, const char *label, Int_t sigDigits,
1216  const char* cutSpec, const char* cutRange, const RooCmdArg* formatCmd)
1217 {
1218  Bool_t showLabel= (label != 0 && strlen(label) > 0);
1219 
1220  TString whatStr(what) ;
1221  whatStr.ToUpper() ;
1222  Bool_t showN = whatStr.Contains("N") ;
1223  Bool_t showR = whatStr.Contains("R") ;
1224  Bool_t showM = whatStr.Contains("M") ;
1225  Int_t nPar= 0;
1226  if (showN) nPar++ ;
1227  if (showR) nPar++ ;
1228  if (showM) nPar++ ;
1229 
1230  // calculate the box's size
1231  Double_t dy(0.06), ymin(ymax-nPar*dy);
1232  if(showLabel) ymin-= dy;
1233 
1234  // create the box and set its options
1235  TPaveText *box= new TPaveText(xmin,ymax,xmax,ymin,"BRNDC");
1236  if(!box) return 0;
1237  box->SetName(Form("%s_statBox",GetName())) ;
1238  box->SetFillColor(0);
1239  box->SetBorderSize(1);
1240  box->SetTextAlign(12);
1241  box->SetTextSize(0.04F);
1242  box->SetFillStyle(1001);
1243 
1244  // add formatted text for each statistic
1245  RooRealVar N("N","Number of Events",sumEntries(cutSpec,cutRange));
1246  N.setPlotLabel("Entries") ;
1247  RooRealVar *meanv= meanVar(*(RooRealVar*)frame->getPlotVar(),cutSpec,cutRange);
1248  meanv->setPlotLabel("Mean") ;
1249  RooRealVar *rms= rmsVar(*(RooRealVar*)frame->getPlotVar(),cutSpec,cutRange);
1250  rms->setPlotLabel("RMS") ;
1251  TString *rmsText, *meanText, *NText ;
1252  if (options) {
1253  rmsText= rms->format(sigDigits,options);
1254  meanText= meanv->format(sigDigits,options);
1255  NText= N.format(sigDigits,options);
1256  } else {
1257  rmsText= rms->format(*formatCmd);
1258  meanText= meanv->format(*formatCmd);
1259  NText= N.format(*formatCmd);
1260  }
1261  if (showR) box->AddText(rmsText->Data());
1262  if (showM) box->AddText(meanText->Data());
1263  if (showN) box->AddText(NText->Data());
1264 
1265  // cleanup heap memory
1266  delete NText;
1267  delete meanText;
1268  delete rmsText;
1269  delete meanv;
1270  delete rms;
1271 
1272  // add the optional label if specified
1273  if(showLabel) box->AddText(label);
1274 
1275  frame->addObject(box) ;
1276  return frame ;
1277 }
1278 
1279 ////////////////////////////////////////////////////////////////////////////////
1280 /// Loop over columns of our tree data and fill the input histogram. Returns a pointer to the
1281 /// input histogram, or zero in case of an error. The input histogram can be any TH1 subclass, and
1282 /// therefore of arbitrary dimension. Variables are matched with the (x,y,...) dimensions of the input
1283 /// histogram according to the order in which they appear in the input plotVars list.
1284 
1285 TH1 *RooAbsData::fillHistogram(TH1 *hist, const RooArgList &plotVars, const char *cuts, const char* cutRange) const
1286 {
1287  // Do we have a valid histogram to use?
1288  if(0 == hist) {
1289  coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: no valid histogram to fill" << endl;
1290  return 0;
1291  }
1292 
1293  // Check that the number of plotVars matches the input histogram's dimension
1294  Int_t hdim= hist->GetDimension();
1295  if(hdim != plotVars.getSize()) {
1296  coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: plotVars has the wrong dimension" << endl;
1297  return 0;
1298  }
1299 
1300  // Check that the plot variables are all actually RooAbsReal's and print a warning if we do not
1301  // explicitly depend on one of them. Clone any variables that we do not contain directly and
1302  // redirect them to use our event data.
1303  RooArgSet plotClones,localVars;
1304  for(Int_t index= 0; index < plotVars.getSize(); index++) {
1305  const RooAbsArg *var= plotVars.at(index);
1306  const RooAbsReal *realVar= dynamic_cast<const RooAbsReal*>(var);
1307  if(0 == realVar) {
1308  coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: cannot plot variable \"" << var->GetName()
1309  << "\" of type " << var->ClassName() << endl;
1310  return 0;
1311  }
1312  RooAbsArg *found= _vars.find(realVar->GetName());
1313  if(!found) {
1314  RooAbsArg *clone= plotClones.addClone(*realVar,kTRUE); // do not complain about duplicates
1315  assert(0 != clone);
1316  if(!clone->dependsOn(_vars)) {
1317  coutE(InputArguments) << ClassName() << "::" << GetName()
1318  << ":fillHistogram: Data does not contain the variable '" << realVar->GetName() << "'." << endl;
1319  return nullptr;
1320  }
1321  else {
1323  }
1324  localVars.add(*clone);
1325  }
1326  else {
1327  localVars.add(*found);
1328  }
1329  }
1330 
1331  // Create selection formula if selection cuts are specified
1332  std::unique_ptr<RooFormula> select;
1333  if (cuts != nullptr && strlen(cuts) > 0) {
1334  select.reset(new RooFormula(cuts, cuts, _vars, false));
1335  if (!select || !select->ok()) {
1336  coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: invalid cuts \"" << cuts << "\"" << endl;
1337  return 0 ;
1338  }
1339  }
1340 
1341  // Lookup each of the variables we are binning in our tree variables
1342  const RooAbsReal *xvar = 0;
1343  const RooAbsReal *yvar = 0;
1344  const RooAbsReal *zvar = 0;
1345  switch(hdim) {
1346  case 3:
1347  zvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(2)->GetName()));
1348  assert(0 != zvar);
1349  // fall through to next case...
1350  case 2:
1351  yvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(1)->GetName()));
1352  assert(0 != yvar);
1353  // fall through to next case...
1354  case 1:
1355  xvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(0)->GetName()));
1356  assert(0 != xvar);
1357  break;
1358  default:
1359  coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: cannot fill histogram with "
1360  << hdim << " dimensions" << endl;
1361  break;
1362  }
1363 
1364  // Parse cutRange specification
1365  const auto cutVec = RooHelpers::tokenise(cutRange ? cutRange : "", ",", /*returnEmptyToken =*/ false);
1366 
1367  // Loop over events and fill the histogram
1368  if (hist->GetSumw2()->fN==0) {
1369  hist->Sumw2() ;
1370  }
1371  Int_t nevent= numEntries() ; //(Int_t)_tree->GetEntries();
1372  for(Int_t i=0; i < nevent; ++i) {
1373 
1374  //Int_t entryNumber= _tree->GetEntryNumber(i);
1375  //if (entryNumber<0) break;
1376  get(i);
1377 
1378  // Apply expression based selection criteria
1379  if (select && select->eval()==0) {
1380  continue ;
1381  }
1382 
1383 
1384  // Apply range based selection criteria
1385  Bool_t selectByRange = kTRUE ;
1386  if (cutRange) {
1387  for (const auto arg : _vars) {
1388  Bool_t selectThisArg = kFALSE ;
1389  UInt_t icut ;
1390  for (icut=0 ; icut<cutVec.size() ; icut++) {
1391  if (arg->inRange(cutVec[icut].c_str())) {
1392  selectThisArg = kTRUE ;
1393  break ;
1394  }
1395  }
1396  if (!selectThisArg) {
1397  selectByRange = kFALSE ;
1398  break ;
1399  }
1400  }
1401  }
1402 
1403  if (!selectByRange) {
1404  // Go to next event in loop over events
1405  continue ;
1406  }
1407 
1408  Int_t bin(0);
1409  switch(hdim) {
1410  case 1:
1411  bin= hist->FindBin(xvar->getVal());
1412  hist->Fill(xvar->getVal(),weight()) ;
1413  break;
1414  case 2:
1415  bin= hist->FindBin(xvar->getVal(),yvar->getVal());
1416  static_cast<TH2*>(hist)->Fill(xvar->getVal(),yvar->getVal(),weight()) ;
1417  break;
1418  case 3:
1419  bin= hist->FindBin(xvar->getVal(),yvar->getVal(),zvar->getVal());
1420  static_cast<TH3*>(hist)->Fill(xvar->getVal(),yvar->getVal(),zvar->getVal(),weight()) ;
1421  break;
1422  default:
1423  assert(hdim < 3);
1424  break;
1425  }
1426 
1427 
1428  Double_t error2 = TMath::Power(hist->GetBinError(bin),2)-TMath::Power(weight(),2) ;
1430  if (we==0) we = weight() ;
1431  error2 += TMath::Power(we,2) ;
1432 
1433 
1434 // Double_t we = weightError(RooAbsData::SumW2) ;
1435 // Double_t error2(0) ;
1436 // if (we==0) {
1437 // we = weight() ; //sqrt(weight()) ;
1438 // error2 = TMath::Power(hist->GetBinError(bin),2)-TMath::Power(weight(),2) + TMath::Power(we,2) ;
1439 // } else {
1440 // error2 = TMath::Power(hist->GetBinError(bin),2)-TMath::Power(weight(),2) + TMath::Power(we,2) ;
1441 // }
1442  //hist->AddBinContent(bin,weight());
1443  hist->SetBinError(bin,sqrt(error2)) ;
1444 
1445  //cout << "RooTreeData::fillHistogram() bin = " << bin << " weight() = " << weight() << " we = " << we << endl ;
1446 
1447  }
1448 
1449  return hist;
1450 }
1451 
1452 ////////////////////////////////////////////////////////////////////////////////
1453 /// Split dataset into subsets based on states of given splitCat in this dataset.
1454 /// A TList of RooDataSets is returned in which each RooDataSet is named
1455 /// after the state name of splitCat of which it contains the dataset subset.
1456 /// The observables splitCat itself is no longer present in the sub datasets.
1457 /// If createEmptyDataSets is kFALSE (default) this method only creates datasets for states
1458 /// which have at least one entry The caller takes ownership of the returned list and its contents
1459 
1460 TList* RooAbsData::split(const RooAbsCategory& splitCat, Bool_t createEmptyDataSets) const
1461 {
1462  // Sanity check
1463  if (!splitCat.dependsOn(*get())) {
1464  coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") ERROR category " << splitCat.GetName()
1465  << " doesn't depend on any variable in this dataset" << endl ;
1466  return 0 ;
1467  }
1468 
1469  // Clone splitting category and attach to self
1470  RooAbsCategory* cloneCat =0;
1471  RooArgSet* cloneSet = 0;
1472  if (splitCat.isDerived()) {
1473  cloneSet = (RooArgSet*) RooArgSet(splitCat).snapshot(kTRUE) ;
1474  if (!cloneSet) {
1475  coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") Couldn't deep-clone splitting category, abort." << endl ;
1476  return 0 ;
1477  }
1478  cloneCat = (RooAbsCategory*) cloneSet->find(splitCat.GetName()) ;
1479  cloneCat->attachDataSet(*this) ;
1480  } else {
1481  cloneCat = dynamic_cast<RooAbsCategory*>(get()->find(splitCat.GetName())) ;
1482  if (!cloneCat) {
1483  coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") ERROR category " << splitCat.GetName()
1484  << " is fundamental and does not appear in this dataset" << endl ;
1485  return 0 ;
1486  }
1487  }
1488 
1489  // Split a dataset in a series of subsets, each corresponding
1490  // to a state of splitCat
1491  TList* dsetList = new TList ;
1492 
1493  // Construct set of variables to be included in split sets = full set - split category
1494  RooArgSet subsetVars(*get()) ;
1495  if (splitCat.isDerived()) {
1496  RooArgSet* vars = splitCat.getVariables() ;
1497  subsetVars.remove(*vars,kTRUE,kTRUE) ;
1498  delete vars ;
1499  } else {
1500  subsetVars.remove(splitCat,kTRUE,kTRUE) ;
1501  }
1502 
1503  // Add weight variable explicitly if dataset has weights, but no top-level weight
1504  // variable exists (can happen with composite datastores)
1505  Bool_t addWV(kFALSE) ;
1506  RooRealVar newweight("weight","weight",-1e9,1e9) ;
1507  if (isWeighted() && !IsA()->InheritsFrom(RooDataHist::Class())) {
1508  subsetVars.add(newweight) ;
1509  addWV = kTRUE ;
1510  }
1511 
1512  // If createEmptyDataSets is true, prepopulate with empty sets corresponding to all states
1513  if (createEmptyDataSets) {
1514  for (const auto& nameIdx : *cloneCat) {
1515  RooAbsData* subset = emptyClone(nameIdx.first.c_str(), nameIdx.first.c_str(), &subsetVars,(addWV?"weight":0)) ;
1516  dsetList->Add((RooAbsArg*)subset) ;
1517  }
1518  }
1519 
1520 
1521  // Loop over dataset and copy event to matching subset
1522  const bool propWeightSquared = isWeighted();
1523  for (Int_t i = 0; i < numEntries(); ++i) {
1524  const RooArgSet* row = get(i);
1525  RooAbsData* subset = (RooAbsData*) dsetList->FindObject(cloneCat->getCurrentLabel());
1526  if (!subset) {
1527  subset = emptyClone(cloneCat->getCurrentLabel(),cloneCat->getCurrentLabel(),&subsetVars,(addWV?"weight":0));
1528  dsetList->Add((RooAbsArg*)subset);
1529  }
1530  if (!propWeightSquared) {
1531  subset->add(*row, weight());
1532  } else {
1533  subset->add(*row, weight(), weightSquared());
1534  }
1535  }
1536 
1537  delete cloneSet;
1538  return dsetList;
1539 }
1540 
1541 ////////////////////////////////////////////////////////////////////////////////
1542 /// Plot dataset on specified frame.
1543 ///
1544 /// By default:
1545 /// - An unbinned dataset will use the default binning of the target frame.
1546 /// - A binned dataset will retain its intrinsic binning.
1547 ///
1548 /// The following optional named arguments can be used to modify the behaviour:
1549 ///
1550 /// <table>
1551 /// <tr><th> <th> Data representation options
1552 /// <tr><td> `Asymmetry(const RooCategory& c)` <td> Show the asymmetry of the data in given two-state category [F(+)-F(-)] / [F(+)+F(-)].
1553 /// Category must have two states with indices -1 and +1 or three states with indices -1,0 and +1.
1554 /// <tr><td> `Efficiency(const RooCategory& c)` <td> Show the efficiency F(acc)/[F(acc)+F(rej)]. Category must have two states with indices 0 and 1
1555 /// <tr><td> `DataError(RooAbsData::EType)` <td> Select the type of error drawn:
1556 /// - `Auto(default)` results in Poisson for unweighted data and SumW2 for weighted data
1557 /// - `Poisson` draws asymmetric Poisson confidence intervals.
1558 /// - `SumW2` draws symmetric sum-of-weights error ( \f$ \left( \sum w \right)^2 / \sum\left(w^2\right) \f$ )
1559 /// - `None` draws no error bars
1560 /// <tr><td> `Binning(int nbins, double xlo, double xhi)` <td> Use specified binning to draw dataset
1561 /// <tr><td> `Binning(const RooAbsBinning&)` <td> Use specified binning to draw dataset
1562 /// <tr><td> `Binning(const char* name)` <td> Use binning with specified name to draw dataset
1563 /// <tr><td> `RefreshNorm(Bool_t flag)` <td> Force refreshing for PDF normalization information in frame.
1564 /// If set, any subsequent PDF will normalize to this dataset, even if it is
1565 /// not the first one added to the frame. By default only the 1st dataset
1566 /// added to a frame will update the normalization information
1567 /// <tr><td> `Rescale(Double_t f)` <td> Rescale drawn histogram by given factor.
1568 /// <tr><td> `Cut(const char*)` <td> Only plot entries that pass the given cut.
1569 /// Apart from cutting in continuous variables `Cut("x>5")`, this can also be used to plot a specific
1570 /// category state. Use something like `Cut("myCategory == myCategory::stateA")`, where
1571 /// `myCategory` resolves to the state number for a given entry and
1572 /// `myCategory::stateA` resolves to the state number of the state named "stateA".
1573 ///
1574 /// <tr><td> `CutRange(const char*)` <td> Only plot data from given range. Separate multiple ranges with ",".
1575 /// \note This often requires passing the normalisation when plotting the PDF because RooFit does not save
1576 /// how many events were being plotted (it will only work for cutting slices out of uniformly distributed variables).
1577 /// ```
1578 /// data->plotOn(frame01, CutRange("SB1"));
1579 /// const double nData = data->sumEntries("", "SB1");
1580 /// // Make clear that the target normalisation is nData. The enumerator NumEvent
1581 /// // is needed to switch between relative and absolute scaling.
1582 /// model.plotOn(frame01, Normalization(nData, RooAbsReal::NumEvent),
1583 /// ProjectionRange("SB1"));
1584 /// ```
1585 ///
1586 /// <tr><th> <th> Histogram drawing options
1587 /// <tr><td> `DrawOption(const char* opt)` <td> Select ROOT draw option for resulting TGraph object
1588 /// <tr><td> `LineStyle(Int_t style)` <td> Select line style by ROOT line style code, default is solid
1589 /// <tr><td> `LineColor(Int_t color)` <td> Select line color by ROOT color code, default is black
1590 /// <tr><td> `LineWidth(Int_t width)` <td> Select line with in pixels, default is 3
1591 /// <tr><td> `MarkerStyle(Int_t style)` <td> Select the ROOT marker style, default is 21
1592 /// <tr><td> `MarkerColor(Int_t color)` <td> Select the ROOT marker color, default is black
1593 /// <tr><td> `MarkerSize(Double_t size)` <td> Select the ROOT marker size
1594 /// <tr><td> `FillStyle(Int_t style)` <td> Select fill style, default is filled.
1595 /// <tr><td> `FillColor(Int_t color)` <td> Select fill color by ROOT color code
1596 /// <tr><td> `XErrorSize(Double_t frac)` <td> Select size of X error bar as fraction of the bin width, default is 1
1597 ///
1598 ///
1599 /// <tr><th> <th> Misc. other options
1600 /// <tr><td> `Name(const chat* name)` <td> Give curve specified name in frame. Useful if curve is to be referenced later
1601 /// <tr><td> `Invisible()` <td> Add curve to frame, but do not display. Useful in combination AddTo()
1602 /// <tr><td> `AddTo(const char* name, double_t wgtSelf, double_t wgtOther)` <td> Add constructed histogram to already existing histogram with given name and relative weight factors
1603 /// </table>
1604 
1605 RooPlot* RooAbsData::plotOn(RooPlot* frame, const RooLinkedList& argList) const
1606 {
1607  // New experimental plotOn() with varargs...
1608 
1609  // Define configuration for this method
1610  RooCmdConfig pc(Form("RooAbsData::plotOn(%s)",GetName())) ;
1611  pc.defineString("drawOption","DrawOption",0,"P") ;
1612  pc.defineString("cutRange","CutRange",0,"",kTRUE) ;
1613  pc.defineString("cutString","CutSpec",0,"") ;
1614  pc.defineString("histName","Name",0,"") ;
1615  pc.defineObject("cutVar","CutVar",0) ;
1616  pc.defineObject("binning","Binning",0) ;
1617  pc.defineString("binningName","BinningName",0,"") ;
1618  pc.defineInt("nbins","BinningSpec",0,100) ;
1619  pc.defineDouble("xlo","BinningSpec",0,0) ;
1620  pc.defineDouble("xhi","BinningSpec",1,1) ;
1621  pc.defineObject("asymCat","Asymmetry",0) ;
1622  pc.defineObject("effCat","Efficiency",0) ;
1623  pc.defineInt("lineColor","LineColor",0,-999) ;
1624  pc.defineInt("lineStyle","LineStyle",0,-999) ;
1625  pc.defineInt("lineWidth","LineWidth",0,-999) ;
1626  pc.defineInt("markerColor","MarkerColor",0,-999) ;
1627  pc.defineInt("markerStyle","MarkerStyle",0,-999) ;
1628  pc.defineDouble("markerSize","MarkerSize",0,-999) ;
1629  pc.defineInt("fillColor","FillColor",0,-999) ;
1630  pc.defineInt("fillStyle","FillStyle",0,-999) ;
1631  pc.defineInt("errorType","DataError",0,(Int_t)RooAbsData::Auto) ;
1632  pc.defineInt("histInvisible","Invisible",0,0) ;
1633  pc.defineInt("refreshFrameNorm","RefreshNorm",0,1) ;
1634  pc.defineString("addToHistName","AddTo",0,"") ;
1635  pc.defineDouble("addToWgtSelf","AddTo",0,1.) ;
1636  pc.defineDouble("addToWgtOther","AddTo",1,1.) ;
1637  pc.defineDouble("xErrorSize","XErrorSize",0,1.) ;
1638  pc.defineDouble("scaleFactor","Rescale",0,1.) ;
1639  pc.defineMutex("DataError","Asymmetry","Efficiency") ;
1640  pc.defineMutex("Binning","BinningName","BinningSpec") ;
1641 
1642  // Process & check varargs
1643  pc.process(argList) ;
1644  if (!pc.ok(kTRUE)) {
1645  return frame ;
1646  }
1647 
1648  PlotOpt o ;
1649 
1650  // Extract values from named arguments
1651  o.drawOptions = pc.getString("drawOption") ;
1652  o.cuts = pc.getString("cutString") ;
1653  if (pc.hasProcessed("Binning")) {
1654  o.bins = (RooAbsBinning*) pc.getObject("binning") ;
1655  } else if (pc.hasProcessed("BinningName")) {
1656  o.bins = &frame->getPlotVar()->getBinning(pc.getString("binningName")) ;
1657  } else if (pc.hasProcessed("BinningSpec")) {
1658  Double_t xlo = pc.getDouble("xlo") ;
1659  Double_t xhi = pc.getDouble("xhi") ;
1660  o.bins = new RooUniformBinning((xlo==xhi)?frame->getPlotVar()->getMin():xlo,
1661  (xlo==xhi)?frame->getPlotVar()->getMax():xhi,pc.getInt("nbins")) ;
1662  }
1663  const RooAbsCategoryLValue* asymCat = (const RooAbsCategoryLValue*) pc.getObject("asymCat") ;
1664  const RooAbsCategoryLValue* effCat = (const RooAbsCategoryLValue*) pc.getObject("effCat") ;
1665  o.etype = (RooAbsData::ErrorType) pc.getInt("errorType") ;
1666  o.histInvisible = pc.getInt("histInvisible") ;
1667  o.xErrorSize = pc.getDouble("xErrorSize") ;
1668  o.cutRange = pc.getString("cutRange",0,kTRUE) ;
1669  o.histName = pc.getString("histName",0,kTRUE) ;
1670  o.addToHistName = pc.getString("addToHistName",0,kTRUE) ;
1671  o.addToWgtSelf = pc.getDouble("addToWgtSelf") ;
1672  o.addToWgtOther = pc.getDouble("addToWgtOther") ;
1673  o.refreshFrameNorm = pc.getInt("refreshFrameNorm") ;
1674  o.scaleFactor = pc.getDouble("scaleFactor") ;
1675 
1676  // Map auto error type to actual type
1677  if (o.etype == Auto) {
1679  if (o.etype == SumW2) {
1680  coutI(InputArguments) << "RooAbsData::plotOn(" << GetName()
1681  << ") INFO: dataset has non-integer weights, auto-selecting SumW2 errors instead of Poisson errors" << endl ;
1682  }
1683  }
1684 
1685  if (o.addToHistName && !frame->findObject(o.addToHistName,RooHist::Class())) {
1686  coutE(InputArguments) << "RooAbsData::plotOn(" << GetName() << ") cannot find existing histogram " << o.addToHistName
1687  << " to add to in RooPlot" << endl ;
1688  return frame ;
1689  }
1690 
1691  RooPlot* ret ;
1692  if (!asymCat && !effCat) {
1693  ret = plotOn(frame,o) ;
1694  } else if (asymCat) {
1695  ret = plotAsymOn(frame,*asymCat,o) ;
1696  } else {
1697  ret = plotEffOn(frame,*effCat,o) ;
1698  }
1699 
1700  Int_t lineColor = pc.getInt("lineColor") ;
1701  Int_t lineStyle = pc.getInt("lineStyle") ;
1702  Int_t lineWidth = pc.getInt("lineWidth") ;
1703  Int_t markerColor = pc.getInt("markerColor") ;
1704  Int_t markerStyle = pc.getInt("markerStyle") ;
1705  Size_t markerSize = pc.getDouble("markerSize") ;
1706  Int_t fillColor = pc.getInt("fillColor") ;
1707  Int_t fillStyle = pc.getInt("fillStyle") ;
1708  if (lineColor!=-999) ret->getAttLine()->SetLineColor(lineColor) ;
1709  if (lineStyle!=-999) ret->getAttLine()->SetLineStyle(lineStyle) ;
1710  if (lineWidth!=-999) ret->getAttLine()->SetLineWidth(lineWidth) ;
1711  if (markerColor!=-999) ret->getAttMarker()->SetMarkerColor(markerColor) ;
1712  if (markerStyle!=-999) ret->getAttMarker()->SetMarkerStyle(markerStyle) ;
1713  if (markerSize!=-999) ret->getAttMarker()->SetMarkerSize(markerSize) ;
1714  if (fillColor!=-999) ret->getAttFill()->SetFillColor(fillColor) ;
1715  if (fillStyle!=-999) ret->getAttFill()->SetFillStyle(fillStyle) ;
1716 
1717  if (pc.hasProcessed("BinningSpec")) {
1718  delete o.bins ;
1719  }
1720 
1721  return ret ;
1722 }
1723 
1724 ////////////////////////////////////////////////////////////////////////////////
1725 /// Create and fill a histogram of the frame's variable and append it to the frame.
1726 /// The frame variable must be one of the data sets dimensions.
1727 ///
1728 /// The plot range and the number of plot bins is determined by the parameters
1729 /// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins()).
1730 ///
1731 /// The optional cut string expression can be used to select the events to be plotted.
1732 /// The cut specification may refer to any variable contained in the data set.
1733 ///
1734 /// The drawOptions are passed to the TH1::Draw() method.
1735 /// \see RooAbsData::plotOn(RooPlot*,const RooLinkedList&) const
1737 {
1738  if(0 == frame) {
1739  coutE(Plotting) << ClassName() << "::" << GetName() << ":plotOn: frame is null" << endl;
1740  return 0;
1741  }
1742  RooAbsRealLValue *var= (RooAbsRealLValue*) frame->getPlotVar();
1743  if(0 == var) {
1744  coutE(Plotting) << ClassName() << "::" << GetName()
1745  << ":plotOn: frame does not specify a plot variable" << endl;
1746  return 0;
1747  }
1748 
1749  // create and fill a temporary histogram of this variable
1750  TString histName(GetName());
1751  histName.Append("_plot");
1752  std::unique_ptr<TH1> hist;
1753  if (o.bins) {
1754  hist.reset( var->createHistogram(histName.Data(), RooFit::AxisLabel("Events"), RooFit::Binning(*o.bins)) );
1755  } else if (!frame->getPlotVar()->getBinning().isUniform()) {
1756  hist.reset( var->createHistogram(histName.Data(), RooFit::AxisLabel("Events"),
1757  RooFit::Binning(frame->getPlotVar()->getBinning())) );
1758  } else {
1759  hist.reset( var->createHistogram(histName.Data(), "Events",
1760  frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(), frame->GetNbinsX()) );
1761  }
1762 
1763  // Keep track of sum-of-weights error
1764  hist->Sumw2() ;
1765 
1766  if(0 == fillHistogram(hist.get(), RooArgList(*var),o.cuts,o.cutRange)) {
1767  coutE(Plotting) << ClassName() << "::" << GetName()
1768  << ":plotOn: fillHistogram() failed" << endl;
1769  return 0;
1770  }
1771 
1772  // If frame has no predefined bin width (event density) it will be adjusted to
1773  // our histograms bin width so we should force that bin width here
1774  Double_t nomBinWidth ;
1775  if (frame->getFitRangeNEvt()==0 && o.bins) {
1776  nomBinWidth = o.bins->averageBinWidth() ;
1777  } else {
1778  nomBinWidth = o.bins ? frame->getFitRangeBinW() : 0 ;
1779  }
1780 
1781  // convert this histogram to a RooHist object on the heap
1782  RooHist *graph= new RooHist(*hist,nomBinWidth,1,o.etype,o.xErrorSize,o.correctForBinWidth,o.scaleFactor);
1783  if(0 == graph) {
1784  coutE(Plotting) << ClassName() << "::" << GetName()
1785  << ":plotOn: unable to create a RooHist object" << endl;
1786  return 0;
1787  }
1788 
1789  // If the dataset variable has a wide range than the plot variable,
1790  // calculate the number of entries in the dataset in the plot variable fit range
1791  RooAbsRealLValue* dataVar = (RooAbsRealLValue*) _vars.find(var->GetName()) ;
1792  Double_t nEnt(sumEntries()) ;
1793  if (dataVar->getMin()<var->getMin() || dataVar->getMax()>var->getMax()) {
1794  RooAbsData* tmp = ((RooAbsData*)this)->reduce(*var) ;
1795  nEnt = tmp->sumEntries() ;
1796  delete tmp ;
1797  }
1798 
1799  // Store the number of entries before the cut, if any was made
1800  if ((o.cuts && strlen(o.cuts)) || o.cutRange) {
1801  coutI(Plotting) << "RooTreeData::plotOn: plotting " << hist->GetSumOfWeights() << " events out of " << nEnt << " total events" << endl ;
1802  graph->setRawEntries(nEnt) ;
1803  }
1804 
1805  // Add self to other hist if requested
1806  if (o.addToHistName) {
1807  RooHist* otherGraph = static_cast<RooHist*>(frame->findObject(o.addToHistName,RooHist::Class())) ;
1808 
1809  if (!graph->hasIdenticalBinning(*otherGraph)) {
1810  coutE(Plotting) << "RooTreeData::plotOn: ERROR Histogram to be added to, '" << o.addToHistName << "',has different binning" << endl ;
1811  delete graph ;
1812  return frame ;
1813  }
1814 
1815  RooHist* sumGraph = new RooHist(*graph,*otherGraph,o.addToWgtSelf,o.addToWgtOther,o.etype) ;
1816  delete graph ;
1817  graph = sumGraph ;
1818  }
1819 
1820  // Rename graph if requested
1821  if (o.histName) {
1822  graph->SetName(o.histName) ;
1823  } else {
1824  TString hname(Form("h_%s",GetName())) ;
1825  if (o.cutRange && strlen(o.cutRange)>0) {
1826  hname.Append(Form("_CutRange[%s]",o.cutRange)) ;
1827  }
1828  if (o.cuts && strlen(o.cuts)>0) {
1829  hname.Append(Form("_Cut[%s]",o.cuts)) ;
1830  }
1831  graph->SetName(hname.Data()) ;
1832  }
1833 
1834  // initialize the frame's normalization setup, if necessary
1835  frame->updateNormVars(_vars);
1836 
1837 
1838  // add the RooHist to the specified plot
1840 
1841  return frame;
1842 }
1843 
1844 ////////////////////////////////////////////////////////////////////////////////
1845 /// Create and fill a histogram with the asymmetry N[+] - N[-] / ( N[+] + N[-] ),
1846 /// where N(+/-) is the number of data points with asymCat=+1 and asymCat=-1
1847 /// as function of the frames variable. The asymmetry category 'asymCat' must
1848 /// have exactly 2 (or 3) states defined with index values +1,-1 (and 0)
1849 ///
1850 /// The plot range and the number of plot bins is determined by the parameters
1851 /// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins())
1852 ///
1853 /// The optional cut string expression can be used to select the events to be plotted.
1854 /// The cut specification may refer to any variable contained in the data set
1855 ///
1856 /// The drawOptions are passed to the TH1::Draw() method
1857 
1859 {
1860  if(0 == frame) {
1861  coutE(Plotting) << ClassName() << "::" << GetName() << ":plotAsymOn: frame is null" << endl;
1862  return 0;
1863  }
1864  RooAbsRealLValue *var= (RooAbsRealLValue*) frame->getPlotVar();
1865  if(0 == var) {
1866  coutE(Plotting) << ClassName() << "::" << GetName()
1867  << ":plotAsymOn: frame does not specify a plot variable" << endl;
1868  return 0;
1869  }
1870 
1871  // create and fill temporary histograms of this variable for each state
1872  TString hist1Name(GetName()),hist2Name(GetName());
1873  hist1Name.Append("_plot1");
1874  std::unique_ptr<TH1> hist1, hist2;
1875  hist2Name.Append("_plot2");
1876 
1877  if (o.bins) {
1878  hist1.reset( var->createHistogram(hist1Name.Data(), "Events", *o.bins) );
1879  hist2.reset( var->createHistogram(hist2Name.Data(), "Events", *o.bins) );
1880  } else {
1881  hist1.reset( var->createHistogram(hist1Name.Data(), "Events",
1882  frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
1883  frame->GetNbinsX()) );
1884  hist2.reset( var->createHistogram(hist2Name.Data(), "Events",
1885  frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
1886  frame->GetNbinsX()) );
1887  }
1888 
1889  assert(hist1 && hist2);
1890 
1891  TString cuts1,cuts2 ;
1892  if (o.cuts && strlen(o.cuts)) {
1893  cuts1 = Form("(%s)&&(%s>0)",o.cuts,asymCat.GetName());
1894  cuts2 = Form("(%s)&&(%s<0)",o.cuts,asymCat.GetName());
1895  } else {
1896  cuts1 = Form("(%s>0)",asymCat.GetName());
1897  cuts2 = Form("(%s<0)",asymCat.GetName());
1898  }
1899 
1900  if(! fillHistogram(hist1.get(), RooArgList(*var),cuts1.Data(),o.cutRange) ||
1901  ! fillHistogram(hist2.get(), RooArgList(*var),cuts2.Data(),o.cutRange)) {
1902  coutE(Plotting) << ClassName() << "::" << GetName()
1903  << ":plotAsymOn: createHistogram() failed" << endl;
1904  return 0;
1905  }
1906 
1907  // convert this histogram to a RooHist object on the heap
1908  RooHist *graph= new RooHist(*hist1,*hist2,0,1,o.etype,o.xErrorSize,kFALSE,o.scaleFactor);
1909  graph->setYAxisLabel(Form("Asymmetry in %s",asymCat.GetName())) ;
1910 
1911  // initialize the frame's normalization setup, if necessary
1912  frame->updateNormVars(_vars);
1913 
1914  // Rename graph if requested
1915  if (o.histName) {
1916  graph->SetName(o.histName) ;
1917  } else {
1918  TString hname(Form("h_%s_Asym[%s]",GetName(),asymCat.GetName())) ;
1919  if (o.cutRange && strlen(o.cutRange)>0) {
1920  hname.Append(Form("_CutRange[%s]",o.cutRange)) ;
1921  }
1922  if (o.cuts && strlen(o.cuts)>0) {
1923  hname.Append(Form("_Cut[%s]",o.cuts)) ;
1924  }
1925  graph->SetName(hname.Data()) ;
1926  }
1927 
1928  // add the RooHist to the specified plot
1930 
1931  return frame;
1932 }
1933 
1934 ////////////////////////////////////////////////////////////////////////////////
1935 /// Create and fill a histogram with the efficiency N[1] / ( N[1] + N[0] ),
1936 /// where N(1/0) is the number of data points with effCat=1 and effCat=0
1937 /// as function of the frames variable. The efficiency category 'effCat' must
1938 /// have exactly 2 +1 and 0.
1939 ///
1940 /// The plot range and the number of plot bins is determined by the parameters
1941 /// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins())
1942 ///
1943 /// The optional cut string expression can be used to select the events to be plotted.
1944 /// The cut specification may refer to any variable contained in the data set
1945 ///
1946 /// The drawOptions are passed to the TH1::Draw() method
1947 
1949 {
1950  if(0 == frame) {
1951  coutE(Plotting) << ClassName() << "::" << GetName() << ":plotEffOn: frame is null" << endl;
1952  return 0;
1953  }
1954  RooAbsRealLValue *var= (RooAbsRealLValue*) frame->getPlotVar();
1955  if(0 == var) {
1956  coutE(Plotting) << ClassName() << "::" << GetName()
1957  << ":plotEffOn: frame does not specify a plot variable" << endl;
1958  return 0;
1959  }
1960 
1961  // create and fill temporary histograms of this variable for each state
1962  TString hist1Name(GetName()),hist2Name(GetName());
1963  hist1Name.Append("_plot1");
1964  std::unique_ptr<TH1> hist1, hist2;
1965  hist2Name.Append("_plot2");
1966 
1967  if (o.bins) {
1968  hist1.reset( var->createHistogram(hist1Name.Data(), "Events", *o.bins) );
1969  hist2.reset( var->createHistogram(hist2Name.Data(), "Events", *o.bins) );
1970  } else {
1971  hist1.reset( var->createHistogram(hist1Name.Data(), "Events",
1972  frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
1973  frame->GetNbinsX()) );
1974  hist2.reset( var->createHistogram(hist2Name.Data(), "Events",
1975  frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
1976  frame->GetNbinsX()) );
1977  }
1978 
1979  assert(hist1 && hist2);
1980 
1981  TString cuts1,cuts2 ;
1982  if (o.cuts && strlen(o.cuts)) {
1983  cuts1 = Form("(%s)&&(%s==1)",o.cuts,effCat.GetName());
1984  cuts2 = Form("(%s)&&(%s==0)",o.cuts,effCat.GetName());
1985  } else {
1986  cuts1 = Form("(%s==1)",effCat.GetName());
1987  cuts2 = Form("(%s==0)",effCat.GetName());
1988  }
1989 
1990  if(! fillHistogram(hist1.get(), RooArgList(*var),cuts1.Data(),o.cutRange) ||
1991  ! fillHistogram(hist2.get(), RooArgList(*var),cuts2.Data(),o.cutRange)) {
1992  coutE(Plotting) << ClassName() << "::" << GetName()
1993  << ":plotEffOn: createHistogram() failed" << endl;
1994  return 0;
1995  }
1996 
1997  // convert this histogram to a RooHist object on the heap
1998  RooHist *graph= new RooHist(*hist1,*hist2,0,1,o.etype,o.xErrorSize,kTRUE);
1999  graph->setYAxisLabel(Form("Efficiency of %s=%s", effCat.GetName(), effCat.lookupName(1).c_str()));
2000 
2001  // initialize the frame's normalization setup, if necessary
2002  frame->updateNormVars(_vars);
2003 
2004  // Rename graph if requested
2005  if (o.histName) {
2006  graph->SetName(o.histName) ;
2007  } else {
2008  TString hname(Form("h_%s_Eff[%s]",GetName(),effCat.GetName())) ;
2009  if (o.cutRange && strlen(o.cutRange)>0) {
2010  hname.Append(Form("_CutRange[%s]",o.cutRange)) ;
2011  }
2012  if (o.cuts && strlen(o.cuts)>0) {
2013  hname.Append(Form("_Cut[%s]",o.cuts)) ;
2014  }
2015  graph->SetName(hname.Data()) ;
2016  }
2017 
2018  // add the RooHist to the specified plot
2020 
2021  return frame;
2022 }
2023 
2024 ////////////////////////////////////////////////////////////////////////////////
2025 /// Create and fill a 1-dimensional table for given category column
2026 /// This functions is the equivalent of plotOn() for category dimensions.
2027 ///
2028 /// The optional cut string expression can be used to select the events to be tabulated
2029 /// The cut specification may refer to any variable contained in the data set
2030 ///
2031 /// The option string is currently not used
2032 
2033 Roo1DTable* RooAbsData::table(const RooAbsCategory& cat, const char* cuts, const char* /*opts*/) const
2034 {
2035  // First see if var is in data set
2036  RooAbsCategory* tableVar = (RooAbsCategory*) _vars.find(cat.GetName()) ;
2037  RooArgSet *tableSet = 0;
2038  Bool_t ownPlotVar(kFALSE) ;
2039  if (!tableVar) {
2040  if (!cat.dependsOn(_vars)) {
2041  coutE(Plotting) << "RooTreeData::Table(" << GetName() << "): Argument " << cat.GetName()
2042  << " is not in dataset and is also not dependent on data set" << endl ;
2043  return 0 ;
2044  }
2045 
2046  // Clone derived variable
2047  tableSet = (RooArgSet*) RooArgSet(cat).snapshot(kTRUE) ;
2048  if (!tableSet) {
2049  coutE(Plotting) << "RooTreeData::table(" << GetName() << ") Couldn't deep-clone table category, abort." << endl ;
2050  return 0 ;
2051  }
2052  tableVar = (RooAbsCategory*) tableSet->find(cat.GetName()) ;
2053  ownPlotVar = kTRUE ;
2054 
2055  //Redirect servers of derived clone to internal ArgSet representing the data in this set
2056  tableVar->recursiveRedirectServers(_vars) ;
2057  }
2058 
2059  TString tableName(GetName()) ;
2060  if (cuts && strlen(cuts)) {
2061  tableName.Append("(") ;
2062  tableName.Append(cuts) ;
2063  tableName.Append(")") ;
2064  }
2065  Roo1DTable* table2 = tableVar->createTable(tableName) ;
2066 
2067  // Make cut selector if cut is specified
2068  RooFormulaVar* cutVar = 0;
2069  if (cuts && strlen(cuts)) {
2070  cutVar = new RooFormulaVar("cutVar",cuts,_vars) ;
2071  }
2072 
2073  // Dump contents
2074  Int_t nevent= numEntries() ;
2075  for(Int_t i=0; i < nevent; ++i) {
2076  get(i);
2077 
2078  if (cutVar && cutVar->getVal()==0) continue ;
2079 
2080  table2->fill(*tableVar,weight()) ;
2081  }
2082 
2083  if (ownPlotVar) delete tableSet ;
2084  if (cutVar) delete cutVar ;
2085 
2086  return table2 ;
2087 }
2088 
2089 ////////////////////////////////////////////////////////////////////////////////
2090 /// Fill Doubles 'lowest' and 'highest' with the lowest and highest value of
2091 /// observable 'var' in this dataset. If the return value is kTRUE and error
2092 /// occurred
2093 
2094 Bool_t RooAbsData::getRange(const RooAbsRealLValue& var, Double_t& lowest, Double_t& highest, Double_t marginFrac, Bool_t symMode) const
2095 {
2096  // Lookup variable in dataset
2097  const auto arg = _vars.find(var.GetName());
2098  if (!arg) {
2099  coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") ERROR: unknown variable: " << var.GetName() << endl ;
2100  return kTRUE;
2101  }
2102 
2103  auto varPtr = dynamic_cast<const RooRealVar*>(arg);
2104  // Check if found variable is of type RooRealVar
2105  if (!varPtr) {
2106  coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") ERROR: variable " << var.GetName() << " is not of type RooRealVar" << endl ;
2107  return kTRUE;
2108  }
2109 
2110  // Check if dataset is not empty
2111  if(sumEntries() == 0.) {
2112  coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") WARNING: empty dataset" << endl ;
2113  return kTRUE;
2114  }
2115 
2116  // Look for highest and lowest value
2117  lowest = RooNumber::infinity() ;
2118  highest = -RooNumber::infinity() ;
2119  for (Int_t i=0 ; i<numEntries() ; i++) {
2120  get(i) ;
2121  if (varPtr->getVal()<lowest) {
2122  lowest = varPtr->getVal() ;
2123  }
2124  if (varPtr->getVal()>highest) {
2125  highest = varPtr->getVal() ;
2126  }
2127  }
2128 
2129  if (marginFrac>0) {
2130  if (symMode==kFALSE) {
2131 
2132  Double_t margin = marginFrac*(highest-lowest) ;
2133  lowest -= margin ;
2134  highest += margin ;
2135  if (lowest<var.getMin()) lowest = var.getMin() ;
2136  if (highest>var.getMax()) highest = var.getMax() ;
2137 
2138  } else {
2139 
2140  Double_t mom1 = moment(*varPtr,1) ;
2141  Double_t delta = ((highest-mom1)>(mom1-lowest)?(highest-mom1):(mom1-lowest))*(1+marginFrac) ;
2142  lowest = mom1-delta ;
2143  highest = mom1+delta ;
2144  if (lowest<var.getMin()) lowest = var.getMin() ;
2145  if (highest>var.getMax()) highest = var.getMax() ;
2146 
2147  }
2148  }
2149 
2150  return kFALSE ;
2151 }
2152 
2153 ////////////////////////////////////////////////////////////////////////////////
2154 /// Prepare dataset for use with cached constant terms listed in
2155 /// 'cacheList' of expression 'arg'. Deactivate tree branches
2156 /// for any dataset observable that is either not used at all,
2157 /// or is used exclusively by cached branch nodes.
2158 
2159 void RooAbsData::optimizeReadingWithCaching(RooAbsArg& arg, const RooArgSet& cacheList, const RooArgSet& keepObsList)
2160 {
2161  RooArgSet pruneSet ;
2162 
2163  // Add unused observables in this dataset to pruneSet
2164  pruneSet.add(*get()) ;
2165  RooArgSet* usedObs = arg.getObservables(*this) ;
2166  pruneSet.remove(*usedObs,kTRUE,kTRUE) ;
2167 
2168  // Add observables exclusively used to calculate cached observables to pruneSet
2169  TIterator* vIter = get()->createIterator() ;
2170  RooAbsArg *var ;
2171  while ((var=(RooAbsArg*) vIter->Next())) {
2172  if (allClientsCached(var,cacheList)) {
2173  pruneSet.add(*var) ;
2174  }
2175  }
2176  delete vIter ;
2177 
2178 
2179  if (pruneSet.getSize()!=0) {
2180 
2181  // Go over all used observables and check if any of them have parameterized
2182  // ranges in terms of pruned observables. If so, remove those observable
2183  // from the pruning list
2184  TIterator* uIter = usedObs->createIterator() ;
2185  RooAbsArg* obs ;
2186  while((obs=(RooAbsArg*)uIter->Next())) {
2187  RooRealVar* rrv = dynamic_cast<RooRealVar*>(obs) ;
2188  if (rrv && !rrv->getBinning().isShareable()) {
2189  RooArgSet depObs ;
2190  RooAbsReal* loFunc = rrv->getBinning().lowBoundFunc() ;
2191  RooAbsReal* hiFunc = rrv->getBinning().highBoundFunc() ;
2192  if (loFunc) {
2193  loFunc->leafNodeServerList(&depObs,0,kTRUE) ;
2194  }
2195  if (hiFunc) {
2196  hiFunc->leafNodeServerList(&depObs,0,kTRUE) ;
2197  }
2198  if (depObs.getSize()>0) {
2199  pruneSet.remove(depObs,kTRUE,kTRUE) ;
2200  }
2201  }
2202  }
2203  delete uIter ;
2204  }
2205 
2206 
2207  // Remove all observables in keep list from prune list
2208  pruneSet.remove(keepObsList,kTRUE,kTRUE) ;
2209 
2210  if (pruneSet.getSize()!=0) {
2211 
2212  // Deactivate tree branches here
2213  cxcoutI(Optimization) << "RooTreeData::optimizeReadingForTestStatistic(" << GetName() << "): Observables " << pruneSet
2214  << " in dataset are either not used at all, orserving exclusively p.d.f nodes that are now cached, disabling reading of these observables for TTree" << endl ;
2215  setArgStatus(pruneSet,kFALSE) ;
2216  }
2217 
2218  delete usedObs ;
2219 
2220 }
2221 
2222 ////////////////////////////////////////////////////////////////////////////////
2223 /// Utility function that determines if all clients of object 'var'
2224 /// appear in given list of cached nodes.
2225 
2227 {
2228  Bool_t ret(kTRUE), anyClient(kFALSE) ;
2229 
2230  for (const auto client : var->valueClients()) {
2231  anyClient = kTRUE ;
2232  if (!cacheList.find(client->GetName())) {
2233  // If client is not cached recurse
2234  ret &= allClientsCached(client,cacheList) ;
2235  }
2236  }
2237 
2238  return anyClient?ret:kFALSE ;
2239 }
2240 
2241 ////////////////////////////////////////////////////////////////////////////////
2242 
2244 {
2245  _dstore->attachBuffers(extObs) ;
2246 }
2247 
2248 ////////////////////////////////////////////////////////////////////////////////
2249 
2251 {
2252  _dstore->resetBuffers() ;
2253 }
2254 
2255 ////////////////////////////////////////////////////////////////////////////////
2256 
2258 {
2259  if (_ownedComponents.size()>0) {
2260  return kTRUE ;
2261  }
2262  return kFALSE ;
2263 }
2264 
2265 ////////////////////////////////////////////////////////////////////////////////
2266 
2268 {
2269  map<string,RooAbsData*>::iterator i = _ownedComponents.find(name) ;
2270  if (i==_ownedComponents.end()) return 0 ;
2271  return i->second ;
2272 }
2273 
2274 ////////////////////////////////////////////////////////////////////////////////
2275 
2276 void RooAbsData::addOwnedComponent(const char* idxlabel, RooAbsData& data)
2277 {
2278  _ownedComponents[idxlabel]= &data ;
2279 }
2280 
2281 ////////////////////////////////////////////////////////////////////////////////
2282 /// Stream an object of class RooAbsData.
2283 
2284 void RooAbsData::Streamer(TBuffer &R__b)
2285 {
2286  if (R__b.IsReading()) {
2287  R__b.ReadClassBuffer(RooAbsData::Class(),this);
2288 
2289  // Convert on the fly to vector storage if that the current working default
2292  }
2293 
2294  } else {
2295  R__b.WriteClassBuffer(RooAbsData::Class(),this);
2296  }
2297 }
2298 
2299 ////////////////////////////////////////////////////////////////////////////////
2300 
2302 {
2303  _dstore->checkInit() ;
2304 }
2305 
2306 ////////////////////////////////////////////////////////////////////////////////
2307 /// Forward draw command to data store
2308 
2310 {
2311  if (_dstore) _dstore->Draw(option) ;
2312 }
2313 
2314 ////////////////////////////////////////////////////////////////////////////////
2315 
2317 {
2318  return _dstore->hasFilledCache() ;
2319 }
2320 
2321 ////////////////////////////////////////////////////////////////////////////////
2322 /// Return a pointer to the TTree which stores the data. Returns a nullpointer
2323 /// if vector-based storage is used. The RooAbsData remains owner of the tree.
2324 /// GetClonedTree() can be used to get a tree even if the internal storage does not use one.
2325 
2326 const TTree *RooAbsData::tree() const
2327 {
2328  if (storageType == RooAbsData::Tree) {
2329  return _dstore->tree();
2330  } else {
2331  coutW(InputArguments) << "RooAbsData::tree(" << GetName() << ") WARNING: is not of StorageType::Tree. "
2332  << "Use GetClonedTree() instead or convert to tree storage." << endl;
2333  return (TTree *)nullptr;
2334  }
2335 }
2336 
2337 ////////////////////////////////////////////////////////////////////////////////
2338 /// Return a clone of the TTree which stores the data or create such a tree
2339 /// if vector storage is used. The user is responsible for deleting the tree
2340 
2342 {
2343  if (storageType == RooAbsData::Tree) {
2344  auto tmp = const_cast<TTree *>(_dstore->tree());
2345  return tmp->CloneTree();
2346  } else {
2347  RooTreeDataStore buffer(GetName(), GetTitle(), *get(), *_dstore);
2348  return buffer.tree().CloneTree();
2349  }
2350 }
2351 
2352 ////////////////////////////////////////////////////////////////////////////////
2353 /// Convert vector-based storage to tree-based storage
2354 
2356 {
2357  if (storageType != RooAbsData::Tree) {
2359  delete _dstore;
2360  _dstore = newStore;
2362  }
2363 }
2364 
2365 ////////////////////////////////////////////////////////////////////////////////
2366 /// If one of the TObject we have a referenced to is deleted, remove the
2367 /// reference.
2368 
2370 {
2371  for(auto &iter : _ownedComponents) {
2372  if (iter.second == obj) {
2373  iter.second = nullptr;
2374  }
2375  }
2376 }
RooFormulaVar.h
RooAbsData::PlotOpt::scaleFactor
Double_t scaleFactor
Definition: RooAbsData.h:169
l
auto * l
Definition: textangle.C:4
RooAbsData::Tree
@ Tree
Definition: RooAbsData.h:248
TAttMarker::SetMarkerSize
virtual void SetMarkerSize(Size_t msize=1)
Set the marker size.
Definition: TAttMarker.h:41
Util.h
RooFormula
RooFormula internally uses ROOT's TFormula to compute user-defined expressions of RooAbsArgs.
Definition: RooFormula.h:34
RooLinkedList::find
TObject * find(const char *name) const
Return pointer to object with given name in collection.
Definition: RooLinkedList.cxx:603
RooCmdArg
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:27
RooHelpers.h
RooAbsData::printName
virtual void printName(std::ostream &os) const
Print name of dataset.
Definition: RooAbsData.cxx:775
RooRealVar::setVal
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:257
ymax
float ymax
Definition: THbookFile.cxx:95
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
RooCmdConfig.h
RooAbsData::PlotOpt
Definition: RooAbsData.h:152
RooMsgService.h
RooUniformBinning.h
RooAbsData::plotEffOn
virtual RooPlot * plotEffOn(RooPlot *frame, const RooAbsCategoryLValue &effCat, PlotOpt o) const
Create and fill a histogram with the efficiency N[1] / ( N[1] + N[0] ), where N(1/0) is the number of...
Definition: RooAbsData.cxx:1948
RooPlot::updateNormVars
void updateNormVars(const RooArgSet &vars)
Install the given set of observables are reference normalization variables for this frame.
Definition: RooPlot.cxx:380
TAttFill::SetFillStyle
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition: TAttFill.h:39
RooPlot::findObject
TObject * findObject(const char *name, const TClass *clas=0) const
Find the named object in our list of items and return a pointer to it.
Definition: RooPlot.cxx:982
RooCmdArg::getInt
Int_t getInt(Int_t idx) const
Definition: RooCmdArg.h:80
TNamed::SetName
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
Option_t
const char Option_t
Definition: RtypesCore.h:66
RooAbsData::isNonPoissonWeighted
virtual Bool_t isNonPoissonWeighted() const
Definition: RooAbsData.h:130
RooAbsData::PlotOpt::cutRange
const char * cutRange
Definition: RooAbsData.h:160
RooAbsData
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:49
RooAbsData::cacheArgs
virtual void cacheArgs(const RooAbsArg *owner, RooArgSet &varSet, const RooArgSet *nset=0, Bool_t skipZeroWeights=kFALSE)
Internal method – Cache given set of functions with data.
Definition: RooAbsData.cxx:327
RooAbsBinning::highBoundFunc
virtual RooAbsReal * highBoundFunc() const
Return pointer to RooAbsReal parameterized upper bound, if any.
Definition: RooAbsBinning.h:88
RooAbsRealLValue::getMax
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Definition: RooAbsRealLValue.h:89
RooAbsData::SumW2
@ SumW2
Definition: RooAbsData.h:99
Roo1DTable
Roo1DTable implements a one-dimensional table.
Definition: Roo1DTable.h:23
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:61
RooAbsData::weightError
virtual Double_t weightError(ErrorType etype=Poisson) const
Return error on current weight (dummy implementation returning zero)
Definition: RooAbsData.cxx:530
RooAbsBinning::lowBoundFunc
virtual RooAbsReal * lowBoundFunc() const
Return pointer to RooAbsReal parameterized lower bound, if any.
Definition: RooAbsBinning.h:84
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:798
TList::FindObject
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:578
RooTreeDataStore::tree
TTree & tree()
Definition: RooTreeDataStore.h:106
RooRealVar::format
TString * format(const RooCmdArg &formatArg) const
Format contents of RooRealVar for pretty printing on RooPlot parameter boxes.
Definition: RooRealVar.cxx:857
TString::Data
const char * Data() const
Definition: TString.h:369
TNamed::operator=
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
RooAbsData::isWeighted
virtual Bool_t isWeighted() const
Definition: RooAbsData.h:126
RooAbsData::_ownedComponents
std::map< std::string, RooAbsData * > _ownedComponents
Definition: RooAbsData.h:296
RooAbsData::meanVar
RooRealVar * meanVar(const RooRealVar &var, const char *cutSpec=0, const char *cutRange=0) const
Create a RooRealVar containing the mean of observable 'var' in this dataset.
Definition: RooAbsData.cxx:1068
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
RooPlot::getAttFill
TAttFill * getAttFill(const char *name=0) const
Return a pointer to the fill attributes of the named object in this plot, or zero if the named object...
Definition: RooPlot.cxx:876
RooAbsData::PlotOpt::refreshFrameNorm
Bool_t refreshFrameNorm
Definition: RooAbsData.h:167
RooAbsCategory::lookupName
const std::string & lookupName(value_type index) const
Get the name corresponding to the given index.
Definition: RooAbsCategory.cxx:201
RooAbsData::canSplitFast
Bool_t canSplitFast() const
Definition: RooAbsData.cxx:2257
TBuffer::ReadClassBuffer
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
RooAbsData::GetClonedTree
TTree * GetClonedTree() const
Return a clone of the TTree which stores the data or create such a tree if vector storage is used.
Definition: RooAbsData.cxx:2341
xmax
float xmax
Definition: THbookFile.cxx:95
RooAbsData::_dstore
RooAbsDataStore * _dstore
External variables cached with this data set.
Definition: RooAbsData.h:294
RooAbsData::fill
virtual void fill()
Definition: RooAbsData.cxx:297
sum
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345
RooPlot::getFitRangeNEvt
Double_t getFitRangeNEvt() const
Return the number of events in the fit range.
Definition: RooPlot.h:141
coutE
#define coutE(a)
Definition: RooMsgService.h:33
RooTrace::create
static void create(const TObject *obj)
Register creation of object 'obj'.
Definition: RooTrace.cxx:68
coutW
#define coutW(a)
Definition: RooMsgService.h:32
RooAbsDataStore::attachCache
virtual void attachCache(const RooAbsArg *newOwner, const RooArgSet &cachedVars)=0
RooAbsReal::getPlotLabel
const char * getPlotLabel() const
Get the label associated with the variable.
Definition: RooAbsReal.cxx:441
TTree::CloneTree
virtual TTree * CloneTree(Long64_t nentries=-1, Option_t *option="")
Create a clone of this tree and copy nentries.
Definition: TTree.cxx:3101
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
RooAbsData::Composite
@ Composite
Definition: RooAbsData.h:248
RooAbsData::Auto
@ Auto
Definition: RooAbsData.h:99
RooAbsData::sigma
Double_t sigma(const RooRealVar &var, const char *cutSpec=0, const char *cutRange=0) const
Definition: RooAbsData.h:212
RooAbsDataStore::resetCache
virtual void resetCache()=0
RooAbsData::weight
virtual Double_t weight() const =0
RooAbsData::RecursiveRemove
virtual void RecursiveRemove(TObject *obj)
If one of the TObject we have a referenced to is deleted, remove the reference.
Definition: RooAbsData.cxx:2369
RooAbsData::resetBuffers
void resetBuffers()
Definition: RooAbsData.cxx:2250
RooAbsCategory::createTable
Roo1DTable * createTable(const char *label) const
Create a table matching the shape of this category.
Definition: RooAbsCategory.cxx:363
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooFit::Optimization
@ Optimization
Definition: RooGlobalFunc.h:61
operator=
Binding & operator=(OUT(*fun)(void))
Definition: TRInterface_Binding.h:15
RooAbsArg::SetName
void SetName(const char *name)
Set the name of the TNamed.
Definition: RooAbsArg.cxx:2325
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:162
RooAbsBinning::averageBinWidth
virtual Double_t averageBinWidth() const =0
Int_t
int Int_t
Definition: RtypesCore.h:45
RooAbsData::getDefaultStorageType
static StorageType getDefaultStorageType()
Definition: RooAbsData.cxx:90
TH1::GetBinError
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.
Definition: TH1.cxx:8897
RooAbsData::defaultStorageType
static StorageType defaultStorageType
Definition: RooAbsData.h:256
N
#define N
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:585
RooAbsData::PlotOpt::correctForBinWidth
Bool_t correctForBinWidth
Definition: RooAbsData.h:168
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:810
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
RooCurve.h
RooAbsData::claimVars
static void claimVars(RooAbsData *)
Definition: RooAbsData.cxx:97
x
Double_t x[n]
Definition: legend1.C:17
RooAbsData::store
RooAbsDataStore * store()
Definition: RooAbsData.h:68
RooAbsReal::setPlotLabel
void setPlotLabel(const char *label)
Set the label associated with this variable.
Definition: RooAbsReal.cxx:451
RooAbsDataStore::clone
virtual RooAbsDataStore * clone(const char *newname=0) const =0
coutI
#define coutI(a)
Definition: RooMsgService.h:30
TClass.h
indent
static void indent(ostringstream &buf, int indent_level)
Definition: TClingCallFunc.cxx:87
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
TMatrixTSym< Double_t >
RooAbsData::fillHistogram
virtual TH1 * fillHistogram(TH1 *hist, const RooArgList &plotVars, const char *cuts="", const char *cutRange=0) const
Loop over columns of our tree data and fill the input histogram.
Definition: RooAbsData.cxx:1285
RooArgList::at
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:72
TBuffer
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
TArray::fN
Int_t fN
Definition: TArray.h:38
RooAbsData::releaseVars
static Bool_t releaseVars(RooAbsData *)
If return value is true variables can be deleted.
Definition: RooAbsData.cxx:106
RooAbsData::checkInit
void checkInit() const
Definition: RooAbsData.cxx:2301
TAttLine::SetLineColor
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
RooAbsArg::recursiveRedirectServers
Bool_t recursiveRedirectServers(const RooAbsCollection &newServerList, Bool_t mustReplaceAll=kFALSE, Bool_t nameChange=kFALSE, Bool_t recurseInNewSet=kTRUE)
Recursively replace all servers with the new servers in newSet.
Definition: RooAbsArg.cxx:1159
RooCmdArg::getDouble
Double_t getDouble(Int_t idx) const
Definition: RooCmdArg.h:84
RooFit::YVar
RooCmdArg YVar(const RooAbsRealLValue &var, const RooCmdArg &arg=RooCmdArg::none())
Definition: RooGlobalFunc.cxx:279
RooAbsData::printClassName
virtual void printClassName(std::ostream &os) const
Print class name of dataset.
Definition: RooAbsData.cxx:791
TTree.h
RooAbsData::PlotOpt::cuts
const char * cuts
Definition: RooAbsData.h:155
RooFit::Binning
RooCmdArg Binning(const RooAbsBinning &binning)
Definition: RooGlobalFunc.cxx:125
TString
Basic string class.
Definition: TString.h:136
RooLinkedList::FindObject
TObject * FindObject(const char *name) const
Return pointer to obejct with given name.
Definition: RooLinkedList.cxx:531
RooPrintable
RooPlotable is a 'mix-in' base class that define the standard RooFit plotting and printing methods.
Definition: RooPrintable.h:25
RooAbsBinning::isShareable
virtual Bool_t isShareable() const
If true (default), the range definition can be shared across clones of a RooRealVar.
Definition: RooAbsBinning.h:92
RooAbsDataStore::tree
virtual const TTree * tree() const
Definition: RooAbsDataStore.h:120
RooAbsData::ErrorType
ErrorType
Definition: RooAbsData.h:99
RooAbsCategory::getCurrentLabel
virtual const char * getCurrentLabel() const
Return label string of current state.
Definition: RooAbsCategory.cxx:130
TObject::InheritsFrom
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:445
RooTreeDataStore.h
RooAbsArg::valueClients
const RefCountList_t & valueClients() const
List of all value clients of this object. Value clients receive value updates.
Definition: RooAbsArg.h:189
RooAbsData::sumEntries
virtual Double_t sumEntries() const =0
Return effective number of entries in dataset, i.e., sum all weights.
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
TPaveText.h
TMatrixDSym.h
RooCompositeDataStore
RooCompositeDataStore combines several disjunct datasets into one.
Definition: RooCompositeDataStore.h:34
RooLinkedList::Replace
Bool_t Replace(const TObject *oldArg, const TObject *newArg)
Replace object 'oldArg' in collection with new object 'newArg'.
Definition: RooLinkedList.cxx:507
TH1::GetDimension
virtual Int_t GetDimension() const
Definition: TH1.h:282
RooAbsCategory
RooAbsCategory is the base class for objects that represent a discrete value with a finite number of ...
Definition: RooAbsCategory.h:37
RooCmdConfig
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
Definition: RooCmdConfig.h:27
ROOT::Math::Cephes::C
static double C[]
Definition: SpecFuncCephes.cxx:187
RooAbsData::tree
const TTree * tree() const
Return a pointer to the TTree which stores the data.
Definition: RooAbsData.cxx:2326
TString::ToUpper
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1163
RooPlot::addObject
void addObject(TObject *obj, Option_t *drawOptions="", Bool_t invisible=kFALSE)
Add a generic object to this plot.
Definition: RooPlot.cxx:422
RooAbsArg::isDerived
virtual Bool_t isDerived() const
Does value or shape of this arg depend on any other arg?
Definition: RooAbsArg.h:92
RooArgSet::snapshot
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:115
RooAbsData::standMoment
Double_t standMoment(const RooRealVar &var, Double_t order, const char *cutSpec=0, const char *cutRange=0) const
Calculate standardized moment.
Definition: RooAbsData.cxx:820
RooPlot::getFitRangeBinW
Double_t getFitRangeBinW() const
Return the bin width that is being used to normalise the PDF.
Definition: RooPlot.h:144
RooTrace.h
RooAbsDataStore::resetBuffers
virtual void resetBuffers()=0
hi
float type_of_call hi(const int &, const int &)
RooAbsData::allClientsCached
Bool_t allClientsCached(RooAbsArg *, const RooArgSet &)
Utility function that determines if all clients of object 'var' appear in given list of cached nodes.
Definition: RooAbsData.cxx:2226
RooAbsData::getRange
Bool_t getRange(const RooAbsRealLValue &var, Double_t &lowest, Double_t &highest, Double_t marginFrac=0, Bool_t symMode=kFALSE) const
Fill Doubles 'lowest' and 'highest' with the lowest and highest value of observable 'var' in this dat...
Definition: RooAbsData.cxx:2094
RooAbsData::getSimData
RooAbsData * getSimData(const char *idxstate)
Definition: RooAbsData.cxx:2267
RooAbsData::split
virtual TList * split(const RooAbsCategory &splitCat, Bool_t createEmptyDataSets=kFALSE) const
Split dataset into subsets based on states of given splitCat in this dataset.
Definition: RooAbsData.cxx:1460
RooAbsData::_vars
RooArgSet _vars
Definition: RooAbsData.h:291
RooAbsData::table
virtual Roo1DTable * table(const RooArgSet &catSet, const char *cuts="", const char *opts="") const
Construct table for product of categories in catSet.
Definition: RooAbsData.cxx:745
RooPrintable::kName
@ kName
Definition: RooPrintable.h:33
RooAbsData::PlotOpt::addToWgtSelf
Double_t addToWgtSelf
Definition: RooAbsData.h:164
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
RooFormulaVar
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooAbsData::moment
Double_t moment(const RooRealVar &var, Double_t order, const char *cutSpec=0, const char *cutRange=0) const
Calculate moment of requested order.
Definition: RooAbsData.cxx:839
RooAbsData::reduceEng
virtual RooAbsData * reduceEng(const RooArgSet &varSubset, const RooFormulaVar *cutVar, const char *cutRange=0, std::size_t nStart=0, std::size_t=std::numeric_limits< std::size_t >::max(), Bool_t copyCache=kTRUE)=0
TBuffer.h
RooFit::AxisLabel
RooCmdArg AxisLabel(const char *name)
Definition: RooGlobalFunc.cxx:281
TAxis::GetXmin
Double_t GetXmin() const
Definition: TAxis.h:133
RooAbsBinning::isUniform
virtual Bool_t isUniform() const
Definition: RooAbsBinning.h:48
RooPlot::getAttLine
TAttLine * getAttLine(const char *name=0) const
Return a pointer to the line attributes of the named object in this plot, or zero if the named object...
Definition: RooPlot.cxx:866
xmin
float xmin
Definition: THbookFile.cxx:95
TH1::GetSumw2
virtual TArrayD * GetSumw2()
Definition: TH1.h:312
RooMultiCategory
RooMultiCategory connects several RooAbsCategory objects into a single category.
Definition: RooMultiCategory.h:28
TH3
The 3-D histogram classes derived from the 1-D histogram classes.
Definition: TH3.h:31
RooMultiCategory.h
RooAbsData::Vector
@ Vector
Definition: RooAbsData.h:248
Double_t
RooAbsData::optimizeReadingWithCaching
virtual void optimizeReadingWithCaching(RooAbsArg &arg, const RooArgSet &cacheList, const RooArgSet &keepObsList)
Prepare dataset for use with cached constant terms listed in 'cacheList' of expression 'arg'.
Definition: RooAbsData.cxx:2159
TAttMarker::SetMarkerColor
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
Definition: TAttMarker.h:38
RooAbsData::plotOn
virtual RooPlot * plotOn(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
Definition: RooAbsData.cxx:544
RooAbsCollection::addClone
virtual RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE)
Add a clone of the specified argument to list.
Definition: RooAbsCollection.cxx:428
RooFit::ZVar
RooCmdArg ZVar(const RooAbsRealLValue &var, const RooCmdArg &arg=RooCmdArg::none())
Definition: RooGlobalFunc.cxx:280
TNamed
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
TBuffer::WriteClassBuffer
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
RooAbsData::setDirtyProp
void setDirtyProp(Bool_t flag)
Control propagation of dirty flags from observables in dataset.
Definition: RooAbsData.cxx:359
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
RooRealVar::setError
void setError(Double_t value)
Definition: RooRealVar.h:62
RooDataHist.h
RooAbsData::PlotOpt::xErrorSize
Double_t xErrorSize
Definition: RooAbsData.h:166
RooAbsData::get
virtual const RooArgSet * get() const
Definition: RooAbsData.h:92
RooAbsData::storageType
StorageType storageType
Definition: RooAbsData.h:258
TNamed::Clone
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition: TNamed.cxx:74
RooAbsData::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:622
RooFit::Plotting
@ Plotting
Definition: RooGlobalFunc.h:60
RooLinkedList
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:37
RooPlot.h
RooAbsCollection::createIterator
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:178
TMath::Power
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:735
RooAbsDataStore::changeObservableName
virtual Bool_t changeObservableName(const char *from, const char *to)=0
RooAbsBinning
RooAbsBinning is the abstract base class for RooRealVar binning definitions.
Definition: RooAbsBinning.h:26
RooAbsData::PlotOpt::histInvisible
Bool_t histInvisible
Definition: RooAbsData.h:162
RooAbsData::setDefaultStorageType
static void setDefaultStorageType(StorageType s)
Definition: RooAbsData.cxx:79
RooUniformBinning
RooUniformBinning is an implementation of RooAbsBinning that provides a uniform binning in 'n' bins b...
Definition: RooUniformBinning.h:23
RooAbsData::rmsVar
RooRealVar * rmsVar(const RooRealVar &var, const char *cutSpec=0, const char *cutRange=0) const
Create a RooRealVar containing the RMS of observable 'var' in this dataset.
Definition: RooAbsData.cxx:1104
TH1::Fill
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3350
RooPrintable::kArgs
@ kArgs
Definition: RooPrintable.h:33
TH2
Service class for 2-Dim histogram classes.
Definition: TH2.h:30
RooPlot
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:44
what
static const char * what
Definition: stlLoader.cc:6
RooAbsData::convertToVectorStore
void convertToVectorStore()
Convert tree-based storage to vector-based storage.
Definition: RooAbsData.cxx:272
RooAbsCollection::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:455
Roo1DTable::fill
virtual void fill(RooAbsCategory &cat, Double_t weight=1.0)
Increment the counter of the table slot with the name corresponding to that of the current category s...
Definition: Roo1DTable.cxx:103
RooAbsData::numEntries
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Definition: RooAbsData.cxx:304
RooCategory.h
RooLinkedList::Add
virtual void Add(TObject *arg)
Definition: RooLinkedList.h:64
y
Double_t y[n]
Definition: legend1.C:17
RooPrintable::kValue
@ kValue
Definition: RooPrintable.h:33
sqrt
double sqrt(double)
RooAbsArg::attachDataSet
void attachDataSet(const RooAbsData &set)
Replace server nodes with names matching the dataset variable names with those data set variables,...
Definition: RooAbsArg.cxx:1553
RooAbsData::corrcovMatrix
TMatrixDSym * corrcovMatrix(const RooArgList &vars, const char *cutSpec, const char *cutRange, Bool_t corr) const
Return covariance matrix from data for given list of observables.
Definition: RooAbsData.cxx:972
RooAbsData::operator=
RooAbsData & operator=(const RooAbsData &other)
Definition: RooAbsData.cxx:208
RooRealVar.h
RooAbsDataStore::get
virtual const RooArgSet * get(Int_t index) const =0
TNamed::SetTitle
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
RooHist
A RooHist is a graphical representation of binned data based on the TGraphAsymmErrors class.
Definition: RooHist.h:27
TH2.h
RooAbsCollection::releaseOwnership
void releaseOwnership()
Definition: RooAbsCollection.h:299
RooTreeDataStore
RooTreeDataStore is a TTree-backed data storage.
Definition: RooTreeDataStore.h:32
RooAbsDataStore.h
RooAbsData::Poisson
@ Poisson
Definition: RooAbsData.h:99
TIterator::Next
virtual TObject * Next()=0
RooAbsData::hasFilledCache
Bool_t hasFilledCache() const
Definition: RooAbsData.cxx:2316
RooAbsData::PlotOpt::bins
RooAbsBinning * bins
Definition: RooAbsData.h:158
unsigned int
RooHist.h
TH3.h
ymin
float ymin
Definition: THbookFile.cxx:95
RooAbsData::addOwnedComponent
void addOwnedComponent(const char *idxlabel, RooAbsData &data)
Definition: RooAbsData.cxx:2276
RooAbsData::dataRealVar
RooRealVar * dataRealVar(const char *methodname, const RooRealVar &extVar) const
Internal method to check if given RooRealVar maps to a RooRealVar in this dataset.
Definition: RooAbsData.cxx:897
_dcc
static std::map< RooAbsData *, int > _dcc
Definition: RooAbsData.cxx:71
RooGlobalFunc.h
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:67
RooRealVar::getBinning
const RooAbsBinning & getBinning(const char *name=0, Bool_t verbose=kTRUE, Bool_t createOnTheFly=kFALSE) const
Return binning definition with name.
Definition: RooRealVar.cxx:318
TIterator::Reset
virtual void Reset()=0
RooAbsDataStore::numEntries
virtual Int_t numEntries() const =0
TAttLine::SetLineWidth
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition: TAttLine.h:43
TBuffer::IsReading
Bool_t IsReading() const
Definition: TBuffer.h:86
RooAbsData::emptyClone
virtual RooAbsData * emptyClone(const char *newName=0, const char *newTitle=0, const RooArgSet *vars=0, const char *wgtVarName=0) const =0
RooNumber::infinity
static Double_t infinity()
Return internal infinity representation.
Definition: RooNumber.cxx:49
TObject::Draw
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition: TObject.cxx:197
RooAbsDataStore::checkInit
virtual void checkInit() const
Definition: RooAbsDataStore.h:116
RooPrintable::kClassName
@ kClassName
Definition: RooPrintable.h:33
box
Definition: box.py:1
RooAbsData.h
RooAbsData::PlotOpt::addToHistName
const char * addToHistName
Definition: RooAbsData.h:163
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:295
RooAbsDataStore::hasFilledCache
virtual Bool_t hasFilledCache() const
Definition: RooAbsDataStore.h:118
RooAbsArg::leafNodeServerList
void leafNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=0, Bool_t recurseNonDerived=kFALSE) const
Fill supplied list with all leaf nodes of the arg tree, starting with ourself as top node.
Definition: RooAbsArg.cxx:509
RooAbsData::Draw
virtual void Draw(Option_t *option="")
Forward draw command to data store.
Definition: RooAbsData.cxx:2309
RooAbsData::plotAsymOn
virtual RooPlot * plotAsymOn(RooPlot *frame, const RooAbsCategoryLValue &asymCat, PlotOpt o) const
Create and fill a histogram with the asymmetry N[+] - N[-] / ( N[+] + N[-] ), where N(+/-) is the num...
Definition: RooAbsData.cxx:1858
RooFit::AutoBinning
RooCmdArg AutoBinning(Int_t nbins=100, Double_t marginFactor=0.1)
Definition: RooGlobalFunc.cxx:289
TH1::SetBinError
virtual void SetBinError(Int_t bin, Double_t error)
Set the bin Error Note that this resets the bin eror option to be of Normal Type and for the non-empt...
Definition: TH1.cxx:9040
Roo1DTable.h
RooCategory
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:27
RooAbsData::~RooAbsData
virtual ~RooAbsData()
Destructor.
Definition: RooAbsData.cxx:250
TH1::FindBin
virtual Int_t FindBin(Double_t x, Double_t y=0, Double_t z=0)
Return Global bin number corresponding to x,y,z.
Definition: TH1.cxx:3680
RooAbsData::reduce
RooAbsData * reduce(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg())
Create a reduced copy of this dataset.
Definition: RooAbsData.cxx:379
RooAbsData::weightSquared
virtual Double_t weightSquared() const =0
TAttFill::SetFillColor
virtual void SetFillColor(Color_t fcolor)
Set the fill area color.
Definition: TAttFill.h:37
graph
Definition: graph.py:1
RooVectorDataStore
RooVectorDataStore uses std::vectors to store data columns.
Definition: RooVectorDataStore.h:37
TList::Add
virtual void Add(TObject *obj)
Definition: TList.h:87
RooAbsData::PlotOpt::addToWgtOther
Double_t addToWgtOther
Definition: RooAbsData.h:165
RooAbsData::changeObservableName
virtual Bool_t changeObservableName(const char *from, const char *to)
Definition: RooAbsData.cxx:284
TH1::Sumw2
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition: TH1.cxx:8854
TObject
Mother of all ROOT objects.
Definition: TObject.h:37
RooAbsData::resetCache
virtual void resetCache()
Internal method – Remove cached function values.
Definition: RooAbsData.cxx:335
TH1
TH1 is the base class of all histogram classes in ROOT.
Definition: TH1.h:58
RooAbsData::RooAbsData
RooAbsData()
Default constructor.
Definition: RooAbsData.cxx:119
RooTrace::destroy
static void destroy(const TObject *obj)
Register deletion of object 'obj'.
Definition: RooTrace.cxx:81
name
char name[80]
Definition: TGX11.cxx:110
RooAbsDataStore::printMultiline
void printMultiline(std::ostream &os, Int_t content, Bool_t verbose, TString indent) const
Detailed printing interface.
Definition: RooAbsDataStore.cxx:190
RooAbsRealLValue::setConstant
void setConstant(Bool_t value=kTRUE)
Definition: RooAbsRealLValue.h:115
RooAbsDataStore::setDirtyProp
virtual void setDirtyProp(Bool_t flag)
Definition: RooAbsDataStore.h:113
RooAbsDataStore::reset
virtual void reset()=0
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:799
RooAbsData::StorageType
StorageType
Definition: RooAbsData.h:248
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
RooAbsData::PlotOpt::etype
RooAbsData::ErrorType etype
Definition: RooAbsData.h:159
RooPlot::getAttMarker
TAttMarker * getAttMarker(const char *name=0) const
Return a pointer to the marker attributes of the named object in this plot, or zero if the named obje...
Definition: RooPlot.cxx:886
TPaveText
A Pave (see TPave) with text, lines or/and boxes inside.
Definition: TPaveText.h:21
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooAbsData::attachBuffers
void attachBuffers(const RooArgSet &extObs)
Definition: RooAbsData.cxx:2243
RooAbsDataStore::cacheArgs
virtual void cacheArgs(const RooAbsArg *cacheOwner, RooArgSet &varSet, const RooArgSet *nset=0, Bool_t skipZeroWeights=kFALSE)=0
TAttMarker::SetMarkerStyle
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition: TAttMarker.h:40
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooAbsData::convertToTreeStore
virtual void convertToTreeStore()
Convert vector-based storage to tree-based storage.
Definition: RooAbsData.cxx:2355
RooAbsArg::getVariables
RooArgSet * getVariables(Bool_t stripDisconnected=kTRUE) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2008
RooVectorDataStore.h
RooAbsData::add
virtual void add(const RooArgSet &row, Double_t weight=1, Double_t weightError=0)=0
TAxis::GetXmax
Double_t GetXmax() const
Definition: TAxis.h:134
RooAbsDataStore
RooAbsDataStore is the abstract base class for data collection that use a TTree as internal storage m...
Definition: RooAbsDataStore.h:34
RooPlot::GetNbinsX
Int_t GetNbinsX() const
Definition: RooPlot.cxx:1274
RooAbsDataStore::fill
virtual Int_t fill()=0
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:344
RooAbsCategoryLValue
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
Definition: RooAbsCategoryLValue.h:25
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:1366
Class
void Class()
Definition: Class.C:29
TMatrixDSym
TMatrixTSym< Double_t > TMatrixDSym
Definition: TMatrixDSymfwd.h:22
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
ROOT::Math::KahanSum< double >
RooAbsData::statOn
virtual RooPlot * statOn(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())
Add a box with statistics information to the specified frame.
Definition: RooAbsData.cxx:1161
RooAbsData::corrcov
Double_t corrcov(const RooRealVar &x, const RooRealVar &y, const char *cutSpec, const char *cutRange, Bool_t corr) const
Internal method to calculate single correlation and covariance elements.
Definition: RooAbsData.cxx:916
RooAbsCollection::removeAll
virtual void removeAll()
Remove all arguments from our set, deleting them if we own them.
Definition: RooAbsCollection.cxx:644
TAttLine::SetLineStyle
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition: TAttLine.h:42
RooFormula::eval
Double_t eval(const RooArgSet *nset=0) const
Evalute all parameters/observables, and then evaluate formula.
Definition: RooFormula.cxx:342
RooAbsDataStore::setArgStatus
virtual void setArgStatus(const RooArgSet &set, Bool_t active)=0
RooAbsRealLValue
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Definition: RooAbsRealLValue.h:31
TH1.h
TObject::ClassName
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:130
cxcoutI
#define cxcoutI(a)
Definition: RooMsgService.h:85
RooAbsRealLValue.h
RooPlot::getPlotVar
RooAbsRealLValue * getPlotVar() const
Definition: RooPlot.h:139
RooAbsData::PlotOpt::histName
const char * histName
Definition: RooAbsData.h:161
TObject::Clear
virtual void Clear(Option_t *="")
Definition: TObject.h:115
RooCompositeDataStore.h
RooAbsRealLValue::getBinning
virtual const RooAbsBinning & getBinning(const char *name=0, Bool_t verbose=kTRUE, Bool_t createOnTheFly=kFALSE) const =0
Retrive binning configuration with given name or default binning.
RooLinkedList::Delete
void Delete(Option_t *o=0)
Remove all elements in collection and delete all elements NB: Collection does not own elements,...
Definition: RooLinkedList.cxx:575
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:231
RooAbsData::defaultPrintContents
virtual Int_t defaultPrintContents(Option_t *opt) const
Define default print options, for a given print style.
Definition: RooAbsData.cxx:806
TList
A doubly linked list.
Definition: TList.h:44
RooAbsData::_cachedVars
RooArgSet _cachedVars
Definition: RooAbsData.h:292
TMath.h
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:33
fw3dlego::xbins
const double xbins[xbins_n]
Definition: collection_proxies.C:48
RooPlot::GetXaxis
TAxis * GetXaxis() const
Definition: RooPlot.cxx:1270
RooAbsData::PlotOpt::drawOptions
Option_t * drawOptions
Definition: RooAbsData.h:157
int
RooAbsData::setArgStatus
virtual void setArgStatus(const RooArgSet &set, Bool_t active)
Definition: RooAbsData.cxx:351
Size_t
float Size_t
Definition: RtypesCore.h:96
RooAbsData::reset
virtual void reset()
Definition: RooAbsData.cxx:311
RooAbsData::printTitle
virtual void printTitle(std::ostream &os) const
Print title of dataset.
Definition: RooAbsData.cxx:783
RooAbsRealLValue::createHistogram
TH1 * createHistogram(const char *name, 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
Definition: RooAbsRealLValue.cxx:524
RooAbsDataStore::attachBuffers
virtual void attachBuffers(const RooArgSet &extObs)=0