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