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