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