Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooHistPdf.h
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * File: $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#ifndef ROO_HIST_PDF
17#define ROO_HIST_PDF
18
19#include "RooAbsPdf.h"
20#include "RooRealProxy.h"
21#include "RooSetProxy.h"
22#include "RooAICRegistry.h"
23#include "RooDataHist.h"
24
25#include <list>
26
27class RooRealVar;
28class RooAbsReal;
29
30class RooHistPdf : public RooAbsPdf {
31public:
33 RooHistPdf(const char *name, const char *title, const RooArgSet& vars, const RooDataHist& dhist, Int_t intOrder=0);
34 RooHistPdf(const char *name, const char *title, const RooArgList& pdfObs, const RooArgList& histObs, const RooDataHist& dhist, Int_t intOrder=0);
35 RooHistPdf(const char *name, const char *title, const RooArgSet& vars,
36 std::unique_ptr<RooDataHist> dhist, int intOrder=0);
37 RooHistPdf(const char *name, const char *title, const RooArgList& pdfObs, const RooArgList& histObs,
38 std::unique_ptr<RooDataHist> dhist, int intOrder=0);
39 RooHistPdf(const RooHistPdf& other, const char* name=nullptr);
40 TObject* clone(const char* newname) const override { return new RooHistPdf(*this,newname); }
41
43 // Return RooDataHist that is represented
44 return *_dataHist ;
45 }
46 const RooDataHist& dataHist() const {
47 // Return RooDataHist that is represented
48 return *_dataHist ;
49 }
50
51 /// Replaces underlying RooDataHist with a clone, which is now owned, and returns the clone.
52 /// If the underlying RooDataHist is already owned, then that is returned instead of being cloned.
53 RooDataHist* cloneAndOwnDataHist(const char* newname="");
54
56 // Set histogram interpolation order
57 _intOrder = order ;
58 }
60 // Return histogram interpolation order
61 return _intOrder ;
62 }
63
64 Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=nullptr) const override ;
65 double analyticalIntegral(Int_t code, const char* rangeName=nullptr) const override ;
66
67 bool forceAnalyticalInt(const RooAbsArg& dep) const override;
68
69 void setCdfBoundaries(bool flag) {
70 // Set use of special boundary conditions for c.d.f.s
71 _cdfBoundaries = flag ;
72 }
73 bool getCdfBoundaries() const {
74 // If true, special boundary conditions for c.d.f.s are used
75 return _cdfBoundaries ;
76 }
77
78 void setUnitNorm(bool flag) {
79 // Declare contents to have unit normalization
80 _unitNorm = flag ;
81 }
82 bool haveUnitNorm() const {
83 // Return true if contents is declared to be unit normalized
84 return _unitNorm ;
85 }
86
87 bool selfNormalized() const override { return _unitNorm ; }
88
89 Int_t getMaxVal(const RooArgSet& vars) const override ;
90 double maxVal(Int_t code) const override ;
91
92 std::list<double>* plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const override ;
93 std::list<double>* binBoundaries(RooAbsRealLValue& /*obs*/, double /*xlo*/, double /*xhi*/) const override ;
94 bool isBinnedDistribution(const RooArgSet&) const override { return _intOrder==0 ; }
95
96 void doEval(RooFit::EvalContext &) const override;
97
98 RooArgSet const &variables() const { return _pdfObsList; }
99
100protected:
101 bool areIdentical(const RooDataHist& dh1, const RooDataHist& dh2) ;
102
103 bool importWorkspaceHook(RooWorkspace& ws) override ;
104
105 double evaluate() const override;
106 double totalVolume() const ;
107 friend class RooAbsCachedPdf ;
108 double totVolume() const ;
109
110 RooArgSet _histObsList; ///< List of observables defining dimensions of histogram
111 RooSetProxy _pdfObsList; ///< List of observables mapped onto histogram observables
112 RooDataHist* _dataHist = nullptr; ///< Unowned pointer to underlying histogram
113 std::unique_ptr<RooDataHist> _ownedDataHist; ///<! Owned pointer to underlying histogram
114 mutable RooAICRegistry _codeReg ; ///<! Auxiliary class keeping tracking of analytical integration code
115 Int_t _intOrder = 0; ///< Interpolation order
116 bool _cdfBoundaries = false; ///< Use boundary conditions for CDFs.
117 mutable double _totVolume = 0.0; ///<! Total volume of space (product of ranges of observables)
118 bool _unitNorm = false; ///< Assume contents is unit normalized (for use as pdf cache)
119
120private:
121
122 friend class RooHistFunc;
123
124 static bool forceAnalyticalInt(RooArgSet const& pdfObsList, RooAbsArg const& dep);
125
126 static Int_t getAnalyticalIntegral(RooArgSet& allVars,
127 RooArgSet& analVars,
128 const char* rangeName,
129 RooArgSet const& histObsList,
130 RooArgSet const& pdfObsList,
131 Int_t intOrder) ;
132
133 static double analyticalIntegral(Int_t code,
134 const char* rangeName,
135 RooArgSet const& histObsList,
136 RooArgSet const& pdfObsList,
138 bool histFuncMode) ;
139
140 static std::list<double>* plotSamplingHint(RooDataHist const& dataHist,
141 RooArgSet const& pdfObsList,
142 RooArgSet const& histObsList,
143 int intOrder,
144 RooAbsRealLValue& obs,
145 double xlo,
146 double xhi);
147
148 inline void initializeOwnedDataHist(std::unique_ptr<RooDataHist> &&dataHist)
149 {
150 _ownedDataHist = std::move(dataHist);
151 }
152
153 ClassDefOverride(RooHistPdf,4) // Histogram based PDF
154};
155
156#endif
#define ClassDefOverride(name, id)
Definition Rtypes.h:346
char name[80]
Definition TGX11.cxx:110
Utility class for operator p.d.f classes that keeps track of analytical integration codes and associa...
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:77
Abstract base class for p.d.f.s that need or want to cache their evaluate() output in a RooHistPdf de...
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 ...
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
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
Container class to hold N-dimensional binned data.
Definition RooDataHist.h:40
A real-valued function sampled from a multidimensional histogram.
Definition RooHistFunc.h:31
A propability density function sampled from a multidimensional histogram.
Definition RooHistPdf.h:30
Int_t getInterpolationOrder() const
Definition RooHistPdf.h:59
void setUnitNorm(bool flag)
Definition RooHistPdf.h:78
RooArgSet _histObsList
List of observables defining dimensions of histogram.
Definition RooHistPdf.h:110
Int_t _intOrder
Interpolation order.
Definition RooHistPdf.h:115
bool forceAnalyticalInt(const RooAbsArg &dep) const override
bool areIdentical(const RooDataHist &dh1, const RooDataHist &dh2)
RooDataHist * _dataHist
Unowned pointer to underlying histogram.
Definition RooHistPdf.h:112
bool _cdfBoundaries
Use boundary conditions for CDFs.
Definition RooHistPdf.h:116
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Return sampling hint for making curves of (projections) of this function as the recursive division st...
double totVolume() const
Return the total volume spanned by the observables of the RooHistPdf.
void initializeOwnedDataHist(std::unique_ptr< RooDataHist > &&dataHist)
Definition RooHistPdf.h:148
bool importWorkspaceHook(RooWorkspace &ws) override
Check if our datahist is already in the workspace.
RooSetProxy _pdfObsList
List of observables mapped onto histogram observables.
Definition RooHistPdf.h:111
bool isBinnedDistribution(const RooArgSet &) const override
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition RooHistPdf.h:94
double maxVal(Int_t code) const override
Return maximum value for set of observables identified by code assigned in getMaxVal.
double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const override
Return integral identified by 'code'.
RooAICRegistry _codeReg
! Auxiliary class keeping tracking of analytical integration code
Definition RooHistPdf.h:114
TObject * clone(const char *newname) const override
Definition RooHistPdf.h:40
void setCdfBoundaries(bool flag)
Definition RooHistPdf.h:69
std::list< double > * plotSamplingHint(RooAbsRealLValue &obs, double xlo, double xhi) const override
Return sampling hint for making curves of (projections) of this function as the recursive division st...
RooDataHist & dataHist()
Definition RooHistPdf.h:42
bool haveUnitNorm() const
Definition RooHistPdf.h:82
Int_t getMaxVal(const RooArgSet &vars) const override
Only handle case of maximum in all variables.
double _totVolume
! Total volume of space (product of ranges of observables)
Definition RooHistPdf.h:117
RooDataHist * cloneAndOwnDataHist(const char *newname="")
Replaces underlying RooDataHist with a clone, which is now owned, and returns the clone.
std::unique_ptr< RooDataHist > _ownedDataHist
! Owned pointer to underlying histogram
Definition RooHistPdf.h:113
double totalVolume() const
void doEval(RooFit::EvalContext &) const override
Base function for computing multiple values of a RooAbsReal.
bool getCdfBoundaries() const
Definition RooHistPdf.h:73
RooArgSet const & variables() const
Definition RooHistPdf.h:98
bool selfNormalized() const override
Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
Definition RooHistPdf.h:87
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Determine integration scenario.
const RooDataHist & dataHist() const
Definition RooHistPdf.h:46
double evaluate() const override
Return the current value: The value of the bin enclosing the current coordinates of the observables,...
void setInterpolationOrder(Int_t order)
Definition RooHistPdf.h:55
bool _unitNorm
Assume contents is unit normalized (for use as pdf cache)
Definition RooHistPdf.h:118
Variable that can be changed from the outside.
Definition RooRealVar.h:37
Persistable container for RooFit projects.
Mother of all ROOT objects.
Definition TObject.h:41