Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooAbsReal.h
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * File: $Id: RooAbsReal.h,v 1.75 2007/07/13 21:50:24 wouter Exp $
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#ifndef ROO_ABS_REAL
17#define ROO_ABS_REAL
18
19#include "RooAbsArg.h"
20#include "RooArgList.h"
21#include "RooArgProxy.h"
22#include "RooArgSet.h"
23#include "RooCmdArg.h"
24#include "RooCurve.h"
26#include "RooFit/EvalContext.h"
27#include "RooGlobalFunc.h"
28
29#include <ROOT/RSpan.hxx>
30
31#include <TList.h>
32#include <TObjString.h>
33
34class RooDataSet ;
35class RooPlot;
36class RooRealVar;
37class RooAbsFunc;
39class RooLinkedList ;
40class RooNumIntConfig ;
41class RooDataHist ;
42class RooFunctor ;
43class RooFitResult ;
44class RooAbsMoment ;
45class RooDerivative ;
47struct TreeReadBuffer; /// A space to attach TBranches
48namespace RooBatchCompute {
49struct RunContext;
50}
51
52class TH1;
53class TH1F;
54class TH2F;
55class TH3F;
56
57#include <iostream>
58#include <list>
59#include <map>
60#include <string>
61#include <sstream>
62
63class RooAbsReal : public RooAbsArg {
64public:
66
67 /// A RooAbsReal::Ref can be constructed from a `RooAbsReal&` or a `double`
68 /// that will be implicitly converted to a RooConstVar&. The RooAbsReal::Ref
69 /// can be used as a replacement for `RooAbsReal&`. With this type
70 /// definition, you can write RooFit interfaces that accept both RooAbsReal,
71 /// or simply a number that will be implicitly converted to a RooConstVar&.
72 class Ref {
73 public:
74 inline Ref(RooAbsReal &ref) : _ref{ref} {}
75 Ref(double val);
76 inline operator RooAbsReal &() const { return _ref; }
77
78 private:
80 };
81
82 // Constructors, assignment etc
83 RooAbsReal() ;
84 RooAbsReal(const char *name, const char *title, const char *unit= "") ;
85 RooAbsReal(const char *name, const char *title, double minVal, double maxVal,
86 const char *unit= "") ;
87 RooAbsReal(const RooAbsReal& other, const char* name=nullptr);
88 ~RooAbsReal() override;
89
90
91
92
93 //////////////////////////////////////////////////////////////////////////////////
94 /// Evaluate object. Returns either cached value or triggers a recalculation.
95 /// The recalculation happens by calling getValV(), which in the end calls the
96 /// virtual evaluate() functions of the respective PDFs.
97 /// \param[in] normalisationSet getValV() reacts differently depending on the value of the normalisation set.
98 /// If the set is `nullptr`, an unnormalised value is returned.
99 /// \note The normalisation is arbitrary, because it is up to the implementation
100 /// of the PDF to e.g. leave out normalisation constants for speed reasons. The range
101 /// of the variables is also ignored.
102 ///
103 /// To normalise the result properly, a RooArgSet has to be passed, which contains
104 /// the variables to normalise over.
105 /// These are integrated over their current ranges to compute the normalisation constant,
106 /// and the unnormalised result is divided by this value.
107 inline double getVal(const RooArgSet* normalisationSet = nullptr) const {
108 // Sometimes, the calling code uses an empty RooArgSet to request evaluation
109 // without normalization set instead of following the `nullptr` convention.
110 // To remove this ambiguity which might not always be correctly handled in
111 // downstream code, we set `normalisationSet` to nullptr if it is pointing
112 // to an empty set.
113 if(normalisationSet && normalisationSet->empty()) {
114 normalisationSet = nullptr;
115 }
116#ifdef ROOFIT_CHECK_CACHED_VALUES
118#else
119
120#ifndef _WIN32
122#else
123 return (_fast && !inhibitDirty()) ? _value : getValV(normalisationSet) ;
124#endif
125
126#endif
127 }
128
129 /// Like getVal(const RooArgSet*), but always requires an argument for normalisation.
130 inline double getVal(const RooArgSet& normalisationSet) const {
131 // Sometimes, the calling code uses an empty RooArgSet to request evaluation
132 // without normalization set instead of following the `nullptr` convention.
133 // To remove this ambiguity which might not always be correctly handled in
134 // downstream code, we set `normalisationSet` to nullptr if it is an empty set.
135 return _fast ? _value : getValV(normalisationSet.empty() ? nullptr : &normalisationSet) ;
136 }
137
138 double getVal(RooArgSet &&) const;
139
140 virtual double getValV(const RooArgSet* normalisationSet = nullptr) const ;
141
142 double getPropagatedError(const RooFitResult &fr, const RooArgSet &nset = {}) const;
143
144 bool operator==(double value) const ;
145 bool operator==(const RooAbsArg& other) const override;
146 bool isIdentical(const RooAbsArg& other, bool assumeSameType=false) const override;
147
148
149 inline const Text_t *getUnit() const {
150 // Return string with unit description
151 return _unit.Data();
152 }
153 inline void setUnit(const char *unit) {
154 // Set unit description to given string
155 _unit= unit;
156 }
157 TString getTitle(bool appendUnit= false) const;
158
159 // Lightweight interface adaptors (caller takes ownership)
160 RooFit::OwningPtr<RooAbsFunc> bindVars(const RooArgSet &vars, const RooArgSet* nset=nullptr, bool clipInvalid=false) const;
161
162 // Create a fundamental-type object that can hold our value.
163 RooFit::OwningPtr<RooAbsArg> createFundamental(const char* newname=nullptr) const override;
164
165 // Analytical integration support
166 virtual Int_t getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName=nullptr) const ;
167 virtual double analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=nullptr) const ;
168 virtual Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=nullptr) const ;
169 virtual double analyticalIntegral(Int_t code, const char* rangeName=nullptr) const ;
170 virtual bool forceAnalyticalInt(const RooAbsArg& /*dep*/) const {
171 // Interface to force RooRealIntegral to offer given observable for internal integration
172 // even if this is deemed unsafe. This default implementation returns always false
173 return false ;
174 }
175 virtual void forceNumInt(bool flag=true) {
176 // If flag is true, all advertised analytical integrals will be ignored
177 // and all integrals are calculated numerically
179 }
180 bool getForceNumInt() const { return _forceNumInt ; }
181
182 // Chi^2 fits to histograms
184 const RooCmdArg& arg3={}, const RooCmdArg& arg4={}, const RooCmdArg& arg5={},
185 const RooCmdArg& arg6={}, const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) ;
187
190 const RooCmdArg& arg3={}, const RooCmdArg& arg4={}, const RooCmdArg& arg5={},
191 const RooCmdArg& arg6={}, const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) ;
192
193 // Chi^2 fits to X-Y datasets
195 const RooCmdArg& arg3={}, const RooCmdArg& arg4={}, const RooCmdArg& arg5={},
196 const RooCmdArg& arg6={}, const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) ;
198
201 const RooCmdArg& arg3={}, const RooCmdArg& arg4={}, const RooCmdArg& arg5={},
202 const RooCmdArg& arg6={}, const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) ;
203
205
206
208 const RooCmdArg& arg3={}, const RooCmdArg& arg4={},
209 const RooCmdArg& arg5={}, const RooCmdArg& arg6={},
210 const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) const ;
211
212 /// Create integral over observables in iset in range named rangeName.
214 return createIntegral(iset,nullptr,nullptr,rangeName) ;
215 }
216 /// Create integral over observables in iset in range named rangeName with integrand normalized over observables in nset
217 RooFit::OwningPtr<RooAbsReal> createIntegral(const RooArgSet& iset, const RooArgSet& nset, const char* rangeName=nullptr) const {
218 return createIntegral(iset,&nset,nullptr,rangeName) ;
219 }
220 /// Create integral over observables in iset in range named rangeName with integrand normalized over observables in nset while
221 /// using specified configuration for any numeric integration.
222 RooFit::OwningPtr<RooAbsReal> createIntegral(const RooArgSet& iset, const RooArgSet& nset, const RooNumIntConfig& cfg, const char* rangeName=nullptr) const {
223 return createIntegral(iset,&nset,&cfg,rangeName) ;
224 }
225 /// Create integral over observables in iset in range named rangeName using specified configuration for any numeric integration.
227 return createIntegral(iset,nullptr,&cfg,rangeName) ;
228 }
229 virtual RooFit::OwningPtr<RooAbsReal> createIntegral(const RooArgSet& iset, const RooArgSet* nset=nullptr, const RooNumIntConfig* cfg=nullptr, const char* rangeName=nullptr) const ;
230
231
233
234 // Create running integrals
237 const RooCmdArg& arg3={}, const RooCmdArg& arg4={},
238 const RooCmdArg& arg5={}, const RooCmdArg& arg6={},
239 const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) ;
242
243
244 // Optimized accept/reject generator support
245 virtual Int_t getMaxVal(const RooArgSet& vars) const ;
246 virtual double maxVal(Int_t code) const ;
247 virtual Int_t minTrialSamples(const RooArgSet& /*arGenObs*/) const { return 0 ; }
248
249
250 // Plotting options
251 void setPlotLabel(const char *label);
252 const char *getPlotLabel() const;
253
254 virtual double defaultErrorLevel() const {
255 // Return default level for MINUIT error analysis
256 return 1.0 ;
257 }
258
259 const RooNumIntConfig* getIntegratorConfig() const ;
264 void setIntegratorConfig() ;
265 void setIntegratorConfig(const RooNumIntConfig& config) ;
266
267 virtual void fixAddCoefNormalization(const RooArgSet& addNormSet=RooArgSet(),bool force=true) ;
268 virtual void fixAddCoefRange(const char* rangeName=nullptr,bool force=true) ;
269
270 virtual void preferredObservableScanOrder(const RooArgSet& obs, RooArgSet& orderedObs) const ;
271
272 // User entry point for plotting
273 virtual RooPlot* plotOn(RooPlot* frame,
274 const RooCmdArg& arg1={}, const RooCmdArg& arg2={},
275 const RooCmdArg& arg3={}, const RooCmdArg& arg4={},
276 const RooCmdArg& arg5={}, const RooCmdArg& arg6={},
277 const RooCmdArg& arg7={}, const RooCmdArg& arg8={},
278 const RooCmdArg& arg9={}, const RooCmdArg& arg10={}
279 ) const ;
280
281
283
284 // Fill an existing histogram
285 TH1 *fillHistogram(TH1 *hist, const RooArgList &plotVars,
286 double scaleFactor= 1, const RooArgSet *projectedVars= nullptr, bool scaling=true,
287 const RooArgSet* condObs=nullptr, bool setError=true) const;
288
289 // Create 1,2, and 3D histograms from and fill it
291 TH1* createHistogram(const char *name, const RooAbsRealLValue& xvar, RooLinkedList& argList) const ;
292 TH1 *createHistogram(const char *name, const RooAbsRealLValue& xvar,
293 const RooCmdArg& arg1={}, const RooCmdArg& arg2={},
294 const RooCmdArg& arg3={}, const RooCmdArg& arg4={},
295 const RooCmdArg& arg5={}, const RooCmdArg& arg6={},
296 const RooCmdArg& arg7={}, const RooCmdArg& arg8={}) const ;
297
298 // Fill a RooDataHist
299 RooDataHist* fillDataHist(RooDataHist *hist, const RooArgSet* nset, double scaleFactor,
300 bool correctForBinVolume=false, bool showProgress=false) const ;
301
302 // I/O streaming interface (machine readable)
303 bool readFromStream(std::istream& is, bool compact, bool verbose=false) override ;
304 void writeToStream(std::ostream& os, bool compact) const override ;
305
306 // Printing interface (human readable)
307 void printValue(std::ostream& os) const override ;
308 void printMultiline(std::ostream& os, Int_t contents, bool verbose=false, TString indent="") const override ;
309
310 inline void setCachedValue(double value, bool notifyClients = true) final;
311
312 // Evaluation error logging
314 public:
316 EvalError(const EvalError& other) : _msg(other._msg), _srvval(other._srvval) { }
317 void setMessage(const char* tmp) { std::string s(tmp); s.swap(_msg); }
318 void setServerValues(const char* tmp) { std::string s(tmp); s.swap(_srvval); }
319 std::string _msg;
320 std::string _srvval;
321 } ;
322
324
325 /// Context to temporarily change the error logging mode as long as the context is alive.
339
342 void logEvalError(const char* message, const char* serverValueString=nullptr) const ;
343 static void logEvalError(const RooAbsReal* originator, const char* origName, const char* message, const char* serverValueString=nullptr) ;
344 static void printEvalErrors(std::ostream&os=std::cout, Int_t maxPerNode=10000000) ;
345 static Int_t numEvalErrors() ;
346 static Int_t numEvalErrorItems();
347 static std::map<const RooAbsArg *, std::pair<std::string, std::list<RooAbsReal::EvalError>>>::iterator evalErrorIter();
348
349 static void clearEvalErrorLog() ;
350
351 /// Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
352 virtual bool isBinnedDistribution(const RooArgSet& /*obs*/) const { return false ; }
353 virtual std::list<double>* binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const;
354 virtual std::list<double>* plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const;
355
356 RooFunctor* functor(const RooArgList& obs, const RooArgList& pars=RooArgList(), const RooArgSet& nset=RooArgSet()) const ;
357 TF1* asTF(const RooArgList& obs, const RooArgList& pars=RooArgList(), const RooArgSet& nset=RooArgSet()) const ;
358
359 RooDerivative* derivative(RooRealVar& obs, Int_t order=1, double eps=0.001) ;
360 RooDerivative* derivative(RooRealVar& obs, const RooArgSet& normSet, Int_t order, double eps=0.001) ;
361
362 RooAbsMoment* moment(RooRealVar& obs, Int_t order, bool central, bool takeRoot) ;
363 RooAbsMoment* moment(RooRealVar& obs, const RooArgSet& normObs, Int_t order, bool central, bool takeRoot, bool intNormObs) ;
364
365 RooAbsMoment* mean(RooRealVar& obs) { return moment(obs,1,false,false) ; }
366 RooAbsMoment* mean(RooRealVar& obs, const RooArgSet& nset) { return moment(obs,nset,1,false,false,true) ; }
367 RooAbsMoment* sigma(RooRealVar& obs) { return moment(obs,2,true,true) ; }
368 RooAbsMoment* sigma(RooRealVar& obs, const RooArgSet& nset) { return moment(obs,nset,2,true,true,true) ; }
369
370 double findRoot(RooRealVar& x, double xmin, double xmax, double yval) ;
371
372
373 virtual bool setData(RooAbsData& /*data*/, bool /*cloneData*/=true) { return true ; }
374
375 virtual void enableOffsetting(bool);
376 virtual bool isOffsetting() const { return false ; }
377 virtual double offset() const { return 0 ; }
378
379 static void setHideOffset(bool flag);
380 static bool hideOffset() ;
381
382 bool isSelectedComp() const ;
383 void selectComp(bool flag) {
384 // If flag is true, only selected component will be included in evaluates of RooAddPdf components
385 _selectComp = flag ;
386 }
387
390 RooArgSet *&cloneSet, const char* rangeName=nullptr, const RooArgSet* condObs=nullptr) const;
391 virtual void doEval(RooFit::EvalContext &) const;
392
393 virtual bool hasGradient() const { return false; }
394 virtual void gradient(double *) const {
395 if(!hasGradient()) throw std::runtime_error("RooAbsReal::gradient(double *) not implemented by this class!");
396 }
397
398 // PlotOn with command list
399 virtual RooPlot* plotOn(RooPlot* frame, RooLinkedList& cmdList) const ;
400
401protected:
403 friend class RooVectorDataStore;
404 friend class RooRealBinding;
405 friend class RooRealSumPdf;
406 friend class RooRealSumFunc;
407 friend class RooAddHelpers;
408 friend class RooAddPdf;
409 friend class RooAddModel;
410 friend class AddCacheElem;
412
413 // Hook for objects with normalization-dependent parameters interpretation
414 virtual void selectNormalization(const RooArgSet* depSet=nullptr, bool force=false) ;
415 virtual void selectNormalizationRange(const char* rangeName=nullptr, bool force=false) ;
416
417 // Helper functions for plotting
418 bool plotSanityChecks(RooPlot* frame) const ;
419 void makeProjectionSet(const RooAbsArg* plotVar, const RooArgSet* allVars,
420 RooArgSet& projectedVars, bool silent) const ;
421
422 TString integralNameSuffix(const RooArgSet& iset, const RooArgSet* nset=nullptr, const char* rangeName=nullptr, bool omitEmpty=false) const ;
423
424 void plotOnCompSelect(RooArgSet* selNodes) const ;
425 RooPlot* plotOnWithErrorBand(RooPlot* frame,const RooFitResult& fr, double Z, const RooArgSet* params, const RooLinkedList& argList, bool method1) const ;
426
427 template<typename... Proxies>
428 bool matchArgs(const RooArgSet& allDeps, RooArgSet& analDeps, const RooArgProxy& a, const Proxies&... proxies) const
429 {
431 // Fold expression to add all proxy names to the list
432 nameList.Add(new TObjString(a.absArg()->GetName()));
433 (nameList.Add(new TObjString(proxies.absArg()->GetName())), ...);
434
436 nameList.Delete(); // Clean up the list contents
437 return result;
438 }
439
441 const RooArgSet& set) const ;
442
443 RooFit::OwningPtr<RooAbsReal> createIntObj(const RooArgSet& iset, const RooArgSet* nset, const RooNumIntConfig* cfg, const char* rangeName) const ;
444 void findInnerMostIntegration(const RooArgSet& allObs, RooArgSet& innerObs, const char* rangeName) const ;
445
446 // Internal consistency checking (needed by RooDataSet)
447 /// Check if current value is valid.
448 bool isValid() const override { return isValidReal(_value); }
449 /// Interface function to check if given value is a valid value for this object. Returns true unless overridden.
450 virtual bool isValidReal(double /*value*/, bool printError = false) const { (void)printError; return true; }
451
452 // Function evaluation and error tracing
453 double traceEval(const RooArgSet* set) const ;
454
455 /// Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
456 virtual double evaluate() const = 0;
457
458 // Hooks for RooDataSet interface
459 void syncCache(const RooArgSet* set=nullptr) override { getVal(set) ; }
460 void copyCache(const RooAbsArg* source, bool valueOnly=false, bool setValDirty=true) override ;
461 void attachToTree(TTree& t, Int_t bufSize=32000) override ;
463 void setTreeBranchStatus(TTree& t, bool active) override ;
464 void fillTreeBranch(TTree& t) override ;
465
466 struct PlotOpt {
468 double scaleFactor = 1.0;
470 const RooAbsData *projData = nullptr;
471 bool binProjData = false;
472 const RooArgSet *projSet = nullptr;
473 double precision = 1e-3;
474 bool shiftToZero = false;
475 const RooArgSet *projDataSet = nullptr;
476 const char *normRangeName = nullptr;
477 double rangeLo = 0.0;
478 double rangeHi = 0.0;
479 bool postRangeFracScale = false;
481 const char *projectionRangeName = nullptr;
482 bool curveInvisible = false;
483 const char *curveName = nullptr;
484 const char *addToCurveName = nullptr;
485 double addToWgtSelf = 1.0;
486 double addToWgtOther = 1.0;
489 const char *curveNameSuffix = "";
491 double eeval = 0.0;
492 bool doeeval = false;
493 bool progress = false;
494 const RooFitResult *errorFR = nullptr;
495 };
496
497 // Plot implementation functions
498 virtual RooPlot *plotOn(RooPlot* frame, PlotOpt o) const;
499
500 virtual RooPlot *plotAsymOn(RooPlot *frame, const RooAbsCategoryLValue& asymCat, PlotOpt o) const;
501
503
505 bool nameChange, bool isRecursiveStep) override;
506
507 static void globalSelectComp(bool flag) ;
508
509 // This struct can be used to flip the global switch to select components.
510 // Doing this with RAII prevents forgetting to reset the state.
525
526
527private:
528
529 /// Debug version of getVal(), which is slow and does error checking.
530 double _DEBUG_getVal(const RooArgSet* normalisationSet) const;
531
532 //--------------------------------------------------------------------
533
534 protected:
535
536 double _plotMin = 0.0; ///< Minimum of plot range
537 double _plotMax = 0.0; ///< Maximum of plot range
538 Int_t _plotBins = 100; ///< Number of plot bins
539 mutable double _value = 0.0; ///< Cache for current value of object
540 TString _unit; ///< Unit for objects value
541 TString _label; ///< Plot label for objects value
542 bool _forceNumInt = false; ///< Force numerical integration if flag set
543 std::unique_ptr<RooNumIntConfig> _specIntegratorConfig; // Numeric integrator configuration specific for this object
544 TreeReadBuffer *_treeReadBuffer = nullptr; //! A buffer for reading values from trees
545 bool _selectComp = true; //! Component selection flag for RooAbsPdf::plotCompOn
547
548 static bool _globalSelectComp; // Global activation switch for component selection
549 static bool _hideOffset; ///< Offset hiding flag
550
551 ClassDefOverride(RooAbsReal,3); // Abstract real-valued variable
552};
553
554
555////////////////////////////////////////////////////////////////////////////////
556/// Overwrite the value stored in this object's cache.
557/// This can be used to fake a computation that resulted in `value`.
558/// \param[in] value Value to write.
559/// \param[in] notifyClients If true, notify users of this object that its value changed.
560/// This is the default.
562 _value = value;
563
564 if (notifyClients) {
566 _valueDirty = false;
567 }
568}
569
570
571#endif
#define a(i)
Definition RSha256.hxx:99
#define e(i)
Definition RSha256.hxx:103
char Text_t
General string (char)
Definition RtypesCore.h:76
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
#define ClassDefOverride(name, id)
Definition Rtypes.h:348
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.
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 result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
char name[80]
Definition TGX11.cxx:110
float xmin
float xmax
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:76
bool _fast
Definition RooAbsArg.h:645
bool _valueDirty
Definition RooAbsArg.h:641
void setValueDirty()
Mark the element dirty. This forces a re-evaluation when a value is requested.
Definition RooAbsArg.h:425
bool inhibitDirty() const
Delete watch flag.
static bool _inhibitDirty
Definition RooAbsArg.h:625
Abstract base class for objects that represent a discrete value that can be set from the outside,...
Abstract container object that can hold multiple RooAbsArg objects.
Abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
Abstract interface for evaluating a real-valued function of one real variable and performing numerica...
Definition RooAbsFunc.h:27
Abstract base class for objects that represent a real value that may appear on the left hand side of ...
Context to temporarily change the error logging mode as long as the context is alive.
Definition RooAbsReal.h:326
EvalErrorContext(ErrorLoggingMode m)
Definition RooAbsReal.h:328
EvalErrorContext & operator=(EvalErrorContext const &)=delete
EvalErrorContext & operator=(EvalErrorContext &&)=delete
EvalErrorContext(EvalErrorContext const &)=delete
EvalErrorContext(EvalErrorContext &&)=delete
void setServerValues(const char *tmp)
Definition RooAbsReal.h:318
void setMessage(const char *tmp)
Definition RooAbsReal.h:317
EvalError(const EvalError &other)
Definition RooAbsReal.h:316
A RooAbsReal::Ref can be constructed from a RooAbsReal& or a double that will be implicitly converted...
Definition RooAbsReal.h:72
Ref(RooAbsReal &ref)
Definition RooAbsReal.h:74
RooAbsReal & _ref
Definition RooAbsReal.h:79
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:63
RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, double scaleFactor, bool correctForBinVolume=false, bool showProgress=false) const
Fill a RooDataHist with values sampled from this function at the bin centers.
virtual void selectNormalizationRange(const char *rangeName=nullptr, bool force=false)
Interface function to force use of a given normalization range to interpret function value.
void plotOnCompSelect(RooArgSet *selNodes) const
Helper function for plotting of composite p.d.fs.
bool isSelectedComp() const
If true, the current pdf is a selected component (for use in plotting)
virtual std::list< double > * binBoundaries(RooAbsRealLValue &obs, double xlo, double xhi) const
Retrieve bin boundaries if this distribution is binned in obs.
void selectComp(bool flag)
Definition RooAbsReal.h:383
TString _label
Plot label for objects value.
Definition RooAbsReal.h:541
bool _selectComp
A buffer for reading values from trees.
Definition RooAbsReal.h:545
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:107
double findRoot(RooRealVar &x, double xmin, double xmax, double yval)
Return value of x (in range xmin,xmax) at which function equals yval.
RooAbsMoment * mean(RooRealVar &obs)
Definition RooAbsReal.h:365
TreeReadBuffer * _treeReadBuffer
Definition RooAbsReal.h:544
virtual bool isOffsetting() const
Definition RooAbsReal.h:376
virtual double getValV(const RooArgSet *normalisationSet=nullptr) const
Return value of object.
static Int_t numEvalErrorItems()
RooAbsReal()
coverity[UNINIT_CTOR] Default constructor
virtual void fixAddCoefNormalization(const RooArgSet &addNormSet=RooArgSet(), bool force=true)
Fix the interpretation of the coefficient of any RooAddPdf component in the expression tree headed by...
virtual bool hasGradient() const
Definition RooAbsReal.h:393
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=nullptr) const
Variant of getAnalyticalIntegral that is also passed the normalization set that should be applied to ...
bool isValid() const override
Check if current value is valid.
Definition RooAbsReal.h:448
bool _forceNumInt
Force numerical integration if flag set.
Definition RooAbsReal.h:542
~RooAbsReal() override
Destructor.
void setParameterizeIntegral(const RooArgSet &paramVars)
bool matchArgsByName(const RooArgSet &allArgs, RooArgSet &matchedArgs, const TList &nameList) const
Check if allArgs contains matching elements for each name in nameList.
static bool hideOffset()
void setTreeBranchStatus(TTree &t, bool active) override
(De)Activate associated tree branch
TH1 * fillHistogram(TH1 *hist, const RooArgList &plotVars, double scaleFactor=1, const RooArgSet *projectedVars=nullptr, bool scaling=true, const RooArgSet *condObs=nullptr, bool setError=true) const
Fill the ROOT histogram 'hist' with values sampled from this function at the bin centers.
RooFit::OwningPtr< RooAbsReal > createScanRI(const RooArgSet &iset, const RooArgSet &nset, Int_t numScanBins, Int_t intOrder)
Utility function for createRunningIntegral that construct an object implementing the numeric scanning...
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const char *rangeName) const
Create integral over observables in iset in range named rangeName.
Definition RooAbsReal.h:213
virtual double defaultErrorLevel() const
Definition RooAbsReal.h:254
double _DEBUG_getVal(const RooArgSet *normalisationSet) const
Debug version of getVal(), which is slow and does error checking.
RooFit::OwningPtr< RooAbsArg > createFundamental(const char *newname=nullptr) const override
Create a RooRealVar fundamental object with our properties.
bool plotSanityChecks(RooPlot *frame) const
Utility function for plotOn(), perform general sanity check on frame to ensure safe plotting operatio...
RooFit::OwningPtr< RooAbsFunc > bindVars(const RooArgSet &vars, const RooArgSet *nset=nullptr, bool clipInvalid=false) const
Create an interface adaptor f(vars) that binds us to the specified variables (in arbitrary order).
virtual void selectNormalization(const RooArgSet *depSet=nullptr, bool force=false)
Interface function to force use of a given set of observables to interpret function value.
friend class AddCacheElem
Definition RooAbsReal.h:410
RooDerivative * derivative(RooRealVar &obs, Int_t order=1, double eps=0.001)
Return function representing first, second or third order derivative of this function.
virtual RooFit::OwningPtr< RooFitResult > chi2FitTo(RooDataHist &data, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={})
Perform a fit to given histogram.
TF1 * asTF(const RooArgList &obs, const RooArgList &pars=RooArgList(), const RooArgSet &nset=RooArgSet()) const
Return a ROOT TF1,2,3 object bound to this RooAbsReal with given definition of observables and parame...
virtual bool setData(RooAbsData &, bool=true)
Definition RooAbsReal.h:373
TString _unit
Unit for objects value.
Definition RooAbsReal.h:540
static RooNumIntConfig * defaultIntegratorConfig()
Returns the default numeric integration configuration for all RooAbsReals.
double getVal(const RooArgSet &normalisationSet) const
Like getVal(const RooArgSet*), but always requires an argument for normalisation.
Definition RooAbsReal.h:130
bool readFromStream(std::istream &is, bool compact, bool verbose=false) override
Read object contents from stream (dummy for now)
void fillTreeBranch(TTree &t) override
Fill the tree branch that associated with this object with its current value.
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooArgSet &nset, const char *rangeName=nullptr) const
Create integral over observables in iset in range named rangeName with integrand normalized over obse...
Definition RooAbsReal.h:217
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Structure printing.
virtual RooPlot * plotAsymOn(RooPlot *frame, const RooAbsCategoryLValue &asymCat, PlotOpt o) const
bool operator==(double value) const
Equality operator comparing to a double.
static ErrorLoggingMode evalErrorLoggingMode()
Return current evaluation error logging mode.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursiveStep) override
Function that is called at the end of redirectServers().
virtual Int_t minTrialSamples(const RooArgSet &) const
Definition RooAbsReal.h:247
virtual bool isValidReal(double, bool printError=false) const
Interface function to check if given value is a valid value for this object. Returns true unless over...
Definition RooAbsReal.h:450
void setIntegratorConfig()
Remove the specialized numeric integration configuration associated with this object.
void syncCache(const RooArgSet *set=nullptr) override
Definition RooAbsReal.h:459
void printValue(std::ostream &os) const override
Print object value.
virtual bool forceAnalyticalInt(const RooAbsArg &) const
Definition RooAbsReal.h:170
bool isIdentical(const RooAbsArg &other, bool assumeSameType=false) const override
void setUnit(const char *unit)
Definition RooAbsReal.h:153
friend class RooAddHelpers
Definition RooAbsReal.h:407
bool getForceNumInt() const
Definition RooAbsReal.h:180
virtual RooFit::OwningPtr< RooAbsReal > createProfile(const RooArgSet &paramsOfInterest)
Create a RooProfileLL object that eliminates all nuisance parameters in the present function.
static bool _hideOffset
Offset hiding flag.
Definition RooAbsReal.h:549
void attachToVStore(RooVectorDataStore &vstore) override
void copyCache(const RooAbsArg *source, bool valueOnly=false, bool setValDirty=true) override
Copy the cached value of another RooAbsArg to our cache.
virtual void gradient(double *) const
Definition RooAbsReal.h:394
TH1 * createHistogram(RooStringView varNameList, Int_t xbins=0, Int_t ybins=0, Int_t zbins=0) const
Create and fill a ROOT histogram TH1, TH2 or TH3 with the values of this function for the variables w...
virtual void fixAddCoefRange(const char *rangeName=nullptr, bool force=true)
Fix the interpretation of the coefficient of any RooAddPdf component in the expression tree headed by...
double _value
Cache for current value of object.
Definition RooAbsReal.h:539
virtual double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
void attachToTree(TTree &t, Int_t bufSize=32000) override
Attach object to a branch of given TTree.
friend class BatchInterfaceAccessor
Definition RooAbsReal.h:402
RooNumIntConfig * specialIntegratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
void writeToStream(std::ostream &os, bool compact) const override
Write object contents to stream (dummy for now)
double traceEval(const RooArgSet *set) const
Calculate current value of object, with error tracing wrapper.
double getPropagatedError(const RooFitResult &fr, const RooArgSet &nset={}) const
Propagates parameter uncertainties to an uncertainty estimate for this RooAbsReal.
static void setHideOffset(bool flag)
static void globalSelectComp(bool flag)
Global switch controlling the activation of the selectComp() functionality.
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooArgSet &nset, const RooNumIntConfig &cfg, const char *rangeName=nullptr) const
Create integral over observables in iset in range named rangeName with integrand normalized over obse...
Definition RooAbsReal.h:222
RooAbsMoment * moment(RooRealVar &obs, Int_t order, bool central, bool takeRoot)
Return function representing moment of function of given order.
RooPlot * plotOnWithErrorBand(RooPlot *frame, const RooFitResult &fr, double Z, const RooArgSet *params, const RooLinkedList &argList, bool method1) const
Plot function or PDF on frame with support for visualization of the uncertainty encoded in the given ...
RooFit::UniqueId< RooArgSet >::Value_t _lastNormSetId
Component selection flag for RooAbsPdf::plotCompOn.
Definition RooAbsReal.h:546
const char * getPlotLabel() const
Get the label associated with the variable.
virtual void forceNumInt(bool flag=true)
Definition RooAbsReal.h:175
RooFit::OwningPtr< RooAbsReal > createRunningIntegral(const RooArgSet &iset, const RooArgSet &nset={})
Calls createRunningIntegral(const RooArgSet&, const RooCmdArg&, const RooCmdArg&, const RooCmdArg&,...
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooNumIntConfig &cfg, const char *rangeName=nullptr) const
Create integral over observables in iset in range named rangeName using specified configuration for a...
Definition RooAbsReal.h:226
std::unique_ptr< RooNumIntConfig > _specIntegratorConfig
Definition RooAbsReal.h:543
virtual Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
static std::map< constRooAbsArg *, std::pair< std::string, std::list< RooAbsReal::EvalError > > >::iterator evalErrorIter()
static Int_t numEvalErrors()
Return the number of logged evaluation errors since the last clearing.
static void setEvalErrorLoggingMode(ErrorLoggingMode m)
Set evaluation error logging mode.
double _plotMax
Maximum of plot range.
Definition RooAbsReal.h:537
virtual void preferredObservableScanOrder(const RooArgSet &obs, RooArgSet &orderedObs) const
Interface method for function objects to indicate their preferred order of observables for scanning t...
virtual double maxVal(Int_t code) const
Return maximum value for set of observables identified by code assigned in getMaxVal.
void findInnerMostIntegration(const RooArgSet &allObs, RooArgSet &innerObs, const char *rangeName) const
Utility function for createIntObj() that aids in the construct of recursive integrals over functions ...
TString integralNameSuffix(const RooArgSet &iset, const RooArgSet *nset=nullptr, const char *rangeName=nullptr, bool omitEmpty=false) const
Construct string with unique suffix name to give to integral object that encodes integrated observabl...
virtual double evaluate() const =0
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
TString getTitle(bool appendUnit=false) const
Return this variable's title string.
void logEvalError(const char *message, const char *serverValueString=nullptr) const
Log evaluation error message.
virtual double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
const Text_t * getUnit() const
Definition RooAbsReal.h:149
const RooNumIntConfig * getIntegratorConfig() const
Return the numeric integration configuration used for this object.
virtual bool isBinnedDistribution(const RooArgSet &) const
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition RooAbsReal.h:352
static void printEvalErrors(std::ostream &os=std::cout, Int_t maxPerNode=10000000)
Print all outstanding logged evaluation error on the given ostream.
RooFit::OwningPtr< RooAbsReal > createIntRI(const RooArgSet &iset, const RooArgSet &nset={})
Utility function for createRunningIntegral.
virtual void enableOffsetting(bool)
static void clearEvalErrorLog()
Clear the stack of evaluation error messages.
RooFit::OwningPtr< RooAbsReal > createIntObj(const RooArgSet &iset, const RooArgSet *nset, const RooNumIntConfig *cfg, const char *rangeName) const
Internal utility function for createIntegral() that creates the actual integral object.
RooAbsMoment * sigma(RooRealVar &obs, const RooArgSet &nset)
Definition RooAbsReal.h:368
RooFunctor * functor(const RooArgList &obs, const RooArgList &pars=RooArgList(), const RooArgSet &nset=RooArgSet()) const
Return a RooFunctor object bound to this RooAbsReal with given definition of observables and paramete...
void setCachedValue(double value, bool notifyClients=true) final
Overwrite the value stored in this object's cache.
Definition RooAbsReal.h:561
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 RooCmdArg &arg9={}, const RooCmdArg &arg10={}) const
Plot (project) PDF on specified frame.
Int_t _plotBins
Number of plot bins.
Definition RooAbsReal.h:538
const RooAbsReal * createPlotProjection(const RooArgSet &depVars, const RooArgSet &projVars, RooArgSet *&cloneSet) const
Utility function for plotOn() that creates a projection of a function or p.d.f to be plotted on a Roo...
virtual std::list< double > * plotSamplingHint(RooAbsRealLValue &obs, double xlo, double xhi) const
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
void setPlotLabel(const char *label)
Set the label associated with this variable.
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, 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 an object that represents the integral of the function over one or more observables listed in ...
virtual void doEval(RooFit::EvalContext &) const
Base function for computing multiple values of a RooAbsReal.
void makeProjectionSet(const RooAbsArg *plotVar, const RooArgSet *allVars, RooArgSet &projectedVars, bool silent) const
Utility function for plotOn() that constructs the set of observables to project when plotting ourselv...
virtual Int_t getMaxVal(const RooArgSet &vars) const
Advertise capability to determine maximum value of function for given set of observables.
double _plotMin
Minimum of plot range.
Definition RooAbsReal.h:536
bool matchArgs(const RooArgSet &allDeps, RooArgSet &analDeps, const RooArgProxy &a, const Proxies &... proxies) const
Definition RooAbsReal.h:428
RooAbsMoment * mean(RooRealVar &obs, const RooArgSet &nset)
Definition RooAbsReal.h:366
virtual RooFit::OwningPtr< RooAbsReal > createChi2(RooDataHist &data, const RooLinkedList &cmdList)
virtual double offset() const
Definition RooAbsReal.h:377
RooAbsMoment * sigma(RooRealVar &obs)
Definition RooAbsReal.h:367
static bool _globalSelectComp
Definition RooAbsReal.h:548
RooAddModel is an efficient implementation of a sum of PDFs of the form.
Definition RooAddModel.h:27
Efficient implementation of a sum of PDFs of the form.
Definition RooAddPdf.h:33
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
Abstract interface for RooAbsArg proxy classes.
Definition RooArgProxy.h:24
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
Named container for two doubles, two integers two object points and three string pointers that can be...
Definition RooCmdArg.h:26
@ Extended
Definition RooCurve.h:39
Container class to hold N-dimensional binned data.
Definition RooDataHist.h:40
Container class to hold unbinned data.
Definition RooDataSet.h:32
Represents the first, second, or third order derivative of any RooAbsReal as calculated (numerically)...
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Lightweight interface adaptor that exports a RooAbsPdf as a functor.
Definition RooFunctor.h:25
Collection class for internal use, storing a collection of RooAbsArg pointers in a doubly linked list...
Holds the configuration parameters of the various numeric integrators used by RooRealIntegral.
Plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
Lightweight interface adaptor that binds a RooAbsReal object to a subset of its servers and present i...
Implements a PDF constructed from a sum of functions:
Variable that can be changed from the outside.
Definition RooRealVar.h:37
The RooStringView is a wrapper around a C-style string that can also be constructed from a std::strin...
Uses std::vector to store data columns.
1-Dim function class
Definition TF1.h:182
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:879
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:109
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:307
3-D histogram with a float per channel (see TH1 documentation)
Definition TH3.h:326
A doubly linked list.
Definition TList.h:38
Collectable string class.
Definition TObjString.h:28
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
A TTree represents a columnar dataset.
Definition TTree.h:89
Double_t x[n]
Definition legend1.C:17
Namespace for dispatching RooFit computations to various backends.
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
RooCurve::WingMode wmode
Definition RooAbsReal.h:480
const char * normRangeName
Definition RooAbsReal.h:476
RooFit::MPSplit interleave
Definition RooAbsReal.h:488
const char * projectionRangeName
Definition RooAbsReal.h:481
const RooArgSet * projDataSet
Definition RooAbsReal.h:475
const char * curveNameSuffix
Definition RooAbsReal.h:489
const char * addToCurveName
Definition RooAbsReal.h:484
const RooFitResult * errorFR
Definition RooAbsReal.h:494
const RooArgSet * projSet
Definition RooAbsReal.h:472
const char * curveName
Definition RooAbsReal.h:483
const RooAbsData * projData
Definition RooAbsReal.h:470
Option_t * drawOptions
Definition RooAbsReal.h:467
A UniqueId can be added as a class member to enhance any class with a unique identifier for each inst...
Definition UniqueId.h:39
unsigned long Value_t
Definition UniqueId.h:41
TMarker m
Definition textangle.C:8