Logo ROOT  
Reference Guide
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 ~RooHistPdf() override ;
42
44 // Return RooDataHist that is represented
45 return *_dataHist ;
46 }
47 const RooDataHist& dataHist() const {
48 // Return RooDataHist that is represented
49 return *_dataHist ;
50 }
51
53 // Set histogram interpolation order
54 _intOrder = order ;
55 }
57 // Return histogram interpolation order
58 return _intOrder ;
59 }
60
62 RooArgSet& analVars,
63 const char* rangeName,
64 RooArgSet const& histObsList,
65 RooSetProxy const& pdfObsList,
66 Int_t intOrder) ;
67
68 static double analyticalIntegral(Int_t code,
69 const char* rangeName,
70 RooArgSet const& histObsList,
71 RooSetProxy const& pdfObsList,
73 bool histFuncMode) ;
74
75 Int_t getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName=nullptr) const override ;
76 double analyticalIntegral(Int_t code, const char* rangeName=nullptr) const override ;
77
78 void setCdfBoundaries(bool flag) {
79 // Set use of special boundary conditions for c.d.f.s
80 _cdfBoundaries = flag ;
81 }
82 bool getCdfBoundaries() const {
83 // If true, special boundary conditions for c.d.f.s are used
84 return _cdfBoundaries ;
85 }
86
87 void setUnitNorm(bool flag) {
88 // Declare contents to have unit normalization
89 _unitNorm = flag ;
90 }
91 bool haveUnitNorm() const {
92 // Return true if contents is declared to be unit normalized
93 return _unitNorm ;
94 }
95
96 bool selfNormalized() const override { return _unitNorm ; }
97
98 Int_t getMaxVal(const RooArgSet& vars) const override ;
99 double maxVal(Int_t code) const override ;
100
101 std::list<double>* plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const override ;
102 std::list<double>* binBoundaries(RooAbsRealLValue& /*obs*/, double /*xlo*/, double /*xhi*/) const override ;
103 bool isBinnedDistribution(const RooArgSet&) const override { return _intOrder==0 ; }
104
105 void computeBatch(cudaStream_t*, double* output, size_t size, RooFit::Detail::DataMap const&) const override;
106
107protected:
108
109 bool areIdentical(const RooDataHist& dh1, const RooDataHist& dh2) ;
110
111 bool importWorkspaceHook(RooWorkspace& ws) override ;
112
113 double evaluate() const override;
114 double totalVolume() const ;
115 friend class RooAbsCachedPdf ;
116 double totVolume() const ;
117
118 RooArgSet _histObsList; ///< List of observables defining dimensions of histogram
119 RooSetProxy _pdfObsList; ///< List of observables mapped onto histogram observables
120 RooDataHist* _dataHist = nullptr; ///< Unowned pointer to underlying histogram
121 std::unique_ptr<RooDataHist> _ownedDataHist; ///<! Owned pointer to underlying histogram
122 mutable RooAICRegistry _codeReg ; ///<! Auxiliary class keeping tracking of analytical integration code
123 Int_t _intOrder = 0; ///< Interpolation order
124 bool _cdfBoundaries = false; ///< Use boundary conditions for CDFs.
125 mutable double _totVolume = 0.0; ///<! Total volume of space (product of ranges of observables)
126 bool _unitNorm = false; ///< Assume contents is unit normalized (for use as pdf cache)
127
128 ClassDefOverride(RooHistPdf,4) // Histogram based PDF
129};
130
131#endif
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
#define ClassDefOverride(name, id)
Definition: Rtypes.h:339
char name[80]
Definition: TGX11.cxx:110
RooAICRegistry is a utility class for operator p.d.f classes that keeps track of analytical integrati...
RooAbsCachedPdf is the abstract base class for p.d.f.s that need or want to cache their evaluate() ou...
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:62
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:56
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:39
RooHistPdf implements a probablity density function sampled from a multidimensional histogram.
Definition: RooHistPdf.h:30
Int_t getInterpolationOrder() const
Definition: RooHistPdf.h:56
void setUnitNorm(bool flag)
Definition: RooHistPdf.h:87
RooArgSet _histObsList
List of observables defining dimensions of histogram.
Definition: RooHistPdf.h:118
Int_t _intOrder
Interpolation order.
Definition: RooHistPdf.h:123
bool areIdentical(const RooDataHist &dh1, const RooDataHist &dh2)
Definition: RooHistPdf.cxx:510
RooDataHist * _dataHist
Unowned pointer to underlying histogram.
Definition: RooHistPdf.h:120
bool _cdfBoundaries
Use boundary conditions for CDFs.
Definition: RooHistPdf.h:124
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...
Definition: RooHistPdf.cxx:443
double totVolume() const
Return the total volume spanned by the observables of the RooHistPdf.
Definition: RooHistPdf.cxx:234
bool importWorkspaceHook(RooWorkspace &ws) override
Check if our datahist is already in the workspace.
Definition: RooHistPdf.cxx:527
RooSetProxy _pdfObsList
List of observables mapped onto histogram observables.
Definition: RooHistPdf.h:119
static double analyticalIntegral(Int_t code, const char *rangeName, RooArgSet const &histObsList, RooSetProxy const &pdfObsList, RooDataHist &dataHist, bool histFuncMode)
Definition: RooHistPdf.cxx:317
bool isBinnedDistribution(const RooArgSet &) const override
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition: RooHistPdf.h:103
double maxVal(Int_t code) const override
Return maximum value for set of observables identified by code assigned in getMaxVal.
Definition: RooHistPdf.cxx:491
RooAICRegistry _codeReg
! Auxiliary class keeping tracking of analytical integration code
Definition: RooHistPdf.h:122
TObject * clone(const char *newname) const override
Definition: RooHistPdf.h:40
void setCdfBoundaries(bool flag)
Definition: RooHistPdf.h:78
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...
Definition: RooHistPdf.cxx:385
void computeBatch(cudaStream_t *, double *output, size_t size, RooFit::Detail::DataMap const &) const override
Base function for computing multiple values of a RooAbsReal.
Definition: RooHistPdf.cxx:191
RooDataHist & dataHist()
Definition: RooHistPdf.h:43
bool haveUnitNorm() const
Definition: RooHistPdf.h:91
Int_t getMaxVal(const RooArgSet &vars) const override
Only handle case of maximum in all variables.
Definition: RooHistPdf.cxx:479
double _totVolume
! Total volume of space (product of ranges of observables)
Definition: RooHistPdf.h:125
std::unique_ptr< RooDataHist > _ownedDataHist
! Owned pointer to underlying histogram
Definition: RooHistPdf.h:121
double totalVolume() const
static Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName, RooArgSet const &histObsList, RooSetProxy const &pdfObsList, Int_t intOrder)
Definition: RooHistPdf.cxx:276
bool getCdfBoundaries() const
Definition: RooHistPdf.h:82
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:96
const RooDataHist & dataHist() const
Definition: RooHistPdf.h:47
~RooHistPdf() override
Destructor.
Definition: RooHistPdf.cxx:185
double evaluate() const override
Return the current value: The value of the bin enclosing the current coordinates of the observables,...
Definition: RooHistPdf.cxx:209
void setInterpolationOrder(Int_t order)
Definition: RooHistPdf.h:52
bool _unitNorm
Assume contents is unit normalized (for use as pdf cache)
Definition: RooHistPdf.h:126
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:43
Mother of all ROOT objects.
Definition: TObject.h:41
static void output()