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