Logo ROOT  
Reference Guide
RooPolynomial.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/** \class RooPolynomial
18 \ingroup Roofit
19
20RooPolynomial implements a polynomial p.d.f of the form
21\f[ f(x) = \mathcal{N} \cdot \sum_{i} a_{i} * x^i \f]
22By default, the coefficient \f$ a_0 \f$ is chosen to be 1, as polynomial
23probability density functions have one degree of freedom
24less than polynomial functions due to the normalisation condition. \f$ \mathcal{N} \f$
25is a normalisation constant that is automatically calculated when the polynomial is used
26in computations.
27
28The sum can be truncated at the low end. See the main constructor
29RooPolynomial::RooPolynomial(const char*, const char*, RooAbsReal&, const RooArgList&, Int_t)
30**/
31
32#include "RooPolynomial.h"
33#include "RooArgList.h"
34#include "RooMsgService.h"
35#include "RooBatchCompute.h"
36
37#include "TError.h"
38#include <vector>
39using namespace std;
40
42
43////////////////////////////////////////////////////////////////////////////////
44/// coverity[UNINIT_CTOR]
45
47{
48}
49
50////////////////////////////////////////////////////////////////////////////////
51/// Create a polynomial in the variable `x`.
52/// \param[in] name Name of the PDF
53/// \param[in] title Title for plotting the PDF
54/// \param[in] x The variable of the polynomial
55/// \param[in] coefList The coefficients \f$ a_i \f$
56/// \param[in] lowestOrder [optional] Truncate the sum such that it skips the lower orders:
57/// \f[
58/// 1. + \sum_{i=0}^{\mathrm{coefList.size()}} a_{i} * x^{(i + \mathrm{lowestOrder})}
59/// \f]
60///
61/// This means that
62/// \code{.cpp}
63/// RooPolynomial pol("pol", "pol", x, RooArgList(a, b), lowestOrder = 2)
64/// \endcode
65/// computes
66/// \f[
67/// \mathrm{pol}(x) = 1 * x^0 + (0 * x^{\ldots}) + a * x^2 + b * x^3.
68/// \f]
69
70
71RooPolynomial::RooPolynomial(const char* name, const char* title,
72 RooAbsReal& x, const RooArgList& coefList, Int_t lowestOrder) :
73 RooAbsPdf(name, title),
74 _x("x", "Dependent", this, x),
75 _coefList("coefList","List of coefficients",this),
76 _lowestOrder(lowestOrder)
77{
78 // Check lowest order
79 if (_lowestOrder<0) {
80 coutE(InputArguments) << "RooPolynomial::ctor(" << GetName()
81 << ") WARNING: lowestOrder must be >=0, setting value to 0" << endl ;
82 _lowestOrder=0 ;
83 }
84
85 for (auto *coef : coefList) {
86 if (!dynamic_cast<RooAbsReal*>(coef)) {
87 coutE(InputArguments) << "RooPolynomial::ctor(" << GetName() << ") ERROR: coefficient " << coef->GetName()
88 << " is not of type RooAbsReal" << endl ;
89 R__ASSERT(0) ;
90 }
91 _coefList.add(*coef) ;
92 }
93}
94
95////////////////////////////////////////////////////////////////////////////////
96
97RooPolynomial::RooPolynomial(const char* name, const char* title,
98 RooAbsReal& x) :
99 RooAbsPdf(name, title),
100 _x("x", "Dependent", this, x),
101 _coefList("coefList","List of coefficients",this),
102 _lowestOrder(1)
103{ }
104
105////////////////////////////////////////////////////////////////////////////////
106/// Copy constructor
107
109 RooAbsPdf(other, name),
110 _x("x", this, other._x),
111 _coefList("coefList",this,other._coefList),
112 _lowestOrder(other._lowestOrder)
113{ }
114
115////////////////////////////////////////////////////////////////////////////////
116/// Destructor
117
119{ }
120
121////////////////////////////////////////////////////////////////////////////////
122
124{
125 // Calculate and return value of polynomial
126
127 const unsigned sz = _coefList.getSize();
128 const int lowestOrder = _lowestOrder;
129 if (!sz) return lowestOrder ? 1. : 0.;
130 _wksp.clear();
131 _wksp.reserve(sz);
132 {
133 const RooArgSet* nset = _coefList.nset();
134 for (auto *c : static_range_cast<RooAbsReal *>(_coefList)) {
135 _wksp.push_back(c->getVal(nset));
136 }
137 }
138 const double x = _x;
139 double retVal = _wksp[sz - 1];
140 for (unsigned i = sz - 1; i--; ) retVal = _wksp[i] + x * retVal;
141 return retVal * std::pow(x, lowestOrder) + (lowestOrder ? 1.0 : 0.0);
142}
143
144// The batch mode support for RooPolynomial was commented out, because that
145// implementation can't deal with observables used as polynomial coefficients
146// yet.
147
148//////////////////////////////////////////////////////////////////////////////////
149///// Compute multiple values of Polynomial.
150//void RooPolynomial::computeBatch(cudaStream_t* stream, double* output, size_t nEvents, RooFit::DataMap& dataMap) const
151//{
152 //RooBatchCompute::ArgVector extraArgs;
153 //for (auto* coef:_coefList)
154 //extraArgs.push_back( static_cast<const RooAbsReal*>(coef)->getVal() );
155 //extraArgs.push_back(_lowestOrder);
156 //auto dispatch = stream ? RooBatchCompute::dispatchCUDA : RooBatchCompute::dispatchCPU;
157 //dispatch->compute(stream, RooBatchCompute::Polynomial, output, nEvents, dataMap, {&*_x,&*_norm}, extraArgs);
158//}
159
160////////////////////////////////////////////////////////////////////////////////
161/// Advertise to RooFit that this function can be analytically integrated.
162Int_t RooPolynomial::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
163{
164 if (matchArgs(allVars, analVars, _x)) return 1;
165 return 0;
166}
167
168////////////////////////////////////////////////////////////////////////////////
169/// Do the analytical integral according to the code that was returned by getAnalyticalIntegral().
170double RooPolynomial::analyticalIntegral(Int_t code, const char* rangeName) const
171{
172 R__ASSERT(code==1) ;
173
174 const double xmin = _x.min(rangeName), xmax = _x.max(rangeName);
175 const int lowestOrder = _lowestOrder;
176 const unsigned sz = _coefList.getSize();
177 if (!sz) return xmax - xmin;
178 _wksp.clear();
179 _wksp.reserve(sz);
180 {
181 const RooArgSet* nset = _coefList.nset();
182 unsigned i = 1 + lowestOrder;
183 for (auto *c : static_range_cast<RooAbsReal *>(_coefList)) {
184 _wksp.push_back(c->getVal(nset) / double(i));
185 ++i;
186 }
187 }
188 double min = _wksp[sz - 1], max = _wksp[sz - 1];
189 for (unsigned i = sz - 1; i--; )
190 min = _wksp[i] + xmin * min, max = _wksp[i] + xmax * max;
191 return max * std::pow(xmax, 1 + lowestOrder) - min * std::pow(xmin, 1 + lowestOrder) +
192 (lowestOrder ? (xmax - xmin) : 0.);
193}
#define c(i)
Definition: RSha256.hxx:101
#define coutE(a)
Definition: RooMsgService.h:37
#define ClassImp(name)
Definition: Rtypes.h:375
#define R__ASSERT(e)
Definition: TError.h:118
char name[80]
Definition: TGX11.cxx:110
float xmin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
Int_t getSize() const
Return the number of elements in the collection.
const RooArgSet * nset() const
Definition: RooAbsProxy.h:48
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:62
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
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
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
RooPolynomial implements a polynomial p.d.f of the form.
Definition: RooPolynomial.h:28
double evaluate() const override
do not persist
RooRealProxy _x
Definition: RooPolynomial.h:45
std::vector< double > _wksp
Definition: RooPolynomial.h:49
Int_t _lowestOrder
Definition: RooPolynomial.h:47
RooPolynomial()
coverity[UNINIT_CTOR]
~RooPolynomial() override
Destructor.
RooListProxy _coefList
Definition: RooPolynomial.h:46
double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const override
Do the analytical integral according to the code that was returned by getAnalyticalIntegral().
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Advertise to RooFit that this function can be analytically integrated.
double max(const char *rname=nullptr) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
double min(const char *rname=nullptr) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
RVec< PromoteTypes< T0, T1 > > pow(const T0 &x, const RVec< T1 > &v)
Definition: RVec.hxx:1770
Double_t x[n]
Definition: legend1.C:17
@ InputArguments
Definition: RooGlobalFunc.h:64