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