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 "TMath.h"
84#include "TTree.h"
85
86#include "RooFormulaVar.h"
87#include "RooCmdConfig.h"
88#include "RooAbsRealLValue.h"
89#include "RooMsgService.h"
90#include "RooMultiCategory.h"
91#include "Roo1DTable.h"
92#include "RooAbsDataStore.h"
93#include "RooVectorDataStore.h"
94#include "RooTreeDataStore.h"
95#include "RooDataHist.h"
97#include "RooCategory.h"
98#include "RooTrace.h"
99#include "RooUniformBinning.h"
100#include "RooSimultaneous.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) ;
176
177 RooTrace::create(this) ;
178}
179
180////////////////////////////////////////////////////////////////////////////////
181/// Constructor from a set of variables. Only fundamental elements of vars
182/// (RooRealVar,RooCategory etc) are stored as part of the dataset
183
185 TNamed(name,title),
186 _vars("Dataset Variables"),
187 _cachedVars("Cached Variables"),
188 _dstore(dstore)
189{
190 if (dynamic_cast<RooTreeDataStore *>(dstore)) {
192 } else if (dynamic_cast<RooVectorDataStore *>(dstore)) {
194 } else {
196 }
197 // cout << "created dataset " << this << endl ;
198 claimVars(this);
199
200 // clone the fundamentals of the given data set into internal buffer
201 for (const auto var : vars) {
202 if (!var->isFundamental()) {
203 coutE(InputArguments) << "RooAbsDataStore::initialize(" << GetName()
204 << "): Data set cannot contain non-fundamental types, ignoring " << var->GetName()
205 << endl;
206 throw std::invalid_argument(std::string("Only fundamental variables can be placed into datasets. This is violated for ") + var->GetName());
207 } else {
208 _vars.addClone(*var);
209 }
210 }
211
212 // reconnect any parameterized ranges to internal dataset observables
213 for (auto var : _vars) {
214 var->attachArgs(_vars);
215 }
216
218
219 RooTrace::create(this);
220}
221
222////////////////////////////////////////////////////////////////////////////////
223/// Copy constructor
224
225RooAbsData::RooAbsData(const RooAbsData& other, const char* newname) :
226 TNamed(newname ? newname : other.GetName(),other.GetTitle()),
227 RooPrintable(other), _vars(),
228 _cachedVars("Cached Variables"),
229 _namePtr(newname ? RooNameReg::instance().constPtr(newname) : other._namePtr)
230{
231 //cout << "created dataset " << this << endl ;
232 claimVars(this) ;
233 _vars.addClone(other._vars) ;
234
235 // reconnect any parameterized ranges to internal dataset observables
236 for (auto var : _vars) {
237 var->attachArgs(_vars);
238 }
239
240
241 if (other._ownedComponents.size()>0) {
242
243 // copy owned components here
244
245 map<string,RooAbsDataStore*> smap ;
246 for (auto& itero : other._ownedComponents) {
247 RooAbsData* dclone = (RooAbsData*) itero.second->Clone();
248 _ownedComponents[itero.first] = dclone;
249 smap[itero.first] = dclone->store();
250 }
251
253 _dstore = std::make_unique<RooCompositeDataStore>(newname?newname:other.GetName(),other.GetTitle(),_vars,*idx,smap) ;
255
256 } else {
257
258 // Convert to vector store if default is vector
259 _dstore.reset(other._dstore->clone(_vars,newname?newname:other.GetName()));
260 storageType = other.storageType;
261 }
262
264
265 RooTrace::create(this) ;
266}
267
269 TNamed::operator=(other);
271
272 claimVars(this);
273 _vars.Clear();
274 _vars.addClone(other._vars);
275 _namePtr = other._namePtr;
276
277 // reconnect any parameterized ranges to internal dataset observables
278 for (const auto var : _vars) {
279 var->attachDataSet(*this) ;
280 }
281
282
283 if (other._ownedComponents.size()>0) {
284
285 // copy owned components here
286
287 map<string,RooAbsDataStore*> smap ;
288 for (auto& itero : other._ownedComponents) {
289 RooAbsData* dclone = (RooAbsData*) itero.second->Clone();
290 _ownedComponents[itero.first] = dclone;
291 smap[itero.first] = dclone->store();
292 }
293
295 _dstore = std::make_unique<RooCompositeDataStore>(GetName(), GetTitle(), _vars, *idx, smap);
297
298 } else {
299
300 // Convert to vector store if default is vector
301 _dstore.reset(other._dstore->clone(_vars));
302 storageType = other.storageType;
303 }
304
306
307 return *this;
308}
309
310
312 if (other._globalObservables) {
313 if(_globalObservables == nullptr) _globalObservables = std::make_unique<RooArgSet>();
314 else _globalObservables->clear();
315 other._globalObservables->snapshot(*_globalObservables);
316 } else {
317 _globalObservables.reset(nullptr);
318 }
319}
320
321
322////////////////////////////////////////////////////////////////////////////////
323/// Destructor
324
326{
327 if (releaseVars(this)) {
328 // will cause content to be deleted subsequently in dtor
329 } else {
331 }
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}
341////////////////////////////////////////////////////////////////////////////////
342/// Convert tree-based storage to vector-based storage
343
345{
346 if (auto treeStore = dynamic_cast<RooTreeDataStore*>(_dstore.get())) {
347 _dstore = std::make_unique<RooVectorDataStore>(*treeStore, _vars, GetName());
349 }
350}
351
352////////////////////////////////////////////////////////////////////////////////
353
354bool RooAbsData::changeObservableName(const char* from, const char* to)
355{
356 bool ret = _dstore->changeObservableName(from,to) ;
357
358 RooAbsArg* tmp = _vars.find(from) ;
359 if (tmp) {
360 tmp->SetName(to) ;
361 }
362 return ret ;
363}
364
365////////////////////////////////////////////////////////////////////////////////
366
368{
369 _dstore->fill() ;
370}
371
372////////////////////////////////////////////////////////////////////////////////
373
375{
376 return nullptr != _dstore ? _dstore->numEntries() : 0;
377}
378
379////////////////////////////////////////////////////////////////////////////////
380
382{
383 _dstore->reset() ;
384}
385
386////////////////////////////////////////////////////////////////////////////////
387
389{
390 checkInit() ;
391 return _dstore->get(index) ;
392}
393
394////////////////////////////////////////////////////////////////////////////////
395/// Internal method -- Cache given set of functions with data
396
397void RooAbsData::cacheArgs(const RooAbsArg* cacheOwner, RooArgSet& varSet, const RooArgSet* nset, bool skipZeroWeights)
398{
399 _dstore->cacheArgs(cacheOwner,varSet,nset,skipZeroWeights) ;
400}
401
402////////////////////////////////////////////////////////////////////////////////
403/// Internal method -- Remove cached function values
404
406{
407 _dstore->resetCache() ;
409}
410
411////////////////////////////////////////////////////////////////////////////////
412/// Internal method -- Attach dataset copied with cache contents to copied instances of functions
413
414void RooAbsData::attachCache(const RooAbsArg* newOwner, const RooArgSet& cachedVars)
415{
416 _dstore->attachCache(newOwner, cachedVars) ;
417}
418
419////////////////////////////////////////////////////////////////////////////////
420
421void RooAbsData::setArgStatus(const RooArgSet& set, bool active)
422{
423 _dstore->setArgStatus(set,active) ;
424}
425
426////////////////////////////////////////////////////////////////////////////////
427/// Control propagation of dirty flags from observables in dataset
428
430{
431 _dstore->setDirtyProp(flag) ;
432}
433
434////////////////////////////////////////////////////////////////////////////////
435/// Create a reduced copy of this dataset. The caller takes ownership of the returned dataset
436///
437/// The following optional named arguments are accepted
438/// <table>
439/// <tr><td> `SelectVars(const RooArgSet& vars)` <td> Only retain the listed observables in the output dataset
440/// <tr><td> `Cut(const char* expression)` <td> Only retain event surviving the given cut expression
441/// <tr><td> `Cut(const RooFormulaVar& expr)` <td> Only retain event surviving the given cut formula
442/// <tr><td> `CutRange(const char* name)` <td> Only retain events inside range with given name. Multiple CutRange
443/// arguments may be given to select multiple ranges
444/// <tr><td> `EventRange(int lo, int hi)` <td> Only retain events with given sequential event numbers
445/// <tr><td> `Name(const char* name)` <td> Give specified name to output dataset
446/// <tr><td> `Title(const char* name)` <td> Give specified title to output dataset
447/// </table>
448
449RooAbsData* RooAbsData::reduce(const RooCmdArg& arg1,const RooCmdArg& arg2,const RooCmdArg& arg3,const RooCmdArg& arg4,
450 const RooCmdArg& arg5,const RooCmdArg& arg6,const RooCmdArg& arg7,const RooCmdArg& arg8)
451{
452 // Define configuration for this method
453 RooCmdConfig pc(Form("RooAbsData::reduce(%s)",GetName())) ;
454 pc.defineString("name","Name",0,"") ;
455 pc.defineString("title","Title",0,"") ;
456 pc.defineString("cutRange","CutRange",0,"") ;
457 pc.defineString("cutSpec","CutSpec",0,"") ;
458 pc.defineObject("cutVar","CutVar",0,0) ;
459 pc.defineInt("evtStart","EventRange",0,0) ;
460 pc.defineInt("evtStop","EventRange",1,std::numeric_limits<int>::max()) ;
461 pc.defineSet("varSel","SelectVars",0,0) ;
462 pc.defineMutex("CutVar","CutSpec") ;
463
464 // Process & check varargs
465 pc.process(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8) ;
466 if (!pc.ok(true)) {
467 return nullptr;
468 }
469
470 // Extract values from named arguments
471 const char* cutRange = pc.getString("cutRange",0,true) ;
472 const char* cutSpec = pc.getString("cutSpec",0,true) ;
473 RooFormulaVar* cutVar = static_cast<RooFormulaVar*>(pc.getObject("cutVar",0)) ;
474 Int_t nStart = pc.getInt("evtStart",0) ;
475 Int_t nStop = pc.getInt("evtStop",std::numeric_limits<int>::max()) ;
476 RooArgSet* varSet = pc.getSet("varSel");
477 const char* name = pc.getString("name",0,true) ;
478 const char* title = pc.getString("title",0,true) ;
479
480 // Make sure varSubset doesn't contain any variable not in this dataset
481 RooArgSet varSubset ;
482 if (varSet) {
483 varSubset.add(*varSet) ;
484 for (const auto arg : varSubset) {
485 if (!_vars.find(arg->GetName())) {
486 coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
487 << arg->GetName() << " not in dataset, ignored" << endl ;
488 varSubset.remove(*arg) ;
489 }
490 }
491 } else {
492 varSubset.add(*get()) ;
493 }
494
495 RooAbsData* ret = nullptr;
496 if (cutSpec) {
497
498 RooFormulaVar cutVarTmp(cutSpec,cutSpec,*get()) ;
499 ret = reduceEng(varSubset,&cutVarTmp,cutRange,nStart,nStop) ;
500
501 } else if (cutVar) {
502
503 ret = reduceEng(varSubset,cutVar,cutRange,nStart,nStop) ;
504
505 } else {
506
507 ret = reduceEng(varSubset,0,cutRange,nStart,nStop) ;
508
509 }
510
511 if (!ret) return nullptr;
512
513 if (name) ret->SetName(name) ;
514 if (title) ret->SetTitle(title) ;
515
516 ret->copyGlobalObservables(*this);
517 return ret ;
518}
519
520////////////////////////////////////////////////////////////////////////////////
521/// Create a subset of the data set by applying the given cut on the data points.
522/// The cut expression can refer to any variable in the data set. For cuts involving
523/// other variables, such as intermediate formula objects, use the equivalent
524/// reduce method specifying the as a RooFormulVar reference.
525
527{
528 RooFormulaVar cutVar(cut,cut,*get()) ;
529 RooAbsData* ret = reduceEng(*get(),&cutVar,0,0,std::numeric_limits<std::size_t>::max()) ;
530 ret->copyGlobalObservables(*this);
531 return ret;
532}
533
534////////////////////////////////////////////////////////////////////////////////
535/// Create a subset of the data set by applying the given cut on the data points.
536/// The 'cutVar' formula variable is used to select the subset of data points to be
537/// retained in the reduced data collection.
538
540{
541 RooAbsData* ret = reduceEng(*get(),&cutVar,0,0,std::numeric_limits<std::size_t>::max()) ;
542 ret->copyGlobalObservables(*this);
543 return ret;
544}
545
546////////////////////////////////////////////////////////////////////////////////
547/// Create a subset of the data set by applying the given cut on the data points
548/// and reducing the dimensions to the specified set.
549///
550/// The cut expression can refer to any variable in the data set. For cuts involving
551/// other variables, such as intermediate formula objects, use the equivalent
552/// reduce method specifying the as a RooFormulVar reference.
553
554RooAbsData* RooAbsData::reduce(const RooArgSet& varSubset, const char* cut)
555{
556 // Make sure varSubset doesn't contain any variable not in this dataset
557 RooArgSet varSubset2(varSubset) ;
558 for (const auto arg : varSubset) {
559 if (!_vars.find(arg->GetName())) {
560 coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
561 << arg->GetName() << " not in dataset, ignored" << endl ;
562 varSubset2.remove(*arg) ;
563 }
564 }
565
566 RooAbsData* ret = nullptr;
567 if (cut && strlen(cut)>0) {
568 RooFormulaVar cutVar(cut, cut, *get(), false);
569 ret = reduceEng(varSubset2,&cutVar,0,0,std::numeric_limits<std::size_t>::max());
570 } else {
571 ret = reduceEng(varSubset2,0,0,0,std::numeric_limits<std::size_t>::max());
572 }
573 ret->copyGlobalObservables(*this);
574 return ret;
575}
576
577////////////////////////////////////////////////////////////////////////////////
578/// Create a subset of the data set by applying the given cut on the data points
579/// and reducing the dimensions to the specified set.
580///
581/// The 'cutVar' formula variable is used to select the subset of data points to be
582/// retained in the reduced data collection.
583
584RooAbsData* RooAbsData::reduce(const RooArgSet& varSubset, const RooFormulaVar& cutVar)
585{
586 // Make sure varSubset doesn't contain any variable not in this dataset
587 RooArgSet varSubset2(varSubset) ;
588 for(RooAbsArg * arg : varSubset) {
589 if (!_vars.find(arg->GetName())) {
590 coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
591 << arg->GetName() << " not in dataset, ignored" << endl ;
592 varSubset2.remove(*arg) ;
593 }
594 }
595
596 RooAbsData* ret = reduceEng(varSubset2,&cutVar,0,0,std::numeric_limits<std::size_t>::max()) ;
597 ret->copyGlobalObservables(*this);
598 return ret;
599}
600
601
602RooPlot* RooAbsData::plotOn(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
603 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
604 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8) const
605{
607 l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
608 l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
609 l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
610 l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
611 return plotOn(frame,l) ;
612}
613
614////////////////////////////////////////////////////////////////////////////////
615/// Create and fill a ROOT histogram TH1,TH2 or TH3 with the values of this dataset for the variables with given names
616/// The range of each observable that is histogrammed is always automatically calculated from the distribution in
617/// the dataset. The number of bins can be controlled using the [xyz]bins parameters. For a greater degree of control
618/// use the createHistogram() method below with named arguments
619///
620/// The caller takes ownership of the returned histogram
621///
622/// \deprecated This function is deprecated and will be removed in ROOT v6.30.
623/// Use the overload
624/// RooAbsData::createHistogram(const char*, const RooCmdArg&, const RooCmdArg&, const RooCmdArg&)
625/// that takes RooFit command arguments.
626
627TH1 *RooAbsData::createHistogram(const char* varNameList, Int_t xbins, Int_t ybins, Int_t zbins) const
628{
629 coutW(DataHandling) << "'RooAbsData::createHistogram' is deprecated and will be removed in ROOT v6.30: "
630 << "Use the overload of 'RooAbsData::createHistogram' that takes RooFit command arguments."
631 << std::endl;
632
633 // Parse list of variable names
634 const auto varNames = ROOT::Split(varNameList, ",:");
635 RooLinkedList argList;
636 RooRealVar* vars[3] = {nullptr, nullptr, nullptr};
637
638 for (unsigned int i = 0; i < varNames.size(); ++i) {
639 if (i >= 3) {
640 coutW(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): Can only create 3-dimensional histograms. Variable "
641 << i << " " << varNames[i] << " unused." << std::endl;
642 continue;
643 }
644
645 vars[i] = static_cast<RooRealVar*>( get()->find(varNames[i].data()) );
646 if (!vars[i]) {
647 coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << ") ERROR: dataset does not contain an observable named " << varNames[i] << std::endl;
648 return nullptr;
649 }
650 }
651
652 if (!vars[0]) {
653 coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): No variable to be histogrammed in list '" << varNameList << "'" << std::endl;
654 return nullptr;
655 }
656
657 if (xbins<=0 || !vars[0]->hasMax() || !vars[0]->hasMin() ) {
658 argList.Add(RooFit::AutoBinning(xbins==0?vars[0]->numBins():abs(xbins)).Clone()) ;
659 } else {
660 argList.Add(RooFit::Binning(xbins).Clone()) ;
661 }
662
663 if (vars[1]) {
664 if (ybins<=0 || !vars[1]->hasMax() || !vars[1]->hasMin() ) {
665 argList.Add(RooFit::YVar(*vars[1],RooFit::AutoBinning(ybins==0?vars[1]->numBins():abs(ybins))).Clone()) ;
666 } else {
667 argList.Add(RooFit::YVar(*vars[1],RooFit::Binning(ybins)).Clone()) ;
668 }
669 }
670 if (vars[2]) {
671 if (zbins<=0 || !vars[2]->hasMax() || !vars[2]->hasMin() ) {
672 argList.Add(RooFit::ZVar(*vars[2],RooFit::AutoBinning(zbins==0?vars[2]->numBins():abs(zbins))).Clone()) ;
673 } else {
674 argList.Add(RooFit::ZVar(*vars[2],RooFit::Binning(zbins)).Clone()) ;
675 }
676 }
677
678
679 // Call implementation function
680 TH1* result = createHistogram(GetName(), *vars[0], argList);
681
682 // Delete temporary list of RooCmdArgs
683 argList.Delete() ;
684
685 return result ;
686}
687
688
690 const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3, const RooCmdArg& arg4,
691 const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8) const
692{
694 l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
695 l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
696 l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
697 l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
698
699 return createHistogram(name,xvar,l) ;
700}
701
702////////////////////////////////////////////////////////////////////////////////
703/// Create and fill a ROOT histogram TH1,TH2 or TH3 with the values of this
704/// dataset for the variables with given names.
705///
706/// \param[in] varNameList Comma-separated variable names.
707/// \param[in] binArgX Control the binning for the `x` variable.
708/// \param[in] binArgY Control the binning for the `y` variable.
709/// \param[in] binArgZ Control the binning for the `z` variable.
710/// \return Histogram now owned by user.
711///
712/// The possible binning command arguments for each axis are:
713///
714/// <table>
715/// <tr><td> `AutoBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin, set binning to nbins
716/// <tr><td> `AutoSymBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin,
717/// with additional constraint that mean of data is in center of range, set binning to nbins
718/// <tr><td> `Binning(const char* name)` <td> Apply binning with given name to x axis of histogram
719/// <tr><td> `Binning(RooAbsBinning& binning)` <td> Apply specified binning to x axis of histogram
720/// <tr><td> `Binning(int nbins, double lo, double hi)` <td> Apply specified binning to x axis of histogram
721///
722/// <tr><td> `YVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on y axis of ROOT histogram
723/// <tr><td> `ZVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on z axis of ROOT histogram
724/// </table>
725
726TH1 *RooAbsData::createHistogram(const char* varNameList,
727 const RooCmdArg& binArgX,
728 const RooCmdArg& binArgY,
729 const RooCmdArg& binArgZ) const
730{
731 // Parse list of variable names
732 const auto varNames = ROOT::Split(varNameList, ",:");
733 RooRealVar* vars[3] = {nullptr, nullptr, nullptr};
734
735 for (unsigned int i = 0; i < varNames.size(); ++i) {
736 if (i >= 3) {
737 coutW(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): Can only create 3-dimensional histograms. Variable "
738 << i << " " << varNames[i] << " unused." << std::endl;
739 continue;
740 }
741
742 vars[i] = static_cast<RooRealVar*>( get()->find(varNames[i].data()) );
743 if (!vars[i]) {
744 coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << ") ERROR: dataset does not contain an observable named " << varNames[i] << std::endl;
745 return nullptr;
746 }
747 }
748
749 if (!vars[0]) {
750 coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): No variable to be histogrammed in list '" << varNameList << "'" << std::endl;
751 return nullptr;
752 }
753
754 // Fill command argument list
755 RooLinkedList argList;
756 argList.Add(binArgX.Clone());
757 if (vars[1]) {
758 argList.Add(RooFit::YVar(*vars[1],binArgY).Clone());
759 }
760 if (vars[2]) {
761 argList.Add(RooFit::ZVar(*vars[2],binArgZ).Clone());
762 }
763
764 // Call implementation function
765 TH1* result = createHistogram(GetName(), *vars[0], argList);
766
767 // Delete temporary list of RooCmdArgs
768 argList.Delete() ;
769
770 return result ;
771}
772
773////////////////////////////////////////////////////////////////////////////////
774///
775/// This function accepts the following arguments
776///
777/// \param[in] name Name of the ROOT histogram
778/// \param[in] xvar Observable to be mapped on x axis of ROOT histogram
779/// \param[in] argListIn list of input arguments
780/// \return Histogram now owned by user.
781///
782/// <table>
783/// <tr><td> `AutoBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin, set binning to nbins
784/// <tr><td> `AutoSymBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin,
785/// with additional constraint that mean of data is in center of range, set binning to nbins
786/// <tr><td> `Binning(const char* name)` <td> Apply binning with given name to x axis of histogram
787/// <tr><td> `Binning(RooAbsBinning& binning)` <td> Apply specified binning to x axis of histogram
788/// <tr><td> `Binning(int nbins, double lo, double hi)` <td> Apply specified binning to x axis of histogram
789///
790/// <tr><td> `YVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on y axis of ROOT histogram
791/// <tr><td> `ZVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on z axis of ROOT histogram
792/// </table>
793///
794/// 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.
795/// ```
796/// createHistogram("histo",x,Binning(-1,1,20), YVar(y,Binning(-1,1,30)), ZVar(z,Binning("zbinning")))
797/// ```
798///
799/// The caller takes ownership of the returned histogram
800
801TH1 *RooAbsData::createHistogram(const char *name, const RooAbsRealLValue& xvar, const RooLinkedList& argListIn) const
802{
803 RooLinkedList argList(argListIn) ;
804
805 // Define configuration for this method
806 RooCmdConfig pc(Form("RooAbsData::createHistogram(%s)",GetName())) ;
807 pc.defineString("cutRange","CutRange",0,"",true) ;
808 pc.defineString("cutString","CutSpec",0,"") ;
809 pc.defineObject("yvar","YVar",0,0) ;
810 pc.defineObject("zvar","ZVar",0,0) ;
811 pc.allowUndefined() ;
812
813 // Process & check varargs
814 pc.process(argList) ;
815 if (!pc.ok(true)) {
816 return nullptr;
817 }
818
819 const char* cutSpec = pc.getString("cutString",0,true) ;
820 const char* cutRange = pc.getString("cutRange",0,true) ;
821
822 RooArgList vars(xvar) ;
823 RooAbsArg* yvar = static_cast<RooAbsArg*>(pc.getObject("yvar")) ;
824 if (yvar) {
825 vars.add(*yvar) ;
826 }
827 RooAbsArg* zvar = static_cast<RooAbsArg*>(pc.getObject("zvar")) ;
828 if (zvar) {
829 vars.add(*zvar) ;
830 }
831
832 RooCmdConfig::stripCmdList(argList,"CutRange,CutSpec") ;
833
834 // Swap Auto(Sym)RangeData with a Binning command
835 RooLinkedList ownedCmds ;
836 RooCmdArg* autoRD = (RooCmdArg*) argList.find("AutoRangeData") ;
837 if (autoRD) {
838 double xmin,xmax ;
839 if (!getRange((RooRealVar&)xvar,xmin,xmax,autoRD->getDouble(0),autoRD->getInt(0))) {
840 RooCmdArg* bincmd = (RooCmdArg*) RooFit::Binning(autoRD->getInt(1),xmin,xmax).Clone() ;
841 ownedCmds.Add(bincmd) ;
842 argList.Replace(autoRD,bincmd) ;
843 }
844 }
845
846 if (yvar) {
847 RooCmdArg* autoRDY = (RooCmdArg*) ((RooCmdArg*)argList.find("YVar"))->subArgs().find("AutoRangeData") ;
848 if (autoRDY) {
849 double ymin,ymax ;
850 if (!getRange((RooRealVar&)(*yvar),ymin,ymax,autoRDY->getDouble(0),autoRDY->getInt(0))) {
851 RooCmdArg* bincmd = (RooCmdArg*) RooFit::Binning(autoRDY->getInt(1),ymin,ymax).Clone() ;
852 //ownedCmds.Add(bincmd) ;
853 ((RooCmdArg*)argList.find("YVar"))->subArgs().Replace(autoRDY,bincmd) ;
854 }
855 delete autoRDY ;
856 }
857 }
858
859 if (zvar) {
860 RooCmdArg* autoRDZ = (RooCmdArg*) ((RooCmdArg*)argList.find("ZVar"))->subArgs().find("AutoRangeData") ;
861 if (autoRDZ) {
862 double zmin,zmax ;
863 if (!getRange((RooRealVar&)(*zvar),zmin,zmax,autoRDZ->getDouble(0),autoRDZ->getInt(0))) {
864 RooCmdArg* bincmd = (RooCmdArg*) RooFit::Binning(autoRDZ->getInt(1),zmin,zmax).Clone() ;
865 //ownedCmds.Add(bincmd) ;
866 ((RooCmdArg*)argList.find("ZVar"))->subArgs().Replace(autoRDZ,bincmd) ;
867 }
868 delete autoRDZ ;
869 }
870 }
871
872
873 TH1* histo = xvar.createHistogram(name,argList) ;
874 fillHistogram(histo,vars,cutSpec,cutRange) ;
875
876 ownedCmds.Delete() ;
877
878 return histo ;
879}
880
881////////////////////////////////////////////////////////////////////////////////
882/// Construct table for product of categories in catSet
883
884Roo1DTable* RooAbsData::table(const RooArgSet& catSet, const char* cuts, const char* opts) const
885{
886 RooArgSet catSet2 ;
887
888 string prodName("(") ;
889 for(auto * arg : catSet) {
890 if (dynamic_cast<RooAbsCategory*>(arg)) {
891 RooAbsCategory* varsArg = dynamic_cast<RooAbsCategory*>(_vars.find(arg->GetName())) ;
892 if (varsArg != 0) catSet2.add(*varsArg) ;
893 else catSet2.add(*arg) ;
894 if (prodName.length()>1) {
895 prodName += " x " ;
896 }
897 prodName += arg->GetName() ;
898 } else {
899 coutW(InputArguments) << "RooAbsData::table(" << GetName() << ") non-RooAbsCategory input argument " << arg->GetName() << " ignored" << endl ;
900 }
901 }
902 prodName += ")" ;
903
904 RooMultiCategory tmp(prodName.c_str(),prodName.c_str(),catSet2) ;
905 return table(tmp,cuts,opts) ;
906}
907
908////////////////////////////////////////////////////////////////////////////////
909/// Print name of dataset
910
911void RooAbsData::printName(ostream& os) const
912{
913 os << GetName() ;
914}
915
916////////////////////////////////////////////////////////////////////////////////
917/// Print title of dataset
918
919void RooAbsData::printTitle(ostream& os) const
920{
921 os << GetTitle() ;
922}
923
924////////////////////////////////////////////////////////////////////////////////
925/// Print class name of dataset
926
927void RooAbsData::printClassName(ostream& os) const
928{
929 os << ClassName() ;
930}
931
932////////////////////////////////////////////////////////////////////////////////
933
934void RooAbsData::printMultiline(ostream& os, Int_t contents, bool verbose, TString indent) const
935{
936 _dstore->printMultiline(os,contents,verbose,indent) ;
937}
938
939////////////////////////////////////////////////////////////////////////////////
940/// Define default print options, for a given print style
941
943{
945}
946
947////////////////////////////////////////////////////////////////////////////////
948/// Calculate standardized moment.
949///
950/// \param[in] var Variable to be used for calculating the moment.
951/// \param[in] order Order of the moment.
952/// \param[in] cutSpec If specified, the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
953/// \param[in] cutRange If specified, calculate inside the range named 'cutRange' (also applies cut spec)
954/// \return \f$ \frac{\left< \left( X - \left< X \right> \right)^n \right>}{\sigma^n} \f$, where n = order.
955
956double RooAbsData::standMoment(const RooRealVar &var, double order, const char* cutSpec, const char* cutRange) const
957{
958 // Hardwire invariant answer for first and second moment
959 if (order==1) return 0 ;
960 if (order==2) return 1 ;
961
962 return moment(var,order,cutSpec,cutRange) / TMath::Power(sigma(var,cutSpec,cutRange),order) ;
963}
964
965////////////////////////////////////////////////////////////////////////////////
966/// Calculate moment of requested order.
967///
968/// \param[in] var Variable to be used for calculating the moment.
969/// \param[in] order Order of the moment.
970/// \param[in] cutSpec If specified, the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
971/// \param[in] cutRange If specified, calculate inside the range named 'cutRange' (also applies cut spec)
972/// \return \f$ \left< \left( X - \left< X \right> \right)^n \right> \f$ of order \f$n\f$.
973///
974
975double RooAbsData::moment(const RooRealVar& var, double order, const char* cutSpec, const char* cutRange) const
976{
977 double offset = order>1 ? moment(var,1,cutSpec,cutRange) : 0 ;
978 return moment(var,order,offset,cutSpec,cutRange) ;
979
980}
981
982////////////////////////////////////////////////////////////////////////////////
983/// Return the 'order'-ed moment of observable 'var' in this dataset. If offset is non-zero it is subtracted
984/// from the values of 'var' prior to the moment calculation. If cutSpec and/or cutRange are specified
985/// the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
986/// and/or are inside the range named 'cutRange'
987
988double RooAbsData::moment(const RooRealVar& var, double order, double offset, const char* cutSpec, const char* cutRange) const
989{
990 // Lookup variable in dataset
991 auto arg = _vars.find(var.GetName());
992 if (!arg) {
993 coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") ERROR: unknown variable: " << var.GetName() << std::endl;
994 return 0;
995 }
996
997 auto varPtr = dynamic_cast<const RooRealVar*>(arg);
998 // Check if found variable is of type RooRealVar
999 if (!varPtr) {
1000 coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") ERROR: variable " << var.GetName() << " is not of type RooRealVar" << endl ;
1001 return 0;
1002 }
1003
1004 // Check if dataset is not empty
1005 if(sumEntries(cutSpec, cutRange) == 0.) {
1006 coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") WARNING: empty dataset" << endl ;
1007 return 0;
1008 }
1009
1010 // Setup RooFormulaVar for cutSpec if it is present
1011 std::unique_ptr<RooFormula> select;
1012 if (cutSpec) {
1013 select.reset(new RooFormula("select",cutSpec,*get()));
1014 }
1015
1016
1017 // Calculate requested moment
1019 for(Int_t index= 0; index < numEntries(); index++) {
1020 const RooArgSet* vars = get(index) ;
1021 if (select && select->eval()==0) continue ;
1022 if (cutRange && vars->allInRange(cutRange)) continue ;
1023
1024 sum += weight() * TMath::Power(varPtr->getVal() - offset,order);
1025 }
1026
1027 return sum/sumEntries(cutSpec, cutRange);
1028}
1029
1030////////////////////////////////////////////////////////////////////////////////
1031/// Internal method to check if given RooRealVar maps to a RooRealVar in this dataset
1032
1033RooRealVar* RooAbsData::dataRealVar(const char* methodname, const RooRealVar& extVar) const
1034{
1035 // Lookup variable in dataset
1036 RooRealVar *xdata = (RooRealVar*) _vars.find(extVar.GetName());
1037 if(!xdata) {
1038 coutE(InputArguments) << "RooDataSet::" << methodname << "(" << GetName() << ") ERROR: variable : " << extVar.GetName() << " is not in data" << endl ;
1039 return nullptr;
1040 }
1041 // Check if found variable is of type RooRealVar
1042 if (!dynamic_cast<RooRealVar*>(xdata)) {
1043 coutE(InputArguments) << "RooDataSet::" << methodname << "(" << GetName() << ") ERROR: variable : " << extVar.GetName() << " is not of type RooRealVar in data" << endl ;
1044 return nullptr;
1045 }
1046 return xdata;
1047}
1048
1049////////////////////////////////////////////////////////////////////////////////
1050/// Internal method to calculate single correlation and covariance elements
1051
1052double RooAbsData::corrcov(const RooRealVar &x, const RooRealVar &y, const char* cutSpec, const char* cutRange, bool corr) const
1053{
1054 // Lookup variable in dataset
1055 RooRealVar *xdata = dataRealVar(corr?"correlation":"covariance",x) ;
1056 RooRealVar *ydata = dataRealVar(corr?"correlation":"covariance",y) ;
1057 if (!xdata||!ydata) return 0 ;
1058
1059 // Check if dataset is not empty
1060 if(sumEntries(cutSpec, cutRange) == 0.) {
1061 coutW(InputArguments) << "RooDataSet::" << (corr?"correlation":"covariance") << "(" << GetName() << ") WARNING: empty dataset, returning zero" << endl ;
1062 return 0;
1063 }
1064
1065 // Setup RooFormulaVar for cutSpec if it is present
1066 RooFormula* select = cutSpec ? new RooFormula("select",cutSpec,*get()) : 0 ;
1067
1068 // Calculate requested moment
1069 double xysum(0),xsum(0),ysum(0),x2sum(0),y2sum(0);
1070 const RooArgSet* vars ;
1071 for(Int_t index= 0; index < numEntries(); index++) {
1072 vars = get(index) ;
1073 if (select && select->eval()==0) continue ;
1074 if (cutRange && vars->allInRange(cutRange)) continue ;
1075
1076 xysum += weight()*xdata->getVal()*ydata->getVal() ;
1077 xsum += weight()*xdata->getVal() ;
1078 ysum += weight()*ydata->getVal() ;
1079 if (corr) {
1080 x2sum += weight()*xdata->getVal()*xdata->getVal() ;
1081 y2sum += weight()*ydata->getVal()*ydata->getVal() ;
1082 }
1083 }
1084
1085 // Normalize entries
1086 xysum/=sumEntries(cutSpec, cutRange) ;
1087 xsum/=sumEntries(cutSpec, cutRange) ;
1088 ysum/=sumEntries(cutSpec, cutRange) ;
1089 if (corr) {
1090 x2sum/=sumEntries(cutSpec, cutRange) ;
1091 y2sum/=sumEntries(cutSpec, cutRange) ;
1092 }
1093
1094 // Cleanup
1095 if (select) delete select ;
1096
1097 // Return covariance or correlation as requested
1098 if (corr) {
1099 return (xysum-xsum*ysum)/(sqrt(x2sum-(xsum*xsum))*sqrt(y2sum-(ysum*ysum))) ;
1100 } else {
1101 return (xysum-xsum*ysum);
1102 }
1103}
1104
1105////////////////////////////////////////////////////////////////////////////////
1106/// Return covariance matrix from data for given list of observables
1107
1108TMatrixDSym* RooAbsData::corrcovMatrix(const RooArgList& vars, const char* cutSpec, const char* cutRange, bool corr) const
1109{
1110 RooArgList varList ;
1111 for(auto * var : static_range_cast<RooRealVar*>(vars)) {
1112 RooRealVar* datavar = dataRealVar("covarianceMatrix",*var) ;
1113 if (!datavar) {
1114 return nullptr;
1115 }
1116 varList.add(*datavar) ;
1117 }
1118
1119
1120 // Check if dataset is not empty
1121 if(sumEntries(cutSpec, cutRange) == 0.) {
1122 coutW(InputArguments) << "RooDataSet::covariance(" << GetName() << ") WARNING: empty dataset, returning zero" << endl ;
1123 return nullptr;
1124 }
1125
1126 // Setup RooFormulaVar for cutSpec if it is present
1127 std::unique_ptr<RooFormula> select = cutSpec ? std::make_unique<RooFormula>("select",cutSpec,*get()) : nullptr;
1128
1129 TMatrixDSym xysum(varList.size()) ;
1130 std::vector<double> xsum(varList.size()) ;
1131 std::vector<double> x2sum(varList.size()) ;
1132
1133 // Calculate <x_i> and <x_i y_j>
1134 for(Int_t index= 0; index < numEntries(); index++) {
1135 const RooArgSet* dvars = get(index) ;
1136 if (select && select->eval()==0) continue ;
1137 if (cutRange && dvars->allInRange(cutRange)) continue ;
1138
1139 for(std::size_t ix = 0; ix < varList.size(); ++ix) {
1140 auto varx = static_cast<RooRealVar const&>(varList[ix]);
1141 xsum[ix] += weight() * varx.getVal() ;
1142 if (corr) {
1143 x2sum[ix] += weight() * varx.getVal() * varx.getVal();
1144 }
1145
1146 for(std::size_t iy = ix; iy < varList.size(); ++iy) {
1147 auto vary = static_cast<RooRealVar const&>(varList[iy]);
1148 xysum(ix,iy) += weight() * varx.getVal() * vary.getVal();
1149 xysum(iy,ix) = xysum(ix,iy) ;
1150 }
1151 }
1152
1153 }
1154
1155 // Normalize sums
1156 for (std::size_t ix=0 ; ix<varList.size() ; ix++) {
1157 xsum[ix] /= sumEntries(cutSpec, cutRange) ;
1158 if (corr) {
1159 x2sum[ix] /= sumEntries(cutSpec, cutRange) ;
1160 }
1161 for (std::size_t iy=0 ; iy<varList.size() ; iy++) {
1162 xysum(ix,iy) /= sumEntries(cutSpec, cutRange) ;
1163 }
1164 }
1165
1166 // Calculate covariance matrix
1167 TMatrixDSym* C = new TMatrixDSym(varList.size()) ;
1168 for (std::size_t ix=0 ; ix<varList.size() ; ix++) {
1169 for (std::size_t iy=0 ; iy<varList.size() ; iy++) {
1170 (*C)(ix,iy) = xysum(ix,iy)-xsum[ix]*xsum[iy] ;
1171 if (corr) {
1172 (*C)(ix,iy) /= std::sqrt((x2sum[ix]-(xsum[ix]*xsum[ix]))*(x2sum[iy]-(xsum[iy]*xsum[iy]))) ;
1173 }
1174 }
1175 }
1176
1177 return C ;
1178}
1179
1180////////////////////////////////////////////////////////////////////////////////
1181/// Create a RooRealVar containing the mean of observable 'var' in
1182/// this dataset. If cutSpec and/or cutRange are specified the
1183/// moment is calculated on the subset of the data which pass the C++
1184/// cut specification expression 'cutSpec' and/or are inside the
1185/// range named 'cutRange'
1186
1187RooRealVar* RooAbsData::meanVar(const RooRealVar &var, const char* cutSpec, const char* cutRange) const
1188{
1189 // Create a new variable with appropriate strings. The error is calculated as
1190 // RMS/Sqrt(N) which is generally valid.
1191
1192 // Create holder variable for mean
1193 std::string name = std::string{var.GetName()} + "Mean";
1194 std::string title = std::string{"Mean of "} + var.GetTitle();
1195 auto *meanv= new RooRealVar(name.c_str(), title.c_str(), 0) ;
1196 meanv->setConstant(false) ;
1197
1198 // Adjust plot label
1199 std::string label = "<" + std::string{var.getPlotLabel()} + ">";
1200 meanv->setPlotLabel(label.c_str());
1201
1202 // fill in this variable's value and error
1203 double meanVal=moment(var,1,0,cutSpec,cutRange) ;
1204 double N(sumEntries(cutSpec,cutRange)) ;
1205
1206 double rmsVal= sqrt(moment(var,2,meanVal,cutSpec,cutRange)*N/(N-1));
1207 meanv->setVal(meanVal) ;
1208 meanv->setError(N > 0 ? rmsVal/sqrt(N) : 0);
1209
1210 return meanv;
1211}
1212
1213////////////////////////////////////////////////////////////////////////////////
1214/// Create a RooRealVar containing the RMS of observable 'var' in
1215/// this dataset. If cutSpec and/or cutRange are specified the
1216/// moment is calculated on the subset of the data which pass the C++
1217/// cut specification expression 'cutSpec' and/or are inside the
1218/// range named 'cutRange'
1219
1220RooRealVar* RooAbsData::rmsVar(const RooRealVar &var, const char* cutSpec, const char* cutRange) const
1221{
1222 // Create a new variable with appropriate strings. The error is calculated as
1223 // RMS/(2*Sqrt(N)) which is only valid if the variable has a Gaussian distribution.
1224
1225 // Create RMS value holder
1226 std::string name(var.GetName()),title("RMS of ") ;
1227 name += "RMS";
1228 title += var.GetTitle();
1229 auto *rms= new RooRealVar(name.c_str(), title.c_str(), 0) ;
1230 rms->setConstant(false) ;
1231
1232 // Adjust plot label
1233 std::string label(var.getPlotLabel());
1234 label += "_{RMS}";
1235 rms->setPlotLabel(label.c_str());
1236
1237 // Fill in this variable's value and error
1238 double meanVal(moment(var,1,0,cutSpec,cutRange)) ;
1239 double N(sumEntries(cutSpec, cutRange));
1240 double rmsVal= sqrt(moment(var,2,meanVal,cutSpec,cutRange)*N/(N-1));
1241 rms->setVal(rmsVal) ;
1242 rms->setError(rmsVal/sqrt(2*N));
1243
1244 return rms;
1245}
1246
1247////////////////////////////////////////////////////////////////////////////////
1248/// Add a box with statistics information to the specified frame. By default a box with the
1249/// event count, mean and rms of the plotted variable is added.
1250///
1251/// The following optional named arguments are accepted
1252/// <table>
1253/// <tr><td> `What(const char* whatstr)` <td> Controls what is printed: "N" = count, "M" is mean, "R" is RMS.
1254/// <tr><td> `Format(const char* optStr)` <td> \deprecated Classing parameter formatting options, provided for backward compatibility
1255///
1256/// <tr><td> `Format(const char* what,...)` <td> Parameter formatting options.
1257/// <table>
1258/// <tr><td> const char* what <td> Controls what is shown:
1259/// - "N" adds name
1260/// - "E" adds error
1261/// - "A" shows asymmetric error
1262/// - "U" shows unit
1263/// - "H" hides the value
1264/// <tr><td> `FixedPrecision(int n)` <td> Controls precision, set fixed number of digits
1265/// <tr><td> `AutoPrecision(int n)` <td> Controls precision. Number of shown digits is calculated from error + n specified additional digits (1 is sensible default)
1266/// <tr><td> `VerbatimName(bool flag)` <td> Put variable name in a \\verb+ + clause.
1267/// </table>
1268/// <tr><td> `Label(const chat* label)` <td> Add header label to parameter box
1269/// <tr><td> `Layout(double xmin, double xmax, double ymax)` <td> Specify relative position of left,right side of box and top of box. Position of
1270/// bottom of box is calculated automatically from number lines in box
1271/// <tr><td> `Cut(const char* expression)` <td> Apply given cut expression to data when calculating statistics
1272/// <tr><td> `CutRange(const char* rangeName)` <td> Only consider events within given range when calculating statistics. Multiple
1273/// CutRange() argument may be specified to combine ranges.
1274///
1275/// </table>
1276
1277RooPlot* RooAbsData::statOn(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
1278 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
1279 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
1280{
1281 // Stuff all arguments in a list
1282 RooLinkedList cmdList;
1283 cmdList.Add(const_cast<RooCmdArg*>(&arg1)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg2)) ;
1284 cmdList.Add(const_cast<RooCmdArg*>(&arg3)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg4)) ;
1285 cmdList.Add(const_cast<RooCmdArg*>(&arg5)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg6)) ;
1286 cmdList.Add(const_cast<RooCmdArg*>(&arg7)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg8)) ;
1287
1288 // Select the pdf-specific commands
1289 RooCmdConfig pc(Form("RooTreeData::statOn(%s)",GetName())) ;
1290 pc.defineString("what","What",0,"MNR") ;
1291 pc.defineString("label","Label",0,"") ;
1292 pc.defineDouble("xmin","Layout",0,0.65) ;
1293 pc.defineDouble("xmax","Layout",1,0.99) ;
1294 pc.defineInt("ymaxi","Layout",0,Int_t(0.95*10000)) ;
1295 pc.defineString("formatStr","Format",0,"NELU") ;
1296 pc.defineInt("sigDigit","Format",0,2) ;
1297 pc.defineInt("dummy","FormatArgs",0,0) ;
1298 pc.defineString("cutRange","CutRange",0,"",true) ;
1299 pc.defineString("cutString","CutSpec",0,"") ;
1300 pc.defineMutex("Format","FormatArgs") ;
1301
1302 // Process and check varargs
1303 pc.process(cmdList) ;
1304 if (!pc.ok(true)) {
1305 return frame ;
1306 }
1307
1308 const char* label = pc.getString("label") ;
1309 double xmin = pc.getDouble("xmin") ;
1310 double xmax = pc.getDouble("xmax") ;
1311 double ymax = pc.getInt("ymaxi") / 10000. ;
1312 const char* formatStr = pc.getString("formatStr") ;
1313 Int_t sigDigit = pc.getInt("sigDigit") ;
1314 const char* what = pc.getString("what") ;
1315
1316 const char* cutSpec = pc.getString("cutString",0,true) ;
1317 const char* cutRange = pc.getString("cutRange",0,true) ;
1318
1319 if (pc.hasProcessed("FormatArgs")) {
1320 RooCmdArg* formatCmd = static_cast<RooCmdArg*>(cmdList.FindObject("FormatArgs")) ;
1321 return statOn(frame,what,label,0,0,xmin,xmax,ymax,cutSpec,cutRange,formatCmd) ;
1322 } else {
1323 return statOn(frame,what,label,sigDigit,formatStr,xmin,xmax,ymax,cutSpec,cutRange) ;
1324 }
1325}
1326
1327////////////////////////////////////////////////////////////////////////////////
1328/// Implementation back-end of statOn() method with named arguments
1329
1330RooPlot* RooAbsData::statOn(RooPlot* frame, const char* what, const char *label, Int_t sigDigits,
1331 Option_t *options, double xmin, double xmax, double ymax,
1332 const char* cutSpec, const char* cutRange, const RooCmdArg* formatCmd)
1333{
1334 bool showLabel= (label != nullptr && strlen(label) > 0);
1335
1336 std::string whatStr{what};
1337 std::transform(whatStr.begin(), whatStr.end(), whatStr.begin(), [](unsigned char c){ return std::toupper(c); });
1338 bool showN = whatStr.find('N') != std::string::npos;
1339 bool showR = whatStr.find('R') != std::string::npos;
1340 bool showM = whatStr.find('M') != std::string::npos;
1341 Int_t nPar= 0;
1342 if (showN) nPar++ ;
1343 if (showR) nPar++ ;
1344 if (showM) nPar++ ;
1345
1346 // calculate the box's size
1347 double dy(0.06), ymin(ymax-nPar*dy);
1348 if(showLabel) ymin-= dy;
1349
1350 // create the box and set its options
1351 TPaveText *box= new TPaveText(xmin,ymax,xmax,ymin,"BRNDC");
1352 if(!box) return nullptr;
1353 box->SetName(Form("%s_statBox",GetName())) ;
1354 box->SetFillColor(0);
1355 box->SetBorderSize(1);
1356 box->SetTextAlign(12);
1357 box->SetTextSize(0.04F);
1358 box->SetFillStyle(1001);
1359
1360 // add formatted text for each statistic
1361 RooRealVar N("N","Number of Events",sumEntries(cutSpec,cutRange));
1362 N.setPlotLabel("Entries") ;
1363 std::unique_ptr<RooRealVar> meanv{meanVar(*(RooRealVar*)frame->getPlotVar(),cutSpec,cutRange)};
1364 meanv->setPlotLabel("Mean") ;
1365 std::unique_ptr<RooRealVar> rms{rmsVar(*(RooRealVar*)frame->getPlotVar(),cutSpec,cutRange)};
1366 rms->setPlotLabel("RMS") ;
1367 std::unique_ptr<TString> rmsText, meanText, NText;
1368 if (options) {
1369 rmsText.reset(rms->format(sigDigits,options));
1370 meanText.reset(meanv->format(sigDigits,options));
1371 NText.reset(N.format(sigDigits,options));
1372 } else {
1373 rmsText.reset(rms->format(*formatCmd));
1374 meanText.reset(meanv->format(*formatCmd));
1375 NText.reset(N.format(*formatCmd));
1376 }
1377 if (showR) box->AddText(rmsText->Data());
1378 if (showM) box->AddText(meanText->Data());
1379 if (showN) box->AddText(NText->Data());
1380
1381 // add the optional label if specified
1382 if(showLabel) box->AddText(label);
1383
1384 frame->addObject(box) ;
1385 return frame ;
1386}
1387
1388////////////////////////////////////////////////////////////////////////////////
1389/// Loop over columns of our tree data and fill the input histogram. Returns a pointer to the
1390/// input histogram, or zero in case of an error. The input histogram can be any TH1 subclass, and
1391/// therefore of arbitrary dimension. Variables are matched with the (x,y,...) dimensions of the input
1392/// histogram according to the order in which they appear in the input plotVars list.
1393
1394TH1 *RooAbsData::fillHistogram(TH1 *hist, const RooArgList &plotVars, const char *cuts, const char* cutRange) const
1395{
1396 // Do we have a valid histogram to use?
1397 if(nullptr == hist) {
1398 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: no valid histogram to fill" << endl;
1399 return nullptr;
1400 }
1401
1402 // Check that the number of plotVars matches the input histogram's dimension
1403 std::size_t hdim= hist->GetDimension();
1404 if(hdim != plotVars.size()) {
1405 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: plotVars has the wrong dimension" << endl;
1406 return nullptr;
1407 }
1408
1409 // Check that the plot variables are all actually RooAbsReal's and print a warning if we do not
1410 // explicitly depend on one of them. Clone any variables that we do not contain directly and
1411 // redirect them to use our event data.
1412 RooArgSet plotClones,localVars;
1413 for(std::size_t index= 0; index < plotVars.size(); index++) {
1414 const RooAbsArg *var= plotVars.at(index);
1415 const RooAbsReal *realVar= dynamic_cast<const RooAbsReal*>(var);
1416 if(0 == realVar) {
1417 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: cannot plot variable \"" << var->GetName()
1418 << "\" of type " << var->ClassName() << endl;
1419 return nullptr;
1420 }
1421 RooAbsArg *found= _vars.find(realVar->GetName());
1422 if(!found) {
1423 RooAbsArg *clone= plotClones.addClone(*realVar,true); // do not complain about duplicates
1424 assert(0 != clone);
1425 if(!clone->dependsOn(_vars)) {
1426 coutE(InputArguments) << ClassName() << "::" << GetName()
1427 << ":fillHistogram: Data does not contain the variable '" << realVar->GetName() << "'." << endl;
1428 return nullptr;
1429 }
1430 else {
1432 }
1433 localVars.add(*clone);
1434 }
1435 else {
1436 localVars.add(*found);
1437 }
1438 }
1439
1440 // Create selection formula if selection cuts are specified
1441 std::unique_ptr<RooFormula> select;
1442 if (cuts != nullptr && strlen(cuts) > 0) {
1443 select = std::make_unique<RooFormula>(cuts, cuts, _vars, false);
1444 if (!select || !select->ok()) {
1445 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: invalid cuts \"" << cuts << "\"" << endl;
1446 return nullptr;
1447 }
1448 }
1449
1450 // Lookup each of the variables we are binning in our tree variables
1451 const RooAbsReal *xvar = 0;
1452 const RooAbsReal *yvar = 0;
1453 const RooAbsReal *zvar = 0;
1454 switch(hdim) {
1455 case 3:
1456 zvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(2)->GetName()));
1457 assert(0 != zvar);
1458 // fall through to next case...
1459 case 2:
1460 yvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(1)->GetName()));
1461 assert(0 != yvar);
1462 // fall through to next case...
1463 case 1:
1464 xvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(0)->GetName()));
1465 assert(0 != xvar);
1466 break;
1467 default:
1468 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: cannot fill histogram with "
1469 << hdim << " dimensions" << endl;
1470 break;
1471 }
1472
1473 // Parse cutRange specification
1474 const auto cutVec = ROOT::Split(cutRange ? cutRange : "", ",");
1475
1476 // Loop over events and fill the histogram
1477 if (hist->GetSumw2()->fN==0) {
1478 hist->Sumw2() ;
1479 }
1480 Int_t nevent= numEntries() ; //(Int_t)_tree->GetEntries();
1481 for(Int_t i=0; i < nevent; ++i) {
1482
1483 //Int_t entryNumber= _tree->GetEntryNumber(i);
1484 //if (entryNumber<0) break;
1485 get(i);
1486
1487 // Apply expression based selection criteria
1488 if (select && select->eval()==0) {
1489 continue ;
1490 }
1491
1492
1493 // Apply range based selection criteria
1494 bool selectByRange = true ;
1495 if (cutRange) {
1496 for (const auto arg : _vars) {
1497 bool selectThisArg = false ;
1498 for (auto const& cut : cutVec) {
1499 if (!cut.empty() && arg->inRange(cut.c_str())) {
1500 selectThisArg = true ;
1501 break ;
1502 }
1503 }
1504 if (!selectThisArg) {
1505 selectByRange = false ;
1506 break ;
1507 }
1508 }
1509 }
1510
1511 if (!selectByRange) {
1512 // Go to next event in loop over events
1513 continue ;
1514 }
1515
1516 Int_t bin(0);
1517 switch(hdim) {
1518 case 1:
1519 bin= hist->FindBin(xvar->getVal());
1520 hist->Fill(xvar->getVal(),weight()) ;
1521 break;
1522 case 2:
1523 bin= hist->FindBin(xvar->getVal(),yvar->getVal());
1524 static_cast<TH2*>(hist)->Fill(xvar->getVal(),yvar->getVal(),weight()) ;
1525 break;
1526 case 3:
1527 bin= hist->FindBin(xvar->getVal(),yvar->getVal(),zvar->getVal());
1528 static_cast<TH3*>(hist)->Fill(xvar->getVal(),yvar->getVal(),zvar->getVal(),weight()) ;
1529 break;
1530 default:
1531 assert(hdim < 3);
1532 break;
1533 }
1534
1535
1536 double error2 = TMath::Power(hist->GetBinError(bin),2)-TMath::Power(weight(),2) ;
1537 double we = weightError(RooAbsData::SumW2) ;
1538 if (we==0) we = weight() ;
1539 error2 += TMath::Power(we,2) ;
1540
1541
1542// double we = weightError(RooAbsData::SumW2) ;
1543// double error2(0) ;
1544// if (we==0) {
1545// we = weight() ; //sqrt(weight()) ;
1546// error2 = TMath::Power(hist->GetBinError(bin),2)-TMath::Power(weight(),2) + TMath::Power(we,2) ;
1547// } else {
1548// error2 = TMath::Power(hist->GetBinError(bin),2)-TMath::Power(weight(),2) + TMath::Power(we,2) ;
1549// }
1550 //hist->AddBinContent(bin,weight());
1551 hist->SetBinError(bin,sqrt(error2)) ;
1552
1553 //cout << "RooTreeData::fillHistogram() bin = " << bin << " weight() = " << weight() << " we = " << we << endl ;
1554
1555 }
1556
1557 return hist;
1558}
1559
1560////////////////////////////////////////////////////////////////////////////////
1561/// Split dataset into subsets based on states of given splitCat in this dataset.
1562/// A TList of RooDataSets is returned in which each RooDataSet is named
1563/// after the state name of splitCat of which it contains the dataset subset.
1564/// The observables splitCat itself is no longer present in the sub datasets.
1565/// If createEmptyDataSets is false (default) this method only creates datasets for states
1566/// which have at least one entry The caller takes ownership of the returned list and its contents
1567
1568TList* RooAbsData::split(const RooAbsCategory& splitCat, bool createEmptyDataSets) const
1569{
1570 // Sanity check
1571 if (!splitCat.dependsOn(*get())) {
1572 coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") ERROR category " << splitCat.GetName()
1573 << " doesn't depend on any variable in this dataset" << endl ;
1574 return nullptr;
1575 }
1576
1577 // Clone splitting category and attach to self
1578 RooAbsCategory* cloneCat =0;
1579 std::unique_ptr<RooArgSet> cloneSet;
1580 if (splitCat.isDerived()) {
1581 cloneSet.reset(static_cast<RooArgSet*>(RooArgSet(splitCat).snapshot(true)));
1582 if (!cloneSet) {
1583 coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") Couldn't deep-clone splitting category, abort." << endl ;
1584 return nullptr;
1585 }
1586 cloneCat = (RooAbsCategory*) cloneSet->find(splitCat.GetName()) ;
1587 cloneCat->attachDataSet(*this) ;
1588 } else {
1589 cloneCat = dynamic_cast<RooAbsCategory*>(get()->find(splitCat.GetName())) ;
1590 if (!cloneCat) {
1591 coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") ERROR category " << splitCat.GetName()
1592 << " is fundamental and does not appear in this dataset" << endl ;
1593 return nullptr;
1594 }
1595 }
1596
1597 // Split a dataset in a series of subsets, each corresponding
1598 // to a state of splitCat
1599 TList* dsetList = new TList ;
1600
1601 // Construct set of variables to be included in split sets = full set - split category
1602 RooArgSet subsetVars(*get()) ;
1603 if (splitCat.isDerived()) {
1604 std::unique_ptr<RooArgSet> vars{splitCat.getVariables()};
1605 subsetVars.remove(*vars,true,true) ;
1606 } else {
1607 subsetVars.remove(splitCat,true,true) ;
1608 }
1609
1610 // Add weight variable explicitly if dataset has weights, but no top-level weight
1611 // variable exists (can happen with composite datastores)
1612 bool addWV(false) ;
1613 RooRealVar newweight("weight","weight",-1e9,1e9) ;
1615 subsetVars.add(newweight) ;
1616 addWV = true ;
1617 }
1618
1619 // If createEmptyDataSets is true, prepopulate with empty sets corresponding to all states
1620 if (createEmptyDataSets) {
1621 for (const auto& nameIdx : *cloneCat) {
1622 RooAbsData* subset = emptyClone(nameIdx.first.c_str(), nameIdx.first.c_str(), &subsetVars,(addWV?"weight":0)) ;
1623 dsetList->Add((RooAbsArg*)subset) ;
1624 }
1625 }
1626
1627
1628 // Loop over dataset and copy event to matching subset
1629 const bool propWeightSquared = isWeighted();
1630 for (Int_t i = 0; i < numEntries(); ++i) {
1631 const RooArgSet* row = get(i);
1632 RooAbsData* subset = (RooAbsData*) dsetList->FindObject(cloneCat->getCurrentLabel());
1633 if (!subset) {
1634 subset = emptyClone(cloneCat->getCurrentLabel(),cloneCat->getCurrentLabel(),&subsetVars,(addWV?"weight":0));
1635 dsetList->Add((RooAbsArg*)subset);
1636 }
1637 if (!propWeightSquared) {
1638 subset->add(*row, weight());
1639 } else {
1640 subset->add(*row, weight(), weightSquared());
1641 }
1642 }
1643
1644 return dsetList;
1645}
1646
1647////////////////////////////////////////////////////////////////////////////////
1648/// Split dataset into subsets based on the categorisation of the RooSimultaneous
1649/// A TList of RooDataSets is returned in which each RooDataSet is named
1650/// after the state name of splitCat of which it contains the dataset subset.
1651/// The observables splitCat itself is no longer present in the sub datasets, as well as the
1652/// observables of the other categories.
1653/// If createEmptyDataSets is false (default) this method only creates datasets for states
1654/// which have at least one entry The caller takes ownership of the returned list and its contents
1655
1656TList* RooAbsData::split(const RooSimultaneous& simpdf, bool createEmptyDataSets) const
1657{
1658 auto& splitCat = const_cast<RooAbsCategoryLValue&>(simpdf.indexCat());
1659
1660 // Sanity check
1661 if (!splitCat.dependsOn(*get())) {
1662 coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") ERROR category " << splitCat.GetName()
1663 << " doesn't depend on any variable in this dataset" << endl ;
1664 return nullptr;
1665 }
1666
1667 // Clone splitting category and attach to self
1668 RooAbsCategory* cloneCat =0;
1669 std::unique_ptr<RooArgSet> cloneSet;
1670 if (splitCat.isDerived()) {
1671 cloneSet.reset(static_cast<RooArgSet*>(RooArgSet(splitCat).snapshot(true)));
1672 if (!cloneSet) {
1673 coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") Couldn't deep-clone splitting category, abort." << endl ;
1674 return nullptr;
1675 }
1676 cloneCat = (RooAbsCategory*) cloneSet->find(splitCat.GetName()) ;
1677 cloneCat->attachDataSet(*this) ;
1678 } else {
1679 cloneCat = dynamic_cast<RooAbsCategory*>(get()->find(splitCat.GetName())) ;
1680 if (!cloneCat) {
1681 coutE(InputArguments) << "RooTreeData::split(" << GetName() << ") ERROR category " << splitCat.GetName()
1682 << " is fundamental and does not appear in this dataset" << endl ;
1683 return nullptr;
1684 }
1685 }
1686
1687 // Split a dataset in a series of subsets, each corresponding
1688 // to a state of splitCat
1689 TList* dsetList = new TList ;
1690
1691 // Construct set of variables to be included in split sets = full set - split category
1692 RooArgSet subsetVars(*get()) ;
1693 if (splitCat.isDerived()) {
1694 std::unique_ptr<RooArgSet> vars{splitCat.getVariables()};
1695 subsetVars.remove(*vars,true,true) ;
1696 } else {
1697 subsetVars.remove(splitCat,true,true) ;
1698 }
1699
1700 // Add weight variable explicitly if dataset has weights, but no top-level weight
1701 // variable exists (can happen with composite datastores)
1702 bool addWV(false) ;
1703 RooRealVar newweight("weight","weight",-1e9,1e9) ;
1705 subsetVars.add(newweight) ;
1706 addWV = true ;
1707 }
1708
1709 // Get the observables for a given pdf in the RooSimultaneous, or an empty
1710 // RooArgSet if no pdf is set
1711 auto getPdfObservables = [this, &simpdf](const char * label) {
1712 RooArgSet obsSet;
1713 if(RooAbsPdf* catPdf = simpdf.getPdf(label)) {
1714 catPdf->getObservables(this->get(), obsSet);
1715 }
1716 return obsSet;
1717 };
1718
1719 // By default, remove all category observables from the subdatasets
1720 RooArgSet allObservables;
1721 for( const auto& catPair : splitCat) {
1722 allObservables.add(getPdfObservables(catPair.first.c_str()));
1723 }
1724 subsetVars.remove(allObservables, true, true);
1725
1726
1727 // If createEmptyDataSets is true, prepopulate with empty sets corresponding to all states
1728 if (createEmptyDataSets) {
1729 for (const auto& nameIdx : *cloneCat) {
1730 // Add in the subset only the observables corresponding to this category
1731 RooArgSet subsetVarsCat(subsetVars);
1732 subsetVarsCat.add(getPdfObservables(nameIdx.first.c_str()));
1733 RooAbsData* subset = emptyClone(nameIdx.first.c_str(), nameIdx.first.c_str(), &subsetVarsCat,(addWV?"weight":0)) ;
1734 dsetList->Add(subset) ;
1735 }
1736 }
1737
1738
1739 // Loop over dataset and copy event to matching subset
1740 const bool propWeightSquared = isWeighted();
1741 for (Int_t i = 0; i < numEntries(); ++i) {
1742 const RooArgSet* row = get(i);
1743 RooAbsData* subset = (RooAbsData*) dsetList->FindObject(cloneCat->getCurrentLabel());
1744 if (!subset) {
1745 // Add in the subset only the observables corresponding to this category
1746 RooArgSet subsetVarsCat(subsetVars);
1747 subsetVarsCat.add(getPdfObservables(cloneCat->getCurrentLabel()));
1748 subset = emptyClone(cloneCat->getCurrentLabel(),cloneCat->getCurrentLabel(),&subsetVarsCat,(addWV?"weight":0));
1749 dsetList->Add(subset);
1750 }
1751 if (!propWeightSquared) {
1752 subset->add(*row, weight());
1753 } else {
1754 subset->add(*row, weight(), weightSquared());
1755 }
1756 }
1757
1758 return dsetList;
1759}
1760
1761////////////////////////////////////////////////////////////////////////////////
1762/// Plot dataset on specified frame.
1763///
1764/// By default:
1765/// - An unbinned dataset will use the default binning of the target frame.
1766/// - A binned dataset will retain its intrinsic binning.
1767///
1768/// The following optional named arguments can be used to modify the behaviour:
1769/// \note Please follow the function links in the left column to learn about PyROOT specifics for a given option.
1770///
1771/// <table>
1772///
1773/// <tr><th> <th> Data representation options
1774/// <tr><td> RooFit::Asymmetry(const RooCategory& c)
1775/// <td> Show the asymmetry of the data in given two-state category [F(+)-F(-)] / [F(+)+F(-)].
1776/// Category must have two states with indices -1 and +1 or three states with indices -1,0 and +1.
1777/// <tr><td> RooFit::Efficiency(const RooCategory& c)
1778/// <td> Show the efficiency F(acc)/[F(acc)+F(rej)]. Category must have two states with indices 0 and 1
1779/// <tr><td> RooFit::DataError(Int_t)
1780/// <td> Select the type of error drawn:
1781/// - `Auto(default)` results in Poisson for unweighted data and SumW2 for weighted data
1782/// - `Poisson` draws asymmetric Poisson confidence intervals.
1783/// - `SumW2` draws symmetric sum-of-weights error ( \f$ \left( \sum w \right)^2 / \sum\left(w^2\right) \f$ )
1784/// - `None` draws no error bars
1785/// <tr><td> RooFit::Binning(int nbins, double xlo, double xhi)
1786/// <td> Use specified binning to draw dataset
1787/// <tr><td> RooFit::Binning(const RooAbsBinning&)
1788/// <td> Use specified binning to draw dataset
1789/// <tr><td> RooFit::Binning(const char* name)
1790/// <td> Use binning with specified name to draw dataset
1791/// <tr><td> RooFit::RefreshNorm()
1792/// <td> Force refreshing for PDF normalization information in frame.
1793/// If set, any subsequent PDF will normalize to this dataset, even if it is
1794/// not the first one added to the frame. By default only the 1st dataset
1795/// added to a frame will update the normalization information
1796/// <tr><td> RooFit::Rescale(double f)
1797/// <td> Rescale drawn histogram by given factor.
1798/// <tr><td> RooFit::Cut(const char*)
1799/// <td> Only plot entries that pass the given cut.
1800/// Apart from cutting in continuous variables `Cut("x>5")`, this can also be used to plot a specific
1801/// category state. Use something like `Cut("myCategory == myCategory::stateA")`, where
1802/// `myCategory` resolves to the state number for a given entry and
1803/// `myCategory::stateA` resolves to the state number of the state named "stateA".
1804///
1805/// <tr><td> RooFit::CutRange(const char*)
1806/// <td> Only plot data from given range. Separate multiple ranges with ",".
1807/// \note This often requires passing the normalisation when plotting the PDF because RooFit does not save
1808/// how many events were being plotted (it will only work for cutting slices out of uniformly distributed
1809/// variables).
1810/// ```
1811/// data->plotOn(frame01, CutRange("SB1"));
1812/// const double nData = data->sumEntries("", "SB1");
1813/// // Make clear that the target normalisation is nData. The enumerator NumEvent
1814/// // is needed to switch between relative and absolute scaling.
1815/// model.plotOn(frame01, Normalization(nData, RooAbsReal::NumEvent),
1816/// ProjectionRange("SB1"));
1817/// ```
1818///
1819/// <tr><th> <th> Histogram drawing options
1820/// <tr><td> RooFit::DrawOption(const char* opt)
1821/// <td> Select ROOT draw option for resulting TGraph object
1822/// <tr><td> RooFit::LineStyle(Style_t style)
1823/// <td> Select line style by ROOT line style code, default is solid
1824/// <tr><td> RooFit::LineColor(Color_t color)
1825/// <td> Select line color by ROOT color code, default is black
1826/// <tr><td> RooFit::LineWidth(Width_t width)
1827/// <td> Select line with in pixels, default is 3
1828/// <tr><td> RooFit::MarkerStyle(Style_t style)
1829/// <td> Select the ROOT marker style, default is 21
1830/// <tr><td> RooFit::MarkerColor(Color_t color)
1831/// <td> Select the ROOT marker color, default is black
1832/// <tr><td> RooFit::MarkerSize(Size_t size)
1833/// <td> Select the ROOT marker size
1834/// <tr><td> RooFit::FillStyle(Style_t style)
1835/// <td> Select fill style, default is filled.
1836/// <tr><td> RooFit::FillColor(Color_t color)
1837/// <td> Select fill color by ROOT color code
1838/// <tr><td> RooFit::XErrorSize(double frac)
1839/// <td> Select size of X error bar as fraction of the bin width, default is 1
1840///
1841/// <tr><th> <th> Misc. other options
1842/// <tr><td> RooFit::Name(const char* name)
1843/// <td> Give curve specified name in frame. Useful if curve is to be referenced later
1844/// <tr><td> RooFit::Invisible()
1845/// <td> Add curve to frame, but do not display. Useful in combination AddTo()
1846/// <tr><td> RooFit::AddTo(const char* name, double wgtSel, double wgtOther)
1847/// <td> Add constructed histogram to already existing histogram with given name and relative weight factors
1848///
1849/// </table>
1850
1851RooPlot* RooAbsData::plotOn(RooPlot* frame, const RooLinkedList& argList) const
1852{
1853 // New experimental plotOn() with varargs...
1854
1855 // Define configuration for this method
1856 RooCmdConfig pc(Form("RooAbsData::plotOn(%s)",GetName())) ;
1857 pc.defineString("drawOption","DrawOption",0,"P") ;
1858 pc.defineString("cutRange","CutRange",0,"",true) ;
1859 pc.defineString("cutString","CutSpec",0,"") ;
1860 pc.defineString("histName","Name",0,"") ;
1861 pc.defineObject("cutVar","CutVar",0) ;
1862 pc.defineObject("binning","Binning",0) ;
1863 pc.defineString("binningName","BinningName",0,"") ;
1864 pc.defineInt("nbins","BinningSpec",0,100) ;
1865 pc.defineDouble("xlo","BinningSpec",0,0) ;
1866 pc.defineDouble("xhi","BinningSpec",1,1) ;
1867 pc.defineObject("asymCat","Asymmetry",0) ;
1868 pc.defineObject("effCat","Efficiency",0) ;
1869 pc.defineInt("lineColor","LineColor",0,-999) ;
1870 pc.defineInt("lineStyle","LineStyle",0,-999) ;
1871 pc.defineInt("lineWidth","LineWidth",0,-999) ;
1872 pc.defineInt("markerColor","MarkerColor",0,-999) ;
1873 pc.defineInt("markerStyle","MarkerStyle",0,-999) ;
1874 pc.defineDouble("markerSize","MarkerSize",0,-999) ;
1875 pc.defineInt("fillColor","FillColor",0,-999) ;
1876 pc.defineInt("fillStyle","FillStyle",0,-999) ;
1877 pc.defineInt("errorType","DataError",0,(Int_t)RooAbsData::Auto) ;
1878 pc.defineInt("histInvisible","Invisible",0,0) ;
1879 pc.defineInt("refreshFrameNorm","RefreshNorm",0,1) ;
1880 pc.defineString("addToHistName","AddTo",0,"") ;
1881 pc.defineDouble("addToWgtSelf","AddTo",0,1.) ;
1882 pc.defineDouble("addToWgtOther","AddTo",1,1.) ;
1883 pc.defineDouble("xErrorSize","XErrorSize",0,1.) ;
1884 pc.defineDouble("scaleFactor","Rescale",0,1.) ;
1885 pc.defineMutex("DataError","Asymmetry","Efficiency") ;
1886 pc.defineMutex("Binning","BinningName","BinningSpec") ;
1887
1888 // Process & check varargs
1889 pc.process(argList) ;
1890 if (!pc.ok(true)) {
1891 return frame ;
1892 }
1893
1894 PlotOpt o ;
1895
1896 // Extract values from named arguments
1897 o.drawOptions = pc.getString("drawOption") ;
1898 o.cuts = pc.getString("cutString") ;
1899 if (pc.hasProcessed("Binning")) {
1900 o.bins = (RooAbsBinning*) pc.getObject("binning") ;
1901 } else if (pc.hasProcessed("BinningName")) {
1902 o.bins = &frame->getPlotVar()->getBinning(pc.getString("binningName")) ;
1903 } else if (pc.hasProcessed("BinningSpec")) {
1904 double xlo = pc.getDouble("xlo") ;
1905 double xhi = pc.getDouble("xhi") ;
1906 o.bins = new RooUniformBinning((xlo==xhi)?frame->getPlotVar()->getMin():xlo,
1907 (xlo==xhi)?frame->getPlotVar()->getMax():xhi,pc.getInt("nbins")) ;
1908 }
1909 const RooAbsCategoryLValue* asymCat = (const RooAbsCategoryLValue*) pc.getObject("asymCat") ;
1910 const RooAbsCategoryLValue* effCat = (const RooAbsCategoryLValue*) pc.getObject("effCat") ;
1911 o.etype = (RooAbsData::ErrorType) pc.getInt("errorType") ;
1912 o.histInvisible = pc.getInt("histInvisible") ;
1913 o.xErrorSize = pc.getDouble("xErrorSize") ;
1914 o.cutRange = pc.getString("cutRange",0,true) ;
1915 o.histName = pc.getString("histName",0,true) ;
1916 o.addToHistName = pc.getString("addToHistName",0,true) ;
1917 o.addToWgtSelf = pc.getDouble("addToWgtSelf") ;
1918 o.addToWgtOther = pc.getDouble("addToWgtOther") ;
1919 o.refreshFrameNorm = pc.getInt("refreshFrameNorm") ;
1920 o.scaleFactor = pc.getDouble("scaleFactor") ;
1921
1922 // Map auto error type to actual type
1923 if (o.etype == Auto) {
1925 if (o.etype == SumW2) {
1926 coutI(InputArguments) << "RooAbsData::plotOn(" << GetName()
1927 << ") INFO: dataset has non-integer weights, auto-selecting SumW2 errors instead of Poisson errors" << endl ;
1928 }
1929 }
1930
1931 if (o.addToHistName && !frame->findObject(o.addToHistName,RooHist::Class())) {
1932 coutE(InputArguments) << "RooAbsData::plotOn(" << GetName() << ") cannot find existing histogram " << o.addToHistName
1933 << " to add to in RooPlot" << endl ;
1934 return frame ;
1935 }
1936
1937 RooPlot* ret ;
1938 if (!asymCat && !effCat) {
1939 ret = plotOn(frame,o) ;
1940 } else if (asymCat) {
1941 ret = plotAsymOn(frame,*asymCat,o) ;
1942 } else {
1943 ret = plotEffOn(frame,*effCat,o) ;
1944 }
1945
1946 Int_t lineColor = pc.getInt("lineColor") ;
1947 Int_t lineStyle = pc.getInt("lineStyle") ;
1948 Int_t lineWidth = pc.getInt("lineWidth") ;
1949 Int_t markerColor = pc.getInt("markerColor") ;
1950 Int_t markerStyle = pc.getInt("markerStyle") ;
1951 Size_t markerSize = pc.getDouble("markerSize") ;
1952 Int_t fillColor = pc.getInt("fillColor") ;
1953 Int_t fillStyle = pc.getInt("fillStyle") ;
1954 if (lineColor!=-999) ret->getAttLine()->SetLineColor(lineColor) ;
1955 if (lineStyle!=-999) ret->getAttLine()->SetLineStyle(lineStyle) ;
1956 if (lineWidth!=-999) ret->getAttLine()->SetLineWidth(lineWidth) ;
1957 if (markerColor!=-999) ret->getAttMarker()->SetMarkerColor(markerColor) ;
1958 if (markerStyle!=-999) ret->getAttMarker()->SetMarkerStyle(markerStyle) ;
1959 if (markerSize!=-999) ret->getAttMarker()->SetMarkerSize(markerSize) ;
1960 if (fillColor!=-999) ret->getAttFill()->SetFillColor(fillColor) ;
1961 if (fillStyle!=-999) ret->getAttFill()->SetFillStyle(fillStyle) ;
1962
1963 if (pc.hasProcessed("BinningSpec")) {
1964 delete o.bins ;
1965 }
1966
1967 return ret ;
1968}
1969
1970////////////////////////////////////////////////////////////////////////////////
1971/// Create and fill a histogram of the frame's variable and append it to the frame.
1972/// The frame variable must be one of the data sets dimensions.
1973///
1974/// The plot range and the number of plot bins is determined by the parameters
1975/// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins()).
1976///
1977/// The optional cut string expression can be used to select the events to be plotted.
1978/// The cut specification may refer to any variable contained in the data set.
1979///
1980/// The drawOptions are passed to the TH1::Draw() method.
1981/// \see RooAbsData::plotOn(RooPlot*,const RooLinkedList&) const
1983{
1984 if(0 == frame) {
1985 coutE(Plotting) << ClassName() << "::" << GetName() << ":plotOn: frame is null" << endl;
1986 return nullptr;
1987 }
1989 if(0 == var) {
1990 coutE(Plotting) << ClassName() << "::" << GetName()
1991 << ":plotOn: frame does not specify a plot variable" << endl;
1992 return nullptr;
1993 }
1994
1995 // create and fill a temporary histogram of this variable
1996 const std::string histName = std::string{GetName()} + "_plot";
1997 std::unique_ptr<TH1> hist;
1998 if (o.bins) {
1999 hist.reset( var->createHistogram(histName.c_str(), RooFit::AxisLabel("Events"), RooFit::Binning(*o.bins)) );
2000 } else if (!frame->getPlotVar()->getBinning().isUniform()) {
2001 hist.reset( var->createHistogram(histName.c_str(), RooFit::AxisLabel("Events"),
2002 RooFit::Binning(frame->getPlotVar()->getBinning())) );
2003 } else {
2004 hist.reset( var->createHistogram(histName.c_str(), "Events",
2005 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(), frame->GetNbinsX()) );
2006 }
2007
2008 // Keep track of sum-of-weights error
2009 hist->Sumw2() ;
2010
2011 if(0 == fillHistogram(hist.get(), RooArgList(*var),o.cuts,o.cutRange)) {
2012 coutE(Plotting) << ClassName() << "::" << GetName()
2013 << ":plotOn: fillHistogram() failed" << endl;
2014 return nullptr;
2015 }
2016
2017 // If frame has no predefined bin width (event density) it will be adjusted to
2018 // our histograms bin width so we should force that bin width here
2019 double nomBinWidth ;
2020 if (frame->getFitRangeNEvt()==0 && o.bins) {
2021 nomBinWidth = o.bins->averageBinWidth() ;
2022 } else {
2023 nomBinWidth = o.bins ? frame->getFitRangeBinW() : 0 ;
2024 }
2025
2026 // convert this histogram to a RooHist object on the heap
2027 RooHist *graph= new RooHist(*hist,nomBinWidth,1,o.etype,o.xErrorSize,o.correctForBinWidth,o.scaleFactor);
2028 if(0 == graph) {
2029 coutE(Plotting) << ClassName() << "::" << GetName()
2030 << ":plotOn: unable to create a RooHist object" << endl;
2031 return nullptr;
2032 }
2033
2034 // If the dataset variable has a wide range than the plot variable,
2035 // calculate the number of entries in the dataset in the plot variable fit range
2036 RooAbsRealLValue* dataVar = (RooAbsRealLValue*) _vars.find(var->GetName()) ;
2037 double nEnt(sumEntries()) ;
2038 if (dataVar->getMin()<var->getMin() || dataVar->getMax()>var->getMax()) {
2039 std::unique_ptr<RooAbsData> tmp{const_cast<RooAbsData*>(this)->reduce(*var)};
2040 nEnt = tmp->sumEntries() ;
2041 }
2042
2043 // Store the number of entries before the cut, if any was made
2044 if ((o.cuts && strlen(o.cuts)) || o.cutRange) {
2045 coutI(Plotting) << "RooTreeData::plotOn: plotting " << hist->GetSumOfWeights() << " events out of " << nEnt << " total events" << endl ;
2046 graph->setRawEntries(nEnt) ;
2047 }
2048
2049 // Add self to other hist if requested
2050 if (o.addToHistName) {
2051 RooHist* otherGraph = static_cast<RooHist*>(frame->findObject(o.addToHistName,RooHist::Class())) ;
2052
2053 if (!graph->hasIdenticalBinning(*otherGraph)) {
2054 coutE(Plotting) << "RooTreeData::plotOn: ERROR Histogram to be added to, '" << o.addToHistName << "',has different binning" << endl ;
2055 delete graph ;
2056 return frame ;
2057 }
2058
2059 RooHist* sumGraph = new RooHist(*graph,*otherGraph,o.addToWgtSelf,o.addToWgtOther,o.etype) ;
2060 delete graph ;
2061 graph = sumGraph ;
2062 }
2063
2064 // Rename graph if requested
2065 if (o.histName) {
2066 graph->SetName(o.histName) ;
2067 } else {
2068 std::string hname = std::string{"h_"} + GetName();
2069 if (o.cutRange && strlen(o.cutRange)>0) {
2070 hname += std::string{"_CutRange["} + o.cutRange + "]";
2071 }
2072 if (o.cuts && strlen(o.cuts)>0) {
2073 hname += std::string{"_Cut["} + o.cuts + "]";
2074 }
2075 graph->SetName(hname.c_str()) ;
2076 }
2077
2078 // initialize the frame's normalization setup, if necessary
2079 frame->updateNormVars(_vars);
2080
2081
2082 // add the RooHist to the specified plot
2084
2085 return frame;
2086}
2087
2088////////////////////////////////////////////////////////////////////////////////
2089/// Create and fill a histogram with the asymmetry N[+] - N[-] / ( N[+] + N[-] ),
2090/// where N(+/-) is the number of data points with asymCat=+1 and asymCat=-1
2091/// as function of the frames variable. The asymmetry category 'asymCat' must
2092/// have exactly 2 (or 3) states defined with index values +1,-1 (and 0)
2093///
2094/// The plot range and the number of plot bins is determined by the parameters
2095/// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins())
2096///
2097/// The optional cut string expression can be used to select the events to be plotted.
2098/// The cut specification may refer to any variable contained in the data set
2099///
2100/// The drawOptions are passed to the TH1::Draw() method
2101
2103{
2104 if(0 == frame) {
2105 coutE(Plotting) << ClassName() << "::" << GetName() << ":plotAsymOn: frame is null" << endl;
2106 return nullptr;
2107 }
2109 if(0 == var) {
2110 coutE(Plotting) << ClassName() << "::" << GetName()
2111 << ":plotAsymOn: frame does not specify a plot variable" << endl;
2112 return nullptr;
2113 }
2114
2115 // create and fill temporary histograms of this variable for each state
2116 std::string hist1Name(GetName()),hist2Name(GetName());
2117 hist1Name += "_plot1";
2118 std::unique_ptr<TH1> hist1, hist2;
2119 hist2Name += "_plot2";
2120
2121 if (o.bins) {
2122 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events", *o.bins) );
2123 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events", *o.bins) );
2124 } else {
2125 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events",
2126 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
2127 frame->GetNbinsX()) );
2128 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events",
2129 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
2130 frame->GetNbinsX()) );
2131 }
2132
2133 assert(hist1 && hist2);
2134
2135 std::string cuts1,cuts2 ;
2136 if (o.cuts && strlen(o.cuts)) {
2137 cuts1 = Form("(%s)&&(%s>0)",o.cuts,asymCat.GetName());
2138 cuts2 = Form("(%s)&&(%s<0)",o.cuts,asymCat.GetName());
2139 } else {
2140 cuts1 = Form("(%s>0)",asymCat.GetName());
2141 cuts2 = Form("(%s<0)",asymCat.GetName());
2142 }
2143
2144 if(! fillHistogram(hist1.get(), RooArgList(*var),cuts1.c_str(),o.cutRange) ||
2145 ! fillHistogram(hist2.get(), RooArgList(*var),cuts2.c_str(),o.cutRange)) {
2146 coutE(Plotting) << ClassName() << "::" << GetName()
2147 << ":plotAsymOn: createHistogram() failed" << endl;
2148 return nullptr;
2149 }
2150
2151 // convert this histogram to a RooHist object on the heap
2152 RooHist *graph= new RooHist(*hist1,*hist2,0,1,o.etype,o.xErrorSize,false,o.scaleFactor);
2153 graph->setYAxisLabel(Form("Asymmetry in %s",asymCat.GetName())) ;
2154
2155 // initialize the frame's normalization setup, if necessary
2156 frame->updateNormVars(_vars);
2157
2158 // Rename graph if requested
2159 if (o.histName) {
2160 graph->SetName(o.histName) ;
2161 } else {
2162 std::string hname{Form("h_%s_Asym[%s]",GetName(),asymCat.GetName())};
2163 if (o.cutRange && strlen(o.cutRange)>0) {
2164 hname += Form("_CutRange[%s]",o.cutRange);
2165 }
2166 if (o.cuts && strlen(o.cuts)>0) {
2167 hname += Form("_Cut[%s]",o.cuts);
2168 }
2169 graph->SetName(hname.c_str()) ;
2170 }
2171
2172 // add the RooHist to the specified plot
2174
2175 return frame;
2176}
2177
2178////////////////////////////////////////////////////////////////////////////////
2179/// Create and fill a histogram with the efficiency N[1] / ( N[1] + N[0] ),
2180/// where N(1/0) is the number of data points with effCat=1 and effCat=0
2181/// as function of the frames variable. The efficiency category 'effCat' must
2182/// have exactly 2 +1 and 0.
2183///
2184/// The plot range and the number of plot bins is determined by the parameters
2185/// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins())
2186///
2187/// The optional cut string expression can be used to select the events to be plotted.
2188/// The cut specification may refer to any variable contained in the data set
2189///
2190/// The drawOptions are passed to the TH1::Draw() method
2191
2193{
2194 if(0 == frame) {
2195 coutE(Plotting) << ClassName() << "::" << GetName() << ":plotEffOn: frame is null" << endl;
2196 return nullptr;
2197 }
2199 if(0 == var) {
2200 coutE(Plotting) << ClassName() << "::" << GetName()
2201 << ":plotEffOn: frame does not specify a plot variable" << endl;
2202 return nullptr;
2203 }
2204
2205 // create and fill temporary histograms of this variable for each state
2206 std::string hist1Name(GetName()),hist2Name(GetName());
2207 hist1Name += "_plot1";
2208 std::unique_ptr<TH1> hist1, hist2;
2209 hist2Name += "_plot2";
2210
2211 if (o.bins) {
2212 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events", *o.bins) );
2213 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events", *o.bins) );
2214 } else {
2215 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events",
2216 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
2217 frame->GetNbinsX()) );
2218 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events",
2219 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
2220 frame->GetNbinsX()) );
2221 }
2222
2223 assert(hist1 && hist2);
2224
2225 std::string cuts1,cuts2 ;
2226 if (o.cuts && strlen(o.cuts)) {
2227 cuts1 = Form("(%s)&&(%s==1)",o.cuts,effCat.GetName());
2228 cuts2 = Form("(%s)&&(%s==0)",o.cuts,effCat.GetName());
2229 } else {
2230 cuts1 = Form("(%s==1)",effCat.GetName());
2231 cuts2 = Form("(%s==0)",effCat.GetName());
2232 }
2233
2234 if(! fillHistogram(hist1.get(), RooArgList(*var),cuts1.c_str(),o.cutRange) ||
2235 ! fillHistogram(hist2.get(), RooArgList(*var),cuts2.c_str(),o.cutRange)) {
2236 coutE(Plotting) << ClassName() << "::" << GetName()
2237 << ":plotEffOn: createHistogram() failed" << endl;
2238 return nullptr;
2239 }
2240
2241 // convert this histogram to a RooHist object on the heap
2242 RooHist *graph= new RooHist(*hist1,*hist2,0,1,o.etype,o.xErrorSize,true);
2243 graph->setYAxisLabel(Form("Efficiency of %s=%s", effCat.GetName(), effCat.lookupName(1).c_str()));
2244
2245 // initialize the frame's normalization setup, if necessary
2246 frame->updateNormVars(_vars);
2247
2248 // Rename graph if requested
2249 if (o.histName) {
2250 graph->SetName(o.histName) ;
2251 } else {
2252 std::string hname(Form("h_%s_Eff[%s]",GetName(),effCat.GetName())) ;
2253 if (o.cutRange && strlen(o.cutRange)>0) {
2254 hname += Form("_CutRange[%s]",o.cutRange);
2255 }
2256 if (o.cuts && strlen(o.cuts)>0) {
2257 hname += Form("_Cut[%s]",o.cuts);
2258 }
2259 graph->SetName(hname.c_str()) ;
2260 }
2261
2262 // add the RooHist to the specified plot
2264
2265 return frame;
2266}
2267
2268////////////////////////////////////////////////////////////////////////////////
2269/// Create and fill a 1-dimensional table for given category column
2270/// This functions is the equivalent of plotOn() for category dimensions.
2271///
2272/// The optional cut string expression can be used to select the events to be tabulated
2273/// The cut specification may refer to any variable contained in the data set
2274///
2275/// The option string is currently not used
2276
2277Roo1DTable* RooAbsData::table(const RooAbsCategory& cat, const char* cuts, const char* /*opts*/) const
2278{
2279 // First see if var is in data set
2280 RooAbsCategory* tableVar = (RooAbsCategory*) _vars.find(cat.GetName()) ;
2281 std::unique_ptr<RooArgSet> tableSet;
2282 if (!tableVar) {
2283 if (!cat.dependsOn(_vars)) {
2284 coutE(Plotting) << "RooTreeData::Table(" << GetName() << "): Argument " << cat.GetName()
2285 << " is not in dataset and is also not dependent on data set" << endl ;
2286 return nullptr;
2287 }
2288
2289 // Clone derived variable
2290 tableSet.reset(static_cast<RooArgSet*>(RooArgSet(cat).snapshot(true)));
2291 if (!tableSet) {
2292 coutE(Plotting) << "RooTreeData::table(" << GetName() << ") Couldn't deep-clone table category, abort." << std::endl;
2293 return nullptr;
2294 }
2295 tableVar = (RooAbsCategory*) tableSet->find(cat.GetName()) ;
2296
2297 //Redirect servers of derived clone to internal ArgSet representing the data in this set
2298 tableVar->recursiveRedirectServers(_vars) ;
2299 }
2300
2301 std::unique_ptr<RooFormulaVar> cutVar;
2302 std::string tableName{GetName()};
2303 if (cuts && strlen(cuts)) {
2304 tableName += "(";
2305 tableName += cuts;
2306 tableName += ")";
2307 // Make cut selector if cut is specified
2308 cutVar = std::make_unique<RooFormulaVar>("cutVar",cuts,_vars) ;
2309 }
2310 Roo1DTable* table2 = tableVar->createTable(tableName.c_str());
2311
2312 // Dump contents
2313 Int_t nevent= numEntries() ;
2314 for(Int_t i=0; i < nevent; ++i) {
2315 get(i);
2316
2317 if (cutVar && cutVar->getVal()==0) continue ;
2318
2319 table2->fill(*tableVar,weight()) ;
2320 }
2321
2322 return table2 ;
2323}
2324
2325////////////////////////////////////////////////////////////////////////////////
2326/// Fill Doubles 'lowest' and 'highest' with the lowest and highest value of
2327/// observable 'var' in this dataset. If the return value is true and error
2328/// occurred
2329
2330bool RooAbsData::getRange(const RooAbsRealLValue& var, double& lowest, double& highest, double marginFrac, bool symMode) const
2331{
2332 // Lookup variable in dataset
2333 const auto arg = _vars.find(var.GetName());
2334 if (!arg) {
2335 coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") ERROR: unknown variable: " << var.GetName() << endl ;
2336 return true;
2337 }
2338
2339 auto varPtr = dynamic_cast<const RooRealVar*>(arg);
2340 // Check if found variable is of type RooRealVar
2341 if (!varPtr) {
2342 coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") ERROR: variable " << var.GetName() << " is not of type RooRealVar" << endl ;
2343 return true;
2344 }
2345
2346 // Check if dataset is not empty
2347 if(sumEntries() == 0.) {
2348 coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") WARNING: empty dataset" << endl ;
2349 return true;
2350 }
2351
2352 // Look for highest and lowest value
2353 lowest = RooNumber::infinity() ;
2354 highest = -RooNumber::infinity() ;
2355 for (Int_t i=0 ; i<numEntries() ; i++) {
2356 get(i) ;
2357 if (varPtr->getVal()<lowest) {
2358 lowest = varPtr->getVal() ;
2359 }
2360 if (varPtr->getVal()>highest) {
2361 highest = varPtr->getVal() ;
2362 }
2363 }
2364
2365 if (marginFrac>0) {
2366 if (symMode==false) {
2367
2368 double margin = marginFrac*(highest-lowest) ;
2369 lowest -= margin ;
2370 highest += margin ;
2371 if (lowest<var.getMin()) lowest = var.getMin() ;
2372 if (highest>var.getMax()) highest = var.getMax() ;
2373
2374 } else {
2375
2376 double mom1 = moment(*varPtr,1) ;
2377 double delta = ((highest-mom1)>(mom1-lowest)?(highest-mom1):(mom1-lowest))*(1+marginFrac) ;
2378 lowest = mom1-delta ;
2379 highest = mom1+delta ;
2380 if (lowest<var.getMin()) lowest = var.getMin() ;
2381 if (highest>var.getMax()) highest = var.getMax() ;
2382
2383 }
2384 }
2385
2386 return false ;
2387}
2388
2389////////////////////////////////////////////////////////////////////////////////
2390/// Prepare dataset for use with cached constant terms listed in
2391/// 'cacheList' of expression 'arg'. Deactivate tree branches
2392/// for any dataset observable that is either not used at all,
2393/// or is used exclusively by cached branch nodes.
2394
2395void RooAbsData::optimizeReadingWithCaching(RooAbsArg& arg, const RooArgSet& cacheList, const RooArgSet& keepObsList)
2396{
2397 RooArgSet pruneSet ;
2398
2399 // Add unused observables in this dataset to pruneSet
2400 pruneSet.add(*get()) ;
2401 std::unique_ptr<RooArgSet> usedObs{arg.getObservables(*this)};
2402 pruneSet.remove(*usedObs,true,true) ;
2403
2404 // Add observables exclusively used to calculate cached observables to pruneSet
2405 for(auto * var : *get()) {
2406 if (allClientsCached(var,cacheList)) {
2407 pruneSet.add(*var) ;
2408 }
2409 }
2410
2411
2412 if (!pruneSet.empty()) {
2413
2414 // Go over all used observables and check if any of them have parameterized
2415 // ranges in terms of pruned observables. If so, remove those observable
2416 // from the pruning list
2417 for(auto const* rrv : dynamic_range_cast<RooRealVar*>(*usedObs)) {
2418 if (rrv && !rrv->getBinning().isShareable()) {
2419 RooArgSet depObs ;
2420 RooAbsReal* loFunc = rrv->getBinning().lowBoundFunc() ;
2421 RooAbsReal* hiFunc = rrv->getBinning().highBoundFunc() ;
2422 if (loFunc) {
2423 loFunc->leafNodeServerList(&depObs,0,true) ;
2424 }
2425 if (hiFunc) {
2426 hiFunc->leafNodeServerList(&depObs,0,true) ;
2427 }
2428 if (!depObs.empty()) {
2429 pruneSet.remove(depObs,true,true) ;
2430 }
2431 }
2432 }
2433 }
2434
2435
2436 // Remove all observables in keep list from prune list
2437 pruneSet.remove(keepObsList,true,true) ;
2438
2439 if (!pruneSet.empty()) {
2440
2441 // Deactivate tree branches here
2442 cxcoutI(Optimization) << "RooTreeData::optimizeReadingForTestStatistic(" << GetName() << "): Observables " << pruneSet
2443 << " 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 ;
2444 setArgStatus(pruneSet,false) ;
2445 }
2446}
2447
2448////////////////////////////////////////////////////////////////////////////////
2449/// Utility function that determines if all clients of object 'var'
2450/// appear in given list of cached nodes.
2451
2453{
2454 bool ret(true), anyClient(false) ;
2455
2456 for (const auto client : var->valueClients()) {
2457 anyClient = true ;
2458 if (!cacheList.find(client->GetName())) {
2459 // If client is not cached recurse
2460 ret &= allClientsCached(client,cacheList) ;
2461 }
2462 }
2463
2464 return anyClient?ret:false ;
2465}
2466
2467////////////////////////////////////////////////////////////////////////////////
2468
2470{
2471 _dstore->attachBuffers(extObs) ;
2472}
2473
2474////////////////////////////////////////////////////////////////////////////////
2475
2477{
2478 _dstore->resetBuffers() ;
2479}
2480
2481////////////////////////////////////////////////////////////////////////////////
2482
2484{
2485 return !_ownedComponents.empty();
2486}
2487
2488////////////////////////////////////////////////////////////////////////////////
2489
2491{
2492 auto i = _ownedComponents.find(name);
2493 return i==_ownedComponents.end() ? nullptr : i->second;
2494}
2495
2496////////////////////////////////////////////////////////////////////////////////
2497
2499{
2500 _ownedComponents[idxlabel]= &data ;
2501}
2502
2503////////////////////////////////////////////////////////////////////////////////
2504/// Stream an object of class RooAbsData.
2505
2507{
2508 if (R__b.IsReading()) {
2511
2512 // Convert on the fly to vector storage if that the current working default
2515 }
2516
2517 } else {
2519 }
2520}
2521
2522////////////////////////////////////////////////////////////////////////////////
2523
2525{
2526 _dstore->checkInit() ;
2527}
2528
2529////////////////////////////////////////////////////////////////////////////////
2530/// Forward draw command to data store
2531
2533{
2534 if (_dstore) _dstore->Draw(option) ;
2535}
2536
2537////////////////////////////////////////////////////////////////////////////////
2538
2540{
2541 return _dstore->hasFilledCache() ;
2542}
2543
2544////////////////////////////////////////////////////////////////////////////////
2545/// Return a pointer to the TTree which stores the data. Returns a nullpointer
2546/// if vector-based storage is used. The RooAbsData remains owner of the tree.
2547/// GetClonedTree() can be used to get a tree even if the internal storage does not use one.
2548
2550{
2552 return _dstore->tree();
2553 } else {
2554 coutW(InputArguments) << "RooAbsData::tree(" << GetName() << ") WARNING: is not of StorageType::Tree. "
2555 << "Use GetClonedTree() instead or convert to tree storage." << endl;
2556 return nullptr;
2557 }
2558}
2559
2560////////////////////////////////////////////////////////////////////////////////
2561/// Return a clone of the TTree which stores the data or create such a tree
2562/// if vector storage is used. The user is responsible for deleting the tree
2563
2565{
2567 auto tmp = const_cast<TTree *>(_dstore->tree());
2568 return tmp->CloneTree();
2569 } else {
2570 RooTreeDataStore buffer(GetName(), GetTitle(), *get(), *_dstore);
2571 return buffer.tree().CloneTree();
2572 }
2573}
2574
2575////////////////////////////////////////////////////////////////////////////////
2576/// Convert vector-based storage to tree-based storage
2577
2579{
2581 _dstore = std::make_unique<RooTreeDataStore>(GetName(), GetTitle(), _vars, *_dstore);
2583 }
2584}
2585
2586////////////////////////////////////////////////////////////////////////////////
2587/// If one of the TObject we have a referenced to is deleted, remove the
2588/// reference.
2589
2591{
2592 for(auto &iter : _ownedComponents) {
2593 if (iter.second == obj) {
2594 iter.second = nullptr;
2595 }
2596 }
2597}
2598
2599
2600////////////////////////////////////////////////////////////////////////////////
2601/// Sets the global observables stored in this data. A snapshot of the
2602/// observables will be saved.
2603/// \param[in] globalObservables The set of global observables to take a snapshot of.
2604
2605void RooAbsData::setGlobalObservables(RooArgSet const& globalObservables) {
2606 if(_globalObservables == nullptr) _globalObservables = std::make_unique<RooArgSet>();
2607 else _globalObservables->clear();
2608 globalObservables.snapshot(*_globalObservables);
2609 for(auto * arg : *_globalObservables) {
2610 arg->setAttribute("global",true);
2611 // Global observables are also always constant in fits
2612 if(auto lval = dynamic_cast<RooAbsRealLValue*>(arg)) lval->setConstant(true);
2613 if(auto lval = dynamic_cast<RooAbsCategoryLValue*>(arg)) lval->setConstant(true);
2614 }
2615}
2616
2617
2618////////////////////////////////////////////////////////////////////////////////
2619
2620void RooAbsData::SetName(const char* name)
2621{
2623 auto newPtr = RooNameReg::instance().constPtr(GetName()) ;
2624 if (newPtr != _namePtr) {
2625 //cout << "Rename '" << _namePtr->GetName() << "' to '" << name << "' (set flag in new name)" << endl;
2626 _namePtr = newPtr;
2629 }
2630}
2631
2632
2633
2634
2635////////////////////////////////////////////////////////////////////////////////
2636
2637void RooAbsData::SetNameTitle(const char *name, const char *title)
2638{
2639 TNamed::SetTitle(title) ;
2640 SetName(name);
2641}
2642
2643
2644
2645////////////////////////////////////////////////////////////////////////////////
2646/// Return sum of squared weights of this data.
2647
2649 const RooSpan<const double> eventWeights = getWeightBatch(0, numEntries(), /*sumW2=*/true);
2650 if (eventWeights.empty()) {
2651 return numEntries() * weightSquared();
2652 }
2653
2655 for (std::size_t i = 0; i < eventWeights.size(); ++i) {
2656 kahanWeight.AddIndexed(eventWeights[i], i);
2657 }
2658 return kahanWeight.Sum();
2659}
2660
2661
2662////////////////////////////////////////////////////////////////////////////////
2663/// Write information to retrieve data columns into `evalData.spans`.
2664/// All spans belonging to variables of this dataset are overwritten. Spans to other
2665/// variables remain intact.
2666/// \param[out] evalData Store references to all data batches in this struct's `spans`.
2667/// The key to retrieve an item is the pointer of the variable that owns the data.
2668/// \param begin Index of first event that ends up in the batch.
2669/// \param len Number of events in each batch.
2670RooAbsData::RealSpans RooAbsData::getBatches(std::size_t begin, std::size_t len) const {
2671 return store()->getBatches(begin, len);
2672}
2673
2674
2676 return store()->getCategoryBatches(first, len);
2677}
2678
2679////////////////////////////////////////////////////////////////////////////////
2680/// Create a TH2F histogram of the distribution of the specified variable
2681/// using this dataset. Apply any cuts to select which events are used.
2682/// The variable being plotted can either be contained directly in this
2683/// dataset, or else be a function of the variables in this dataset.
2684/// The histogram will be created using RooAbsReal::createHistogram() with
2685/// the name provided (with our dataset name prepended).
2686
2687TH2F *RooAbsData::createHistogram(const RooAbsRealLValue &var1, const RooAbsRealLValue &var2, const char *cuts,
2688 const char *name) const
2689{
2690 checkInit();
2691 return createHistogram(var1, var2, var1.getBins(), var2.getBins(), cuts, name);
2692}
2693
2694////////////////////////////////////////////////////////////////////////////////
2695/// Create a TH2F histogram of the distribution of the specified variable
2696/// using this dataset. Apply any cuts to select which events are used.
2697/// The variable being plotted can either be contained directly in this
2698/// dataset, or else be a function of the variables in this dataset.
2699/// The histogram will be created using RooAbsReal::createHistogram() with
2700/// the name provided (with our dataset name prepended).
2701
2702TH2F *RooAbsData::createHistogram(const RooAbsRealLValue &var1, const RooAbsRealLValue &var2, int nx, int ny,
2703 const char *cuts, const char *name) const
2704{
2705 checkInit();
2706 static int counter(0);
2707
2708 std::unique_ptr<RooAbsReal> ownedPlotVarX;
2709 // Is this variable in our dataset?
2710 auto *plotVarX = static_cast<RooAbsReal *>(_vars.find(var1.GetName()));
2711 if (plotVarX == nullptr) {
2712 // Is this variable a client of our dataset?
2713 if (!var1.dependsOn(_vars)) {
2714 coutE(InputArguments) << GetName() << "::createHistogram: Argument " << var1.GetName()
2715 << " is not in dataset and is also not dependent on data set" << std::endl;
2716 return nullptr;
2717 }
2718
2719 // Clone derived variable
2720 ownedPlotVarX.reset(static_cast<RooAbsReal *>(var1.Clone()));
2721 plotVarX = ownedPlotVarX.get();
2722
2723 // Redirect servers of derived clone to internal ArgSet representing the data in this set
2724 plotVarX->redirectServers(const_cast<RooArgSet &>(_vars));
2725 }
2726
2727 std::unique_ptr<RooAbsReal> ownedPlotVarY;
2728 // Is this variable in our dataset?
2729 RooAbsReal *plotVarY = (RooAbsReal *)_vars.find(var2.GetName());
2730 if (plotVarY == nullptr) {
2731 // Is this variable a client of our dataset?
2732 if (!var2.dependsOn(_vars)) {
2733 coutE(InputArguments) << GetName() << "::createHistogram: Argument " << var2.GetName()
2734 << " is not in dataset and is also not dependent on data set" << std::endl;
2735 return nullptr;
2736 }
2737
2738 // Clone derived variable
2739 ownedPlotVarY.reset(static_cast<RooAbsReal *>(var2.Clone()));
2740 plotVarY = ownedPlotVarY.get();
2741
2742 // Redirect servers of derived clone to internal ArgSet representing the data in this set
2743 plotVarY->redirectServers(const_cast<RooArgSet &>(_vars));
2744 }
2745
2746 // Create selection formula if selection cuts are specified
2747 std::unique_ptr<RooFormula> select;
2748 if (0 != cuts && strlen(cuts)) {
2749 select = std::make_unique<RooFormula>(cuts, cuts, _vars);
2750 if (!select->ok()) {
2751 return nullptr;
2752 }
2753 }
2754
2755 const std::string histName = std::string{GetName()} + "_" + name + "_" + Form("%08x", counter++);
2756
2757 // create the histogram
2758 auto *histogram =
2759 new TH2F(histName.c_str(), "Events", nx, var1.getMin(), var1.getMax(), ny, var2.getMin(), var2.getMax());
2760 if (!histogram) {
2761 coutE(DataHandling) << GetName() << "::createHistogram: unable to create a new histogram" << endl;
2762 return nullptr;
2763 }
2764
2765 // Dump contents
2766 Int_t nevent = numEntries();
2767 for (Int_t i = 0; i < nevent; ++i) {
2768 get(i);
2769
2770 if (select && select->eval() == 0)
2771 continue;
2772 histogram->Fill(plotVarX->getVal(), plotVarY->getVal(), weight());
2773 }
2774
2775 return histogram;
2776}
#define c(i)
Definition: RSha256.hxx:101
static std::map< RooAbsData *, int > _dcc
Definition: RooAbsData.cxx:126
#define coutI(a)
Definition: RooMsgService.h:34
#define cxcoutI(a)
Definition: RooMsgService.h:89
#define coutW(a)
Definition: RooMsgService.h:36
#define coutE(a)
Definition: RooMsgService.h:37
int Int_t
Definition: RtypesCore.h:45
float Size_t
Definition: RtypesCore.h:96
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
#define N
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h offset
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t UChar_t len
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
TMatrixTSym< Double_t > TMatrixDSym
Binding & operator=(OUT(*fun)(void))
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2452
T Sum() const
Definition: Util.h:240
void AddIndexed(T input, std::size_t index)
Add input to the sum.
Definition: Util.h:231
Roo1DTable implements a one-dimensional table.
Definition: Roo1DTable.h:23
void fill(RooAbsCategory &cat, double weight=1.0) override
Increment the counter of the table slot with the name corresponding to that of the current category s...
Definition: Roo1DTable.cxx:102
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:71
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=nullptr, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:805
bool recursiveRedirectServers(const RooAbsCollection &newServerList, bool mustReplaceAll=false, bool nameChange=false, bool recurseInNewSet=true)
Recursively replace all servers with the new servers in newSet.
Definition: RooAbsArg.cxx:1165
void attachArgs(const RooAbsCollection &set)
Bind this node to objects in set.
Definition: RooAbsArg.cxx:1609
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:293
void SetName(const char *name) override
Set the name of the TNamed.
Definition: RooAbsArg.cxx:2314
const RefCountList_t & valueClients() const
List of all value clients of this object. Value clients receive value updates.
Definition: RooAbsArg.h:188
bool redirectServers(const RooAbsCollection &newServerList, bool mustReplaceAll=false, bool nameChange=false, bool isRecursionStep=false)
Replace all direct servers of this object with the new servers in newServerList.
Definition: RooAbsArg.cxx:999
RooArgSet * getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
Definition: RooAbsArg.cxx:2057
virtual bool isDerived() const
Does value or shape of this arg depend on any other arg?
Definition: RooAbsArg.h:91
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition: RooAbsArg.h:83
void leafNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=nullptr, bool recurseNonDerived=false) const
Fill supplied list with all leaf nodes of the arg tree, starting with ourself as top node.
Definition: RooAbsArg.cxx:472
void attachDataSet(const RooAbsData &set)
Replace server nodes with names matching the dataset variable names with those data set variables,...
Definition: RooAbsArg.cxx:1624
RooAbsBinning is the abstract base class for RooRealVar binning definitions.
Definition: RooAbsBinning.h:25
virtual double averageBinWidth() const =0
virtual bool isUniform() const
Definition: RooAbsBinning.h:70
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
A space to attach TBranches.
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.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
bool allInRange(const char *rangeSpec) const
Return true if all contained object report to have their value inside the specified range.
bool empty() const
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Storage_t::size_type size() const
virtual RooAbsArg * addClone(const RooAbsArg &var, bool silent=false)
Add a clone of the specified argument to list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
RooAbsDataStore is the abstract base class for data collection that use a TTree as internal storage m...
virtual RooAbsData::CategorySpans getCategoryBatches(std::size_t, std::size_t) const
virtual RooAbsData::RealSpans getBatches(std::size_t first, std::size_t len) const =0
Retrieve batches for all observables in this data store.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:62
virtual double weight() const =0
virtual double sumEntries() const =0
Return effective number of entries in dataset, i.e., sum all weights.
virtual const RooArgSet * get() const
Definition: RooAbsData.h:106
RooRealVar * meanVar(const RooRealVar &var, const char *cutSpec=nullptr, const char *cutRange=nullptr) const
Create a RooRealVar containing the mean of observable 'var' in this dataset.
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Interface for detailed printing of object.
Definition: RooAbsData.cxx:934
const TNamed * _namePtr
! De-duplicated name pointer. This will be equal for all objects with the same name.
Definition: RooAbsData.h:371
RooAbsData()
Default constructor.
Definition: RooAbsData.cxx:172
static void setDefaultStorageType(StorageType s)
Definition: RooAbsData.cxx:132
void SetName(const char *name) override
Set the name of the TNamed.
CategorySpans getCategoryBatches(std::size_t first=0, std::size_t len=std::numeric_limits< std::size_t >::max()) const
RooRealVar * dataRealVar(const char *methodname, const RooRealVar &extVar) const
Internal method to check if given RooRealVar maps to a RooRealVar in this dataset.
virtual Roo1DTable * table(const RooArgSet &catSet, const char *cuts="", const char *opts="") const
Construct table for product of categories in catSet.
Definition: RooAbsData.cxx:884
void setGlobalObservables(RooArgSet const &globalObservables)
Sets the global observables stored in this data.
RooAbsDataStore * store()
Definition: RooAbsData.h:82
virtual RooAbsData * emptyClone(const char *newName=nullptr, const char *newTitle=nullptr, const RooArgSet *vars=nullptr, const char *wgtVarName=nullptr) const =0
std::map< RooFit::Detail::DataKey, RooSpan< const double > > RealSpans
Definition: RooAbsData.h:136
void printClassName(std::ostream &os) const override
Print class name of dataset.
Definition: RooAbsData.cxx:927
virtual void reset()
Definition: RooAbsData.cxx:381
RooRealVar * rmsVar(const RooRealVar &var, const char *cutSpec=nullptr, const char *cutRange=nullptr) const
Create a RooRealVar containing the RMS of observable 'var' in this dataset.
double standMoment(const RooRealVar &var, double order, const char *cutSpec=nullptr, const char *cutRange=nullptr) const
Calculate standardized moment.
Definition: RooAbsData.cxx:956
void Draw(Option_t *option="") override
Forward draw command to data store.
virtual bool changeObservableName(const char *from, const char *to)
Definition: RooAbsData.cxx:354
void printTitle(std::ostream &os) const override
Print title of dataset.
Definition: RooAbsData.cxx:919
void RecursiveRemove(TObject *obj) override
If one of the TObject we have a referenced to is deleted, remove the reference.
virtual double weightError(ErrorType=Poisson) const
Return the symmetric error on the current weight.
Definition: RooAbsData.h:117
void setDirtyProp(bool flag)
Control propagation of dirty flags from observables in dataset.
Definition: RooAbsData.cxx:429
virtual void add(const RooArgSet &row, double weight=1, double weightError=0.0)=0
virtual TH1 * fillHistogram(TH1 *hist, const RooArgList &plotVars, const char *cuts="", const char *cutRange=nullptr) 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.
TClass * IsA() const override
Definition: RooAbsData.h:378
void checkInit() const
virtual void setArgStatus(const RooArgSet &set, bool active)
Definition: RooAbsData.cxx:421
virtual void cacheArgs(const RooAbsArg *owner, RooArgSet &varSet, const RooArgSet *nset=nullptr, bool skipZeroWeights=false)
Internal method – Cache given set of functions with data.
Definition: RooAbsData.cxx:397
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...
RealSpans getBatches(std::size_t first=0, std::size_t len=std::numeric_limits< std::size_t >::max()) const
Write information to retrieve data columns into evalData.spans.
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
static StorageType defaultStorageType
Definition: RooAbsData.h:328
double corrcov(const RooRealVar &x, const RooRealVar &y, const char *cutSpec, const char *cutRange, bool corr) const
Internal method to calculate single correlation and covariance elements.
bool allClientsCached(RooAbsArg *, const RooArgSet &)
Utility function that determines if all clients of object 'var' appear in given list of cached nodes.
std::unique_ptr< RooAbsDataStore > _dstore
Data storage implementation.
Definition: RooAbsData.h:365
double sigma(const RooRealVar &var, const char *cutSpec=nullptr, const char *cutRange=nullptr) const
Definition: RooAbsData.h:258
static TClass * Class()
void addOwnedComponent(const char *idxlabel, RooAbsData &data)
virtual void fill()
Definition: RooAbsData.cxx:367
RooArgSet _vars
Dimensions of this data set.
Definition: RooAbsData.h:362
bool canSplitFast() const
static bool releaseVars(RooAbsData *)
If return value is true variables can be deleted.
Definition: RooAbsData.cxx:159
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)
void copyGlobalObservables(const RooAbsData &other)
Definition: RooAbsData.cxx:311
virtual bool isNonPoissonWeighted() const
Definition: RooAbsData.h:161
bool hasFilledCache() const
double sumEntriesW2() const
Return sum of squared weights of this data.
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:414
std::map< RooFit::Detail::DataKey, RooSpan< const RooAbsCategory::value_type > > CategorySpans
Definition: RooAbsData.h:137
void convertToVectorStore()
Convert tree-based storage to vector-based storage.
Definition: RooAbsData.cxx:344
virtual RooSpan< const double > getWeightBatch(std::size_t first, std::size_t len, bool sumW2=false) const =0
Return event weights of all events in range [first, first+len).
bool getRange(const RooAbsRealLValue &var, double &lowest, double &highest, double marginFrac=0.0, bool symMode=false) const
Fill Doubles 'lowest' and 'highest' with the lowest and highest value of observable 'var' in this dat...
RooArgSet _cachedVars
! External variables cached with this data set
Definition: RooAbsData.h:363
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Definition: RooAbsData.cxx:374
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:602
double moment(const RooRealVar &var, double order, const char *cutSpec=nullptr, const char *cutRange=nullptr) const
Calculate moment of requested order.
Definition: RooAbsData.cxx:975
StorageType storageType
Definition: RooAbsData.h:330
RooAbsData & operator=(const RooAbsData &other)
Definition: RooAbsData.cxx:268
void SetNameTitle(const char *name, const char *title) override
Set all the TNamed parameters (name and title).
virtual void resetCache()
Internal method – Remove cached function values.
Definition: RooAbsData.cxx:405
virtual TList * split(const RooAbsCategory &splitCat, bool createEmptyDataSets=false) const
Split dataset into subsets based on states of given splitCat in this dataset.
Int_t defaultPrintContents(Option_t *opt) const override
Define default print options, for a given print style.
Definition: RooAbsData.cxx:942
std::unique_ptr< RooArgSet > _globalObservables
Snapshot of global observables.
Definition: RooAbsData.h:369
virtual double weightSquared() const =0
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:449
void attachBuffers(const RooArgSet &extObs)
virtual RooAbsData * reduceEng(const RooArgSet &varSubset, const RooFormulaVar *cutVar, const char *cutRange=nullptr, std::size_t nStart=0, std::size_t=std::numeric_limits< std::size_t >::max())=0
std::map< std::string, RooAbsData * > _ownedComponents
Owned external components.
Definition: RooAbsData.h:367
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:689
static StorageType getDefaultStorageType()
Definition: RooAbsData.cxx:143
virtual bool isWeighted() const
Definition: RooAbsData.h:157
void Streamer(TBuffer &) override
Stream an object of class RooAbsData.
void resetBuffers()
TMatrixDSym * corrcovMatrix(const RooArgList &vars, const char *cutSpec, const char *cutRange, bool corr) const
Return covariance matrix from data for given list of observables.
void printName(std::ostream &os) const override
Print name of dataset.
Definition: RooAbsData.cxx:911
const TTree * tree() const
Return a pointer to the TTree which stores the data.
~RooAbsData() override
Destructor.
Definition: RooAbsData.cxx:325
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
virtual Int_t getBins(const char *name=nullptr) const
Get number of bins of currently defined range.
virtual const RooAbsBinning & getBinning(const char *name=nullptr, bool verbose=true, bool createOnTheFly=false) const =0
Retrive binning configuration with given name or default binning.
virtual double getMax(const char *name=nullptr) const
Get maximum of currently defined range.
virtual double getMin(const char *name=nullptr) const
Get minimum 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:62
double 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:415
void setPlotLabel(const char *label)
Set the label associated with this variable.
Definition: RooAbsReal.cxx:425
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:56
RooArgSet * snapshot(bool deepCopy=true) const
Use RooAbsCollection::snapshot(), but return as RooArgSet.
Definition: RooArgSet.h:179
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:28
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:26
double getDouble(Int_t idx) const
Return double stored in slot idx.
Definition: RooCmdArg.h:91
Int_t getInt(Int_t idx) const
Definition: RooCmdArg.h:86
TObject * Clone(const char *newName=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition: RooCmdArg.h:57
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
Definition: RooCmdConfig.h:31
static void stripCmdList(RooLinkedList &cmdList, const char *cmdsToPurge)
Utility function that strips command names listed (comma separated) in cmdsToPurge from cmdList.
RooCompositeDataStore combines several disjunct datasets into one.
static TClass * Class()
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:33
double eval(const RooArgSet *nset=nullptr) const
Evalute all parameters/observables, and then evaluate formula.
Definition: RooFormula.cxx:431
A RooHist is a graphical representation of binned data based on the TGraphAsymmErrors class.
Definition: RooHist.h:29
static TClass * Class()
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:38
bool Replace(const TObject *oldArg, const TObject *newArg)
Replace object 'oldArg' in collection with new object 'newArg'.
void Delete(Option_t *o=nullptr) override
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
TObject * FindObject(const char *name) const override
Return pointer to obejct with given name.
RooMultiCategory connects several RooAbsCategory objects into a single category.
RooNameReg is a registry for const char* names.
Definition: RooNameReg.h:25
const TNamed * constPtr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:60
static RooNameReg & instance()
Return reference to singleton instance.
Definition: RooNameReg.cxx:50
@ kRenamedArg
TNamed flag to indicate that some RooAbsArg has been renamed (flag set in new name)
Definition: RooNameReg.h:44
static void incrementRenameCounter()
The renaming counter has to be incremented every time a RooAbsArg is renamed.
Definition: RooNameReg.cxx:106
static double infinity()
Return internal infinity representation.
Definition: RooNumber.cxx:48
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:43
TObject * findObject(const char *name, const TClass *clas=nullptr) const
Find the named object in our list of items and return a pointer to it.
Definition: RooPlot.cxx:984
void addObject(TObject *obj, Option_t *drawOptions="", bool invisible=false)
Add a generic object to this plot.
Definition: RooPlot.cxx:410
double getFitRangeNEvt() const
Return the number of events in the fit range.
Definition: RooPlot.h:142
TAttLine * getAttLine(const char *name=nullptr) const
Return a pointer to the line attributes of the named object in this plot, or zero if the named object...
Definition: RooPlot.cxx:849
TAttFill * getAttFill(const char *name=nullptr) const
Return a pointer to the fill attributes of the named object in this plot, or zero if the named object...
Definition: RooPlot.cxx:859
RooAbsRealLValue * getPlotVar() const
Definition: RooPlot.h:140
TAttMarker * getAttMarker(const char *name=nullptr) const
Return a pointer to the marker attributes of the named object in this plot, or zero if the named obje...
Definition: RooPlot.cxx:869
TAxis * GetXaxis() const
Definition: RooPlot.cxx:1283
void updateNormVars(const RooArgSet &vars)
Install the given set of observables are reference normalization variables for this frame.
Definition: RooPlot.cxx:368
Int_t GetNbinsX() const
Definition: RooPlot.cxx:1287
void addPlotable(RooPlotable *plotable, Option_t *drawOptions="", bool invisible=false, bool refreshNorm=false)
Add the specified plotable object to our plot.
Definition: RooPlot.cxx:558
double getFitRangeBinW() const
Return the bin width that is being used to normalise the PDF.
Definition: RooPlot.h:145
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:40
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
const RooAbsCategoryLValue & indexCat() const
RooAbsPdf * getPdf(const char *catName) const
Return the p.d.f associated with the given index category name.
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
constexpr std::span< T >::index_type size() const noexcept
Definition: RooSpan.h:121
constexpr bool empty() const noexcept
Definition: RooSpan.h:125
The RooStringView is a wrapper around a C-syle string that can also be constructed from a std::string...
Definition: RooStringView.h:27
static void destroy(const TObject *obj)
Register deletion of object 'obj'.
Definition: RooTrace.cxx:79
static void create(const TObject *obj)
Register creation of object 'obj'.
Definition: RooTrace.cxx:66
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:135
Double_t GetXmin() const
Definition: TAxis.h:134
Buffer base class used for serializing objects.
Definition: TBuffer.h:43
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=nullptr)=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:8930
virtual Int_t GetDimension() const
Definition: TH1.h:281
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:9073
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3338
virtual TArrayD * GetSumw2()
Definition: TH1.h:311
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:3668
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition: TH1.cxx:8887
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:257
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
A doubly linked list.
Definition: TList.h:38
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition: TList.cxx:578
void Add(TObject *obj) override
Definition: TList.h:81
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
TObject * Clone(const char *newname="") const override
Make a clone of an object using the Streamer facility.
Definition: TNamed.cxx:74
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition: TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition: TNamed.h:48
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition: TNamed.cxx:140
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
Mother of all ROOT objects.
Definition: TObject.h:41
virtual void Clear(Option_t *="")
Definition: TObject.h:119
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:200
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:768
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:519
A Pave (see TPave) with text, lines or/and boxes inside.
Definition: TPaveText.h:21
Basic string class.
Definition: TString.h:136
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 AutoBinning(Int_t nbins=100, double marginFactor=0.1)
RooCmdArg AxisLabel(const char *name)
RooCmdArg Binning(const RooAbsBinning &binning)
RVec< PromoteType< T > > abs(const RVec< T > &v)
Definition: RVec.hxx:1756
Double_t y[n]
Definition: legend1.C:17
Double_t x[n]
Definition: legend1.C:17
static double C[]
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
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
static Roo_reg_AGKInteg1D instance
@ Optimization
Definition: RooGlobalFunc.h:62
@ DataHandling
Definition: RooGlobalFunc.h:63
@ InputArguments
Definition: RooGlobalFunc.h:62
static constexpr double s
static constexpr double pc
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Definition: TMath.h:719
Definition: first.py:1
const double xbins[xbins_n]
Definition: graph.py:1
static const char * what
Definition: stlLoader.cc:6
const char * cuts
Definition: RooAbsData.h:187
const char * cutRange
Definition: RooAbsData.h:191
const char * histName
Definition: RooAbsData.h:192
const char * addToHistName
Definition: RooAbsData.h:194
RooAbsData::ErrorType etype
Definition: RooAbsData.h:190
RooAbsBinning * bins
Definition: RooAbsData.h:189
Option_t * drawOptions
Definition: RooAbsData.h:188
TLine l
Definition: textangle.C:4
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345