Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
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 "RooCmdArg.h"
21 #include "RooCurve.h"
22 #include "RooArgSet.h"
23 #include "RooArgList.h"
24 #include "RooGlobalFunc.h"
25 #include "RooSpan.h"
26 #include "BatchData.h"
27 
28 class RooArgList ;
29 class RooDataSet ;
30 class RooPlot;
31 class RooRealVar;
32 class RooAbsFunc;
33 class RooAbsCategoryLValue ;
34 class RooCategory ;
35 class RooLinkedList ;
36 class RooNumIntConfig ;
37 class RooDataHist ;
38 class RooFunctor ;
39 class RooGenFunction ;
40 class RooMultiGenFunction ;
41 class RooFitResult ;
42 class RooAbsMoment ;
43 class RooDerivative ;
44 class RooVectorDataStore ;
45 namespace RooHelpers {
46 class BatchInterfaceAccessor;
47 }
48 
49 class TH1;
50 class TH1F;
51 class TH2F;
52 class TH3F;
53 
54 #include <list>
55 #include <string>
56 #include <iostream>
57 #include <sstream>
58 
59 class RooAbsReal : public RooAbsArg {
60 public:
61  // Constructors, assignment etc
62  RooAbsReal() ;
63  RooAbsReal(const char *name, const char *title, const char *unit= "") ;
64  RooAbsReal(const char *name, const char *title, Double_t minVal, Double_t maxVal,
65  const char *unit= "") ;
66  RooAbsReal(const RooAbsReal& other, const char* name=0);
67  RooAbsReal& operator=(const RooAbsReal& other);
68  virtual ~RooAbsReal();
69 
70 
71 
72 
73  //////////////////////////////////////////////////////////////////////////////////
74  /// Evaluate object. Returns either cached value or triggers a recalculation.
75  /// The recalculation happens by calling getValV(), which in the end calls the
76  /// virtual evaluate() functions of the respective PDFs.
77  /// \param[in] normalisationSet getValV() reacts differently depending on the value of the normalisation set.
78  /// If the set is `nullptr`, an unnormalised value is returned.
79  /// \note The normalisation is arbitrary, because it is up to the implementation
80  /// of the PDF to e.g. leave out normalisation constants for speed reasons. The range
81  /// of the variables is also ignored.
82  ///
83  /// To normalise the result properly, a RooArgSet has to be passed, which contains
84  /// the variables to normalise over.
85  /// These are integrated over their current ranges to compute the normalisation constant,
86  /// and the unnormalised result is divided by this value.
87  inline Double_t getVal(const RooArgSet* normalisationSet = nullptr) const {
88 #ifdef ROOFIT_CHECK_CACHED_VALUES
89  return _DEBUG_getVal(normalisationSet);
90 #else
91 
92 #ifndef _WIN32
93  return (_fast && !_inhibitDirty) ? _value : getValV(normalisationSet) ;
94 #else
95  return (_fast && !inhibitDirty()) ? _value : getValV(normalisationSet) ;
96 #endif
97 
98 #endif
99  }
100 
101  /// Like getVal(const RooArgSet*), but always requires an argument for normalisation.
102  inline Double_t getVal(const RooArgSet& normalisationSet) const { return _fast ? _value : getValV(&normalisationSet) ; }
103 
104  virtual Double_t getValV(const RooArgSet* normalisationSet = nullptr) const ;
105 
106  virtual RooSpan<const double> getValBatch(std::size_t begin, std::size_t maxSize, const RooArgSet* normSet = nullptr) const;
107 
108  Double_t getPropagatedError(const RooFitResult &fr, const RooArgSet &nset = RooArgSet()) const;
109 
110  Bool_t operator==(Double_t value) const ;
111  virtual Bool_t operator==(const RooAbsArg& other) ;
112  virtual Bool_t isIdentical(const RooAbsArg& other, Bool_t assumeSameType=kFALSE) ;
113 
114 
115  inline const Text_t *getUnit() const {
116  // Return string with unit description
117  return _unit.Data();
118  }
119  inline void setUnit(const char *unit) {
120  // Set unit description to given string
121  _unit= unit;
122  }
123  TString getTitle(Bool_t appendUnit= kFALSE) const;
124 
125  // Lightweight interface adaptors (caller takes ownership)
126  RooAbsFunc *bindVars(const RooArgSet &vars, const RooArgSet* nset=0, Bool_t clipInvalid=kFALSE) const;
127 
128  // Create a fundamental-type object that can hold our value.
129  RooAbsArg *createFundamental(const char* newname=0) const;
130 
131  // Analytical integration support
132  virtual Int_t getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName=0) const ;
133  virtual Double_t analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=0) const ;
134  virtual Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=0) const ;
135  virtual Double_t analyticalIntegral(Int_t code, const char* rangeName=0) const ;
136  virtual Bool_t forceAnalyticalInt(const RooAbsArg& /*dep*/) const {
137  // Interface to force RooRealIntegral to offer given observable for internal integration
138  // even if this is deemed unsafe. This default implementation returns always flase
139  return kFALSE ;
140  }
141  virtual void forceNumInt(Bool_t flag=kTRUE) {
142  // If flag is true, all advertised analytical integrals will be ignored
143  // and all integrals are calculated numerically
144  _forceNumInt = flag ;
145  }
146  Bool_t getForceNumInt() const { return _forceNumInt ; }
147 
148  // Chi^2 fits to histograms
149  virtual RooFitResult* chi2FitTo(RooDataHist& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
150  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
151  const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
152  virtual RooFitResult* chi2FitTo(RooDataHist& data, const RooLinkedList& cmdList) ;
153 
154  virtual RooAbsReal* createChi2(RooDataHist& data, const RooLinkedList& cmdList) ;
155  virtual RooAbsReal* createChi2(RooDataHist& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
156  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
157  const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
158 
159  // Chi^2 fits to X-Y datasets
160  virtual RooFitResult* chi2FitTo(RooDataSet& xydata, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
161  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
162  const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
163  virtual RooFitResult* chi2FitTo(RooDataSet& xydata, const RooLinkedList& cmdList) ;
164 
165  virtual RooAbsReal* createChi2(RooDataSet& data, const RooLinkedList& cmdList) ;
166  virtual RooAbsReal* createChi2(RooDataSet& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
167  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
168  const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
169 
170 
171  virtual RooAbsReal* createProfile(const RooArgSet& paramsOfInterest) ;
172 
173 
174  RooAbsReal* createIntegral(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2=RooCmdArg::none(),
175  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
176  const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
177  const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) const ;
178 
179  /// Create integral over observables in iset in range named rangeName.
180  RooAbsReal* createIntegral(const RooArgSet& iset, const char* rangeName) const {
181  return createIntegral(iset,0,0,rangeName) ;
182  }
183  /// Create integral over observables in iset in range named rangeName with integrand normalized over observables in nset
184  RooAbsReal* createIntegral(const RooArgSet& iset, const RooArgSet& nset, const char* rangeName=0) const {
185  return createIntegral(iset,&nset,0,rangeName) ;
186  }
187  /// Create integral over observables in iset in range named rangeName with integrand normalized over observables in nset while
188  /// using specified configuration for any numeric integration.
189  RooAbsReal* createIntegral(const RooArgSet& iset, const RooArgSet& nset, const RooNumIntConfig& cfg, const char* rangeName=0) const {
190  return createIntegral(iset,&nset,&cfg,rangeName) ;
191  }
192  /// Create integral over observables in iset in range named rangeName using specified configuration for any numeric integration.
193  RooAbsReal* createIntegral(const RooArgSet& iset, const RooNumIntConfig& cfg, const char* rangeName=0) const {
194  return createIntegral(iset,0,&cfg,rangeName) ;
195  }
196  virtual RooAbsReal* createIntegral(const RooArgSet& iset, const RooArgSet* nset=0, const RooNumIntConfig* cfg=0, const char* rangeName=0) const ;
197 
198 
199  void setParameterizeIntegral(const RooArgSet& paramVars) ;
200 
201  // Create running integrals
202  RooAbsReal* createRunningIntegral(const RooArgSet& iset, const RooArgSet& nset=RooArgSet()) ;
203  RooAbsReal* createRunningIntegral(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2=RooCmdArg::none(),
204  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
205  const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
206  const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
207  RooAbsReal* createIntRI(const RooArgSet& iset, const RooArgSet& nset=RooArgSet()) ;
208  RooAbsReal* createScanRI(const RooArgSet& iset, const RooArgSet& nset, Int_t numScanBins, Int_t intOrder) ;
209 
210 
211  // Optimized accept/reject generator support
212  virtual Int_t getMaxVal(const RooArgSet& vars) const ;
213  virtual Double_t maxVal(Int_t code) const ;
214  virtual Int_t minTrialSamples(const RooArgSet& /*arGenObs*/) const { return 0 ; }
215 
216 
217  // Plotting options
218  void setPlotLabel(const char *label);
219  const char *getPlotLabel() const;
220 
221  virtual Double_t defaultErrorLevel() const {
222  // Return default level for MINUIT error analysis
223  return 1.0 ;
224  }
225 
226  const RooNumIntConfig* getIntegratorConfig() const ;
227  RooNumIntConfig* getIntegratorConfig() ;
228  static RooNumIntConfig* defaultIntegratorConfig() ;
229  RooNumIntConfig* specialIntegratorConfig() const ;
230  RooNumIntConfig* specialIntegratorConfig(Bool_t createOnTheFly) ;
231  void setIntegratorConfig() ;
232  void setIntegratorConfig(const RooNumIntConfig& config) ;
233 
234  virtual void fixAddCoefNormalization(const RooArgSet& addNormSet=RooArgSet(),Bool_t force=kTRUE) ;
235  virtual void fixAddCoefRange(const char* rangeName=0,Bool_t force=kTRUE) ;
236 
237  virtual void preferredObservableScanOrder(const RooArgSet& obs, RooArgSet& orderedObs) const ;
238 
239  // User entry point for plotting
240  virtual RooPlot* plotOn(RooPlot* frame,
241  const RooCmdArg& arg1=RooCmdArg(), const RooCmdArg& arg2=RooCmdArg(),
242  const RooCmdArg& arg3=RooCmdArg(), const RooCmdArg& arg4=RooCmdArg(),
243  const RooCmdArg& arg5=RooCmdArg(), const RooCmdArg& arg6=RooCmdArg(),
244  const RooCmdArg& arg7=RooCmdArg(), const RooCmdArg& arg8=RooCmdArg(),
245  const RooCmdArg& arg9=RooCmdArg(), const RooCmdArg& arg10=RooCmdArg()
246  ) const ;
247 
248 
249  enum ScaleType { Raw, Relative, NumEvent, RelativeExpected } ;
250 
251  // Forwarder function for backward compatibility
252  virtual RooPlot *plotSliceOn(RooPlot *frame, const RooArgSet& sliceSet, Option_t* drawOptions="L",
253  Double_t scaleFactor=1.0, ScaleType stype=Relative, const RooAbsData* projData=0) const;
254 
255  // Fill an existing histogram
256  TH1 *fillHistogram(TH1 *hist, const RooArgList &plotVars,
257  Double_t scaleFactor= 1, const RooArgSet *projectedVars= 0, Bool_t scaling=kTRUE,
258  const RooArgSet* condObs=0, Bool_t setError=kTRUE) const;
259 
260  // Create 1,2, and 3D histograms from and fill it
261  TH1 *createHistogram(const char* varNameList, Int_t xbins=0, Int_t ybins=0, Int_t zbins=0) const ;
262  TH1* createHistogram(const char *name, const RooAbsRealLValue& xvar, RooLinkedList& argList) const ;
263  TH1 *createHistogram(const char *name, const RooAbsRealLValue& xvar,
264  const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
265  const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
266  const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
267  const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) const ;
268 
269  // Fill a RooDataHist
270  RooDataHist* fillDataHist(RooDataHist *hist, const RooArgSet* nset, Double_t scaleFactor,
271  Bool_t correctForBinVolume=kFALSE, Bool_t showProgress=kFALSE) const ;
272 
273  // I/O streaming interface (machine readable)
274  virtual Bool_t readFromStream(std::istream& is, Bool_t compact, Bool_t verbose=kFALSE) ;
275  virtual void writeToStream(std::ostream& os, Bool_t compact) const ;
276 
277  // Printing interface (human readable)
278  virtual void printValue(std::ostream& os) const ;
279  virtual void printMultiline(std::ostream& os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const ;
280 
281  static void setCacheCheck(Bool_t flag) ;
282 
283  // Evaluation error logging
284  class EvalError {
285  public:
286  EvalError() { }
287  EvalError(const EvalError& other) : _msg(other._msg), _srvval(other._srvval) { }
288  void setMessage(const char* tmp) { std::string s(tmp); s.swap(_msg); }
289  void setServerValues(const char* tmp) { std::string s(tmp); s.swap(_srvval); }
290  std::string _msg;
291  std::string _srvval;
292  } ;
293 
294  enum ErrorLoggingMode { PrintErrors, CollectErrors, CountErrors, Ignore } ;
295  static ErrorLoggingMode evalErrorLoggingMode() ;
296  static void setEvalErrorLoggingMode(ErrorLoggingMode m) ;
297  void logEvalError(const char* message, const char* serverValueString=0) const ;
298  static void logEvalError(const RooAbsReal* originator, const char* origName, const char* message, const char* serverValueString=0) ;
299  static void printEvalErrors(std::ostream&os=std::cout, Int_t maxPerNode=10000000) ;
300  static Int_t numEvalErrors() ;
301  static Int_t numEvalErrorItems() ;
302 
303 
304  typedef std::map<const RooAbsArg*,std::pair<std::string,std::list<EvalError> > >::const_iterator EvalErrorIter ;
305  static EvalErrorIter evalErrorIter() ;
306 
307  static void clearEvalErrorLog() ;
308 
309  virtual Bool_t isBinnedDistribution(const RooArgSet& /*obs*/) const { return kFALSE ; }
310  virtual std::list<Double_t>* binBoundaries(RooAbsRealLValue& /*obs*/, Double_t /*xlo*/, Double_t /*xhi*/) const { return 0 ; }
311  virtual std::list<Double_t>* plotSamplingHint(RooAbsRealLValue& /*obs*/, Double_t /*xlo*/, Double_t /*xhi*/) const {
312  // Interface for returning an optional hint for initial sampling points when constructing a curve
313  // projected on observable.
314  return 0 ;
315  }
316 
317  RooGenFunction* iGenFunction(RooRealVar& x, const RooArgSet& nset=RooArgSet()) ;
318  RooMultiGenFunction* iGenFunction(const RooArgSet& observables, const RooArgSet& nset=RooArgSet()) ;
319 
320  RooFunctor* functor(const RooArgList& obs, const RooArgList& pars=RooArgList(), const RooArgSet& nset=RooArgSet()) const ;
321  TF1* asTF(const RooArgList& obs, const RooArgList& pars=RooArgList(), const RooArgSet& nset=RooArgSet()) const ;
322 
323  RooDerivative* derivative(RooRealVar& obs, Int_t order=1, Double_t eps=0.001) ;
324  RooDerivative* derivative(RooRealVar& obs, const RooArgSet& normSet, Int_t order, Double_t eps=0.001) ;
325 
326  RooAbsMoment* moment(RooRealVar& obs, Int_t order, Bool_t central, Bool_t takeRoot) ;
327  RooAbsMoment* moment(RooRealVar& obs, const RooArgSet& normObs, Int_t order, Bool_t central, Bool_t takeRoot, Bool_t intNormObs) ;
328 
329  RooAbsMoment* mean(RooRealVar& obs) { return moment(obs,1,kFALSE,kFALSE) ; }
330  RooAbsMoment* mean(RooRealVar& obs, const RooArgSet& nset) { return moment(obs,nset,1,kFALSE,kFALSE,kTRUE) ; }
331  RooAbsMoment* sigma(RooRealVar& obs) { return moment(obs,2,kTRUE,kTRUE) ; }
332  RooAbsMoment* sigma(RooRealVar& obs, const RooArgSet& nset) { return moment(obs,nset,2,kTRUE,kTRUE,kTRUE) ; }
333 
334  Double_t findRoot(RooRealVar& x, Double_t xmin, Double_t xmax, Double_t yval) ;
335 
336 
337  virtual Bool_t setData(RooAbsData& /*data*/, Bool_t /*cloneData*/=kTRUE) { return kTRUE ; }
338 
339  virtual void enableOffsetting(Bool_t) {} ;
340  virtual Bool_t isOffsetting() const { return kFALSE ; }
341  virtual Double_t offset() const { return 0 ; }
342 
343  static void setHideOffset(Bool_t flag);
344  static Bool_t hideOffset() ;
345 
346 protected:
347  // Hook for objects with normalization-dependent parameters interperetation
348  virtual void selectNormalization(const RooArgSet* depSet=0, Bool_t force=kFALSE) ;
349  virtual void selectNormalizationRange(const char* rangeName=0, Bool_t force=kFALSE) ;
350 
351  // Helper functions for plotting
352  Bool_t plotSanityChecks(RooPlot* frame) const ;
353  void makeProjectionSet(const RooAbsArg* plotVar, const RooArgSet* allVars,
354  RooArgSet& projectedVars, Bool_t silent) const ;
355 
356  TString integralNameSuffix(const RooArgSet& iset, const RooArgSet* nset=0, const char* rangeName=0, Bool_t omitEmpty=kFALSE) const ;
357 
358 
359  Bool_t isSelectedComp() const ;
360 
361 
362  public:
363  const RooAbsReal* createPlotProjection(const RooArgSet& depVars, const RooArgSet& projVars) const ;
364  const RooAbsReal* createPlotProjection(const RooArgSet& depVars, const RooArgSet& projVars, RooArgSet*& cloneSet) const ;
365  const RooAbsReal *createPlotProjection(const RooArgSet &dependentVars, const RooArgSet *projectedVars,
366  RooArgSet *&cloneSet, const char* rangeName=0, const RooArgSet* condObs=0) const;
367  protected:
368 
369  RooFitResult* chi2FitDriver(RooAbsReal& fcn, RooLinkedList& cmdList) ;
370 
371  void plotOnCompSelect(RooArgSet* selNodes) const ;
372  RooPlot* plotOnWithErrorBand(RooPlot* frame,const RooFitResult& fr, Double_t Z, const RooArgSet* params, const RooLinkedList& argList, Bool_t method1) const ;
373 
374  // Support interface for subclasses to advertise their analytic integration
375  // and generator capabilities in their analticalIntegral() and generateEvent()
376  // implementations.
377  Bool_t matchArgs(const RooArgSet& allDeps, RooArgSet& numDeps,
378  const RooArgProxy& a) const ;
379  Bool_t matchArgs(const RooArgSet& allDeps, RooArgSet& numDeps,
380  const RooArgProxy& a, const RooArgProxy& b) const ;
381  Bool_t matchArgs(const RooArgSet& allDeps, RooArgSet& numDeps,
382  const RooArgProxy& a, const RooArgProxy& b, const RooArgProxy& c) const ;
383  Bool_t matchArgs(const RooArgSet& allDeps, RooArgSet& numDeps,
384  const RooArgProxy& a, const RooArgProxy& b,
385  const RooArgProxy& c, const RooArgProxy& d) const ;
386 
387  Bool_t matchArgs(const RooArgSet& allDeps, RooArgSet& numDeps,
388  const RooArgSet& set) const ;
389 
390 
391  RooAbsReal* createIntObj(const RooArgSet& iset, const RooArgSet* nset, const RooNumIntConfig* cfg, const char* rangeName) const ;
392  void findInnerMostIntegration(const RooArgSet& allObs, RooArgSet& innerObs, const char* rangeName) const ;
393 
394 
395  // Internal consistency checking (needed by RooDataSet)
396  virtual Bool_t isValid() const ;
397  virtual Bool_t isValidReal(Double_t value, Bool_t printError=kFALSE) const ;
398 
399  // Function evaluation and error tracing
400  Double_t traceEval(const RooArgSet* set) const ;
401 
402  /// Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
403  virtual Double_t evaluate() const = 0;
404  virtual RooSpan<double> evaluateBatch(std::size_t begin, std::size_t maxSize) const;
405 
406  //---------- Interface to access batch data ---------------------------
407  //
408  friend class RooHelpers::BatchInterfaceAccessor;
409  void clearBatchMemory() {
410  _batchData.clear();
411  for (auto arg : _serverList) {
412  //TODO get rid of this cast?
413  auto absReal = dynamic_cast<RooAbsReal*>(arg);
414  if (absReal)
415  absReal->clearBatchMemory();
416  }
417  }
418 
419  private:
420  void checkBatchComputation(std::size_t evtNo, const RooArgSet* normSet = nullptr, double relAccuracy = 1.E-13) const;
421 
422  const BatchHelpers::BatchData& batchData() const {
423  return _batchData;
424  }
425 
426  /// Debug version of getVal(), which is slow and does error checking.
427  Double_t _DEBUG_getVal(const RooArgSet* normalisationSet) const;
428 
429  //--------------------------------------------------------------------
430 
431  protected:
432  // Hooks for RooDataSet interface
433  friend class RooRealIntegral ;
434  friend class RooVectorDataStore ;
435  virtual void syncCache(const RooArgSet* set=0) { getVal(set) ; }
436  virtual void copyCache(const RooAbsArg* source, Bool_t valueOnly=kFALSE, Bool_t setValDirty=kTRUE) ;
437  virtual void attachToTree(TTree& t, Int_t bufSize=32000) ;
438  virtual void attachToVStore(RooVectorDataStore& vstore) ;
439  virtual void setTreeBranchStatus(TTree& t, Bool_t active) ;
440  virtual void fillTreeBranch(TTree& t) ;
441 
442  friend class RooRealBinding ;
443  Double_t _plotMin ; // Minimum of plot range
444  Double_t _plotMax ; // Maximum of plot range
445  Int_t _plotBins ; // Number of plot bins
446  mutable Double_t _value ; // Cache for current value of object
447  mutable BatchHelpers::BatchData _batchData; //! Value storage for batches of events
448  TString _unit ; // Unit for objects value
449  TString _label ; // Plot label for objects value
450  Bool_t _forceNumInt ; // Force numerical integration if flag set
451 
452  mutable Float_t _floatValue{0.}; //! Transient cache for floating point values from tree branches
453  mutable Int_t _intValue{0}; //! Transient cache for integer values from tree branches
454  mutable Bool_t _boolValue{false}; //! Transient cache for bool values from tree branches
455  mutable UChar_t _byteValue{'\0'}; //! Transient cache for byte values from tree branches
456  mutable Char_t _sbyteValue{'\0'}; //! Transient cache for signed byte values from tree branches
457  mutable UInt_t _uintValue{0u}; //! Transient cache for unsigned integer values from tree branches
458 
459  friend class RooAbsPdf ;
460  friend class RooAbsAnaConvPdf ;
461 
462  RooNumIntConfig* _specIntegratorConfig ; // Numeric integrator configuration specific for this object
463 
464  Bool_t _treeVar ; // !do not persist
465 
466  friend class RooDataProjBinding ;
467  friend class RooAbsOptGoodnessOfFit ;
468 
469  struct PlotOpt {
470  PlotOpt() : drawOptions("L"), scaleFactor(1.0), stype(Relative), projData(0), binProjData(kFALSE), projSet(0), precision(1e-3),
471  shiftToZero(kFALSE),projDataSet(0),normRangeName(0),rangeLo(0),rangeHi(0),postRangeFracScale(kFALSE),wmode(RooCurve::Extended),
472  projectionRangeName(0),curveInvisible(kFALSE), curveName(0),addToCurveName(0),addToWgtSelf(1.),addToWgtOther(1.),
473  numCPU(1),interleave(RooFit::Interleave),curveNameSuffix(""), numee(10), eeval(0), doeeval(kFALSE), progress(kFALSE), errorFR(0) {} ;
474  Option_t* drawOptions ;
475  Double_t scaleFactor ;
476  ScaleType stype ;
477  const RooAbsData* projData ;
478  Bool_t binProjData ;
479  const RooArgSet* projSet ;
480  Double_t precision ;
481  Bool_t shiftToZero ;
482  const RooArgSet* projDataSet ;
483  const char* normRangeName ;
484  Double_t rangeLo ;
485  Double_t rangeHi ;
486  Bool_t postRangeFracScale ;
487  RooCurve::WingMode wmode ;
488  const char* projectionRangeName ;
489  Bool_t curveInvisible ;
490  const char* curveName ;
491  const char* addToCurveName ;
492  Double_t addToWgtSelf ;
493  Double_t addToWgtOther ;
494  Int_t numCPU ;
495  RooFit::MPSplit interleave ;
496  const char* curveNameSuffix ;
497  Int_t numee ;
498  Double_t eeval ;
499  Bool_t doeeval ;
500  Bool_t progress ;
501  const RooFitResult* errorFR ;
502  } ;
503 
504  // Plot implementation functions
505  virtual RooPlot *plotOn(RooPlot* frame, PlotOpt o) const;
506 
507 public:
508  // PlotOn with command list
509  virtual RooPlot* plotOn(RooPlot* frame, RooLinkedList& cmdList) const ;
510 
511  protected:
512  virtual RooPlot *plotAsymOn(RooPlot *frame, const RooAbsCategoryLValue& asymCat, PlotOpt o) const;
513 
514 
515 private:
516 
517  static ErrorLoggingMode _evalErrorMode ;
518  static std::map<const RooAbsArg*,std::pair<std::string,std::list<EvalError> > > _evalErrorList ;
519  static Int_t _evalErrorCount ;
520 
521  Bool_t matchArgsByName(const RooArgSet &allArgs, RooArgSet &matchedArgs, const TList &nameList) const;
522 
523 protected:
524 
525 
526  friend class RooRealSumPdf ;
527  friend class RooRealSumFunc;
528  friend class RooAddPdf ;
529  friend class RooAddModel ;
530  void selectComp(Bool_t flag) {
531  // If flag is true, only selected component will be included in evaluates of RooAddPdf components
532  _selectComp = flag ;
533  }
534  static void globalSelectComp(Bool_t flag) ;
535  Bool_t _selectComp ; //! Component selection flag for RooAbsPdf::plotCompOn
536  static Bool_t _globalSelectComp ; // Global activation switch for component selection
537  // This struct can be used to flip the global switch to select components.
538  // Doing this with RAII prevents forgetting to reset the state.
539  struct GlobalSelectComponentRAII {
540  GlobalSelectComponentRAII(bool state) :
541  _oldState{_globalSelectComp} {
542  if (state != RooAbsReal::_globalSelectComp)
543  RooAbsReal::_globalSelectComp = state;
544  }
545 
546  ~GlobalSelectComponentRAII() {
547  if (RooAbsReal::_globalSelectComp != _oldState)
548  RooAbsReal::_globalSelectComp = _oldState;
549  }
550 
551  bool _oldState;
552  };
553 
554 
555  mutable RooArgSet* _lastNSet ; //!
556  static Bool_t _hideOffset ; // Offset hiding flag
557 
558  ClassDef(RooAbsReal,2) // Abstract real-valued variable
559 };
560 
561 #endif