ROOT   Reference Guide
Searching...
No Matches
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 *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
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 "RooAbsReal.h"
34#include "RooArgList.h"
35#include "RooMsgService.h"
36#include "RooBatchCompute.h"
37
38#include "TError.h"
39#include <vector>
40using namespace std;
41
43
44////////////////////////////////////////////////////////////////////////////////
45/// coverity[UNINIT_CTOR]
46
48{
49}
50
51////////////////////////////////////////////////////////////////////////////////
52/// Create a polynomial in the variable x.
53/// \param[in] name Name of the PDF
54/// \param[in] title Title for plotting the PDF
55/// \param[in] x The variable of the polynomial
56/// \param[in] coefList The coefficients \f$a_i \f$
57/// \param[in] lowestOrder [optional] Truncate the sum such that it skips the lower orders:
58/// \f[
59/// 1. + \sum_{i=0}^{\mathrm{coefList.size()}} a_{i} * x^{(i + \mathrm{lowestOrder})}
60/// \f]
61///
62/// This means that
63/// \code{.cpp}
64/// RooPolynomial pol("pol", "pol", x, RooArgList(a, b), lowestOrder = 2)
65/// \endcode
66/// computes
67/// \f[
68/// \mathrm{pol}(x) = 1 * x^0 + (0 * x^{\ldots}) + a * x^2 + b * x^3.
69/// \f]
70
71
72RooPolynomial::RooPolynomial(const char* name, const char* title,
73 RooAbsReal& x, const RooArgList& coefList, Int_t lowestOrder) :
74 RooAbsPdf(name, title),
75 _x("x", "Dependent", this, x),
76 _coefList("coefList","List of coefficients",this),
77 _lowestOrder(lowestOrder)
78{
79 // Check lowest order
80 if (_lowestOrder<0) {
81 coutE(InputArguments) << "RooPolynomial::ctor(" << GetName()
82 << ") WARNING: lowestOrder must be >=0, setting value to 0" << endl ;
83 _lowestOrder=0 ;
84 }
85
86 RooFIter coefIter = coefList.fwdIterator() ;
87 RooAbsArg* coef ;
88 while((coef = (RooAbsArg*)coefIter.next())) {
89 if (!dynamic_cast<RooAbsReal*>(coef)) {
90 coutE(InputArguments) << "RooPolynomial::ctor(" << GetName() << ") ERROR: coefficient " << coef->GetName()
91 << " is not of type RooAbsReal" << endl ;
92 R__ASSERT(0) ;
93 }
95 }
96}
97
98////////////////////////////////////////////////////////////////////////////////
99
100RooPolynomial::RooPolynomial(const char* name, const char* title,
101 RooAbsReal& x) :
102 RooAbsPdf(name, title),
103 _x("x", "Dependent", this, x),
104 _coefList("coefList","List of coefficients",this),
105 _lowestOrder(1)
106{ }
107
108////////////////////////////////////////////////////////////////////////////////
109/// Copy constructor
110
112 RooAbsPdf(other, name),
113 _x("x", this, other._x),
114 _coefList("coefList",this,other._coefList),
115 _lowestOrder(other._lowestOrder)
116{ }
117
118////////////////////////////////////////////////////////////////////////////////
119/// Destructor
120
122{ }
123
124////////////////////////////////////////////////////////////////////////////////
125
127{
128 // Calculate and return value of polynomial
129
130 const unsigned sz = _coefList.getSize();
131 const int lowestOrder = _lowestOrder;
132 if (!sz) return lowestOrder ? 1. : 0.;
133 _wksp.clear();
134 _wksp.reserve(sz);
135 {
136 const RooArgSet* nset = _coefList.nset();
138 RooAbsReal* c;
139 while ((c = (RooAbsReal*) it.next())) _wksp.push_back(c->getVal(nset));
140 }
141 const Double_t x = _x;
142 Double_t retVal = _wksp[sz - 1];
143 for (unsigned i = sz - 1; i--; ) retVal = _wksp[i] + x * retVal;
144 return retVal * std::pow(x, lowestOrder) + (lowestOrder ? 1.0 : 0.0);
145}
146
147////////////////////////////////////////////////////////////////////////////////
148/// Compute multiple values of Polynomial.
150 RooSpan<const double> xData = _x->getValues(evalData, normSet);
151 int batchSize = xData.size();
152 RooSpan<double> output = evalData.makeBatch(this, batchSize);
153
154 const int nCoef = _coefList.getSize();
155 const RooArgSet* listNormSet = _coefList.nset();
157 for (int i=0; i<nCoef; i++) {
158 auto valBatch = static_cast<RooAbsReal&>(_coefList[i]).getValues(evalData, listNormSet);
159 coefList.emplace_back(valBatch);
160 }
161
162 RooBatchCompute::dispatch->computePolynomial(batchSize, output.data(), xData.data(), _lowestOrder, coefList);
163 return output;
164}
165
166////////////////////////////////////////////////////////////////////////////////
167/// Advertise to RooFit that this function can be analytically integrated.
168Int_t RooPolynomial::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
169{
170 if (matchArgs(allVars, analVars, _x)) return 1;
171 return 0;
172}
173
174////////////////////////////////////////////////////////////////////////////////
175/// Do the analytical integral according to the code that was returned by getAnalyticalIntegral().
176Double_t RooPolynomial::analyticalIntegral(Int_t code, const char* rangeName) const
177{
178 R__ASSERT(code==1) ;
179
180 const Double_t xmin = _x.min(rangeName), xmax = _x.max(rangeName);
181 const int lowestOrder = _lowestOrder;
182 const unsigned sz = _coefList.getSize();
183 if (!sz) return xmax - xmin;
184 _wksp.clear();
185 _wksp.reserve(sz);
186 {
187 const RooArgSet* nset = _coefList.nset();
189 unsigned i = 1 + lowestOrder;
190 RooAbsReal* c;
191 while ((c = (RooAbsReal*) it.next())) {
192 _wksp.push_back(c->getVal(nset) / Double_t(i));
193 ++i;
194 }
195 }
196 Double_t min = _wksp[sz - 1], max = _wksp[sz - 1];
197 for (unsigned i = sz - 1; i--; )
198 min = _wksp[i] + xmin * min, max = _wksp[i] + xmax * max;
199 return max * std::pow(xmax, 1 + lowestOrder) - min * std::pow(xmin, 1 + lowestOrder) +
200 (lowestOrder ? (xmax - xmin) : 0.);
201}
#define c(i)
Definition RSha256.hxx:101
#define coutE(a)
#define ClassImp(name)
Definition Rtypes.h:364
#define R__ASSERT(e)
Definition TError.h:120
char name[80]
Definition TGX11.cxx:110
float xmin
float xmax
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition RooAbsArg.h:72
Int_t getSize() const
RooFIter fwdIterator() const
One-time forward iterator.
RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Compute batch of values for given input data, and normalise by integrating over the observables in no...
const RooArgSet * nset() const
Definition RooAbsProxy.h:45
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:61
virtual RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet=nullptr) const
by this change, please consult the release notes for ROOT 6.24 for guidance on how to make this trans...
Bool_t 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:21
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:29
virtual void computePolynomial(size_t batchSize, double *__restrict output, const double *__restrict const xData, int lowestOrder, std::vector< BracketAdapterWithMask > &coef)=0
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
RooAbsArg * next()
Return next element or nullptr if at end.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
RooPolynomial implements a polynomial p.d.f of the form.
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Advertise to RooFit that this function can be analytically integrated.
RooRealProxy _x
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Do the analytical integral according to the code that was returned by getAnalyticalIntegral().
RooPolynomial()
coverity[UNINIT_CTOR]
Double_t evaluate() const
do not persist
RooListProxy _coefList
RooSpan< double > evaluateSpan(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Compute multiple values of Polynomial.
std::vector< Double_t > _wksp
virtual ~RooPolynomial()
Destructor.
A simple container to hold a batch of data values.
Definition RooSpan.h:34
constexpr std::span< T >::pointer data() const
Definition RooSpan.h:106
constexpr std::span< T >::index_type size() const noexcept
Definition RooSpan.h:121
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47
Double_t x[n]
Definition legend1.C:17
R__EXTERN RooBatchComputeInterface * dispatch
This dispatch pointer points to an implementation of the compute library, provided one has been loade...
This struct enables passing computation data around between elements of a computation graph.
Definition RunContext.h:31
RooSpan< double > makeBatch(const RooAbsReal *owner, std::size_t size)
Create a writable batch.
static void output(int code)
Definition gifencode.c:226