Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
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 RooFormulaVar cutVar(cut,cut,*get()) ;
477 auto ret = reduceEng(*get(),&cutVar,nullptr,0,std::numeric_limits<std::size_t>::max()) ;
478 ret->copyGlobalObservables(*this);
479 return RooFit::makeOwningPtr(std::move(ret));
480}
481
482////////////////////////////////////////////////////////////////////////////////
483/// Create a subset of the data set by applying the given cut on the data points.
484/// The 'cutVar' formula variable is used to select the subset of data points to be
485/// retained in the reduced data collection.
486
488{
489 auto ret = reduceEng(*get(),&cutVar,nullptr,0,std::numeric_limits<std::size_t>::max()) ;
490 ret->copyGlobalObservables(*this);
491 return RooFit::makeOwningPtr(std::move(ret));
492}
493
494////////////////////////////////////////////////////////////////////////////////
495/// Create a subset of the data set by applying the given cut on the data points
496/// and reducing the dimensions to the specified set.
497///
498/// The cut expression can refer to any variable in the data set. For cuts involving
499/// other variables, such as intermediate formula objects, use the equivalent
500/// reduce method specifying the as a RooFormulVar reference.
501
503{
504 // Make sure varSubset doesn't contain any variable not in this dataset
506 for (const auto arg : varSubset) {
507 if (!_vars.find(arg->GetName())) {
508 coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
509 << arg->GetName() << " not in dataset, ignored" << std::endl ;
510 varSubset2.remove(*arg) ;
511 }
512 }
513
514 std::unique_ptr<RooAbsData> ret;
515 if (cut && strlen(cut)>0) {
516 RooFormulaVar cutVar(cut, cut, *get(), false);
517 ret = reduceEng(varSubset2,&cutVar,nullptr,0,std::numeric_limits<std::size_t>::max());
518 } else {
519 ret = reduceEng(varSubset2,nullptr,nullptr,0,std::numeric_limits<std::size_t>::max());
520 }
521 ret->copyGlobalObservables(*this);
522 return RooFit::makeOwningPtr(std::move(ret));
523}
524
525////////////////////////////////////////////////////////////////////////////////
526/// Create a subset of the data set by applying the given cut on the data points
527/// and reducing the dimensions to the specified set.
528///
529/// The 'cutVar' formula variable is used to select the subset of data points to be
530/// retained in the reduced data collection.
531
533{
534 // Make sure varSubset doesn't contain any variable not in this dataset
536 for(RooAbsArg * arg : varSubset) {
537 if (!_vars.find(arg->GetName())) {
538 coutW(InputArguments) << "RooAbsData::reduce(" << GetName() << ") WARNING: variable "
539 << arg->GetName() << " not in dataset, ignored" << std::endl ;
540 varSubset2.remove(*arg) ;
541 }
542 }
543
544 auto ret = reduceEng(varSubset2,&cutVar,nullptr,0,std::numeric_limits<std::size_t>::max()) ;
545 ret->copyGlobalObservables(*this);
546 return RooFit::makeOwningPtr(std::move(ret));
547}
548
549
551 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
552 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8) const
553{
555 l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
556 l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
557 l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
558 l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
559 return plotOn(frame,l) ;
560}
561
562
564 const RooCmdArg& arg1, const RooCmdArg& arg2, const RooCmdArg& arg3, const RooCmdArg& arg4,
565 const RooCmdArg& arg5, const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8) const
566{
568 l.Add((TObject*)&arg1) ; l.Add((TObject*)&arg2) ;
569 l.Add((TObject*)&arg3) ; l.Add((TObject*)&arg4) ;
570 l.Add((TObject*)&arg5) ; l.Add((TObject*)&arg6) ;
571 l.Add((TObject*)&arg7) ; l.Add((TObject*)&arg8) ;
572
573 return createHistogram(name,xvar,l) ;
574}
575
576////////////////////////////////////////////////////////////////////////////////
577/// Create and fill a ROOT histogram TH1,TH2 or TH3 with the values of this
578/// dataset for the variables with given names.
579///
580/// \param[in] varNameList Comma-separated variable names.
581/// \param[in] binArgX Control the binning for the `x` variable.
582/// \param[in] binArgY Control the binning for the `y` variable.
583/// \param[in] binArgZ Control the binning for the `z` variable.
584/// \return Histogram now owned by user.
585///
586/// The possible binning command arguments for each axis are:
587///
588/// <table>
589/// <tr><td> `AutoBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin, set binning to nbins
590/// <tr><td> `AutoSymBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin,
591/// with additional constraint that mean of data is in center of range, set binning to nbins
592/// <tr><td> `Binning(const char* name)` <td> Apply binning with given name to x axis of histogram
593/// <tr><td> `Binning(RooAbsBinning& binning)` <td> Apply specified binning to x axis of histogram
594/// <tr><td> `Binning(int nbins, double lo, double hi)` <td> Apply specified binning to x axis of histogram
595///
596/// <tr><td> `YVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on y axis of ROOT histogram
597/// <tr><td> `ZVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on z axis of ROOT histogram
598/// </table>
599
601 const RooCmdArg& binArgX,
602 const RooCmdArg& binArgY,
603 const RooCmdArg& binArgZ) const
604{
605 // Parse list of variable names
606 const auto varNames = ROOT::Split(varNameList, ",:");
607 RooRealVar* vars[3] = {nullptr, nullptr, nullptr};
608
609 for (unsigned int i = 0; i < varNames.size(); ++i) {
610 if (i >= 3) {
611 coutW(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): Can only create 3-dimensional histograms. Variable "
612 << i << " " << varNames[i] << " unused." << std::endl;
613 continue;
614 }
615
616 vars[i] = static_cast<RooRealVar*>(get()->find(varNames[i].data()) );
617 if (!vars[i]) {
618 coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << ") ERROR: dataset does not contain an observable named " << varNames[i] << std::endl;
619 return nullptr;
620 }
621 }
622
623 if (!vars[0]) {
624 coutE(InputArguments) << "RooAbsData::createHistogram(" << GetName() << "): No variable to be histogrammed in list '" << varNameList << "'" << std::endl;
625 return nullptr;
626 }
627
628 // Fill command argument list
629 RooLinkedList argList;
630 argList.Add(binArgX.Clone());
631 if (vars[1]) {
632 argList.Add(RooFit::YVar(*vars[1],binArgY).Clone());
633 }
634 if (vars[2]) {
635 argList.Add(RooFit::ZVar(*vars[2],binArgZ).Clone());
636 }
637
638 // Call implementation function
639 TH1* result = createHistogram(GetName(), *vars[0], argList);
640
641 // Delete temporary list of RooCmdArgs
642 argList.Delete() ;
643
644 return result ;
645}
646
647////////////////////////////////////////////////////////////////////////////////
648///
649/// This function accepts the following arguments
650///
651/// \param[in] name Name of the ROOT histogram
652/// \param[in] xvar Observable to be mapped on x axis of ROOT histogram
653/// \param[in] argListIn list of input arguments
654/// \return Histogram now owned by user.
655///
656/// <table>
657/// <tr><td> `AutoBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin, set binning to nbins
658/// <tr><td> `AutoSymBinning(Int_t nbins, Double_y margin)` <td> Automatically calculate range with given added fractional margin,
659/// with additional constraint that mean of data is in center of range, set binning to nbins
660/// <tr><td> `Binning(const char* name)` <td> Apply binning with given name to x axis of histogram
661/// <tr><td> `Binning(RooAbsBinning& binning)` <td> Apply specified binning to x axis of histogram
662/// <tr><td> `Binning(int nbins, double lo, double hi)` <td> Apply specified binning to x axis of histogram
663///
664/// <tr><td> `YVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on y axis of ROOT histogram
665/// <tr><td> `ZVar(const RooAbsRealLValue& var,...)` <td> Observable to be mapped on z axis of ROOT histogram
666/// </table>
667///
668/// 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.
669/// ```
670/// createHistogram("histo",x,Binning(-1,1,20), YVar(y,Binning(-1,1,30)), ZVar(z,Binning("zbinning")))
671/// ```
672///
673/// The caller takes ownership of the returned histogram
674
676{
677 RooLinkedList argList(argListIn) ;
678
679 // Define configuration for this method
680 RooCmdConfig pc("RooAbsData::createHistogram(" + std::string(GetName()) + ")");
681 pc.defineString("cutRange","CutRange",0,"",true) ;
682 pc.defineString("cutString","CutSpec",0,"") ;
683 pc.defineObject("yvar","YVar",0,nullptr) ;
684 pc.defineObject("zvar","ZVar",0,nullptr) ;
685 pc.allowUndefined() ;
686
687 // Process & check varargs
688 pc.process(argList) ;
689 if (!pc.ok(true)) {
690 return nullptr;
691 }
692
693 const char* cutSpec = pc.getString("cutString",nullptr,true) ;
694 const char* cutRange = pc.getString("cutRange",nullptr,true) ;
695
696 RooArgList vars(xvar) ;
697 RooAbsArg* yvar = static_cast<RooAbsArg*>(pc.getObject("yvar")) ;
698 if (yvar) {
699 vars.add(*yvar) ;
700 }
701 RooAbsArg* zvar = static_cast<RooAbsArg*>(pc.getObject("zvar")) ;
702 if (zvar) {
703 vars.add(*zvar) ;
704 }
705
706 RooCmdConfig::stripCmdList(argList,"CutRange,CutSpec") ;
707
708 // Swap Auto(Sym)RangeData with a Binning command
710 RooCmdArg* autoRD = static_cast<RooCmdArg*>(argList.find("AutoRangeData")) ;
711 if (autoRD) {
712 double xmin;
713 double xmax;
714 if (!getRange(static_cast<RooRealVar const&>(xvar),xmin,xmax,autoRD->getDouble(0),autoRD->getInt(0))) {
715 RooCmdArg* bincmd = static_cast<RooCmdArg*>(RooFit::Binning(autoRD->getInt(1),xmin,xmax).Clone()) ;
716 ownedCmds.Add(bincmd) ;
717 argList.Replace(autoRD,bincmd) ;
718 }
719 }
720
721 if (yvar) {
722 std::unique_ptr<RooCmdArg> autoRDY{static_cast<RooCmdArg*>((static_cast<RooCmdArg*>(argList.find("YVar")))->subArgs().find("AutoRangeData"))};
723 if (autoRDY) {
724 double ymin;
725 double ymax;
726 if (!getRange(static_cast<RooRealVar &>(*yvar), ymin, ymax, autoRDY->getDouble(0), autoRDY->getInt(0))) {
727 RooCmdArg *bincmd = static_cast<RooCmdArg *>(RooFit::Binning(autoRDY->getInt(1), ymin, ymax).Clone());
728 // ownedCmds.Add(bincmd) ;
729 (static_cast<RooCmdArg *>(argList.find("YVar")))->subArgs().Replace(autoRDY.get(), bincmd);
730 }
731 }
732 }
733
734 if (zvar) {
735 std::unique_ptr<RooCmdArg> autoRDZ{static_cast<RooCmdArg*>((static_cast<RooCmdArg*>(argList.find("ZVar")))->subArgs().find("AutoRangeData"))};
736 if (autoRDZ) {
737 double zmin;
738 double zmax;
739 if (!getRange(static_cast<RooRealVar&>(*zvar),zmin,zmax,autoRDZ->getDouble(0),autoRDZ->getInt(0))) {
740 RooCmdArg* bincmd = static_cast<RooCmdArg*>(RooFit::Binning(autoRDZ->getInt(1),zmin,zmax).Clone()) ;
741 //ownedCmds.Add(bincmd) ;
742 (static_cast<RooCmdArg*>(argList.find("ZVar")))->subArgs().Replace(autoRDZ.get(),bincmd) ;
743 }
744 }
745 }
746
747
748 TH1* histo = xvar.createHistogram(name,argList) ;
749 fillHistogram(histo,vars,cutSpec,cutRange) ;
750
751 ownedCmds.Delete() ;
752
753 return histo ;
754}
755
756////////////////////////////////////////////////////////////////////////////////
757/// Construct table for product of categories in catSet
758
759Roo1DTable* RooAbsData::table(const RooArgSet& catSet, const char* cuts, const char* opts) const
760{
762
763 std::string prodName("(") ;
764 for(auto * arg : catSet) {
765 if (dynamic_cast<RooAbsCategory*>(arg)) {
766 if (auto varsArg = dynamic_cast<RooAbsCategory*>(_vars.find(arg->GetName()))) catSet2.add(*varsArg) ;
767 else catSet2.add(*arg) ;
768 if (prodName.length()>1) {
769 prodName += " x " ;
770 }
771 prodName += arg->GetName() ;
772 } else {
773 coutW(InputArguments) << "RooAbsData::table(" << GetName() << ") non-RooAbsCategory input argument " << arg->GetName() << " ignored" << std::endl ;
774 }
775 }
776 prodName += ")" ;
777
778 RooMultiCategory tmp(prodName.c_str(),prodName.c_str(),catSet2) ;
779 return table(tmp,cuts,opts) ;
780}
781
782////////////////////////////////////////////////////////////////////////////////
783/// Print name of dataset
784
785void RooAbsData::printName(std::ostream& os) const
786{
787 os << GetName() ;
788}
789
790////////////////////////////////////////////////////////////////////////////////
791/// Print title of dataset
792
793void RooAbsData::printTitle(std::ostream& os) const
794{
795 os << GetTitle() ;
796}
797
798////////////////////////////////////////////////////////////////////////////////
799/// Print class name of dataset
800
801void RooAbsData::printClassName(std::ostream& os) const
802{
803 os << ClassName() ;
804}
805
806////////////////////////////////////////////////////////////////////////////////
807
808void RooAbsData::printMultiline(std::ostream& os, Int_t contents, bool verbose, TString indent) const
809{
810 _dstore->printMultiline(os,contents,verbose,indent) ;
811}
812
813////////////////////////////////////////////////////////////////////////////////
814/// Define default print options, for a given print style
815
820
821////////////////////////////////////////////////////////////////////////////////
822/// Calculate standardized moment.
823///
824/// \param[in] var Variable to be used for calculating the moment.
825/// \param[in] order Order of the moment.
826/// \param[in] cutSpec If specified, the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
827/// \param[in] cutRange If specified, calculate inside the range named 'cutRange' (also applies cut spec)
828/// \return \f$ \frac{\left< \left( X - \left< X \right> \right)^n \right>}{\sigma^n} \f$, where n = order.
829
830double RooAbsData::standMoment(const RooRealVar &var, double order, const char* cutSpec, const char* cutRange) const
831{
832 // Hardwire invariant answer for first and second moment
833 if (order==1) return 0 ;
834 if (order==2) return 1 ;
835
836 return moment(var,order,cutSpec,cutRange) / std::pow(sigma(var,cutSpec,cutRange),order) ;
837}
838
839////////////////////////////////////////////////////////////////////////////////
840/// Calculate moment of requested order.
841///
842/// \param[in] var Variable to be used for calculating the moment.
843/// \param[in] order Order of the moment.
844/// \param[in] cutSpec If specified, the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
845/// \param[in] cutRange If specified, calculate inside the range named 'cutRange' (also applies cut spec)
846/// \return \f$ \left< \left( X - \left< X \right> \right)^n \right> \f$ of order \f$n\f$.
847///
848
849double RooAbsData::moment(const RooRealVar& var, double order, const char* cutSpec, const char* cutRange) const
850{
851 double offset = order>1 ? moment(var,1,cutSpec,cutRange) : 0 ;
852 return moment(var,order,offset,cutSpec,cutRange) ;
853
854}
855
856////////////////////////////////////////////////////////////////////////////////
857/// Return the 'order'-ed moment of observable 'var' in this dataset. If offset is non-zero it is subtracted
858/// from the values of 'var' prior to the moment calculation. If cutSpec and/or cutRange are specified
859/// the moment is calculated on the subset of the data which pass the C++ cut specification expression 'cutSpec'
860/// and/or are inside the range named 'cutRange'
861
862double RooAbsData::moment(const RooRealVar& var, double order, double offset, const char* cutSpec, const char* cutRange) const
863{
864 // Lookup variable in dataset
865 auto arg = _vars.find(var.GetName());
866 if (!arg) {
867 coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") ERROR: unknown variable: " << var.GetName() << std::endl;
868 return 0;
869 }
870
871 auto varPtr = dynamic_cast<const RooRealVar*>(arg);
872 // Check if found variable is of type RooRealVar
873 if (!varPtr) {
874 coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") ERROR: variable " << var.GetName() << " is not of type RooRealVar" << std::endl ;
875 return 0;
876 }
877
878 // Check if dataset is not empty
879 if(sumEntries(cutSpec, cutRange) == 0.) {
880 coutE(InputArguments) << "RooDataSet::moment(" << GetName() << ") WARNING: empty dataset" << std::endl ;
881 return 0;
882 }
883
884 // Setup RooFormulaVar for cutSpec if it is present
885 std::unique_ptr<RooFormula> select;
886 if (cutSpec) {
887 select = std::make_unique<RooFormula>("select",cutSpec,*get());
888 }
889
890
891 // Calculate requested moment
893 for(int index= 0; index < numEntries(); index++) {
894 const RooArgSet* vars = get(index) ;
895 if (select && select->eval()==0) continue ;
896 if (cutRange && vars->allInRange(cutRange)) continue ;
897
898 sum += weight() * std::pow(varPtr->getVal() - offset,order);
899 }
900
901 return sum.Sum()/sumEntries(cutSpec, cutRange);
902}
903
904////////////////////////////////////////////////////////////////////////////////
905/// Internal method to check if given RooRealVar maps to a RooRealVar in this dataset
906
908{
909 // Lookup variable in dataset
910 RooRealVar *xdata = static_cast<RooRealVar*>(_vars.find(extVar.GetName()));
911 if(!xdata) {
912 coutE(InputArguments) << "RooDataSet::" << methodname << "(" << GetName() << ") ERROR: variable : " << extVar.GetName() << " is not in data" << std::endl ;
913 return nullptr;
914 }
915 // Check if found variable is of type RooRealVar
916 if (!dynamic_cast<RooRealVar*>(xdata)) {
917 coutE(InputArguments) << "RooDataSet::" << methodname << "(" << GetName() << ") ERROR: variable : " << extVar.GetName() << " is not of type RooRealVar in data" << std::endl ;
918 return nullptr;
919 }
920 return xdata;
921}
922
923////////////////////////////////////////////////////////////////////////////////
924/// Internal method to calculate single correlation and covariance elements
925
926double RooAbsData::corrcov(const RooRealVar &x, const RooRealVar &y, const char* cutSpec, const char* cutRange, bool corr) const
927{
928 // Lookup variable in dataset
929 RooRealVar *xdata = dataRealVar(corr?"correlation":"covariance",x) ;
930 RooRealVar *ydata = dataRealVar(corr?"correlation":"covariance",y) ;
931 if (!xdata||!ydata) return 0 ;
932
933 // Check if dataset is not empty
934 if(sumEntries(cutSpec, cutRange) == 0.) {
935 coutW(InputArguments) << "RooDataSet::" << (corr?"correlation":"covariance") << "(" << GetName() << ") WARNING: empty dataset, returning zero" << std::endl ;
936 return 0;
937 }
938
939 // Setup RooFormulaVar for cutSpec if it is present
940 std::unique_ptr<RooFormula> select;
941 if (cutSpec) select = std::make_unique<RooFormula>("select",cutSpec,*get());
942
943 // Calculate requested moment
944 double xysum(0);
945 double xsum(0);
946 double ysum(0);
947 double x2sum(0);
948 double y2sum(0);
949 const RooArgSet* vars ;
950 for(int index= 0; index < numEntries(); index++) {
951 vars = get(index) ;
952 if (select && select->eval()==0) continue ;
953 if (cutRange && vars->allInRange(cutRange)) continue ;
954
955 xysum += weight()*xdata->getVal()*ydata->getVal() ;
956 xsum += weight()*xdata->getVal() ;
957 ysum += weight()*ydata->getVal() ;
958 if (corr) {
959 x2sum += weight()*xdata->getVal()*xdata->getVal() ;
960 y2sum += weight()*ydata->getVal()*ydata->getVal() ;
961 }
962 }
963
964 // Normalize entries
965 xysum/=sumEntries(cutSpec, cutRange) ;
966 xsum/=sumEntries(cutSpec, cutRange) ;
967 ysum/=sumEntries(cutSpec, cutRange) ;
968 if (corr) {
969 x2sum/=sumEntries(cutSpec, cutRange) ;
970 y2sum/=sumEntries(cutSpec, cutRange) ;
971 }
972
973 // Return covariance or correlation as requested
974 if (corr) {
975 return (xysum-xsum*ysum)/(sqrt(x2sum-(xsum*xsum))*sqrt(y2sum-(ysum*ysum))) ;
976 } else {
977 return (xysum-xsum*ysum);
978 }
979}
980
981////////////////////////////////////////////////////////////////////////////////
982/// Return covariance matrix from data for given list of observables
983
984RooFit::OwningPtr<TMatrixDSym> RooAbsData::corrcovMatrix(const RooArgList& vars, const char* cutSpec, const char* cutRange, bool corr) const
985{
987 for(auto * var : static_range_cast<RooRealVar*>(vars)) {
988 RooRealVar* datavar = dataRealVar("covarianceMatrix",*var) ;
989 if (!datavar) {
990 return nullptr;
991 }
992 varList.add(*datavar) ;
993 }
994
995
996 // Check if dataset is not empty
997 if(sumEntries(cutSpec, cutRange) == 0.) {
998 coutW(InputArguments) << "RooDataSet::covariance(" << GetName() << ") WARNING: empty dataset, returning zero" << std::endl ;
999 return nullptr;
1000 }
1001
1002 // Setup RooFormulaVar for cutSpec if it is present
1003 std::unique_ptr<RooFormula> select = cutSpec ? std::make_unique<RooFormula>("select",cutSpec,*get()) : nullptr;
1004
1005 TMatrixDSym xysum(varList.size()) ;
1006 std::vector<double> xsum(varList.size()) ;
1007 std::vector<double> x2sum(varList.size()) ;
1008
1009 // Calculate <x_i> and <x_i y_j>
1010 for(int index= 0; index < numEntries(); index++) {
1011 const RooArgSet* dvars = get(index) ;
1012 if (select && select->eval()==0) continue ;
1013 if (cutRange && dvars->allInRange(cutRange)) continue ;
1014
1015 for(std::size_t iX = 0; iX < varList.size(); ++iX) {
1016 auto varx = static_cast<RooRealVar const&>(varList[iX]);
1017 xsum[iX] += weight() * varx.getVal() ;
1018 if (corr) {
1019 x2sum[iX] += weight() * varx.getVal() * varx.getVal();
1020 }
1021
1022 for(std::size_t iY = iX; iY < varList.size(); ++iY) {
1023 auto vary = static_cast<RooRealVar const&>(varList[iY]);
1024 xysum(iX,iY) += weight() * varx.getVal() * vary.getVal();
1025 xysum(iY,iX) = xysum(iX,iY) ;
1026 }
1027 }
1028
1029 }
1030
1031 // Normalize sums
1032 for (std::size_t iX=0 ; iX<varList.size() ; iX++) {
1033 xsum[iX] /= sumEntries(cutSpec, cutRange) ;
1034 if (corr) {
1035 x2sum[iX] /= sumEntries(cutSpec, cutRange) ;
1036 }
1037 for (std::size_t iY=0 ; iY<varList.size() ; iY++) {
1038 xysum(iX,iY) /= sumEntries(cutSpec, cutRange) ;
1039 }
1040 }
1041
1042 // Calculate covariance matrix
1043 auto C = std::make_unique<TMatrixDSym>(varList.size()) ;
1044 for (std::size_t iX=0 ; iX<varList.size() ; iX++) {
1045 for (std::size_t iY=0 ; iY<varList.size() ; iY++) {
1046 (*C)(iX,iY) = xysum(iX,iY)-xsum[iX]*xsum[iY] ;
1047 if (corr) {
1048 (*C)(iX,iY) /= std::sqrt((x2sum[iX]-(xsum[iX]*xsum[iX]))*(x2sum[iY]-(xsum[iY]*xsum[iY]))) ;
1049 }
1050 }
1051 }
1052
1053 return RooFit::makeOwningPtr(std::move(C));
1054}
1055
1056////////////////////////////////////////////////////////////////////////////////
1057/// Create a RooRealVar containing the mean of observable 'var' in
1058/// this dataset. If cutSpec and/or cutRange are specified the
1059/// moment is calculated on the subset of the data which pass the C++
1060/// cut specification expression 'cutSpec' and/or are inside the
1061/// range named 'cutRange'
1062
1063RooRealVar* RooAbsData::meanVar(const RooRealVar &var, const char* cutSpec, const char* cutRange) const
1064{
1065 // Create a new variable with appropriate strings. The error is calculated as
1066 // RMS/Sqrt(N) which is generally valid.
1067
1068 // Create holder variable for mean
1069 std::string name = std::string{var.GetName()} + "Mean";
1070 std::string title = std::string{"Mean of "} + var.GetTitle();
1071 auto *meanv= new RooRealVar(name.c_str(), title.c_str(), 0) ;
1072 meanv->setConstant(false) ;
1073
1074 // Adjust plot label
1075 std::string label = "<" + std::string{var.getPlotLabel()} + ">";
1076 meanv->setPlotLabel(label.c_str());
1077
1078 // fill in this variable's value and error
1079 double meanVal=moment(var,1,0,cutSpec,cutRange) ;
1080 double N(sumEntries(cutSpec,cutRange)) ;
1081
1082 double rmsVal= sqrt(moment(var,2,meanVal,cutSpec,cutRange)*N/(N-1));
1083 meanv->setVal(meanVal) ;
1084 meanv->setError(N > 0 ? rmsVal/sqrt(N) : 0);
1085
1086 return meanv;
1087}
1088
1089////////////////////////////////////////////////////////////////////////////////
1090/// Create a RooRealVar containing the RMS of observable 'var' in
1091/// this dataset. If cutSpec and/or cutRange are specified the
1092/// moment is calculated on the subset of the data which pass the C++
1093/// cut specification expression 'cutSpec' and/or are inside the
1094/// range named 'cutRange'
1095
1096RooRealVar* RooAbsData::rmsVar(const RooRealVar &var, const char* cutSpec, const char* cutRange) const
1097{
1098 // Create a new variable with appropriate strings. The error is calculated as
1099 // RMS/(2*Sqrt(N)) which is only valid if the variable has a Gaussian distribution.
1100
1101 // Create RMS value holder
1102 std::string name(var.GetName());
1103 std::string title("RMS of ");
1104 name += "RMS";
1105 title += var.GetTitle();
1106 auto *rms= new RooRealVar(name.c_str(), title.c_str(), 0) ;
1107 rms->setConstant(false) ;
1108
1109 // Adjust plot label
1110 std::string label(var.getPlotLabel());
1111 label += "_{RMS}";
1112 rms->setPlotLabel(label.c_str());
1113
1114 // Fill in this variable's value and error
1115 double meanVal(moment(var,1,0,cutSpec,cutRange)) ;
1116 double N(sumEntries(cutSpec, cutRange));
1117 double rmsVal= sqrt(moment(var,2,meanVal,cutSpec,cutRange)*N/(N-1));
1118 rms->setVal(rmsVal) ;
1119 rms->setError(rmsVal/sqrt(2*N));
1120
1121 return rms;
1122}
1123
1124////////////////////////////////////////////////////////////////////////////////
1125/// Add a box with statistics information to the specified frame. By default a box with the
1126/// event count, mean and rms of the plotted variable is added.
1127///
1128/// The following optional named arguments are accepted
1129/// <table>
1130/// <tr><td> `What(const char* whatstr)` <td> Controls what is printed: "N" = count, "M" is mean, "R" is RMS.
1131/// <tr><td> `Format(const char* optStr)` <td> \deprecated Classing parameter formatting options, provided for backward compatibility
1132///
1133/// <tr><td> `Format(const char* what,...)` <td> Parameter formatting options.
1134/// <table>
1135/// <tr><td> const char* what <td> Controls what is shown:
1136/// - "N" adds name
1137/// - "E" adds error
1138/// - "A" shows asymmetric error
1139/// - "U" shows unit
1140/// - "H" hides the value
1141/// <tr><td> `FixedPrecision(int n)` <td> Controls precision, set fixed number of digits
1142/// <tr><td> `AutoPrecision(int n)` <td> Controls precision. Number of shown digits is calculated from error + n specified additional digits (1 is sensible default)
1143/// <tr><td> `VerbatimName(bool flag)` <td> Put variable name in a \verb+ + clause.
1144/// </table>
1145/// <tr><td> `Label(const chat* label)` <td> Add header label to parameter box
1146/// <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
1147/// bottom of box is calculated automatically from number lines in box
1148/// <tr><td> `Cut(const char* expression)` <td> Apply given cut expression to data when calculating statistics
1149/// <tr><td> `CutRange(const char* rangeName)` <td> Only consider events within given range when calculating statistics. Multiple
1150/// CutRange() argument may be specified to combine ranges.
1151///
1152/// </table>
1153
1155 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
1156 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
1157{
1158 // Stuff all arguments in a list
1160 cmdList.Add(const_cast<RooCmdArg*>(&arg1)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg2)) ;
1161 cmdList.Add(const_cast<RooCmdArg*>(&arg3)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg4)) ;
1162 cmdList.Add(const_cast<RooCmdArg*>(&arg5)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg6)) ;
1163 cmdList.Add(const_cast<RooCmdArg*>(&arg7)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg8)) ;
1164
1165 // Select the pdf-specific commands
1166 RooCmdConfig pc("RooTreeData::statOn(" + std::string(GetName()) + ")");
1167 pc.defineString("what","What",0,"MNR") ;
1168 pc.defineString("label","Label",0,"") ;
1169 pc.defineDouble("xmin","Layout",0,0.65) ;
1170 pc.defineDouble("xmax","Layout",1,0.99) ;
1171 pc.defineInt("ymaxi","Layout",0,int(0.95*10000)) ;
1172 pc.defineString("formatStr","Format",0,"NELU") ;
1173 pc.defineInt("sigDigit","Format",0,2) ;
1174 pc.defineInt("dummy","FormatArgs",0,0) ;
1175 pc.defineString("cutRange","CutRange",0,"",true) ;
1176 pc.defineString("cutString","CutSpec",0,"") ;
1177 pc.defineMutex("Format","FormatArgs") ;
1178
1179 // Process and check varargs
1180 pc.process(cmdList) ;
1181 if (!pc.ok(true)) {
1182 return frame ;
1183 }
1184
1185 const char* label = pc.getString("label") ;
1186 double xmin = pc.getDouble("xmin") ;
1187 double xmax = pc.getDouble("xmax") ;
1188 double ymax = pc.getInt("ymaxi") / 10000. ;
1189 const char* formatStr = pc.getString("formatStr") ;
1190 int sigDigit = pc.getInt("sigDigit") ;
1191 const char* what = pc.getString("what") ;
1192
1193 const char* cutSpec = pc.getString("cutString",nullptr,true) ;
1194 const char* cutRange = pc.getString("cutRange",nullptr,true) ;
1195
1196 if (pc.hasProcessed("FormatArgs")) {
1197 RooCmdArg* formatCmd = static_cast<RooCmdArg*>(cmdList.FindObject("FormatArgs")) ;
1198 return statOn(frame,what,label,0,nullptr,xmin,xmax,ymax,cutSpec,cutRange,formatCmd) ;
1199 } else {
1200 return statOn(frame,what,label,sigDigit,formatStr,xmin,xmax,ymax,cutSpec,cutRange) ;
1201 }
1202}
1203
1204////////////////////////////////////////////////////////////////////////////////
1205/// Implementation back-end of statOn() method with named arguments
1206
1207RooPlot* RooAbsData::statOn(RooPlot* frame, const char* what, const char *label, Int_t sigDigits,
1208 Option_t *options, double xmin, double xmax, double ymax,
1209 const char* cutSpec, const char* cutRange, const RooCmdArg* formatCmd)
1210{
1211 bool showLabel= (label != nullptr && strlen(label) > 0);
1212
1213 std::string whatStr{what};
1214 std::transform(whatStr.begin(), whatStr.end(), whatStr.begin(), [](unsigned char c){ return std::toupper(c); });
1215 bool showN = whatStr.find('N') != std::string::npos;
1216 bool showR = whatStr.find('R') != std::string::npos;
1217 bool showM = whatStr.find('M') != std::string::npos;
1218 int nPar= 0;
1219 if (showN) nPar++ ;
1220 if (showR) nPar++ ;
1221 if (showM) nPar++ ;
1222
1223 // calculate the box's size
1224 double dy(0.06);
1225 double ymin(ymax - nPar * dy);
1226 if(showLabel) ymin-= dy;
1227
1228 // create the box and set its options
1229 TPaveText *box= new TPaveText(xmin,ymax,xmax,ymin,"BRNDC");
1230 if(!box) return nullptr;
1231 box->SetName(Form("%s_statBox",GetName())) ;
1232 box->SetFillColor(0);
1233 box->SetBorderSize(1);
1234 box->SetTextAlign(12);
1235 box->SetTextSize(0.04F);
1236 box->SetFillStyle(1001);
1237
1238 // add formatted text for each statistic
1239 RooRealVar N("N","Number of Events",sumEntries(cutSpec,cutRange));
1240 N.setPlotLabel("Entries") ;
1241 std::unique_ptr<RooRealVar> meanv{meanVar(*static_cast<RooRealVar*>(frame->getPlotVar()),cutSpec,cutRange)};
1242 meanv->setPlotLabel("Mean") ;
1243 std::unique_ptr<RooRealVar> rms{rmsVar(*static_cast<RooRealVar*>(frame->getPlotVar()),cutSpec,cutRange)};
1244 rms->setPlotLabel("RMS") ;
1245 std::string rmsText = options ? rms->format(sigDigits,options) : rms->format(*formatCmd);
1246 std::string meanText = options ? meanv->format(sigDigits,options) : meanv->format(*formatCmd);
1247 std::string NText = options ? N.format(sigDigits,options) : N.format(*formatCmd);
1248 if (showR) box->AddText(rmsText.c_str());
1249 if (showM) box->AddText(meanText.c_str());
1250 if (showN) box->AddText(NText.c_str());
1251
1252 // add the optional label if specified
1253 if(showLabel) box->AddText(label);
1254
1255 frame->addObject(box) ;
1256 return frame ;
1257}
1258
1259////////////////////////////////////////////////////////////////////////////////
1260/// Loop over columns of our tree data and fill the input histogram. Returns a pointer to the
1261/// input histogram, or zero in case of an error. The input histogram can be any TH1 subclass, and
1262/// therefore of arbitrary dimension. Variables are matched with the (x,y,...) dimensions of the input
1263/// histogram according to the order in which they appear in the input plotVars list.
1264
1265TH1 *RooAbsData::fillHistogram(TH1 *hist, const RooArgList &plotVars, const char *cuts, const char* cutRange) const
1266{
1267 // Do we have a valid histogram to use?
1268 if(nullptr == hist) {
1269 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: no valid histogram to fill" << std::endl;
1270 return nullptr;
1271 }
1272
1273 // Check that the number of plotVars matches the input histogram's dimension
1274 std::size_t hdim= hist->GetDimension();
1275 if(hdim != plotVars.size()) {
1276 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: plotVars has the wrong dimension" << std::endl;
1277 return nullptr;
1278 }
1279
1280 // Check that the plot variables are all actually RooAbsReal's and print a warning if we do not
1281 // explicitly depend on one of them. Clone any variables that we do not contain directly and
1282 // redirect them to use our event data.
1285 for(std::size_t index= 0; index < plotVars.size(); index++) {
1286 const RooAbsArg *var= plotVars.at(index);
1287 const RooAbsReal *realVar= dynamic_cast<const RooAbsReal*>(var);
1288 if(realVar == nullptr) {
1289 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: cannot plot variable \"" << var->GetName()
1290 << "\" of type " << var->ClassName() << std::endl;
1291 return nullptr;
1292 }
1293 RooAbsArg *found= _vars.find(realVar->GetName());
1294 if(!found) {
1295 RooAbsArg *clone= plotClones.addClone(*realVar,true); // do not complain about duplicates
1296 assert(nullptr != clone);
1297 if(!clone->dependsOn(_vars)) {
1298 coutE(InputArguments) << ClassName() << "::" << GetName()
1299 << ":fillHistogram: Data does not contain the variable '" << realVar->GetName() << "'." << std::endl;
1300 return nullptr;
1301 }
1302 else {
1304 }
1305 localVars.add(*clone);
1306 }
1307 else {
1308 localVars.add(*found);
1309 }
1310 }
1311
1312 // Create selection formula if selection cuts are specified
1313 std::unique_ptr<RooFormula> select;
1314 if (cuts != nullptr && strlen(cuts) > 0) {
1315 select = std::make_unique<RooFormula>(cuts, cuts, _vars, false);
1316 if (!select || !select->ok()) {
1317 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: invalid cuts \"" << cuts << "\"" << std::endl;
1318 return nullptr;
1319 }
1320 }
1321
1322 // Lookup each of the variables we are binning in our tree variables
1323 const RooAbsReal *xvar = nullptr;
1324 const RooAbsReal *yvar = nullptr;
1325 const RooAbsReal *zvar = nullptr;
1326 switch(hdim) {
1327 case 3:
1328 zvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(2)->GetName()));
1329 assert(nullptr != zvar);
1330 // fall through to next case...
1331 case 2:
1332 yvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(1)->GetName()));
1333 assert(nullptr != yvar);
1334 // fall through to next case...
1335 case 1:
1336 xvar= dynamic_cast<RooAbsReal*>(localVars.find(plotVars.at(0)->GetName()));
1337 assert(nullptr != xvar);
1338 break;
1339 default:
1340 coutE(InputArguments) << ClassName() << "::" << GetName() << ":fillHistogram: cannot fill histogram with "
1341 << hdim << " dimensions" << std::endl;
1342 break;
1343 }
1344
1345 // Parse cutRange specification
1346 const auto cutVec = ROOT::Split(cutRange ? cutRange : "", ",");
1347
1348 // Loop over events and fill the histogram
1349 if (hist->GetSumw2()->fN==0) {
1350 hist->Sumw2() ;
1351 }
1352 int nevent= numEntries() ; //(int)_tree->GetEntries();
1353 for(int i=0; i < nevent; ++i) {
1354
1355 //int entryNumber= _tree->GetEntryNumber(i);
1356 //if (entryNumber<0) break;
1357 get(i);
1358
1359 // Apply expression based selection criteria
1360 if (select && select->eval()==0) {
1361 continue ;
1362 }
1363
1364
1365 // Apply range based selection criteria
1366 bool selectByRange = true ;
1367 if (cutRange) {
1368 for (const auto arg : _vars) {
1369 bool selectThisArg = false ;
1370 for (auto const& cut : cutVec) {
1371 if (!cut.empty() && arg->inRange(cut.c_str())) {
1373 break ;
1374 }
1375 }
1376 if (!selectThisArg) {
1378 break ;
1379 }
1380 }
1381 }
1382
1383 if (!selectByRange) {
1384 // Go to next event in loop over events
1385 continue ;
1386 }
1387
1388 int bin(0);
1389 switch(hdim) {
1390 case 1:
1391 bin= hist->FindBin(xvar->getVal());
1392 hist->Fill(xvar->getVal(),weight()) ;
1393 break;
1394 case 2:
1395 bin= hist->FindBin(xvar->getVal(),yvar->getVal());
1396 static_cast<TH2*>(hist)->Fill(xvar->getVal(),yvar->getVal(),weight()) ;
1397 break;
1398 case 3:
1399 bin= hist->FindBin(xvar->getVal(),yvar->getVal(),zvar->getVal());
1400 static_cast<TH3*>(hist)->Fill(xvar->getVal(),yvar->getVal(),zvar->getVal(),weight()) ;
1401 break;
1402 default:
1403 assert(hdim < 3);
1404 break;
1405 }
1406
1407
1408 double error2 = std::pow(hist->GetBinError(bin),2)-std::pow(weight(),2) ;
1409 double we = weightError(RooAbsData::SumW2) ;
1410 if (we==0) we = weight() ;
1411 error2 += std::pow(we,2) ;
1412
1413
1414// double we = weightError(RooAbsData::SumW2) ;
1415// double error2(0) ;
1416// if (we==0) {
1417// we = weight() ; //sqrt(weight()) ;
1418// error2 = std::pow(hist->GetBinError(bin),2)-std::pow(weight(),2) + std::pow(we,2) ;
1419// } else {
1420// error2 = std::pow(hist->GetBinError(bin),2)-std::pow(weight(),2) + std::pow(we,2) ;
1421// }
1422 //hist->AddBinContent(bin,weight());
1423 hist->SetBinError(bin,sqrt(error2)) ;
1424
1425 //cout << "RooTreeData::fillHistogram() bin = " << bin << " weight() = " << weight() << " we = " << we << std::endl ;
1426
1427 }
1428
1429 return hist;
1430}
1431
1432
1433namespace {
1434
1435struct SplittingSetup {
1436 RooArgSet ownedSet;
1437 RooAbsCategory *cloneCat = nullptr;
1438 RooArgSet subsetVars;
1439 bool addWeightVar = false;
1440};
1441
1442SplittingSetup initSplit(RooAbsData const &data, RooAbsCategory const &splitCat)
1443{
1444 SplittingSetup setup;
1445
1446 // Sanity check
1447 if (!splitCat.dependsOn(*data.get())) {
1448 oocoutE(&data, InputArguments) << "RooTreeData::split(" << data.GetName() << ") ERROR category "
1449 << splitCat.GetName() << " doesn't depend on any variable in this dataset"
1450 << std::endl;
1451 return setup;
1452 }
1453
1454 // Clone splitting category and attach to self
1455 if (splitCat.isDerived()) {
1456 RooArgSet(splitCat).snapshot(setup.ownedSet, true);
1457 setup.cloneCat = static_cast<RooAbsCategory *>(setup.ownedSet.find(splitCat.GetName()));
1458 setup.cloneCat->attachDataSet(data);
1459 } else {
1460 setup.cloneCat = dynamic_cast<RooAbsCategory *>(data.get()->find(splitCat.GetName()));
1461 if (!setup.cloneCat) {
1462 oocoutE(&data, InputArguments) << "RooTreeData::split(" << data.GetName() << ") ERROR category "
1463 << splitCat.GetName() << " is fundamental and does not appear in this dataset"
1464 << std::endl;
1465 return setup;
1466 }
1467 }
1468
1469 // Construct set of variables to be included in split sets = full set - split category
1470 setup.subsetVars.add(*data.get());
1471 if (splitCat.isDerived()) {
1472 std::unique_ptr<RooArgSet> vars{splitCat.getVariables()};
1473 setup.subsetVars.remove(*vars, true, true);
1474 } else {
1475 setup.subsetVars.remove(splitCat, true, true);
1476 }
1477
1478 // Add weight variable explicitly if dataset has weights, but no top-level weight
1479 // variable exists (can happen with composite datastores)
1480 setup.addWeightVar = data.isWeighted();
1481
1482 return setup;
1483}
1484
1486 std::function<std::unique_ptr<RooAbsData>(const char *label)> createEmptyData)
1487{
1488 auto dsetList = std::make_unique<TList>();
1489
1490 // If createEmptyDataSets is true, prepopulate with empty sets corresponding to all states
1491 if (createEmptyDataSets) {
1492 for (const auto &nameIdx : cloneCat) {
1493 dsetList->Add(createEmptyData(nameIdx.first.c_str()).release());
1494 }
1495 }
1496
1497 bool isDataHist = dynamic_cast<RooDataHist const *>(&data);
1498
1499 // Loop over dataset and copy event to matching subset
1500 for (int i = 0; i < data.numEntries(); ++i) {
1501 const RooArgSet *row = data.get(i);
1502 auto subset = static_cast<RooAbsData *>(dsetList->FindObject(cloneCat.getCurrentLabel()));
1503 if (!subset) {
1504 subset = createEmptyData(cloneCat.getCurrentLabel()).release();
1505 dsetList->Add(subset);
1506 }
1507
1508 // For datasets with weight errors or sumW2, the interface to fill
1509 // RooDataHist and RooDataSet is not the same.
1510 if (isDataHist) {
1511 static_cast<RooDataHist *>(subset)->add(*row, data.weight(), data.weightSquared());
1512 } else {
1513 static_cast<RooDataSet *>(subset)->add(*row, data.weight(), data.weightError());
1514 }
1515 }
1516
1517 return RooFit::makeOwningPtr(std::move(dsetList));
1518}
1519
1520} // namespace
1521
1522
1523/**
1524 * \brief Split the dataset into subsets based on states of a categorical variable in this dataset.
1525 *
1526 * Returns a list of sub-datasets, which each dataset named after a given state
1527 * name in the `splitCat`. The observables `splitCat` itself is no longer present
1528 * in the sub-datasets.
1529 *
1530 * \note If you mean to split a dataset into sub-datasets that correspond to
1531 * the individual channels of a RooSimultaneous, it is better to use
1532 * RooAbsData::split(const RooSimultaneous &, bool), because then the
1533 * sub-datasets only contain variables that the pdf for the corresponding
1534 * channel depends on. This is much faster in case of many channels, and the
1535 * resulting sub-datasets don't waste memory for unused columns.
1536 *
1537 * \param splitCat The categorical variable used for splitting the dataset.
1538 * \param createEmptyDataSets Flag indicating whether to create empty datasets
1539 * for missing categories (`false` by default).
1540 *
1541 * \return An owning pointer to a TList of subsets of the dataset.
1542 * Returns `nullptr` if an error occurs.
1543 */
1544
1546{
1547 SplittingSetup setup = initSplit(*this, splitCat);
1548
1549 // Something went wrong
1550 if (!setup.cloneCat)
1551 return nullptr;
1552
1553 auto createEmptyData = [&](const char *label) -> std::unique_ptr<RooAbsData> {
1554 return std::unique_ptr<RooAbsData>{
1555 emptyClone(label, label, &setup.subsetVars, setup.addWeightVar ? "weight" : nullptr)};
1556 };
1557
1558 return splitImpl(*this, *setup.cloneCat, createEmptyDataSets, createEmptyData);
1559}
1560
1561/**
1562 * \brief Split the dataset into subsets based on the channels of a RooSimultaneous.
1563 *
1564 * Returns a list of sub-datasets, which each dataset named after the
1565 * applicable state name of the RooSimultaneous index category. The index
1566 * category itself is no longer present in the sub-datasets. The sub-datasets
1567 * only contain variables that the pdf for the corresponding channel depends
1568 * on.
1569 *
1570 * \param simPdf The simultaneous pdf used for splitting the dataset.
1571 * \param createEmptyDataSets Flag indicating whether to create empty datasets
1572 * for missing categories (`false` by default).
1573 *
1574 * \return An owning pointer to a TList of subsets of the dataset.
1575 * Returns `nullptr` if an error occurs.
1576 */
1578{
1579 auto &splitCat = const_cast<RooAbsCategoryLValue &>(simPdf.indexCat());
1580
1581 SplittingSetup setup = initSplit(*this, splitCat);
1582
1583 // Something went wrong
1584 if (!setup.cloneCat)
1585 return nullptr;
1586
1587 // Get the observables for a given pdf in the RooSimultaneous, or an empty
1588 // RooArgSet if no pdf is set
1589 auto getPdfObservables = [this, &simPdf](const char *label) {
1591 if (RooAbsPdf *catPdf = simPdf.getPdf(label)) {
1592 catPdf->getObservables(this->get(), obsSet);
1593 }
1594 return obsSet;
1595 };
1596
1597 // By default, remove all category observables from the subdatasets
1599 for (const auto &catPair : splitCat) {
1600 allObservables.add(getPdfObservables(catPair.first.c_str()));
1601 }
1602 setup.subsetVars.remove(allObservables, true, true);
1603
1604 auto createEmptyData = [&](const char *label) -> std::unique_ptr<RooAbsData> {
1605 // Add in the subset only the observables corresponding to this category
1606 RooArgSet subsetVarsCat(setup.subsetVars);
1607 subsetVarsCat.add(getPdfObservables(label));
1608 return std::unique_ptr<RooAbsData>{
1609 this->emptyClone(label, label, &subsetVarsCat, setup.addWeightVar ? "weight" : nullptr)};
1610 };
1611
1612 return splitImpl(*this, *setup.cloneCat, createEmptyDataSets, createEmptyData);
1613}
1614
1615////////////////////////////////////////////////////////////////////////////////
1616/// Plot dataset on specified frame.
1617///
1618/// By default:
1619/// - An unbinned dataset will use the default binning of the target frame.
1620/// - A binned dataset will retain its intrinsic binning.
1621///
1622/// The following optional named arguments can be used to modify the behaviour:
1623/// \note Please follow the function links in the left column to learn about PyROOT specifics for a given option.
1624///
1625/// <table>
1626///
1627/// <tr><th> <th> Data representation options
1628/// <tr><td> RooFit::Asymmetry(const RooCategory& c)
1629/// <td> Show the asymmetry of the data in given two-state category [F(+)-F(-)] / [F(+)+F(-)].
1630/// Category must have two states with indices -1 and +1 or three states with indices -1,0 and +1.
1631/// <tr><td> RooFit::Efficiency(const RooCategory& c)
1632/// <td> Show the efficiency F(acc)/[F(acc)+F(rej)]. Category must have two states with indices 0 and 1
1633/// <tr><td> RooFit::DataError(Int_t)
1634/// <td> Select the type of error drawn:
1635/// - `Auto(default)` results in Poisson for unweighted data and SumW2 for weighted data
1636/// - `Poisson` draws asymmetric Poisson confidence intervals.
1637/// - `SumW2` draws symmetric sum-of-weights error ( \f$ \left( \sum w \right)^2 / \sum\left(w^2\right) \f$ )
1638/// - `None` draws no error bars
1639/// <tr><td> RooFit::Binning(int nbins, double xlo, double xhi)
1640/// <td> Use specified binning to draw dataset
1641/// <tr><td> RooFit::Binning(const RooAbsBinning&)
1642/// <td> Use specified binning to draw dataset
1643/// <tr><td> RooFit::Binning(const char* name)
1644/// <td> Use binning with specified name to draw dataset
1645/// <tr><td> RooFit::RefreshNorm()
1646/// <td> Force refreshing for PDF normalization information in frame.
1647/// If set, any subsequent PDF will normalize to this dataset, even if it is
1648/// not the first one added to the frame. By default only the 1st dataset
1649/// added to a frame will update the normalization information
1650/// <tr><td> RooFit::Rescale(double f)
1651/// <td> Rescale drawn histogram by given factor.
1652/// <tr><td> RooFit::Cut(const char*)
1653/// <td> Only plot entries that pass the given cut.
1654/// Apart from cutting in continuous variables `Cut("x>5")`, this can also be used to plot a specific
1655/// category state. Use something like `Cut("myCategory == myCategory::stateA")`, where
1656/// `myCategory` resolves to the state number for a given entry and
1657/// `myCategory::stateA` resolves to the state number of the state named "stateA".
1658///
1659/// <tr><td> RooFit::CutRange(const char*)
1660/// <td> Only plot data from given range. Separate multiple ranges with ",".
1661/// \note This often requires passing the normalisation when plotting the PDF because RooFit does not save
1662/// how many events were being plotted (it will only work for cutting slices out of uniformly distributed
1663/// variables).
1664/// ```
1665/// data->plotOn(frame01, CutRange("SB1"));
1666/// const double nData = data->sumEntries("", "SB1");
1667/// // Make clear that the target normalisation is nData. The enumerator NumEvent
1668/// // is needed to switch between relative and absolute scaling.
1669/// model.plotOn(frame01, Normalization(nData, RooAbsReal::NumEvent),
1670/// ProjectionRange("SB1"));
1671/// ```
1672///
1673/// <tr><th> <th> Histogram drawing options
1674/// <tr><td> RooFit::DrawOption(const char* opt)
1675/// <td> Select ROOT draw option for resulting TGraph object
1676/// <tr><td> RooFit::LineStyle(Style_t style)
1677/// <td> Select line style by ROOT line style code, default is solid
1678/// <tr><td> RooFit::LineColor(Color_t color)
1679/// <td> Select line color by ROOT color code, default is black
1680/// <tr><td> RooFit::LineWidth(Width_t width)
1681/// <td> Select line with in pixels, default is 3
1682/// <tr><td> RooFit::MarkerStyle(Style_t style)
1683/// <td> Select the ROOT marker style, default is 21
1684/// <tr><td> RooFit::MarkerColor(Color_t color)
1685/// <td> Select the ROOT marker color, default is black
1686/// <tr><td> RooFit::MarkerSize(Size_t size)
1687/// <td> Select the ROOT marker size
1688/// <tr><td> RooFit::FillStyle(Style_t style)
1689/// <td> Select fill style, default is filled.
1690/// <tr><td> RooFit::FillColor(Color_t color)
1691/// <td> Select fill color by ROOT color code
1692/// <tr><td> RooFit::XErrorSize(double frac)
1693/// <td> Select size of X error bar as fraction of the bin width, default is 1
1694///
1695/// <tr><th> <th> Misc. other options
1696/// <tr><td> RooFit::Name(const char* name)
1697/// <td> Give curve specified name in frame. Useful if curve is to be referenced later
1698/// <tr><td> RooFit::Invisible()
1699/// <td> Add curve to frame, but do not display. Useful in combination AddTo()
1700/// <tr><td> RooFit::AddTo(const char* name, double wgtSel, double wgtOther)
1701/// <td> Add constructed histogram to already existing histogram with given name and relative weight factors
1702///
1703/// </table>
1704
1705RooPlot* RooAbsData::plotOn(RooPlot* frame, const RooLinkedList& argList) const
1706{
1707 // New experimental plotOn() with varargs...
1708
1709 // Define configuration for this method
1710 RooCmdConfig pc("RooAbsData::plotOn(" + std::string(GetName()) + ")");
1711 pc.defineString("drawOption","DrawOption",0,"P") ;
1712 pc.defineString("cutRange","CutRange",0,"",true) ;
1713 pc.defineString("cutString","CutSpec",0,"") ;
1714 pc.defineString("histName","Name",0,"") ;
1715 pc.defineObject("cutVar","CutVar",0) ;
1716 pc.defineObject("binning","Binning",0) ;
1717 pc.defineString("binningName","BinningName",0,"") ;
1718 pc.defineInt("nbins","BinningSpec",0,100) ;
1719 pc.defineDouble("xlo","BinningSpec",0,0) ;
1720 pc.defineDouble("xhi","BinningSpec",1,1) ;
1721 pc.defineObject("asymCat","Asymmetry",0) ;
1722 pc.defineObject("effCat","Efficiency",0) ;
1723 pc.defineInt("lineColor","LineColor",0,-999) ;
1724 pc.defineInt("lineStyle","LineStyle",0,-999) ;
1725 pc.defineInt("lineWidth","LineWidth",0,-999) ;
1726 pc.defineInt("markerColor","MarkerColor",0,-999) ;
1727 pc.defineInt("markerStyle","MarkerStyle",0,-999) ;
1728 pc.defineDouble("markerSize","MarkerSize",0,-999) ;
1729 pc.defineInt("fillColor","FillColor",0,-999) ;
1730 pc.defineInt("fillStyle","FillStyle",0,-999) ;
1731 pc.defineInt("errorType","DataError",0,(int)RooAbsData::Auto) ;
1732 pc.defineInt("histInvisible","Invisible",0,0) ;
1733 pc.defineInt("refreshFrameNorm","RefreshNorm",0,1) ;
1734 pc.defineString("addToHistName","AddTo",0,"") ;
1735 pc.defineDouble("addToWgtSelf","AddTo",0,1.) ;
1736 pc.defineDouble("addToWgtOther","AddTo",1,1.) ;
1737 pc.defineDouble("xErrorSize","XErrorSize",0,1.) ;
1738 pc.defineDouble("scaleFactor","Rescale",0,1.) ;
1739 pc.defineMutex("DataError","Asymmetry","Efficiency") ;
1740 pc.defineMutex("Binning","BinningName","BinningSpec") ;
1741
1742 // Process & check varargs
1743 pc.process(argList) ;
1744 if (!pc.ok(true)) {
1745 return frame ;
1746 }
1747
1748 PlotOpt o ;
1749
1750 // Extract values from named arguments
1751 o.drawOptions = pc.getString("drawOption") ;
1752 o.cuts = pc.getString("cutString") ;
1753 if (pc.hasProcessed("Binning")) {
1754 o.bins = static_cast<RooAbsBinning*>(pc.getObject("binning")) ;
1755 } else if (pc.hasProcessed("BinningName")) {
1756 o.bins = &frame->getPlotVar()->getBinning(pc.getString("binningName")) ;
1757 } else if (pc.hasProcessed("BinningSpec")) {
1758 double xlo = pc.getDouble("xlo") ;
1759 double xhi = pc.getDouble("xhi") ;
1760 o.bins = new RooUniformBinning((xlo==xhi)?frame->getPlotVar()->getMin():xlo,
1761 (xlo==xhi)?frame->getPlotVar()->getMax():xhi,pc.getInt("nbins")) ;
1762 }
1763 const RooAbsCategoryLValue* asymCat = static_cast<const RooAbsCategoryLValue*>(pc.getObject("asymCat")) ;
1764 const RooAbsCategoryLValue* effCat = static_cast<const RooAbsCategoryLValue*>(pc.getObject("effCat")) ;
1765 o.etype = (RooAbsData::ErrorType) pc.getInt("errorType") ;
1766 o.histInvisible = pc.getInt("histInvisible") ;
1767 o.xErrorSize = pc.getDouble("xErrorSize") ;
1768 o.cutRange = pc.getString("cutRange",nullptr,true) ;
1769 o.histName = pc.getString("histName",nullptr,true) ;
1770 o.addToHistName = pc.getString("addToHistName",nullptr,true) ;
1771 o.addToWgtSelf = pc.getDouble("addToWgtSelf") ;
1772 o.addToWgtOther = pc.getDouble("addToWgtOther") ;
1773 o.refreshFrameNorm = pc.getInt("refreshFrameNorm") ;
1774 o.scaleFactor = pc.getDouble("scaleFactor") ;
1775
1776 // Map auto error type to actual type
1777 if (o.etype == Auto) {
1779 if (o.etype == SumW2) {
1780 coutI(InputArguments) << "RooAbsData::plotOn(" << GetName()
1781 << ") INFO: dataset has non-integer weights, auto-selecting SumW2 errors instead of Poisson errors" << std::endl ;
1782 }
1783 }
1784
1785 if (o.addToHistName && !frame->findObject(o.addToHistName,RooHist::Class())) {
1786 coutE(InputArguments) << "RooAbsData::plotOn(" << GetName() << ") cannot find existing histogram " << o.addToHistName
1787 << " to add to in RooPlot" << std::endl ;
1788 return frame ;
1789 }
1790
1791 RooPlot* ret ;
1792 if (!asymCat && !effCat) {
1793 ret = plotOn(frame,o) ;
1794 } else if (asymCat) {
1795 ret = plotAsymOn(frame,*asymCat,o) ;
1796 } else {
1797 ret = plotEffOn(frame,*effCat,o) ;
1798 }
1799
1800 int lineColor = pc.getInt("lineColor") ;
1801 int lineStyle = pc.getInt("lineStyle") ;
1802 int lineWidth = pc.getInt("lineWidth") ;
1803 int markerColor = pc.getInt("markerColor") ;
1804 int markerStyle = pc.getInt("markerStyle") ;
1805 Size_t markerSize = pc.getDouble("markerSize") ;
1806 int fillColor = pc.getInt("fillColor") ;
1807 int fillStyle = pc.getInt("fillStyle") ;
1808 if (lineColor!=-999) ret->getAttLine()->SetLineColor(lineColor) ;
1809 if (lineStyle!=-999) ret->getAttLine()->SetLineStyle(lineStyle) ;
1810 if (lineWidth!=-999) ret->getAttLine()->SetLineWidth(lineWidth) ;
1811 if (markerColor!=-999) ret->getAttMarker()->SetMarkerColor(markerColor) ;
1812 if (markerStyle!=-999) ret->getAttMarker()->SetMarkerStyle(markerStyle) ;
1813 if (markerSize!=-999) ret->getAttMarker()->SetMarkerSize(markerSize) ;
1814 if (fillColor!=-999) ret->getAttFill()->SetFillColor(fillColor) ;
1815 if (fillStyle!=-999) ret->getAttFill()->SetFillStyle(fillStyle) ;
1816
1817 if (pc.hasProcessed("BinningSpec")) {
1818 delete o.bins ;
1819 }
1820
1821 return ret ;
1822}
1823
1824////////////////////////////////////////////////////////////////////////////////
1825/// Create and fill a histogram of the frame's variable and append it to the frame.
1826/// The frame variable must be one of the data sets dimensions.
1827///
1828/// The plot range and the number of plot bins is determined by the parameters
1829/// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins()).
1830///
1831/// The optional cut string expression can be used to select the events to be plotted.
1832/// The cut specification may refer to any variable contained in the data set.
1833///
1834/// The drawOptions are passed to the TH1::Draw() method.
1835/// \see RooAbsData::plotOn(RooPlot*,const RooLinkedList&) const
1837{
1838 if(nullptr == frame) {
1839 coutE(Plotting) << ClassName() << "::" << GetName() << ":plotOn: frame is null" << std::endl;
1840 return nullptr;
1841 }
1843 if(nullptr == var) {
1844 coutE(Plotting) << ClassName() << "::" << GetName()
1845 << ":plotOn: frame does not specify a plot variable" << std::endl;
1846 return nullptr;
1847 }
1848
1849 // create and fill a temporary histogram of this variable
1850 const std::string histName = std::string{GetName()} + "_plot";
1851 std::unique_ptr<TH1> hist;
1852 if (o.bins) {
1853 hist.reset( var->createHistogram(histName.c_str(), RooFit::AxisLabel("Events"), RooFit::Binning(*o.bins)) );
1854 } else if (!frame->getPlotVar()->getBinning().isUniform()) {
1855 hist.reset( var->createHistogram(histName.c_str(), RooFit::AxisLabel("Events"),
1856 RooFit::Binning(frame->getPlotVar()->getBinning())) );
1857 } else {
1858 hist.reset( var->createHistogram(histName.c_str(), "Events",
1859 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(), frame->GetNbinsX()) );
1860 }
1861
1862 // Keep track of sum-of-weights error
1863 hist->Sumw2() ;
1864
1865 if(nullptr == fillHistogram(hist.get(), RooArgList(*var),o.cuts,o.cutRange)) {
1866 coutE(Plotting) << ClassName() << "::" << GetName()
1867 << ":plotOn: fillHistogram() failed" << std::endl;
1868 return nullptr;
1869 }
1870
1871 // If frame has no predefined bin width (event density) it will be adjusted to
1872 // our histograms bin width so we should force that bin width here
1873 double nomBinWidth ;
1874 if (frame->getFitRangeNEvt()==0 && o.bins) {
1876 } else {
1877 nomBinWidth = o.bins ? frame->getFitRangeBinW() : 0 ;
1878 }
1879
1880 // convert this histogram to a RooHist object on the heap
1882 if(nullptr == graph) {
1883 coutE(Plotting) << ClassName() << "::" << GetName()
1884 << ":plotOn: unable to create a RooHist object" << std::endl;
1885 return nullptr;
1886 }
1887
1888 // If the dataset variable has a wide range than the plot variable,
1889 // calculate the number of entries in the dataset in the plot variable fit range
1890 RooAbsRealLValue* dataVar = static_cast<RooAbsRealLValue*>(_vars.find(var->GetName())) ;
1891 double nEnt(sumEntries()) ;
1892 if (dataVar->getMin()<var->getMin() || dataVar->getMax()>var->getMax()) {
1893 std::unique_ptr<RooAbsData> tmp{const_cast<RooAbsData*>(this)->reduce(RooFit::SelectVars(*var))};
1894 nEnt = tmp->sumEntries() ;
1895 }
1896
1897 // Store the number of entries before the cut, if any was made
1898 if ((o.cuts && strlen(o.cuts)) || o.cutRange) {
1899 coutI(Plotting) << "RooTreeData::plotOn: plotting " << hist->GetSumOfWeights() << " events out of " << nEnt << " total events" << std::endl ;
1900 graph->setRawEntries(nEnt) ;
1901 }
1902
1903 // Add self to other hist if requested
1904 if (o.addToHistName) {
1905 RooHist* otherGraph = static_cast<RooHist*>(frame->findObject(o.addToHistName,RooHist::Class())) ;
1906
1907 if (!graph->hasIdenticalBinning(*otherGraph)) {
1908 coutE(Plotting) << "RooTreeData::plotOn: ERROR Histogram to be added to, '" << o.addToHistName << "',has different binning" << std::endl ;
1909 delete graph ;
1910 return frame ;
1911 }
1912
1914 delete graph ;
1915 graph = sumGraph ;
1916 }
1917
1918 // Rename graph if requested
1919 if (o.histName) {
1920 graph->SetName(o.histName) ;
1921 } else {
1922 std::string hname = std::string{"h_"} + GetName();
1923 if (o.cutRange && strlen(o.cutRange)>0) {
1924 hname += std::string{"_CutRange["} + o.cutRange + "]";
1925 }
1926 if (o.cuts && strlen(o.cuts)>0) {
1927 hname += std::string{"_Cut["} + o.cuts + "]";
1928 }
1929 graph->SetName(hname.c_str()) ;
1930 }
1931
1932 // initialize the frame's normalization setup, if necessary
1933 frame->updateNormVars(_vars);
1934
1935
1936 // add the RooHist to the specified plot
1938
1939 return frame;
1940}
1941
1942////////////////////////////////////////////////////////////////////////////////
1943/// Create and fill a histogram with the asymmetry N[+] - N[-] / ( N[+] + N[-] ),
1944/// where N(+/-) is the number of data points with asymCat=+1 and asymCat=-1
1945/// as function of the frames variable. The asymmetry category 'asymCat' must
1946/// have exactly 2 (or 3) states defined with index values +1,-1 (and 0)
1947///
1948/// The plot range and the number of plot bins is determined by the parameters
1949/// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins())
1950///
1951/// The optional cut string expression can be used to select the events to be plotted.
1952/// The cut specification may refer to any variable contained in the data set
1953///
1954/// The drawOptions are passed to the TH1::Draw() method
1955
1957{
1958 if(nullptr == frame) {
1959 coutE(Plotting) << ClassName() << "::" << GetName() << ":plotAsymOn: frame is null" << std::endl;
1960 return nullptr;
1961 }
1963 if(nullptr == var) {
1964 coutE(Plotting) << ClassName() << "::" << GetName()
1965 << ":plotAsymOn: frame does not specify a plot variable" << std::endl;
1966 return nullptr;
1967 }
1968
1969 // create and fill temporary histograms of this variable for each state
1970 std::string hist1Name(GetName());
1971 std::string hist2Name(GetName());
1972 hist1Name += "_plot1";
1973 std::unique_ptr<TH1> hist1;
1974 std::unique_ptr<TH1> hist2;
1975 hist2Name += "_plot2";
1976
1977 if (o.bins) {
1978 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events", *o.bins) );
1979 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events", *o.bins) );
1980 } else {
1981 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events",
1982 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
1983 frame->GetNbinsX()) );
1984 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events",
1985 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
1986 frame->GetNbinsX()) );
1987 }
1988
1989 assert(hist1 && hist2);
1990
1991 std::string cuts1;
1992 std::string cuts2;
1993 if (o.cuts && strlen(o.cuts)) {
1994 cuts1 = Form("(%s)&&(%s>0)",o.cuts,asymCat.GetName());
1995 cuts2 = Form("(%s)&&(%s<0)",o.cuts,asymCat.GetName());
1996 } else {
1997 cuts1 = Form("(%s>0)",asymCat.GetName());
1998 cuts2 = Form("(%s<0)",asymCat.GetName());
1999 }
2000
2001 if(! fillHistogram(hist1.get(), RooArgList(*var),cuts1.c_str(),o.cutRange) ||
2002 ! fillHistogram(hist2.get(), RooArgList(*var),cuts2.c_str(),o.cutRange)) {
2003 coutE(Plotting) << ClassName() << "::" << GetName()
2004 << ":plotAsymOn: createHistogram() failed" << std::endl;
2005 return nullptr;
2006 }
2007
2008 // convert this histogram to a RooHist object on the heap
2009 RooHist *graph= new RooHist(*hist1,*hist2,0,1,o.etype,o.xErrorSize,false,o.scaleFactor);
2010 graph->setYAxisLabel(Form("Asymmetry in %s",asymCat.GetName())) ;
2011
2012 // initialize the frame's normalization setup, if necessary
2013 frame->updateNormVars(_vars);
2014
2015 // Rename graph if requested
2016 if (o.histName) {
2017 graph->SetName(o.histName) ;
2018 } else {
2019 std::string hname{Form("h_%s_Asym[%s]",GetName(),asymCat.GetName())};
2020 if (o.cutRange && strlen(o.cutRange)>0) {
2021 hname += Form("_CutRange[%s]",o.cutRange);
2022 }
2023 if (o.cuts && strlen(o.cuts)>0) {
2024 hname += Form("_Cut[%s]",o.cuts);
2025 }
2026 graph->SetName(hname.c_str()) ;
2027 }
2028
2029 // add the RooHist to the specified plot
2031
2032 return frame;
2033}
2034
2035////////////////////////////////////////////////////////////////////////////////
2036/// Create and fill a histogram with the efficiency N[1] / ( N[1] + N[0] ),
2037/// where N(1/0) is the number of data points with effCat=1 and effCat=0
2038/// as function of the frames variable. The efficiency category 'effCat' must
2039/// have exactly 2 +1 and 0.
2040///
2041/// The plot range and the number of plot bins is determined by the parameters
2042/// of the plot variable of the frame (RooAbsReal::setPlotRange(), RooAbsReal::setPlotBins())
2043///
2044/// The optional cut string expression can be used to select the events to be plotted.
2045/// The cut specification may refer to any variable contained in the data set
2046///
2047/// The drawOptions are passed to the TH1::Draw() method
2048
2050{
2051 if(nullptr == frame) {
2052 coutE(Plotting) << ClassName() << "::" << GetName() << ":plotEffOn: frame is null" << std::endl;
2053 return nullptr;
2054 }
2056 if(nullptr == var) {
2057 coutE(Plotting) << ClassName() << "::" << GetName()
2058 << ":plotEffOn: frame does not specify a plot variable" << std::endl;
2059 return nullptr;
2060 }
2061
2062 // create and fill temporary histograms of this variable for each state
2063 std::string hist1Name(GetName());
2064 std::string hist2Name(GetName());
2065 hist1Name += "_plot1";
2066 std::unique_ptr<TH1> hist1;
2067 std::unique_ptr<TH1> hist2;
2068 hist2Name += "_plot2";
2069
2070 if (o.bins) {
2071 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events", *o.bins) );
2072 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events", *o.bins) );
2073 } else {
2074 hist1.reset( var->createHistogram(hist1Name.c_str(), "Events",
2075 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
2076 frame->GetNbinsX()) );
2077 hist2.reset( var->createHistogram(hist2Name.c_str(), "Events",
2078 frame->GetXaxis()->GetXmin(), frame->GetXaxis()->GetXmax(),
2079 frame->GetNbinsX()) );
2080 }
2081
2082 assert(hist1 && hist2);
2083
2084 std::string cuts1;
2085 std::string cuts2;
2086 if (o.cuts && strlen(o.cuts)) {
2087 cuts1 = Form("(%s)&&(%s==1)",o.cuts,effCat.GetName());
2088 cuts2 = Form("(%s)&&(%s==0)",o.cuts,effCat.GetName());
2089 } else {
2090 cuts1 = Form("(%s==1)",effCat.GetName());
2091 cuts2 = Form("(%s==0)",effCat.GetName());
2092 }
2093
2094 if(! fillHistogram(hist1.get(), RooArgList(*var),cuts1.c_str(),o.cutRange) ||
2095 ! fillHistogram(hist2.get(), RooArgList(*var),cuts2.c_str(),o.cutRange)) {
2096 coutE(Plotting) << ClassName() << "::" << GetName()
2097 << ":plotEffOn: createHistogram() failed" << std::endl;
2098 return nullptr;
2099 }
2100
2101 // convert this histogram to a RooHist object on the heap
2102 RooHist *graph= new RooHist(*hist1,*hist2,0,1,o.etype,o.xErrorSize,true);
2103 graph->setYAxisLabel(Form("Efficiency of %s=%s", effCat.GetName(), effCat.lookupName(1).c_str()));
2104
2105 // initialize the frame's normalization setup, if necessary
2106 frame->updateNormVars(_vars);
2107
2108 // Rename graph if requested
2109 if (o.histName) {
2110 graph->SetName(o.histName) ;
2111 } else {
2112 std::string hname(Form("h_%s_Eff[%s]",GetName(),effCat.GetName())) ;
2113 if (o.cutRange && strlen(o.cutRange)>0) {
2114 hname += Form("_CutRange[%s]",o.cutRange);
2115 }
2116 if (o.cuts && strlen(o.cuts)>0) {
2117 hname += Form("_Cut[%s]",o.cuts);
2118 }
2119 graph->SetName(hname.c_str()) ;
2120 }
2121
2122 // add the RooHist to the specified plot
2124
2125 return frame;
2126}
2127
2128////////////////////////////////////////////////////////////////////////////////
2129/// Create and fill a 1-dimensional table for given category column
2130/// This functions is the equivalent of plotOn() for category dimensions.
2131///
2132/// The optional cut string expression can be used to select the events to be tabulated
2133/// The cut specification may refer to any variable contained in the data set
2134///
2135/// The option string is currently not used
2136
2137Roo1DTable* RooAbsData::table(const RooAbsCategory& cat, const char* cuts, const char* /*opts*/) const
2138{
2139 // First see if var is in data set
2140 RooAbsCategory* tableVar = static_cast<RooAbsCategory*>(_vars.find(cat.GetName())) ;
2141 std::unique_ptr<RooArgSet> tableSet;
2142 if (!tableVar) {
2143 if (!cat.dependsOn(_vars)) {
2144 coutE(Plotting) << "RooTreeData::Table(" << GetName() << "): Argument " << cat.GetName()
2145 << " is not in dataset and is also not dependent on data set" << std::endl ;
2146 return nullptr;
2147 }
2148
2149 // Clone derived variable
2150 tableSet = std::make_unique<RooArgSet>();
2151 if (RooArgSet(cat).snapshot(*tableSet, true)) {
2152 coutE(Plotting) << "RooTreeData::table(" << GetName() << ") Couldn't deep-clone table category, abort." << std::endl;
2153 return nullptr;
2154 }
2155 tableVar = static_cast<RooAbsCategory*>(tableSet->find(cat.GetName())) ;
2156
2157 //Redirect servers of derived clone to internal ArgSet representing the data in this set
2158 tableVar->recursiveRedirectServers(_vars) ;
2159 }
2160
2161 std::unique_ptr<RooFormulaVar> cutVar;
2162 std::string tableName{GetName()};
2163 if (cuts && strlen(cuts)) {
2164 tableName += "(";
2165 tableName += cuts;
2166 tableName += ")";
2167 // Make cut selector if cut is specified
2168 cutVar = std::make_unique<RooFormulaVar>("cutVar",cuts,_vars) ;
2169 }
2170 Roo1DTable* table2 = tableVar->createTable(tableName.c_str());
2171
2172 // Dump contents
2173 int nevent= numEntries() ;
2174 for(int i=0; i < nevent; ++i) {
2175 get(i);
2176
2177 if (cutVar && cutVar->getVal()==0) continue ;
2178
2179 table2->fill(*tableVar,weight()) ;
2180 }
2181
2182 return table2 ;
2183}
2184
2185////////////////////////////////////////////////////////////////////////////////
2186/// Fill Doubles 'lowest' and 'highest' with the lowest and highest value of
2187/// observable 'var' in this dataset. If the return value is true and error
2188/// occurred
2189
2190bool RooAbsData::getRange(const RooAbsRealLValue& var, double& lowest, double& highest, double marginFrac, bool symMode) const
2191{
2192 // Lookup variable in dataset
2193 const auto arg = _vars.find(var.GetName());
2194 if (!arg) {
2195 coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") ERROR: unknown variable: " << var.GetName() << std::endl ;
2196 return true;
2197 }
2198
2199 auto varPtr = dynamic_cast<const RooRealVar*>(arg);
2200 // Check if found variable is of type RooRealVar
2201 if (!varPtr) {
2202 coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") ERROR: variable " << var.GetName() << " is not of type RooRealVar" << std::endl ;
2203 return true;
2204 }
2205
2206 // Check if dataset is not empty
2207 if(sumEntries() == 0.) {
2208 coutE(InputArguments) << "RooDataSet::getRange(" << GetName() << ") WARNING: empty dataset" << std::endl ;
2209 return true;
2210 }
2211
2212 // Look for highest and lowest value
2215 for (int i=0 ; i<numEntries() ; i++) {
2216 get(i) ;
2217 if (varPtr->getVal()<lowest) {
2218 lowest = varPtr->getVal() ;
2219 }
2220 if (varPtr->getVal()>highest) {
2221 highest = varPtr->getVal() ;
2222 }
2223 }
2224
2225 if (marginFrac>0) {
2226 if (symMode==false) {
2227
2228 double margin = marginFrac*(highest-lowest) ;
2229 lowest -= margin ;
2230 highest += margin ;
2231 if (lowest<var.getMin()) lowest = var.getMin() ;
2232 if (highest>var.getMax()) highest = var.getMax() ;
2233
2234 } else {
2235
2236 double mom1 = moment(*varPtr,1) ;
2237 double delta = ((highest-mom1)>(mom1-lowest)?(highest-mom1):(mom1-lowest))*(1+marginFrac) ;
2238 lowest = mom1-delta ;
2239 highest = mom1+delta ;
2240 if (lowest<var.getMin()) lowest = var.getMin() ;
2241 if (highest>var.getMax()) highest = var.getMax() ;
2242
2243 }
2244 }
2245
2246 return false ;
2247}
2248
2249////////////////////////////////////////////////////////////////////////////////
2250/// Prepare dataset for use with cached constant terms listed in
2251/// 'cacheList' of expression 'arg'. Deactivate tree branches
2252/// for any dataset observable that is either not used at all,
2253/// or is used exclusively by cached branch nodes.
2254
2256{
2258
2259 // Add unused observables in this dataset to pruneSet
2260 pruneSet.add(*get()) ;
2261 std::unique_ptr<RooArgSet> usedObs{arg.getObservables(*this)};
2262 pruneSet.remove(*usedObs,true,true) ;
2263
2264 // Add observables exclusively used to calculate cached observables to pruneSet
2265 for(auto * var : *get()) {
2266 if (allClientsCached(var,cacheList)) {
2267 pruneSet.add(*var) ;
2268 }
2269 }
2270
2271
2272 if (!pruneSet.empty()) {
2273
2274 // Go over all used observables and check if any of them have parameterized
2275 // ranges in terms of pruned observables. If so, remove those observable
2276 // from the pruning list
2277 for(auto const* rrv : dynamic_range_cast<RooRealVar*>(*usedObs)) {
2278 if (rrv && !rrv->getBinning().isShareable()) {
2280 RooAbsReal* loFunc = rrv->getBinning().lowBoundFunc() ;
2281 RooAbsReal* hiFunc = rrv->getBinning().highBoundFunc() ;
2282 if (loFunc) {
2283 loFunc->leafNodeServerList(&depObs,nullptr,true) ;
2284 }
2285 if (hiFunc) {
2286 hiFunc->leafNodeServerList(&depObs,nullptr,true) ;
2287 }
2288 if (!depObs.empty()) {
2289 pruneSet.remove(depObs,true,true) ;
2290 }
2291 }
2292 }
2293 }
2294
2295
2296 // Remove all observables in keep list from prune list
2297 pruneSet.remove(keepObsList,true,true) ;
2298
2299 if (!pruneSet.empty()) {
2300
2301 // Deactivate tree branches here
2302 cxcoutI(Optimization) << "RooTreeData::optimizeReadingForTestStatistic(" << GetName() << "): Observables " << pruneSet
2303 << " 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 ;
2304 setArgStatus(pruneSet,false) ;
2305 }
2306}
2307
2308////////////////////////////////////////////////////////////////////////////////
2309/// Utility function that determines if all clients of object 'var'
2310/// appear in given list of cached nodes.
2311
2313{
2314 bool ret(true);
2315 bool anyClient(false);
2316
2317 for (const auto client : var->valueClients()) {
2318 anyClient = true ;
2319 if (!cacheList.find(client->GetName())) {
2320 // If client is not cached recurse
2321 ret &= allClientsCached(client,cacheList) ;
2322 }
2323 }
2324
2325 return anyClient?ret:false ;
2326}
2327
2328////////////////////////////////////////////////////////////////////////////////
2329
2331{
2332 _dstore->attachBuffers(extObs) ;
2333}
2334
2335////////////////////////////////////////////////////////////////////////////////
2336
2338{
2339 _dstore->resetBuffers() ;
2340}
2341
2342////////////////////////////////////////////////////////////////////////////////
2343
2345{
2346 return !_ownedComponents.empty();
2347}
2348
2349////////////////////////////////////////////////////////////////////////////////
2350
2352{
2353 auto i = _ownedComponents.find(name);
2354 return i==_ownedComponents.end() ? nullptr : i->second;
2355}
2356
2357////////////////////////////////////////////////////////////////////////////////
2358
2363
2364////////////////////////////////////////////////////////////////////////////////
2365/// Stream an object of class RooAbsData.
2366
2368{
2369 if (R__b.IsReading()) {
2370 R__b.ReadClassBuffer(RooAbsData::Class(),this);
2371 _namePtr = RooNameReg::instance().constPtr(GetName()) ;
2372
2373 // Convert on the fly to vector storage if that the current working default
2376 }
2377
2378 } else {
2379 R__b.WriteClassBuffer(RooAbsData::Class(),this);
2380 }
2381}
2382
2383////////////////////////////////////////////////////////////////////////////////
2384
2386{
2387 _dstore->checkInit() ;
2388}
2389
2390////////////////////////////////////////////////////////////////////////////////
2391/// Forward draw command to data store
2392
2394{
2395 if (_dstore) _dstore->Draw(option) ;
2396}
2397
2398////////////////////////////////////////////////////////////////////////////////
2399
2401{
2402 return _dstore->hasFilledCache() ;
2403}
2404
2405////////////////////////////////////////////////////////////////////////////////
2406/// Return a pointer to the TTree which stores the data. Returns a nullpointer
2407/// if vector-based storage is used. The RooAbsData remains owner of the tree.
2408/// GetClonedTree() can be used to get a tree even if the internal storage does not use one.
2409
2411{
2413 return static_cast<RooTreeDataStore&>(*_dstore).tree();
2414 } else {
2415 coutW(InputArguments) << "RooAbsData::tree(" << GetName() << ") WARNING: is not of StorageType::Tree. "
2416 << "Use GetClonedTree() instead or convert to tree storage." << std::endl;
2417 return nullptr;
2418 }
2419}
2420
2421////////////////////////////////////////////////////////////////////////////////
2422/// Return a clone of the TTree which stores the data or create such a tree
2423/// if vector storage is used. The user is responsible for deleting the tree
2424
2426{
2428 return static_cast<RooTreeDataStore&>(*_dstore).tree()->CloneTree();
2429 } else {
2430 RooTreeDataStore buffer(GetName(), GetTitle(), *get(), *_dstore);
2431 return buffer.tree()->CloneTree();
2432 }
2433}
2434
2435////////////////////////////////////////////////////////////////////////////////
2436/// Convert vector-based storage to tree-based storage
2437
2439{
2441 _dstore = std::make_unique<RooTreeDataStore>(GetName(), GetTitle(), _vars, *_dstore);
2443 }
2444}
2445
2446////////////////////////////////////////////////////////////////////////////////
2447/// If one of the TObject we have a referenced to is deleted, remove the
2448/// reference.
2449
2451{
2452 for(auto &iter : _ownedComponents) {
2453 if (iter.second == obj) {
2454 iter.second = nullptr;
2455 }
2456 }
2457}
2458
2459
2460////////////////////////////////////////////////////////////////////////////////
2461/// Sets the global observables stored in this data. A snapshot of the
2462/// observables will be saved.
2463/// \param[in] globalObservables The set of global observables to take a snapshot of.
2464
2465void RooAbsData::setGlobalObservables(RooArgSet const& globalObservables) {
2466 if(_globalObservables == nullptr) _globalObservables = std::make_unique<RooArgSet>();
2467 else _globalObservables->clear();
2468 globalObservables.snapshot(*_globalObservables);
2469 for(auto * arg : *_globalObservables) {
2470 arg->setAttribute("global",true);
2471 // Global observables are also always constant in fits
2472 if(auto lval = dynamic_cast<RooAbsRealLValue*>(arg)) lval->setConstant(true);
2473 if(auto lval = dynamic_cast<RooAbsCategoryLValue*>(arg)) lval->setConstant(true);
2474 }
2475}
2476
2477
2478////////////////////////////////////////////////////////////////////////////////
2479
2480void RooAbsData::SetName(const char* name)
2481{
2483 auto newPtr = RooNameReg::instance().constPtr(GetName()) ;
2484 if (newPtr != _namePtr) {
2485 //cout << "Rename '" << _namePtr->GetName() << "' to '" << name << "' (set flag in new name)" << std::endl;
2486 _namePtr = newPtr;
2489 }
2490}
2491
2492
2493
2494
2495////////////////////////////////////////////////////////////////////////////////
2496
2497void RooAbsData::SetNameTitle(const char *name, const char *title)
2498{
2499 TNamed::SetTitle(title) ;
2500 SetName(name);
2501}
2502
2503
2504
2505////////////////////////////////////////////////////////////////////////////////
2506/// Return sum of squared weights of this data.
2507
2509 const std::span<const double> eventWeights = getWeightBatch(0, numEntries(), /*sumW2=*/true);
2510 if (eventWeights.empty()) {
2511 return numEntries() * weightSquared();
2512 }
2513
2515 for (std::size_t i = 0; i < eventWeights.size(); ++i) {
2516 kahanWeight.AddIndexed(eventWeights[i], i);
2517 }
2518 return kahanWeight.Sum();
2519}
2520
2521
2522////////////////////////////////////////////////////////////////////////////////
2523/// Write information to retrieve data columns into `evalData.spans`.
2524/// All spans belonging to variables of this dataset are overwritten. Spans to other
2525/// variables remain intact.
2526/// \param begin Index of first event that ends up in the batch.
2527/// \param len Number of events in each batch.
2528RooAbsData::RealSpans RooAbsData::getBatches(std::size_t begin, std::size_t len) const {
2529 return store()->getBatches(begin, len);
2530}
2531
2532
2533RooAbsData::CategorySpans RooAbsData::getCategoryBatches(std::size_t first, std::size_t len) const {
2534 return store()->getCategoryBatches(first, len);
2535}
2536
2537////////////////////////////////////////////////////////////////////////////////
2538/// Create a TH2F histogram of the distribution of the specified variable
2539/// using this dataset. Apply any cuts to select which events are used.
2540/// The variable being plotted can either be contained directly in this
2541/// dataset, or else be a function of the variables in this dataset.
2542/// The histogram will be created using RooAbsReal::createHistogram() with
2543/// the name provided (with our dataset name prepended).
2544
2546 const char *name) const
2547{
2548 checkInit();
2549 return createHistogram(var1, var2, var1.getBins(), var2.getBins(), cuts, name);
2550}
2551
2552////////////////////////////////////////////////////////////////////////////////
2553/// Create a TH2F histogram of the distribution of the specified variable
2554/// using this dataset. Apply any cuts to select which events are used.
2555/// The variable being plotted can either be contained directly in this
2556/// dataset, or else be a function of the variables in this dataset.
2557/// The histogram will be created using RooAbsReal::createHistogram() with
2558/// the name provided (with our dataset name prepended).
2559
2561 const char *cuts, const char *name) const
2562{
2563 checkInit();
2564 static int counter(0);
2565
2566 std::unique_ptr<RooAbsReal> ownedPlotVarX;
2567 // Is this variable in our dataset?
2568 auto *plotVarX = static_cast<RooAbsReal *>(_vars.find(var1.GetName()));
2569 if (plotVarX == nullptr) {
2570 // Is this variable a client of our dataset?
2571 if (!var1.dependsOn(_vars)) {
2572 coutE(InputArguments) << GetName() << "::createHistogram: Argument " << var1.GetName()
2573 << " is not in dataset and is also not dependent on data set" << std::endl;
2574 return nullptr;
2575 }
2576
2577 // Clone derived variable
2578 ownedPlotVarX.reset(static_cast<RooAbsReal *>(var1.Clone()));
2579 plotVarX = ownedPlotVarX.get();
2580
2581 // Redirect servers of derived clone to internal ArgSet representing the data in this set
2582 plotVarX->redirectServers(const_cast<RooArgSet &>(_vars));
2583 }
2584
2585 std::unique_ptr<RooAbsReal> ownedPlotVarY;
2586 // Is this variable in our dataset?
2587 RooAbsReal *plotVarY = static_cast<RooAbsReal *>(_vars.find(var2.GetName()));
2588 if (plotVarY == nullptr) {
2589 // Is this variable a client of our dataset?
2590 if (!var2.dependsOn(_vars)) {
2591 coutE(InputArguments) << GetName() << "::createHistogram: Argument " << var2.GetName()
2592 << " is not in dataset and is also not dependent on data set" << std::endl;
2593 return nullptr;
2594 }
2595
2596 // Clone derived variable
2597 ownedPlotVarY.reset(static_cast<RooAbsReal *>(var2.Clone()));
2598 plotVarY = ownedPlotVarY.get();
2599
2600 // Redirect servers of derived clone to internal ArgSet representing the data in this set
2601 plotVarY->redirectServers(const_cast<RooArgSet &>(_vars));
2602 }
2603
2604 // Create selection formula if selection cuts are specified
2605 std::unique_ptr<RooFormula> select;
2606 if (nullptr != cuts && strlen(cuts)) {
2607 select = std::make_unique<RooFormula>(cuts, cuts, _vars);
2608 if (!select->ok()) {
2609 return nullptr;
2610 }
2611 }
2612
2613 const std::string histName = std::string{GetName()} + "_" + name + "_" + Form("%08x", counter++);
2614
2615 // create the histogram
2616 auto *histogram =
2617 new TH2F(histName.c_str(), "Events", nx, var1.getMin(), var1.getMax(), ny, var2.getMin(), var2.getMax());
2618 if (!histogram) {
2619 coutE(DataHandling) << GetName() << "::createHistogram: unable to create a new histogram" << std::endl;
2620 return nullptr;
2621 }
2622
2623 // Dump contents
2624 int nevent = numEntries();
2625 for (int i = 0; i < nevent; ++i) {
2626 get(i);
2627
2628 if (select && select->eval() == 0)
2629 continue;
2630 histogram->Fill(plotVarX->getVal(), plotVarY->getVal(), weight());
2631 }
2632
2633 return histogram;
2634}
2635
2636////////////////////////////////////////////////////////////////////////////////
2637/// Convert a string to the value of the RooAbsData::ErrorType enum with the
2638/// same name.
2640{
2641 using Map = std::unordered_map<std::string, RooAbsData::ErrorType>;
2642 static Map enumMap{{"Poisson", RooAbsData::Poisson},
2643 {"SumW2", RooAbsData::SumW2},
2644 {"None", RooAbsData::None},
2645 {"Auto", RooAbsData::Auto},
2646 {"Expected", RooAbsData::Expected}};
2647 auto found = enumMap.find(name);
2648 if (found == enumMap.end()) {
2649 std::stringstream msg;
2650 msg << "Unsupported error type type passed to DataError(). "
2651 "Supported decay types are : \"Poisson\", \"SumW2\", \"Auto\", \"Expected\", and None.";
2652 throw std::invalid_argument(msg.str());
2653 }
2654 return found->second;
2655}
#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
Definition RtypesCore.h:89
const char Option_t
Definition RtypesCore.h:66
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:2489
The Kahan summation is a compensated summation algorithm, which significantly reduces numerical error...
Definition Util.h:122
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:77
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 RooFit::OwningPtr< TList > split(const RooAbsCategory &splitCat, bool createEmptyDataSets=false) const
Split the dataset into subsets based on states of a categorical variable in this dataset.
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
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:40
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:59
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:34
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:2386
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:59
virtual Double_t GetBinError(Int_t bin) const
Return value of error associated to bin number bin.
Definition TH1.cxx:9065
virtual Int_t GetDimension() const
Definition TH1.h:299
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:9208
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3315
virtual TArrayD * GetSumw2()
Definition TH1.h:329
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:3649
virtual void Sumw2(Bool_t flag=kTRUE)
Create structure to store sum of squares of weights.
Definition TH1.cxx:9022
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:31
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
TObject * Clone(const char *newname="") const override
Make a clone of an object using the Streamer facility.
Definition TNamed.cxx:74
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
virtual void SetName(const char *name)
Set the name of the TNamed.
Definition TNamed.cxx:140
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition TNamed.cxx:51
Mother of all ROOT objects.
Definition TObject.h:41
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:225
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:798
A Pave (see TPave) with text, lines or/and boxes inside.
Definition TPaveText.h:21
Basic string class.
Definition TString.h:139
A TTree represents a columnar dataset.
Definition TTree.h:79
virtual TTree * CloneTree(Long64_t nentries=-1, Option_t *option="")
Create a clone of this tree and copy nentries.
Definition TTree.cxx:3132
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:2345