Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooLognormal.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * @(#)root/roofit:$Id$ *
4 * *
5 * RooFit Lognormal PDF *
6 * *
7 * Author: Gregory Schott and Stefan Schmitz *
8 * *
9 *****************************************************************************/
10
11/** \class RooLognormal
12 \ingroup Roofit
13
14RooFit Lognormal PDF. The two parameters are:
15 - `m0`: the median of the distribution
16 - `k = exp(sigma)`: sigma is called the shape parameter in the TMath parameterization
17
18\f[
19 \mathrm{RooLognormal}(x \, | \, m_0, k) = \frac{1}{\sqrt{2\pi \cdot \ln(k) \cdot x}} \cdot \exp\left(
20\frac{-\ln^2(\frac{x}{m_0})}{2 \ln^2(k)} \right) \f]
21
22The parameterization here is physics driven and differs from the ROOT::Math::lognormal_pdf() in `x,m,s,x0` with:
23 - `m = log(m0)`
24 - `s = log(k)`
25 - `x0 = 0`
26**/
27
28#include "RooLognormal.h"
29#include "RooRandom.h"
30#include "RooMath.h"
31#include "RooHelpers.h"
32#include "RooBatchCompute.h"
33
35
37
38////////////////////////////////////////////////////////////////////////////////
39
40RooLognormal::RooLognormal(const char *name, const char *title, RooAbsReal &_x, RooAbsReal &_m0, RooAbsReal &_k,
41 bool useStandardParametrization)
42 : RooAbsPdf{name, title},
43 x{"x", "Observable", this, _x},
44 m0{"m0", "m0", this, _m0},
45 k{"k", "k", this, _k},
46 _useStandardParametrization{useStandardParametrization}
47{
48 RooHelpers::checkRangeOfParameters(this, {&_x, &_m0, &_k}, 0.);
49
50 auto par = dynamic_cast<const RooAbsRealLValue *>(&_k);
51 const double unsafeValue = useStandardParametrization ? 0.0 : 1.0;
52 if (par && par->getMin() <= unsafeValue && par->getMax() >= unsafeValue) {
53 coutE(InputArguments) << "The parameter '" << par->GetName() << "' with range [" << par->getMin("") << ", "
54 << par->getMax() << "] of the " << this->ClassName() << " '" << this->GetName()
55 << "' can reach the unsafe value " << unsafeValue << " "
56 << ". Advise to limit its range." << std::endl;
57 }
58}
59
60////////////////////////////////////////////////////////////////////////////////
61
63 : RooAbsPdf(other, name),
64 x("x", this, other.x),
65 m0("m0", this, other.m0),
66 k{"k", this, other.k},
67 _useStandardParametrization{other._useStandardParametrization}
68{
69}
70
71////////////////////////////////////////////////////////////////////////////////
72/// ln(k)<1 would correspond to sigma < 0 in the parameterization
73/// resulting by transforming a normal random variable in its
74/// standard parameterization to a lognormal random variable
75/// => treat ln(k) as -ln(k) for k<1
76
78{
79 const double ln_k = std::abs(_useStandardParametrization ? k : std::log(k));
80 const double ln_m0 = _useStandardParametrization ? m0 : std::log(m0);
81
82 return ROOT::Math::lognormal_pdf(x, ln_m0, ln_k);
83}
84
85////////////////////////////////////////////////////////////////////////////////
86/// Compute multiple values of Lognormal distribution.
88{
90 RooBatchCompute::compute(ctx.config(this), computer, ctx.output(),
91 {ctx.at(x), ctx.at(m0), ctx.at(k)});
92}
93
94////////////////////////////////////////////////////////////////////////////////
95
96Int_t RooLognormal::getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char * /*rangeName*/) const
97{
98 return matchArgs(allVars, analVars, x) ? 1 : 0;
99}
100
101////////////////////////////////////////////////////////////////////////////////
102
103double RooLognormal::analyticalIntegral(Int_t /*code*/, const char *rangeName) const
104{
105 static const double root2 = std::sqrt(2.);
106
107 double ln_k = std::abs(_useStandardParametrization ? k : std::log(k));
108 double scaledMin = _useStandardParametrization ? std::log(x.min(rangeName)) - m0 : std::log(x.min(rangeName) / m0);
109 double scaledMax = _useStandardParametrization ? std::log(x.max(rangeName)) - m0 : std::log(x.max(rangeName) / m0);
110 return 0.5 * (RooMath::erf(scaledMax / (root2 * ln_k)) - RooMath::erf(scaledMin / (root2 * ln_k)));
111}
112
113////////////////////////////////////////////////////////////////////////////////
114
115Int_t RooLognormal::getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool /*staticInitOK*/) const
116{
117 return matchArgs(directVars, generateVars, x) ? 1 : 0;
118}
119
120////////////////////////////////////////////////////////////////////////////////
121
123{
124 const double ln_k = std::abs(_useStandardParametrization ? k : std::log(k));
125 const double ln_m0 = _useStandardParametrization ? m0 : std::log(m0);
126
127 double xgen;
128 while (true) {
129 xgen = std::exp(RooRandom::randomGenerator()->Gaus(ln_m0, ln_k));
130 if (xgen <= x.max() && xgen >= x.min()) {
131 x = xgen;
132 break;
133 }
134 }
135
136 return;
137}
#define coutE(a)
#define ClassImp(name)
Definition Rtypes.h:382
char name[80]
Definition TGX11.cxx:110
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
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
std::span< double > output()
RooBatchCompute::Config config(RooAbsArg const *arg) const
RooFit Lognormal PDF.
void doEval(RooFit::EvalContext &) const override
Compute multiple values of Lognormal distribution.
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool staticInitOK=true) const override
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
RooRealProxy k
the shape parameter, exp(sigma)
void generateEvent(Int_t code) override
Interface for generation of an event using the algorithm corresponding to the specified code.
bool _useStandardParametrization
double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const override
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
bool useStandardParametrization() const
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const override
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
RooRealProxy x
the variable
RooRealProxy m0
the median, exp(mu)
double evaluate() const override
ln(k)<1 would correspond to sigma < 0 in the parameterization resulting by transforming a normal rand...
static std::complex< double > erf(const std::complex< double > z)
complex erf function
Definition RooMath.cxx:59
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition RooRandom.cxx:48
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
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:225
double lognormal_pdf(double x, double m, double s, double x0=0)
Probability density function of the lognormal distribution.
Double_t x[n]
Definition legend1.C:17
void compute(Config cfg, Computer comp, std::span< double > output, VarSpan vars, ArgSpan extraArgs={})
void checkRangeOfParameters(const RooAbsReal *callingClass, std::initializer_list< const RooAbsReal * > pars, double min=-std::numeric_limits< double >::max(), double max=std::numeric_limits< double >::max(), bool limitsInAllowedRange=false, std::string const &extraMessage="")
Check if the parameters have a range, and warn if the range extends below / above the set limits.