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